From f1702e507a41caccf447d8c07bbdeb3fa70ef446 Mon Sep 17 00:00:00 2001 From: Manuel Lera-Ramirez Date: Sun, 10 Nov 2024 08:37:37 +0000 Subject: [PATCH] fix icon link + 2 small bug fixes (#320) * fix icon link * closes #279 again (for circular seqs) (#324) * move requests optional dependency into function body (#325) * closes #262 (#326) --- README.md | 2 +- src/pydna/amplify.py | 8 +++++--- src/pydna/tm.py | 3 ++- tests/test_module_amplify.py | 33 +++++++++++++++++++++++++++++++++ 4 files changed, 41 insertions(+), 5 deletions(-) diff --git a/README.md b/README.md index 4872628c..b884a43b 100755 --- a/README.md +++ b/README.md @@ -1,4 +1,4 @@ -# ![icon](https://raw.githubusercontent.com/bjornFJohansson/pydna/master/docs/pics/pydna.resized.png) pydna +# ![icon](docs/_static/pydna.resized.png) pydna | [![Tests & Coverage](https://github.com/BjornFJohansson/pydna/actions/workflows/pydna_test_and_coverage_workflow.yml/badge.svg?branch=dev_bjorn)](https://github.com/BjornFJohansson/pydna/actions/workflows/pydna_test_and_coverage_workflow.yml) | [![codecov](https://codecov.io/gh/BjornFJohansson/pydna/branch/master/graph/badge.svg)](https://codecov.io/gh/BjornFJohansson/pydna/branch/master) | [![PyPI version](https://badge.fury.io/py/pydna.svg)](https://badge.fury.io/py/pydna) | [![Google group : pydna](https://img.shields.io/badge/Google%20Group-pydna-blue.svg)](https://groups.google.com/g/pydna) | | -------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------- | -------------------------------------------------------------------------------------------------------------------------------------------------- | -------------------------------------------------------------------------------------------------------------------------------------- | ------------------------------------------------------------------------------------------------------------------------------------- | diff --git a/src/pydna/amplify.py b/src/pydna/amplify.py index 3c12854a..907a5b20 100644 --- a/src/pydna/amplify.py +++ b/src/pydna/amplify.py @@ -344,18 +344,20 @@ def products(self): if self.template.circular: shift = fp.position - fp._fp tpl = self.template.shifted(shift) # shift template so that it starts where the fp starts anneling - feats = tpl[: rp.position + rp._fp].features fp.position = fp._fp # New position of fp becomes the footprint length rp.position = (rp.position - shift) % len(self.template) # Shift the rp position as well + feats = tpl[: rp.position + rp._fp].features elif fp.position <= rp.position: # pcr products only formed if fp anneals forward of rp feats = self.template[ fp.position - fp._fp : rp.position + rp._fp ].features # Save features covered by primers - shift_amount = len(fp.tail) - feats = [_shift_feature(f, shift_amount, None) for f in feats] tpl = self.template else: continue + # Shift features to the right if there was a tail + shift_amount = len(fp.tail) + feats = [_shift_feature(f, shift_amount, None) for f in feats] + if tpl.circular and fp.position == rp.position: prd = _Dseqrecord(fp) + _Dseqrecord(rp).reverse_complement() else: diff --git a/src/pydna/tm.py b/src/pydna/tm.py index e7938041..8c4e3423 100644 --- a/src/pydna/tm.py +++ b/src/pydna/tm.py @@ -14,7 +14,6 @@ import textwrap as _textwrap from pydna._pretty import pretty_str as _pretty_str -import requests # See the documentation for Bio.SeqUtils.MeltingTemp for more details # The 10X Taq Buffer with (NH4)2SO4 is commercialized by companies like @@ -331,6 +330,8 @@ def tmbresluc(primer: str, *args, primerc=500.0, saltc=50, **kwargs): def tm_neb(primer, conc=0.5, prodcode="q5-0"): + import requests + """Calculates a single primers melting temp from NEB. Parameters diff --git a/tests/test_module_amplify.py b/tests/test_module_amplify.py index eeb4d895..0fa15014 100644 --- a/tests/test_module_amplify.py +++ b/tests/test_module_amplify.py @@ -788,6 +788,7 @@ def test_annotation(): """ from pydna.amplify import pcr from pydna.dseqrecord import Dseqrecord + from pydna.primer import Primer dsr = Dseqrecord("ATGCAAACAGTAATGATGGATGACATTCAAAGCACTGATTCTATTGCTGAAAAAGATAAT") dsr.add_feature(x=0, y=60, type="gene", label="my_gene") # We add a feature to highlight the sequence as a gene @@ -799,6 +800,38 @@ def test_annotation(): assert pcr_product.features[0].location.extract(pcr_product).seq == dsr.seq + # Also works in circular sequences + dsr_circ = dsr.looped() + pcr_product_circ = pcr(forward_primer, reverse_primer, dsr_circ) + assert str(pcr_product_circ.features[0].location.extract(pcr_product_circ).seq) == str(dsr_circ.seq) + + # Check that annotations are transmitted properly if the PCR product spans + # the origin in a circular sequence + + vector = Dseqrecord( + "ATGCAAACAGTAATGATGGATGACACCAGCTTCATGAAATGGAACAGTGCCAGAAAAAACTTGAAGATGTTCAAAGCACTGATTCTATTGCTGAAAAAGATAAT", + circular=True, + ) + + vector.add_feature(17, 40, type_="test", label=["left"]) + vector.add_feature(41, 90, type_="test", label=["right"]) + vector.add_feature(17, 90, type_="test", label=["all"]) + vector.add_feature(30, 60, type_="test", label=["middle"]) + + feature_seqs = set(str(f.location.extract(vector).seq) for f in vector.features) + + for shift in range(len(vector)): + shifted_vector = vector.shifted(shift) + + primer_f = Primer("acgtGGATGACACCAGCTTCAT") + primer_r = Primer("attacCAATAGAATCAGTGCTTTGAACA") + + product = pcr(primer_f, primer_r, shifted_vector) + + product_seqs = set(str(f.location.extract(product).seq) for f in product.features if f.type == "test") + + assert product_seqs == feature_seqs, f"Shift {shift}" + if __name__ == "__main__": pytest.main([__file__, "-vv", "-s"])