From 5586c9f84d9ce7f09d42edbd0acfd20ee108450c Mon Sep 17 00:00:00 2001 From: David Doty Date: Wed, 30 Nov 2022 11:34:13 -0800 Subject: [PATCH 01/14] bumped version --- scadnano/scadnano.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/scadnano/scadnano.py b/scadnano/scadnano.py index 47d0f0fd..a07f806d 100644 --- a/scadnano/scadnano.py +++ b/scadnano/scadnano.py @@ -54,7 +54,7 @@ # commented out for now to support Py3.6, which does not support this feature # from __future__ import annotations -__version__ = "0.17.6" # version line; WARNING: do not remove or change this line or comment +__version__ = "0.17.7" # version line; WARNING: do not remove or change this line or comment import dataclasses from abc import abstractmethod, ABC, ABCMeta From 64638658728fc964ec33336a6a800acbbfd41f53 Mon Sep 17 00:00:00 2001 From: David Doty Date: Fri, 2 Dec 2022 09:19:34 -0800 Subject: [PATCH 02/14] closes #252: add method `Design.base_pairs()` --- .github/workflows/run_unit_tests.yml | 2 +- scadnano/scadnano.py | 164 +++++++++++++++++----- tests/scadnano_tests.py | 196 +++++++++++++++++++++++---- 3 files changed, 299 insertions(+), 63 deletions(-) diff --git a/.github/workflows/run_unit_tests.yml b/.github/workflows/run_unit_tests.yml index 0cae72c6..a0baddf4 100644 --- a/.github/workflows/run_unit_tests.yml +++ b/.github/workflows/run_unit_tests.yml @@ -8,7 +8,7 @@ jobs: runs-on: ubuntu-latest strategy: matrix: - python-version: [3.7, 3.8, 3.9, 3.10, 3.11] + python-version: [3.7, 3.8, 3.9, "3.10", "3.11"] steps: - uses: actions/checkout@v2 diff --git a/scadnano/scadnano.py b/scadnano/scadnano.py index 7c5917a5..de7ce3fd 100644 --- a/scadnano/scadnano.py +++ b/scadnano/scadnano.py @@ -56,6 +56,7 @@ __version__ = "0.17.6" # version line; WARNING: do not remove or change this line or comment +import collections import dataclasses from abc import abstractmethod, ABC, ABCMeta import json @@ -65,7 +66,7 @@ from builtins import ValueError from dataclasses import dataclass, field, InitVar, replace from typing import Iterator, Tuple, List, Sequence, Iterable, Set, Dict, Union, Optional, Type, cast, Any, \ - TypeVar, Generic, Callable, AbstractSet + TypeVar, Generic, Callable, AbstractSet, Deque from collections import defaultdict, OrderedDict, Counter import sys import os.path @@ -1949,8 +1950,7 @@ def overlaps(self, other: 'Domain') -> bool: has nonempty intersection with those of `other`, and they appear on the same helix, and they point in opposite directions.""" # noqa (suppress PEP warning) - return (self.helix == other.helix and - self.forward == (not other.forward) and + return (self.forward == (not other.forward) and self.compute_overlap(other)[0] >= 0) # remove quotes when Py3.6 support dropped @@ -1965,8 +1965,7 @@ def overlaps_illegally(self, other: 'Domain') -> bool: has nonempty intersection with those of `other`, and they appear on the same helix, and they point in the same direction.""" # noqa (suppress PEP warning) - return (self.helix == other.helix and - self.forward == other.forward and + return (self.forward == other.forward and self.compute_overlap(other)[0] >= 0) # remove quotes when Py3.6 support dropped @@ -1976,6 +1975,8 @@ def compute_overlap(self, other: 'Domain') -> Tuple[int, int]: Return ``(-1,-1)`` if they do not overlap (different helices, or non-overlapping regions of the same helix).""" + if self.helix != other.helix: + return -1, -1 overlap_start = max(self.start, other.start) overlap_end = min(self.end, other.end) if overlap_start >= overlap_end: # overlap is empty @@ -4671,6 +4672,72 @@ def _check_type_is_one_of(obj: Any, expected_types: Iterable) -> None: f'but instead it is of type {type(obj)}') +def _popleft(queue: Deque[T]) -> Optional[T]: + # return left element of queue if queue is nonempty, otherwise return None + try: + elt = queue.popleft() + return elt + except IndexError: + return None + + +def find_overlapping_domains_on_helix(helix: Helix) -> List[Tuple[Domain, Domain]]: + # compute list of pairs of domains that overlap on Helix `helix` + # assumes that `helix.domains` has been populated by calling `Design._build_domains_on_helix_lists()` + forward_domains = [] + reverse_domains = [] + for domain in helix.domains: + if domain.forward: + forward_domains.append(domain) + else: + reverse_domains.append(domain) + + forward_domains.sort(key=lambda domain: domain.start) + reverse_domains.sort(key=lambda domain: domain.start) + + # need to be efficient to remove the front element repeatedly + reverse_domains = collections.deque(reverse_domains) + + overlapping_domains = [] + + for forward_domain in forward_domains: + reverse_domain = reverse_domains[0] + # remove each reverse_domain that strictly precedes forward domain + # they cannot overlap forward_domain nor any domain following it in the list forward_domains + while reverse_domain.end <= forward_domain.start and len(reverse_domains) > 0: + reverse_domains.popleft() + if len(reverse_domains) > 0: + reverse_domain = reverse_domains[0] + else: # if all reverse domains are gone, we're done + return overlapping_domains + + # otherwise we may have found an overlapping reverse_domain, OR forward_domain could precede it + # if forward_domain precedes reverse_domain, next inner loop is skipped, + # and we go to next forward_domain in the outer loop + + # add each reverse_domain that overlaps forward_domain + while forward_domain.overlaps(reverse_domain): + overlapping_domains.append((forward_domain, reverse_domain)) + + if reverse_domain.end <= forward_domain.end: + # [-----f_dom--->[--next_f_dom--> + # [--r_dom---> + # reverse_domain can't overlap *next* forward_domain, so safe to remove + reverse_domains.popleft() + if len(reverse_domains) == 0: + break + else: + reverse_domain = reverse_domains[0] + else: + # [---f_dom---> [---next_f_dom--> + # [----r_dom->[--next_r_dom----> + # reverse_domain possibly overlaps next forward_domain, so keep it in queue + # but this is last reverse_domain overlapping current forward_domain, so safe to break loop + break + + return overlapping_domains + + @dataclass class Design(_JSONSerializable, Generic[StrandLabel, DomainLabel]): """Object representing the entire design of the DNA structure.""" @@ -4707,7 +4774,7 @@ class Design(_JSONSerializable, Generic[StrandLabel, DomainLabel]): def __init__(self, *, helices: Optional[Union[List[Helix], Dict[int, Helix]]] = None, groups: Optional[Dict[str, HelixGroup]] = None, - strands: List[Strand] = None, + strands: Iterable[Strand] = None, grid: Grid = Grid.none, helices_view_order: List[int] = None, geometry: Geometry = None) -> None: @@ -5257,6 +5324,24 @@ def idx_of(helix: Helix, order: int) -> int: return helices + def base_pairs(self) -> List[Tuple[int, int]]: + """ + List of "addresses" of base pairs in this design. Each address is a pair (`helix_idx`, `offset`) + of the address of the *forward* strand in the base pair. + + :return: + all base pairs (`helix_idx`, `offset`) in this :any:`Design` + """ + base_pairs = [] + for idx, helix in self.helices.items(): + overlapping_domains = find_overlapping_domains_on_helix(helix) + for dom1, dom2 in overlapping_domains: + start, end = dom1.compute_overlap(dom2) + for offset in range(start, end): + base_pairs.append((idx, offset)) + + return base_pairs + @staticmethod def assign_modifications_to_strands(strands: List[Strand], strand_jsons: List[dict], all_mods: Dict[str, Modification]) -> None: @@ -6954,7 +7039,8 @@ def _write_plates_default(self, directory: str, filename: Optional[str], strands workbook.save(filename_plate) - def to_oxview_format(self, warn_duplicate_strand_names: bool = True, use_strand_colors: bool = True) -> dict: + def to_oxview_format(self, warn_duplicate_strand_names: bool = True, + use_strand_colors: bool = True) -> dict: """ Exports to oxView format. @@ -6977,23 +7063,23 @@ def to_oxview_format(self, warn_duplicate_strand_names: bool = True, use_strand_ strand_count += 1 oxvnucs: List[Dict[str, Any]] = [] strand_nuc_start.append(nuc_count) - oxvstrand = {'id': strand_count, - 'class': 'NucleicAcidStrand', - 'end5': nuc_count, - 'end3': nuc_count+len(oxdna_strand.nucleotides), + oxvstrand = {'id': strand_count, + 'class': 'NucleicAcidStrand', + 'end5': nuc_count, + 'end3': nuc_count + len(oxdna_strand.nucleotides), 'monomers': oxvnucs} if use_strand_colors and (strand1.color is not None): scolor = strand1.color.to_cadnano_v2_int_hex() else: scolor = None - + for index_in_strand, nuc in enumerate(oxdna_strand.nucleotides): - oxvnuc = {'id': nuc_count, + oxvnuc = {'id': nuc_count, 'p': [nuc.r.x, nuc.r.y, nuc.r.z], 'a1': [nuc.b.x, nuc.b.y, nuc.b.z], 'a3': [nuc.n.x, nuc.n.y, nuc.n.z], - 'class': 'DNA', - 'type': nuc.base, + 'class': 'DNA', + 'type': nuc.base, 'cluster': 1} if index_in_strand != 0: oxvnuc['n5'] = nuc_count - 1 @@ -7009,11 +7095,11 @@ def to_oxview_format(self, warn_duplicate_strand_names: bool = True, use_strand_ for si2, strand2 in enumerate(self.strands): if not strand1.overlaps(strand2): continue - s1_nuc_idx = strand_nuc_start[si1+1] + s1_nuc_idx = strand_nuc_start[si1 + 1] for domain1 in strand1.domains: if isinstance(domain1, (Loopout, Extension)): continue - s2_nuc_idx = strand_nuc_start[si2+1] + s2_nuc_idx = strand_nuc_start[si2 + 1] for domain2 in strand2.domains: if isinstance(domain2, (Loopout, Extension)): continue @@ -7028,35 +7114,37 @@ def to_oxview_format(self, warn_duplicate_strand_names: bool = True, use_strand_ d1range = range(s1_left, s1_right) d2range = range(s2_left, s2_right, -1) else: - d1range = range(s1_right+1, s1_left+1) - d2range = range(s2_right-1, s2_left-1, -1) + d1range = range(s1_right + 1, s1_left + 1) + d2range = range(s2_right - 1, s2_left - 1, -1) assert len(d1range) == len(d2range) # Check for mismatches, and do not add a pair if the bases are *known* # to mismatch. (FIXME: this must be changed if scadnano later supports # degenerate base codes.) for d1, d2 in zip(d1range, d2range): - if ((strand1.dna_sequence is not None) and - (strand2.dna_sequence is not None) and - (strand1.dna_sequence[d1] != "?") and - (strand2.dna_sequence[d2] != "?") and - (wc(strand1.dna_sequence[d1]) != strand2.dna_sequence[d2])): + if ((strand1.dna_sequence is not None) and + (strand2.dna_sequence is not None) and + (strand1.dna_sequence[d1] != "?") and + (strand2.dna_sequence[d2] != "?") and + (wc(strand1.dna_sequence[d1]) != strand2.dna_sequence[d2])): continue oxv_strand1['monomers'][d1]['bp'] = s2_nuc_idx + d2 if 'bp' in oxview_strands[si2]['monomers'][d2]: if oxview_strands[si2]['monomers'][d2]['bp'] != s1_nuc_idx + d1: - print (s2_nuc_idx+d2, s1_nuc_idx+d1, oxview_strands[si2]['monomers'][d2]['bp'], domain1, domain2) + print(s2_nuc_idx + d2, s1_nuc_idx + d1, + oxview_strands[si2]['monomers'][d2]['bp'], domain1, domain2) b = system.compute_bounding_box() - oxvsystem = {'box': [b.x, b.y, b.z], - 'date': datetime.datetime.now().isoformat(), - 'systems': [{'id': 0, 'strands': oxview_strands}], + oxvsystem = {'box': [b.x, b.y, b.z], + 'date': datetime.datetime.now().isoformat(), + 'systems': [{'id': 0, 'strands': oxview_strands}], 'forces': [], 'selections': []} return oxvsystem - def write_oxview_file(self, directory: str = '.', filename: Optional[str] = None, warn_duplicate_strand_names: bool = True, use_strand_colors: bool = True) -> None: + def write_oxview_file(self, directory: str = '.', filename: Optional[str] = None, + warn_duplicate_strand_names: bool = True, use_strand_colors: bool = True) -> None: """Writes an oxView file rerpesenting this design. :param directory: @@ -7070,7 +7158,8 @@ def write_oxview_file(self, directory: str = '.', filename: Optional[str] = None if True (default), sets the color of each nucleotide in a strand in oxView to the color of the strand. """ - oxvsystem = self.to_oxview_format(warn_duplicate_strand_names=warn_duplicate_strand_names, use_strand_colors=use_strand_colors) + oxvsystem = self.to_oxview_format(warn_duplicate_strand_names=warn_duplicate_strand_names, + use_strand_colors=use_strand_colors) write_file_same_name_as_running_python_script(json.dumps(oxvsystem), 'oxview', directory, filename) def to_oxdna_format(self, warn_duplicate_strand_names: bool = True) -> Tuple[str, str]: @@ -7121,8 +7210,10 @@ def write_oxdna_files(self, directory: str = '.', filename_no_extension: Optiona """ dat, top = self.to_oxdna_format(warn_duplicate_strand_names) - write_file_same_name_as_running_python_script(dat, 'dat', directory, filename_no_extension, add_extension=True) - write_file_same_name_as_running_python_script(top, 'top', directory, filename_no_extension, add_extension=True) + write_file_same_name_as_running_python_script(dat, 'dat', directory, filename_no_extension, + add_extension=True) + write_file_same_name_as_running_python_script(top, 'top', directory, filename_no_extension, + add_extension=True) # @_docstring_parameter was used to substitute sc in for the filename extension, but it is # incompatible with .. code-block:: and caused a very strange and hard-to-determine error, @@ -7874,7 +7965,8 @@ def _name_of_this_script() -> str: def write_file_same_name_as_running_python_script(contents: str, extension: str, directory: str = '.', - filename: Optional[str] = None, add_extension: bool = False) -> None: + filename: Optional[str] = None, + add_extension: bool = False) -> None: """ Writes a text file with `contents` whose name is (by default) the same as the name of the currently running script, but with extension ``.py`` changed to `extension`. @@ -7888,13 +7980,15 @@ def write_file_same_name_as_running_python_script(contents: str, extension: str, :param filename: filename to use instead of the currently running script """ - relative_filename = _get_filename_same_name_as_running_python_script(directory, extension, filename, add_extension=add_extension) + relative_filename = _get_filename_same_name_as_running_python_script(directory, extension, filename, + add_extension=add_extension) with open(relative_filename, 'w') as out_file: out_file.write(contents) def _get_filename_same_name_as_running_python_script(directory: str, extension: str, - filename: Optional[str], add_extension: bool = False) -> str: + filename: Optional[str], + add_extension: bool = False) -> str: # if filename is not None, assume it has an extension if filename is None: filename = _name_of_this_script() + '.' + extension diff --git a/tests/scadnano_tests.py b/tests/scadnano_tests.py index f5dfed3e..6a0e8b02 100644 --- a/tests/scadnano_tests.py +++ b/tests/scadnano_tests.py @@ -16,6 +16,7 @@ from scadnano.scadnano import _convert_design_to_oxdna_system + def strand_matching(strands: Iterable[sc.Strand], helix: int, forward: bool, start: int, end: int) -> sc.Strand: """ @@ -1221,10 +1222,12 @@ def test_write_idt_plate_excel_file(self) -> None: # add 10 strands in excess of 3 plates for plate_type in [sc.PlateType.wells96, sc.PlateType.wells384]: + num_strands = 3 * plate_type.num_wells_per_plate() + 10 filename = f'test_excel_export_{plate_type.num_wells_per_plate()}.xls' - helices = [sc.Helix(max_offset=100) for _ in range(1)] + max_offset = num_strands * strand_len + helices = [sc.Helix(max_offset=max_offset) for _ in range(1)] design = sc.Design(helices=helices, strands=[], grid=sc.square) - for strand_idx in range(3 * plate_type.num_wells_per_plate() + 10): + for strand_idx in range(num_strands): design.draw_strand(0, strand_len * strand_idx).move(strand_len).with_name(f's{strand_idx}') design.strands[-1].set_dna_sequence('T' * strand_len) @@ -5470,7 +5473,7 @@ class TestCircularStrandEdits(unittest.TestCase): ''' def setUp(self) -> None: - helices = [sc.Helix(max_offset=10) for _ in range(3)] + helices = [sc.Helix(max_offset=50) for _ in range(3)] self.design = sc.Design(helices=helices, strands=[]) self.design.draw_strand(0, 0).move(10).cross(1).move(-10) self.design.draw_strand(0, 15).move(5).cross(1).move(-10).cross(0).move(5) @@ -6906,9 +6909,9 @@ def test_export(self): for i, (oxdna_strand, oxview_strand, oxview_nocolor_strand, des_strand) in enumerate( - zip(oxdna_system.strands, oxv['systems'][0]['strands'], - oxv_no_color['systems'][0]['strands'], - design.strands)): + zip(oxdna_system.strands, oxv['systems'][0]['strands'], + oxv_no_color['systems'][0]['strands'], + design.strands)): self.assertEqual(i + 1, oxview_strand['id']) if des_strand.color: @@ -6933,14 +6936,17 @@ def test_export(self): self.assertEqual(oxdna_nt.base, oxview_nt['type']) def test_bp(self): - des = sc.Design() + des = sc.Design() des.set_grid(sc.Grid.square) - des.helices = {i :sc.Helix(max_offset=20, idx=i, grid_position=(0,i)) for i in range(3)} - des.draw_strand(0, 0).to(6).with_deletions(4).to(15).cross(1, 9).to(20).with_insertions((15, 2)).cross(0).to(9) + des.helices = {i: sc.Helix(max_offset=20, idx=i, grid_position=(0, i)) for i in range(3)} + des.draw_strand(0, 0).to(6).with_deletions(4).to(15).cross(1, 9).to(20).with_insertions( + (15, 2)).cross(0).to(9) des.draw_strand(1, 0).to(9).cross(0).to(0).with_deletions(4) - des.draw_strand(1, 20).to(2).with_insertions((15, 2)).cross(2, 0).to(20).with_sequence('TTTCTCATGGGAAGCAAACTCGGTTTCCGCGTCGGATAGT') + des.draw_strand(1, 20).to(2).with_insertions((15, 2)).cross(2, 0).to(20).with_sequence( + 'TTTCTCATGGGAAGCAAACTCGGTTTCCGCGTCGGATAGT') des.draw_strand(2, 8).to(5).loopout(2, 5, 4).to(0) - des.draw_strand(2, 20).extension_5p(8).to(12).extension_3p(8).with_sequence('ATACTGGAACTACGCGCGTGAATT', assign_complement=False) + des.draw_strand(2, 20).extension_5p(8).to(12).extension_3p(8).with_sequence( + 'ATACTGGAACTACGCGCGTGAATT', assign_complement=False) oxv = des.to_oxview_format() @@ -6948,41 +6954,41 @@ def test_bp(self): # Basic complements with a deletion (wildcard sequences) for i in range(0, 8): - self.assertEqual(strands[0]['monomers'][i]['bp'], strands[1]['monomers'][-i-1]['id']) - self.assertEqual(strands[1]['monomers'][-i-1]['bp'], strands[0]['monomers'][i]['id']) + self.assertEqual(strands[0]['monomers'][i]['bp'], strands[1]['monomers'][-i - 1]['id']) + self.assertEqual(strands[1]['monomers'][-i - 1]['bp'], strands[0]['monomers'][i]['id']) # Self-complementary strand (wildcard sequences) for i in range(8, 14): - self.assertEqual(strands[0]['monomers'][i]['bp'], strands[0]['monomers'][7-i]['id']) + self.assertEqual(strands[0]['monomers'][i]['bp'], strands[0]['monomers'][7 - i]['id']) # Insertion (defined sequences) for i in range(14, 27): - self.assertEqual(strands[0]['monomers'][i]['bp'], strands[2]['monomers'][26-i]['id']) + self.assertEqual(strands[0]['monomers'][i]['bp'], strands[2]['monomers'][26 - i]['id']) # Before, in, and after a loopout (one strand with no sequence, one with defined sequence) for i in range(0, 3): - self.assertEqual(strands[3]['monomers'][i]['bp'], strands[2]['monomers'][27-i]['id']) + self.assertEqual(strands[3]['monomers'][i]['bp'], strands[2]['monomers'][27 - i]['id']) for i in range(3, 8): self.assertNotIn('bp', strands[3]['monomers'][i]) for i in range(8, 12): - self.assertEqual(strands[3]['monomers'][i]['bp'], strands[2]['monomers'][23+8-i]['id']) + self.assertEqual(strands[3]['monomers'][i]['bp'], strands[2]['monomers'][23 + 8 - i]['id']) # Mismatches should not be paired; also, extensions: - for i in range(0, 8): # 5p extension + for i in range(0, 8): # 5p extension self.assertNotIn('bp', strands[4]['monomers'][i]) - for i in range(8, 12): # complementary + for i in range(8, 12): # complementary print(i) - self.assertEqual(strands[4]['monomers'][i]['bp'], strands[2]['monomers'][40+7-i]['id']) - for i in range(12, 14): # two mismatches + self.assertEqual(strands[4]['monomers'][i]['bp'], strands[2]['monomers'][40 + 7 - i]['id']) + for i in range(12, 14): # two mismatches self.assertNotIn('bp', strands[4]['monomers'][i]) - self.assertNotIn('bp', strands[2]['monomers'][32+15-i]) - for i in range(14, 16): # complementary again - self.assertEqual(strands[4]['monomers'][i]['bp'], strands[2]['monomers'][32+15-i]['id']) - for i in range(16, len(strands[4]['monomers'])): # 3p extension + self.assertNotIn('bp', strands[2]['monomers'][32 + 15 - i]) + for i in range(14, 16): # complementary again + self.assertEqual(strands[4]['monomers'][i]['bp'], strands[2]['monomers'][32 + 15 - i]['id']) + for i in range(16, len(strands[4]['monomers'])): # 3p extension self.assertNotIn('bp', strands[4]['monomers'][i]) - + # Unbound region for i in range(28, 32): self.assertNotIn('bp', strands[2]['monomers'][i]) @@ -7009,6 +7015,7 @@ def test_export_file(self): self.assertEqual(oxv, oxv2) + class TestOxdnaExport(unittest.TestCase): def setUp(self) -> None: self.OX_UNITS_TO_NM = 0.8518 @@ -7704,10 +7711,12 @@ def test_file_output(self) -> None: self.assertEqual(dat, open(tmpdir + '/' + scriptname + '.dat').read()) # Now, write, to a specific filename without extensions - design.write_oxdna_files(directory=tmpdir, filename_no_extension='oxdna-Export with spaces in name') + design.write_oxdna_files(directory=tmpdir, + filename_no_extension='oxdna-Export with spaces in name') self.assertEqual(top, open(tmpdir + '/oxdna-Export with spaces in name.top').read()) self.assertEqual(dat, open(tmpdir + '/oxdna-Export with spaces in name.dat').read()) + class TestPlateMaps(unittest.TestCase): def setUp(self) -> None: @@ -7785,3 +7794,136 @@ def test_to_json_serializable__label_key_contains_non_default_name(self) -> None ext = sc.Extension(5, label="ext1") result = ext.to_json_serializable(False) self.assertEqual(result[sc.domain_label_key], "ext1") + + +class TestBasePairs(unittest.TestCase): + def setUp(self) -> None: + ''' + 111111111122222222223333333333 + 0123456789012345678901234567890123456789 + 0 [-->[--> [--> [--> [--> + <] <--------] <-] <-] <] <] + + 111111111122222222223333333333 + 0123456789012345678901234567890123456789 + 1 [-----------> [--> + <--] <-----------] + ''' + helices = [sc.Helix(max_offset=40) for _ in range(2)] + self.design = sc.Design(helices=helices) + # helix 0 forward + self.design.draw_strand(0, 0).move(4) + self.design.draw_strand(0, 4).move(4) + self.design.draw_strand(0, 12).move(4) + self.design.draw_strand(0, 20).move(4) + self.design.draw_strand(0, 28).move(4) + # helix 0 reverse + self.design.draw_strand(0, 3).move(-2) + self.design.draw_strand(0, 14).to(4) + self.design.draw_strand(0, 20).to(17) + self.design.draw_strand(0, 29).to(26) + self.design.draw_strand(0, 36).to(34) + self.design.draw_strand(0, 39).to(37) + # helix 1 forward + self.design.draw_strand(1, 4).to(17) + self.design.draw_strand(1, 20).to(24) + # helix 1 reverse + self.design.draw_strand(1, 12).to(8) + self.design.draw_strand(1, 26).to(13) + + def test_find_overlapping_domains(self) -> None: + d01f = self.design.strands[0].domains[0] + d02f = self.design.strands[1].domains[0] + d03f = self.design.strands[2].domains[0] + d04f = self.design.strands[3].domains[0] + d05f = self.design.strands[4].domains[0] + + d01r = self.design.strands[5].domains[0] + d02r = self.design.strands[6].domains[0] + d03r = self.design.strands[7].domains[0] + d04r = self.design.strands[8].domains[0] + d05r = self.design.strands[9].domains[0] + d06r = self.design.strands[10].domains[0] + + d11f = self.design.strands[11].domains[0] + d12f = self.design.strands[12].domains[0] + + d11r = self.design.strands[13].domains[0] + d12r = self.design.strands[14].domains[0] + + overlapping_domains_h0 = sc.find_overlapping_domains_on_helix(self.design.helices[0]) + + overlapping_domains_h1 = sc.find_overlapping_domains_on_helix(self.design.helices[1]) + + self.assertEqual(len(overlapping_domains_h0), 4) + self.assertEqual(len(overlapping_domains_h1), 3) + + self.assertIn((d01f, d01r), overlapping_domains_h0) + self.assertIn((d02f, d02r), overlapping_domains_h0) + self.assertIn((d03f, d02r), overlapping_domains_h0) + self.assertIn((d05f, d04r), overlapping_domains_h0) + + self.assertIn((d11f, d11r), overlapping_domains_h1) + self.assertIn((d11f, d12r), overlapping_domains_h1) + self.assertIn((d12f, d12r), overlapping_domains_h1) + ''' + 111111111122222222223333333333 + 0123456789012345678901234567890123456789 + 0 [-->[--> [--> [--> [--> + <] <--------] <-] <-] <] <] + + 111111111122222222223333333333 + 0123456789012345678901234567890123456789 + 1 [-----------> [--> + <--] <-----------] + ''' + + def test_design_base_pairs(self) -> None: + base_pairs = self.design.base_pairs() + self.assertEqual(len(base_pairs), 21) + + # d01f, d01r + self.assertIn((0, 1), base_pairs) + self.assertIn((0, 2), base_pairs) + + # d02f, d02r + self.assertIn((0, 4), base_pairs) + self.assertIn((0, 5), base_pairs) + self.assertIn((0, 6), base_pairs) + self.assertIn((0, 7), base_pairs) + + # d03f, d02r + self.assertIn((0, 12), base_pairs) + self.assertIn((0, 13), base_pairs) + + # d05f, d04r + self.assertIn((0, 28), base_pairs) + + # d11f, d11r + self.assertIn((1, 8), base_pairs) + self.assertIn((1, 9), base_pairs) + self.assertIn((1, 10), base_pairs) + self.assertIn((1, 11), base_pairs) + + # d11f, d12r + self.assertIn((1, 13), base_pairs) + self.assertIn((1, 14), base_pairs) + self.assertIn((1, 15), base_pairs) + self.assertIn((1, 16), base_pairs) + + # d12f, d12r + self.assertIn((1, 20), base_pairs) + self.assertIn((1, 21), base_pairs) + self.assertIn((1, 22), base_pairs) + self.assertIn((1, 23), base_pairs) + ''' + 111111111122222222223333333333 + 0123456789012345678901234567890123456789 + 0 [-->[--> [--> [--> [--> + <] <--------] <-] <-] <] <] + + 111111111122222222223333333333 + 0123456789012345678901234567890123456789 + 1 [-----------> [--> + <--] <-----------] + ''' From fec6f43bd7817de33b484e17f75f8a3de632c215 Mon Sep 17 00:00:00 2001 From: David Doty Date: Fri, 2 Dec 2022 10:27:36 -0800 Subject: [PATCH 03/14] changed return type of `Design.base_pairs()` to dict and added `allow_mismatches` parameter --- scadnano/scadnano.py | 59 ++++++++++++++++++++++++++++++++++++++------ 1 file changed, 52 insertions(+), 7 deletions(-) diff --git a/scadnano/scadnano.py b/scadnano/scadnano.py index de7ce3fd..2ed8f194 100644 --- a/scadnano/scadnano.py +++ b/scadnano/scadnano.py @@ -2780,7 +2780,7 @@ def with_color(self, color: Color) -> 'StrandBuilder[StrandLabel, DomainLabel]': return self # remove quotes when Py3.6 support dropped - def with_sequence(self, sequence: str, assign_complement: bool = True) \ + def with_sequence(self, sequence: str, assign_complement: bool = False) \ -> 'StrandBuilder[StrandLabel, DomainLabel]': """ Assigns `sequence` as DNA sequence of the :any:`Strand` being built. @@ -2802,7 +2802,7 @@ def with_sequence(self, sequence: str, assign_complement: bool = True) \ return self # remove quotes when Py3.6 support dropped - def with_domain_sequence(self, sequence: str, assign_complement: bool = True) \ + def with_domain_sequence(self, sequence: str, assign_complement: bool = False) \ -> 'StrandBuilder[StrandLabel, DomainLabel]': """ Assigns `sequence` as DNA sequence of the most recently created :any:`Domain` in @@ -4738,6 +4738,43 @@ def find_overlapping_domains_on_helix(helix: Helix) -> List[Tuple[Domain, Domain return overlapping_domains +def bases_complementary(base1: str, base2: str) -> bool: + """ + Indicates if `base1` and `base2` are complementary DNA bases. + + :param base1: + first DNA base + :param base2: + second DNA base + :return: + whether `base1` and `base2` are complementary DNA bases + """ + if len(base1) != 1 or len(base2) != 1: + raise ValueError(f'base1 and base2 must each be a single character: ' + f'base1 = {base1}, base2 = {base2}') + base1 = base1.upper() + base2 = base2.upper() + return {base1, base2} == {'A', 'T'} or {base1, base2} == {'C', 'G'} + + +def reverse_complementary(seq1: str, seq2: str) -> bool: + """ + Indicates if `seq1` and `seq2` are reverse complementary DNA sequences. + + :param seq1: + first DNA sequence + :param seq1: + second DNA sequence + :return: + whether `seq1` and `seq2` are reverse complementary DNA sequences + """ + if len(seq1) != len(seq2): + return False + for b1, b2 in zip(seq1, seq2[::]): + if not bases_complementary(b1, b2): + return False + return True + @dataclass class Design(_JSONSerializable, Generic[StrandLabel, DomainLabel]): """Object representing the entire design of the DNA structure.""" @@ -5324,21 +5361,29 @@ def idx_of(helix: Helix, order: int) -> int: return helices - def base_pairs(self) -> List[Tuple[int, int]]: + def base_pairs(self, allow_mismatches: bool = False) -> Dict[int, List[int]]: """ - List of "addresses" of base pairs in this design. Each address is a pair (`helix_idx`, `offset`) - of the address of the *forward* strand in the base pair. + Base pairs in this design, represented as a dict mapping a :data:`Helix.idx` to a list of offsets + on that helix where two strands are. + :param allow_mismatches: + if True, then all offsets on a :any:`Helix` where there is both a forward and reverse + :any:`Domain` will be included. Otherwise, only offsets where the :any:`Domain`'s have + complementary bases will be included. :return: all base pairs (`helix_idx`, `offset`) in this :any:`Design` """ - base_pairs = [] + base_pairs = {} for idx, helix in self.helices.items(): + offsets = base_pairs[idx] = [] overlapping_domains = find_overlapping_domains_on_helix(helix) for dom1, dom2 in overlapping_domains: start, end = dom1.compute_overlap(dom2) for offset in range(start, end): - base_pairs.append((idx, offset)) + base1 = dom1.dna_sequence_in(offset, offset) + base2 = dom2.dna_sequence_in(offset, offset) + if allow_mismatches or bases_complementary(base1, base2): + offsets.append(offset) return base_pairs From e79562700f6c4a2661e130d2409f15cb27454b06 Mon Sep 17 00:00:00 2001 From: David Doty Date: Fri, 2 Dec 2022 10:28:03 -0800 Subject: [PATCH 04/14] fixed unit tests after updating `with_sequence` and `with_domain_sequence` to default to not assigning complement --- tests/scadnano_tests.py | 157 +++++++++++++++++++++++++++------------- 1 file changed, 108 insertions(+), 49 deletions(-) diff --git a/tests/scadnano_tests.py b/tests/scadnano_tests.py index 6a0e8b02..5035144c 100644 --- a/tests/scadnano_tests.py +++ b/tests/scadnano_tests.py @@ -6638,7 +6638,7 @@ def test_method_chaining_with_domain_sequence(self) -> None: None: 098 765 432 109 876 543 210 """ sb = design.draw_strand(0, 9).to(12) - sb.with_domain_sequence('AAC') + sb.with_domain_sequence('AAC', assign_complement=True) sb.cross(0, offset=3) sb.to(6) - sb.with_domain_sequence('TTC') + sb.with_domain_sequence('TTC', assign_complement=True) self.assertEqual('AAC TTC'.replace(' ', ''), design.strands[-1].dna_sequence) self.assertEqual('??? ??? GCA GTT ??? GAA ???'.replace(' ', ''), design.strands[0].dna_sequence) @@ -6682,12 +6682,12 @@ def test_method_chaining_with_domain_sequence(self) -> None: 098 765 432 109 876 543 210 """ design.draw_strand(0, 6).to(9) \ - .with_domain_sequence('GGA') \ + .with_domain_sequence('GGA', assign_complement=True) \ .cross(0, offset=15) \ .to(18) \ - .with_domain_sequence('TTG') \ + .with_domain_sequence('TTG', assign_complement=True) \ .to(21) \ - .with_domain_sequence('GCA') + .with_domain_sequence('GCA', assign_complement=True) self.assertEqual('GGA TTG GCA'.replace(' ', ''), design.strands[-1].dna_sequence) self.assertEqual('TGC CAA GCA GTT TCC GAA ???'.replace(' ', ''), design.strands[0].dna_sequence) @@ -7004,10 +7004,12 @@ def test_export_file(self): oxv = design.to_oxview_format(use_strand_colors=True) - with tempfile.NamedTemporaryFile(mode='w', suffix='.json') as f: + with tempfile.NamedTemporaryFile(mode='w', suffix='.json', delete=False) as f: design.write_oxview_file(filename=f.name) - with open(f.name, 'r') as f2: - oxv2 = json.load(f2) + filename = f.name + with open(filename, 'r') as f2: + oxv2 = json.load(f2) + os.unlink(filename) # The dates won't be equal, so delete them del oxv['date'] @@ -7799,37 +7801,39 @@ def test_to_json_serializable__label_key_contains_non_default_name(self) -> None class TestBasePairs(unittest.TestCase): def setUp(self) -> None: ''' + X shows position of mismatches 111111111122222222223333333333 0123456789012345678901234567890123456789 0 [-->[--> [--> [--> [--> <] <--------] <-] <-] <] <] - + X 111111111122222222223333333333 0123456789012345678901234567890123456789 1 [-----------> [--> <--] <-----------] + X ''' helices = [sc.Helix(max_offset=40) for _ in range(2)] self.design = sc.Design(helices=helices) # helix 0 forward - self.design.draw_strand(0, 0).move(4) - self.design.draw_strand(0, 4).move(4) - self.design.draw_strand(0, 12).move(4) - self.design.draw_strand(0, 20).move(4) - self.design.draw_strand(0, 28).move(4) + self.design.draw_strand(0, 0).move(4).with_sequence('AAAA') + self.design.draw_strand(0, 4).move(4).with_sequence('AAAA') + self.design.draw_strand(0, 12).move(4).with_sequence('AAAA') + self.design.draw_strand(0, 20).move(4).with_sequence('AAAA') + self.design.draw_strand(0, 28).move(4).with_sequence('AAAA') # helix 0 reverse - self.design.draw_strand(0, 3).move(-2) - self.design.draw_strand(0, 14).to(4) - self.design.draw_strand(0, 20).to(17) - self.design.draw_strand(0, 29).to(26) - self.design.draw_strand(0, 36).to(34) - self.design.draw_strand(0, 39).to(37) + self.design.draw_strand(0, 3).move(-2).with_sequence('TT') + self.design.draw_strand(0, 14).to(4).with_sequence('TTTTTTTTCT') + self.design.draw_strand(0, 20).to(17).with_sequence('TTT') + self.design.draw_strand(0, 29).to(26).with_sequence('TTT') + self.design.draw_strand(0, 36).to(34).with_sequence('TT') + self.design.draw_strand(0, 39).to(37).with_sequence('TT') # helix 1 forward - self.design.draw_strand(1, 4).to(17) - self.design.draw_strand(1, 20).to(24) + self.design.draw_strand(1, 4).to(17).with_sequence('A'*13) + self.design.draw_strand(1, 20).to(24).with_sequence('A'*4) # helix 1 reverse - self.design.draw_strand(1, 12).to(8) - self.design.draw_strand(1, 26).to(13) + self.design.draw_strand(1, 12).to(8).with_sequence('TGTT') + self.design.draw_strand(1, 26).to(13).with_sequence('T'*13) def test_find_overlapping_domains(self) -> None: d01f = self.design.strands[0].domains[0] @@ -7878,44 +7882,46 @@ def test_find_overlapping_domains(self) -> None: <--] <-----------] ''' - def test_design_base_pairs(self) -> None: - base_pairs = self.design.base_pairs() - self.assertEqual(len(base_pairs), 21) + def test_design_base_pairs_mismatches(self) -> None: + base_pairs = self.design.base_pairs(allow_mismatches=True) + self.assertEqual(len(base_pairs), 2) + self.assertEqual(len(base_pairs[0]), 9) + self.assertEqual(len(base_pairs[1]), 12) # d01f, d01r - self.assertIn((0, 1), base_pairs) - self.assertIn((0, 2), base_pairs) + self.assertIn(1, base_pairs[0]) + self.assertIn(2, base_pairs[0]) # d02f, d02r - self.assertIn((0, 4), base_pairs) - self.assertIn((0, 5), base_pairs) - self.assertIn((0, 6), base_pairs) - self.assertIn((0, 7), base_pairs) + self.assertIn(4, base_pairs[0]) + self.assertIn(5, base_pairs[0]) + self.assertIn(6, base_pairs[0]) + self.assertIn(7, base_pairs[0]) # d03f, d02r - self.assertIn((0, 12), base_pairs) - self.assertIn((0, 13), base_pairs) + self.assertIn(12, base_pairs[0]) + self.assertIn(13, base_pairs[0]) # d05f, d04r - self.assertIn((0, 28), base_pairs) + self.assertIn(28, base_pairs[0]) # d11f, d11r - self.assertIn((1, 8), base_pairs) - self.assertIn((1, 9), base_pairs) - self.assertIn((1, 10), base_pairs) - self.assertIn((1, 11), base_pairs) + self.assertIn(8, base_pairs[1]) + self.assertIn(9, base_pairs[1]) + self.assertIn(10, base_pairs[1]) + self.assertIn(11, base_pairs[1]) # d11f, d12r - self.assertIn((1, 13), base_pairs) - self.assertIn((1, 14), base_pairs) - self.assertIn((1, 15), base_pairs) - self.assertIn((1, 16), base_pairs) + self.assertIn(13, base_pairs[1]) + self.assertIn(14, base_pairs[1]) + self.assertIn(15, base_pairs[1]) + self.assertIn(16, base_pairs[1]) # d12f, d12r - self.assertIn((1, 20), base_pairs) - self.assertIn((1, 21), base_pairs) - self.assertIn((1, 22), base_pairs) - self.assertIn((1, 23), base_pairs) + self.assertIn(20, base_pairs[1]) + self.assertIn(21, base_pairs[1]) + self.assertIn(22, base_pairs[1]) + self.assertIn(23, base_pairs[1]) ''' 111111111122222222223333333333 0123456789012345678901234567890123456789 @@ -7927,3 +7933,56 @@ def test_design_base_pairs(self) -> None: 1 [-----------> [--> <--] <-----------] ''' + def test_design_base_pairs_no_mismatches(self) -> None: + base_pairs = self.design.base_pairs(allow_mismatches=False) + self.assertEqual(len(base_pairs), 2) + self.assertEqual(len(base_pairs[0]), 8) + self.assertEqual(len(base_pairs[1]), 11) + + # d01f, d01r + self.assertIn(1, base_pairs[0]) + self.assertIn(2, base_pairs[0]) + + # d02f, d02r + self.assertIn(4, base_pairs[0]) + # self.assertIn(5, base_pairs[0]) # mismatch + self.assertIn(6, base_pairs[0]) + self.assertIn(7, base_pairs[0]) + + # d03f, d02r + self.assertIn(12, base_pairs[0]) + self.assertIn(13, base_pairs[0]) + + # d05f, d04r + self.assertIn(28, base_pairs[0]) + + # d11f, d11r + self.assertIn(8, base_pairs[1]) + self.assertIn(9, base_pairs[1]) + # self.assertIn(10, base_pairs[1]) # mismatch + self.assertIn(11, base_pairs[1]) + + # d11f, d12r + self.assertIn(13, base_pairs[1]) + self.assertIn(14, base_pairs[1]) + self.assertIn(15, base_pairs[1]) + self.assertIn(16, base_pairs[1]) + + # d12f, d12r + self.assertIn(20, base_pairs[1]) + self.assertIn(21, base_pairs[1]) + self.assertIn(22, base_pairs[1]) + self.assertIn(23, base_pairs[1]) + ''' + X shows position of mismatches + 111111111122222222223333333333 + 0123456789012345678901234567890123456789 + 0 [-->[--> [--> [--> [--> + <] <--------] <-] <-] <] <] + X + 111111111122222222223333333333 + 0123456789012345678901234567890123456789 + 1 [-----------> [--> + <--] <-----------] + X + ''' From 6b6001a030859a72052400a7d537d5bd013de1b9 Mon Sep 17 00:00:00 2001 From: David Doty Date: Fri, 2 Dec 2022 10:29:14 -0800 Subject: [PATCH 05/14] removed unused _popleft() function --- scadnano/scadnano.py | 10 +--------- 1 file changed, 1 insertion(+), 9 deletions(-) diff --git a/scadnano/scadnano.py b/scadnano/scadnano.py index 2ed8f194..2cdb2e14 100644 --- a/scadnano/scadnano.py +++ b/scadnano/scadnano.py @@ -4672,15 +4672,6 @@ def _check_type_is_one_of(obj: Any, expected_types: Iterable) -> None: f'but instead it is of type {type(obj)}') -def _popleft(queue: Deque[T]) -> Optional[T]: - # return left element of queue if queue is nonempty, otherwise return None - try: - elt = queue.popleft() - return elt - except IndexError: - return None - - def find_overlapping_domains_on_helix(helix: Helix) -> List[Tuple[Domain, Domain]]: # compute list of pairs of domains that overlap on Helix `helix` # assumes that `helix.domains` has been populated by calling `Design._build_domains_on_helix_lists()` @@ -4775,6 +4766,7 @@ def reverse_complementary(seq1: str, seq2: str) -> bool: return False return True + @dataclass class Design(_JSONSerializable, Generic[StrandLabel, DomainLabel]): """Object representing the entire design of the DNA structure.""" From e55609da734df57222c37c419b76e1db38917c17 Mon Sep 17 00:00:00 2001 From: David Doty Date: Fri, 2 Dec 2022 11:17:26 -0800 Subject: [PATCH 06/14] Update scadnano.py --- scadnano/scadnano.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/scadnano/scadnano.py b/scadnano/scadnano.py index 4e579713..0c509ab8 100644 --- a/scadnano/scadnano.py +++ b/scadnano/scadnano.py @@ -4715,10 +4715,10 @@ def find_overlapping_domains_on_helix(helix: Helix) -> List[Tuple[Domain, Domain # [--r_dom---> # reverse_domain can't overlap *next* forward_domain, so safe to remove reverse_domains.popleft() - if len(reverse_domains) == 0: - break - else: + if len(reverse_domains) > 0: reverse_domain = reverse_domains[0] + else: + break else: # [---f_dom---> [---next_f_dom--> # [----r_dom->[--next_r_dom----> From afc44b8c8a3a627a0b70da48595710d0b11bd2ea Mon Sep 17 00:00:00 2001 From: David Doty Date: Fri, 2 Dec 2022 12:16:24 -0800 Subject: [PATCH 07/14] don't put Helix idx in base_pairs dict if Helix has no base pairs --- scadnano/scadnano.py | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/scadnano/scadnano.py b/scadnano/scadnano.py index 0c509ab8..7019e685 100644 --- a/scadnano/scadnano.py +++ b/scadnano/scadnano.py @@ -5358,6 +5358,8 @@ def base_pairs(self, allow_mismatches: bool = False) -> Dict[int, List[int]]: Base pairs in this design, represented as a dict mapping a :data:`Helix.idx` to a list of offsets on that helix where two strands are. + If a :any:`Helix` has no base pairs, then its :data:`Helix.idx` is not a key in the returned dict. + :param allow_mismatches: if True, then all offsets on a :any:`Helix` where there is both a forward and reverse :any:`Domain` will be included. Otherwise, only offsets where the :any:`Domain`'s have @@ -5367,7 +5369,7 @@ def base_pairs(self, allow_mismatches: bool = False) -> Dict[int, List[int]]: """ base_pairs = {} for idx, helix in self.helices.items(): - offsets = base_pairs[idx] = [] + offsets = [] overlapping_domains = find_overlapping_domains_on_helix(helix) for dom1, dom2 in overlapping_domains: start, end = dom1.compute_overlap(dom2) @@ -5376,6 +5378,8 @@ def base_pairs(self, allow_mismatches: bool = False) -> Dict[int, List[int]]: base2 = dom2.dna_sequence_in(offset, offset) if allow_mismatches or bases_complementary(base1, base2): offsets.append(offset) + if len(offsets) > 0: + base_pairs[idx] = offsets return base_pairs From 1049cdc45e3502874c5a1affdfa073f5f54bc7fe Mon Sep 17 00:00:00 2001 From: David Doty Date: Fri, 2 Dec 2022 13:26:05 -0800 Subject: [PATCH 08/14] fixed base pair calculation to allow for deletions/insertions/wildcards/None if DNA sequence is not assigned --- scadnano/scadnano.py | 40 ++++++++++++++--- tests/scadnano_tests.py | 97 +++++++++++++++++++++++++++++++++++++++-- 2 files changed, 128 insertions(+), 9 deletions(-) diff --git a/scadnano/scadnano.py b/scadnano/scadnano.py index 7019e685..02b8e997 100644 --- a/scadnano/scadnano.py +++ b/scadnano/scadnano.py @@ -4729,7 +4729,8 @@ def find_overlapping_domains_on_helix(helix: Helix) -> List[Tuple[Domain, Domain return overlapping_domains -def bases_complementary(base1: str, base2: str) -> bool: +def bases_complementary(base1: str, base2: str, allow_wildcard: bool = False, + allow_none: bool = False) -> bool: """ Indicates if `base1` and `base2` are complementary DNA bases. @@ -4737,9 +4738,21 @@ def bases_complementary(base1: str, base2: str) -> bool: first DNA base :param base2: second DNA base + :param allow_wildcard: + if true a "wildcard" (the symbol '?') is considered to be complementary to anything + :param allow_none: + if true the object None is considered to be complementary to anything :return: whether `base1` and `base2` are complementary DNA bases """ + if allow_none and (base1 is None or base2 is None): + return True + elif not allow_none and (base1 is None or base2 is None): + return False + + if allow_wildcard and (base1 == DNA_base_wildcard or base2 == DNA_base_wildcard): + return True + if len(base1) != 1 or len(base2) != 1: raise ValueError(f'base1 and base2 must each be a single character: ' f'base1 = {base1}, base2 = {base2}') @@ -4748,7 +4761,8 @@ def bases_complementary(base1: str, base2: str) -> bool: return {base1, base2} == {'A', 'T'} or {base1, base2} == {'C', 'G'} -def reverse_complementary(seq1: str, seq2: str) -> bool: +def reverse_complementary(seq1: str, seq2: str, allow_wildcard: bool = False, + allow_none: bool = False) -> bool: """ Indicates if `seq1` and `seq2` are reverse complementary DNA sequences. @@ -4756,13 +4770,22 @@ def reverse_complementary(seq1: str, seq2: str) -> bool: first DNA sequence :param seq1: second DNA sequence + :param allow_wildcard: + if true a "wildcard" (the symbol '?') is considered to be complementary to anything + :param allow_none: + if true the object None is considered to be complementary to anything :return: whether `seq1` and `seq2` are reverse complementary DNA sequences """ + if allow_none and (seq1 is None or seq2 is None): + return True + elif not allow_none and (seq1 is None or seq2 is None): + return False + if len(seq1) != len(seq2): return False for b1, b2 in zip(seq1, seq2[::]): - if not bases_complementary(b1, b2): + if not bases_complementary(b1, b2, allow_wildcard, allow_none): return False return True @@ -5374,9 +5397,14 @@ def base_pairs(self, allow_mismatches: bool = False) -> Dict[int, List[int]]: for dom1, dom2 in overlapping_domains: start, end = dom1.compute_overlap(dom2) for offset in range(start, end): - base1 = dom1.dna_sequence_in(offset, offset) - base2 = dom2.dna_sequence_in(offset, offset) - if allow_mismatches or bases_complementary(base1, base2): + if offset in dom1.deletions or offset in dom2.deletions: + continue + seq1 = dom1.dna_sequence_in(offset, offset) + seq2 = dom2.dna_sequence_in(offset, offset) + # we use reverse_complementary instead of base_complementary here to allow for insertions + # that may give a larger DNA sequence than length 1 at a given offset + if allow_mismatches or reverse_complementary(seq1, seq2, + allow_wildcard=True, allow_none=True): offsets.append(offset) if len(offsets) > 0: base_pairs[idx] = offsets diff --git a/tests/scadnano_tests.py b/tests/scadnano_tests.py index 5035144c..29b0dea7 100644 --- a/tests/scadnano_tests.py +++ b/tests/scadnano_tests.py @@ -7829,11 +7829,11 @@ def setUp(self) -> None: self.design.draw_strand(0, 36).to(34).with_sequence('TT') self.design.draw_strand(0, 39).to(37).with_sequence('TT') # helix 1 forward - self.design.draw_strand(1, 4).to(17).with_sequence('A'*13) - self.design.draw_strand(1, 20).to(24).with_sequence('A'*4) + self.design.draw_strand(1, 4).to(17).with_sequence('A' * 13) + self.design.draw_strand(1, 20).to(24).with_sequence('A' * 4) # helix 1 reverse self.design.draw_strand(1, 12).to(8).with_sequence('TGTT') - self.design.draw_strand(1, 26).to(13).with_sequence('T'*13) + self.design.draw_strand(1, 26).to(13).with_sequence('T' * 13) def test_find_overlapping_domains(self) -> None: d01f = self.design.strands[0].domains[0] @@ -7933,6 +7933,7 @@ def test_design_base_pairs_mismatches(self) -> None: 1 [-----------> [--> <--] <-----------] ''' + def test_design_base_pairs_no_mismatches(self) -> None: base_pairs = self.design.base_pairs(allow_mismatches=False) self.assertEqual(len(base_pairs), 2) @@ -7986,3 +7987,93 @@ def test_design_base_pairs_no_mismatches(self) -> None: <--] <-----------] X ''' + + def test_design_base_pairs_no_dna(self) -> None: + ''' + 0123456789 + 0 [--------> + <---]<---] + ''' + design = sc.Design(helices=[sc.Helix(max_offset=100)]) + design.draw_strand(0, 0).move(10) + design.draw_strand(0, 5).move(-5) + design.draw_strand(0, 10).move(-5) + + base_pairs = design.base_pairs() + self.assertEqual(len(base_pairs), 1) + self.assertEqual(len(base_pairs[0]), 10) + + for offset in range(10): + self.assertIn(offset, base_pairs[0]) + + def test_design_base_pairs_dna_on_some_strands_and_mismatches(self) -> None: + ''' + 0123456789 + AAAAAAAAAA + 0 [--------> + <---]<---] + TTCTT + ''' + design = sc.Design(helices=[sc.Helix(max_offset=100)]) + design.draw_strand(0, 0).move(10).with_sequence('A'*10) + design.draw_strand(0, 5).move(-5).with_sequence('TTCTT') + design.draw_strand(0, 10).move(-5) + + base_pairs = design.base_pairs(allow_mismatches=False) + self.assertEqual(len(base_pairs), 1) + self.assertEqual(len(base_pairs[0]), 9) + + for offset in range(10): + if offset != 2: + self.assertIn(offset, base_pairs[0]) + + def test_design_base_pairs_deletions_insertions(self) -> None: + ''' + 0123456789 + AA + A AAAAAAA + 0 [XX---II-> + <-XX]<-II] + TT TTTTTT + TT + ''' + design = sc.Design(helices=[sc.Helix(max_offset=100)]) + design.draw_strand(0, 0).move(10).with_deletions([1, 2]).with_insertions([(6, 1), (7, 1)]) \ + .with_sequence('A'*10) + design.draw_strand(0, 5).move(-5).with_deletions([2, 3]).with_sequence('TTT') + design.draw_strand(0, 10).move(-5).with_insertions([(7, 1), (8, 1)]).with_sequence('TTTTTTT') + + base_pairs = design.base_pairs(allow_mismatches=False) + self.assertEqual(len(base_pairs), 1) + self.assertEqual(len(base_pairs[0]), 5) + + self.assertIn(0, base_pairs[0]) + self.assertIn(4, base_pairs[0]) + self.assertIn(5, base_pairs[0]) + self.assertIn(7, base_pairs[0]) + self.assertIn(9, base_pairs[0]) + + def test_design_base_pairs_deletions_insertions_mismatch_in_insertion(self) -> None: + ''' + 0123456789 + AA + A AAAAAAA + 0 [XX---II-> + <-XX]<-II] + TT TTTTTT + CT + ''' + design = sc.Design(helices=[sc.Helix(max_offset=100)]) + design.draw_strand(0, 0).move(10).with_deletions([1, 2]).with_insertions([(6, 1), (7, 1)]) \ + .with_sequence('A'*10) + design.draw_strand(0, 5).move(-5).with_deletions([2, 3]).with_sequence('TTT') + design.draw_strand(0, 10).move(-5).with_insertions([(7, 1), (8, 1)]).with_sequence('TTTCTTT') + + base_pairs = design.base_pairs(allow_mismatches=False) + self.assertEqual(len(base_pairs), 1) + self.assertEqual(len(base_pairs[0]), 4) + + self.assertIn(0, base_pairs[0]) + self.assertIn(4, base_pairs[0]) + self.assertIn(5, base_pairs[0]) + self.assertIn(9, base_pairs[0]) From e6c07379223a7fc59d63fc113ee156d34b5241b3 Mon Sep 17 00:00:00 2001 From: David Doty Date: Fri, 2 Dec 2022 13:58:28 -0800 Subject: [PATCH 09/14] added unit test for no base pairs --- tests/scadnano_tests.py | 14 ++++++++++++++ 1 file changed, 14 insertions(+) diff --git a/tests/scadnano_tests.py b/tests/scadnano_tests.py index 29b0dea7..02c5baed 100644 --- a/tests/scadnano_tests.py +++ b/tests/scadnano_tests.py @@ -8077,3 +8077,17 @@ def test_design_base_pairs_deletions_insertions_mismatch_in_insertion(self) -> N self.assertIn(4, base_pairs[0]) self.assertIn(5, base_pairs[0]) self.assertIn(9, base_pairs[0]) + + def test_no_base_pairs(self) -> None: + ''' + 0123456789 + [--> + <--] + ''' + design = sc.Design(helices=[sc.Helix(max_offset=100)]) + design.draw_strand(0, 0).move(4) + design.draw_strand(0, 9).move(-4) + + base_pairs = design.base_pairs(allow_mismatches=False) + self.assertEqual(len(base_pairs), 1) + self.assertEqual(len(base_pairs[0]), 0) From f81b8ff82a218a92cb2899bff2b2e7857ff65e01 Mon Sep 17 00:00:00 2001 From: David Doty Date: Fri, 2 Dec 2022 13:58:59 -0800 Subject: [PATCH 10/14] Update scadnano_tests.py --- tests/scadnano_tests.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/tests/scadnano_tests.py b/tests/scadnano_tests.py index 02c5baed..58edd294 100644 --- a/tests/scadnano_tests.py +++ b/tests/scadnano_tests.py @@ -8089,5 +8089,4 @@ def test_no_base_pairs(self) -> None: design.draw_strand(0, 9).move(-4) base_pairs = design.base_pairs(allow_mismatches=False) - self.assertEqual(len(base_pairs), 1) - self.assertEqual(len(base_pairs[0]), 0) + self.assertEqual(len(base_pairs), 0) From 7b83a43bfd1da6c96a1095f0c1e2eb966eac9a0a Mon Sep 17 00:00:00 2001 From: David Doty Date: Fri, 2 Dec 2022 14:04:28 -0800 Subject: [PATCH 11/14] fixed bug in `Design.base_pairs()` when no reverse domains on helix --- scadnano/scadnano.py | 3 +++ tests/scadnano_tests.py | 11 +++++++++++ 2 files changed, 14 insertions(+) diff --git a/scadnano/scadnano.py b/scadnano/scadnano.py index 02b8e997..86e9f114 100644 --- a/scadnano/scadnano.py +++ b/scadnano/scadnano.py @@ -4686,6 +4686,9 @@ def find_overlapping_domains_on_helix(helix: Helix) -> List[Tuple[Domain, Domain forward_domains.sort(key=lambda domain: domain.start) reverse_domains.sort(key=lambda domain: domain.start) + if len(forward_domains) == 0 or len(reverse_domains) == 0: + return [] + # need to be efficient to remove the front element repeatedly reverse_domains = collections.deque(reverse_domains) diff --git a/tests/scadnano_tests.py b/tests/scadnano_tests.py index 58edd294..9f8bff38 100644 --- a/tests/scadnano_tests.py +++ b/tests/scadnano_tests.py @@ -8090,3 +8090,14 @@ def test_no_base_pairs(self) -> None: base_pairs = design.base_pairs(allow_mismatches=False) self.assertEqual(len(base_pairs), 0) + + def test_no_base_pairs_only_forward_strand(self) -> None: + ''' + 0123456789 + [--> + ''' + design = sc.Design(helices=[sc.Helix(max_offset=100)]) + design.draw_strand(0, 0).move(4) + + base_pairs = design.base_pairs(allow_mismatches=False) + self.assertEqual(len(base_pairs), 0) From 5312a59eddfa0d19535213e98147a45c4118ab6d Mon Sep 17 00:00:00 2001 From: David Doty Date: Fri, 2 Dec 2022 16:23:12 -0800 Subject: [PATCH 12/14] Update scadnano_tests.py --- tests/scadnano_tests.py | 9 ++++----- 1 file changed, 4 insertions(+), 5 deletions(-) diff --git a/tests/scadnano_tests.py b/tests/scadnano_tests.py index 9f8bff38..7f4973da 100644 --- a/tests/scadnano_tests.py +++ b/tests/scadnano_tests.py @@ -7856,7 +7856,6 @@ def test_find_overlapping_domains(self) -> None: d12r = self.design.strands[14].domains[0] overlapping_domains_h0 = sc.find_overlapping_domains_on_helix(self.design.helices[0]) - overlapping_domains_h1 = sc.find_overlapping_domains_on_helix(self.design.helices[1]) self.assertEqual(len(overlapping_domains_h0), 4) @@ -8015,7 +8014,7 @@ def test_design_base_pairs_dna_on_some_strands_and_mismatches(self) -> None: TTCTT ''' design = sc.Design(helices=[sc.Helix(max_offset=100)]) - design.draw_strand(0, 0).move(10).with_sequence('A'*10) + design.draw_strand(0, 0).move(10).with_sequence('A' * 10) design.draw_strand(0, 5).move(-5).with_sequence('TTCTT') design.draw_strand(0, 10).move(-5) @@ -8039,9 +8038,9 @@ def test_design_base_pairs_deletions_insertions(self) -> None: ''' design = sc.Design(helices=[sc.Helix(max_offset=100)]) design.draw_strand(0, 0).move(10).with_deletions([1, 2]).with_insertions([(6, 1), (7, 1)]) \ - .with_sequence('A'*10) + .with_sequence('A' * 10) design.draw_strand(0, 5).move(-5).with_deletions([2, 3]).with_sequence('TTT') - design.draw_strand(0, 10).move(-5).with_insertions([(7, 1), (8, 1)]).with_sequence('TTTTTTT') + design.draw_strand(0, 10).move(-5).with_insertions([(7, 1), (8, 1)]).with_sequence('T' * 7) base_pairs = design.base_pairs(allow_mismatches=False) self.assertEqual(len(base_pairs), 1) @@ -8065,7 +8064,7 @@ def test_design_base_pairs_deletions_insertions_mismatch_in_insertion(self) -> N ''' design = sc.Design(helices=[sc.Helix(max_offset=100)]) design.draw_strand(0, 0).move(10).with_deletions([1, 2]).with_insertions([(6, 1), (7, 1)]) \ - .with_sequence('A'*10) + .with_sequence('A' * 10) design.draw_strand(0, 5).move(-5).with_deletions([2, 3]).with_sequence('TTT') design.draw_strand(0, 10).move(-5).with_insertions([(7, 1), (8, 1)]).with_sequence('TTTCTTT') From 1b7e605026ebadee59eed6394394d889ed041460 Mon Sep 17 00:00:00 2001 From: David Doty Date: Sat, 3 Dec 2022 12:18:04 -0800 Subject: [PATCH 13/14] closes #238: Domain colors --- scadnano/scadnano.py | 83 ++++++++++++++++++++++++++++++++++++----- tests/scadnano_tests.py | 15 +++++--- 2 files changed, 83 insertions(+), 15 deletions(-) diff --git a/scadnano/scadnano.py b/scadnano/scadnano.py index 86e9f114..7128d829 100644 --- a/scadnano/scadnano.py +++ b/scadnano/scadnano.py @@ -197,6 +197,25 @@ def to_json_serializable(self, suppress_indent: bool = True, **kwargs: Any) -> s # return NoIndent(self.__dict__) if suppress_indent else self.__dict__ return f'#{self.r:02x}{self.g:02x}{self.b:02x}' + @staticmethod + def from_json(color_json: Union[str, int, None]) -> Union['Color', None]: + if color_json is None: + return None + + color_str: str + if isinstance(color_json, int): + def decimal_int_to_hex(d: int) -> str: + return "#" + "{0:#08x}".format(d, 8)[2:] # type: ignore + + color_str = decimal_int_to_hex(color_json) + elif isinstance(color_json, str): + color_str = color_json + else: + raise IllegalDesignError(f'color must be a string or int, ' + f'but it is a {type(color_json)}: {color_json}') + color = Color(hex_string=color_str) + return color + def to_cadnano_v2_int_hex(self) -> int: return int(f'{self.r:02x}{self.g:02x}{self.b:02x}', 16) @@ -1679,6 +1698,11 @@ class Domain(_JSONSerializable, Generic[DomainLabel]): """DNA sequence of this Domain, or ``None`` if no DNA sequence has been assigned to this :any:`Domain`'s :any:`Strand`.""" + color: Optional[Color] = None + """ + Color to show this domain in the main view. If specified, overrides the field :data:`Strand.color`. + """ + # not serialized; for efficiency # remove quotes when Py3.6 support dropped _parent_strand: Optional['Strand'] = field(init=False, repr=False, compare=False, default=None) @@ -1701,6 +1725,8 @@ def to_json_serializable(self, suppress_indent: bool = True, dct[insertions_key] = self.insertions if self.label is not None: dct[domain_label_key] = self.label + if self.color is not None: + dct[color_key] = self.color.to_json_serializable(suppress_indent) return NoIndent(dct) if suppress_indent else dct @staticmethod @@ -1714,6 +1740,8 @@ def from_json(json_map: Dict[str, Any]) -> 'Domain': # remove quotes when Py3.6 list(map(tuple, json_map.get(insertions_key, [])))) # type: ignore name = json_map.get(domain_name_key) label = json_map.get(domain_label_key) + color_json = json_map.get(color_key) + color = Color.from_json(color_json) return Domain( helix=helix, forward=forward, @@ -1723,6 +1751,7 @@ def from_json(json_map: Dict[str, Any]) -> 'Domain': # remove quotes when Py3.6 insertions=insertions, name=name, label=label, + color=color, ) def __repr__(self) -> str: @@ -1734,6 +1763,7 @@ def __repr__(self) -> str: f', end={self.end}') + \ (f', deletions={self.deletions}' if len(self.deletions) > 0 else '') + \ (f', insertions={self.insertions}' if len(self.insertions) > 0 else '') + \ + (f', color={self.color}' if self.color is not None else '') + \ ')' return rep @@ -2077,6 +2107,11 @@ class Loopout(_JSONSerializable, Generic[DomainLabel]): dna_sequence: Optional[str] = None """DNA sequence of this :any:`Loopout`, or ``None`` if no DNA sequence has been assigned.""" + color: Optional[Color] = None + """ + Color to show this loopout in the main view. If specified, overrides the field :data:`Strand.color`. + """ + # not serialized; for efficiency # remove quotes when Py3.6 support dropped _parent_strand: Optional['Strand'] = field(init=False, repr=False, compare=False, default=None) @@ -2088,6 +2123,8 @@ def to_json_serializable(self, suppress_indent: bool = True, dct[domain_name_key] = self.name if self.label is not None: dct[domain_label_key] = self.label + if self.color is not None: + dct[color_key] = self.color.to_json_serializable(suppress_indent) return NoIndent(dct) if suppress_indent else dct @staticmethod @@ -2098,7 +2135,9 @@ def from_json(json_map: Dict[str, Any]) -> 'Loopout': # remove quotes when Py3. length = int(length_str) name = json_map.get(domain_name_key) label = json_map.get(domain_label_key) - return Loopout(length=length, name=name, label=label) + color_json = json_map.get(color_key) + color = Color.from_json(color_json) + return Loopout(length=length, name=name, label=label, color=color) def strand(self) -> 'Strand': # remove quotes when Py3.6 support dropped """ @@ -2232,6 +2271,11 @@ class Extension(_JSONSerializable, Generic[DomainLabel]): dna_sequence: Optional[str] = None """DNA sequence of this :any:`Extension`, or ``None`` if no DNA sequence has been assigned.""" + color: Optional[Color] = None + """ + Color to show this extension in the main view. If specified, overrides the field :data:`Strand.color`. + """ + # not serialized; for efficiency # remove quotes when Py3.6 support dropped _parent_strand: Optional['Strand'] = field(init=False, repr=False, compare=False, default=None) @@ -2243,6 +2287,8 @@ def to_json_serializable(self, suppress_indent: bool = True, **kwargs: Any) \ self._add_display_angle_if_not_default(json_map) self._add_name_if_not_default(json_map) self._add_label_if_not_default(json_map) + if self.color is not None: + json_map[color_key] = self.color.to_json_serializable(suppress_indent) return NoIndent(json_map) if suppress_indent else json_map def dna_length(self) -> int: @@ -2269,11 +2315,16 @@ def float_transformer(x): return float(x) json_map, display_angle_key, transformer=float_transformer) name = json_map.get(domain_name_key) label = json_map.get(domain_label_key) + color_json = json_map.get(color_key) + color = Color.from_json(color_json) return Extension( num_bases=num_bases, display_length=display_length, display_angle=display_angle, - label=label, name=name) + label=label, + name=name, + color=color, + ) def _add_display_length_if_not_default(self, json_map) -> None: self._add_key_value_to_json_map_if_not_default( @@ -2839,6 +2890,23 @@ def with_domain_sequence(self, sequence: str, assign_complement: bool = False) \ assign_complement=assign_complement) return self + # remove quotes when Py3.6 support dropped + def with_domain_color(self, color: Color) -> 'StrandBuilder[StrandLabel, DomainLabel]': + """ + Sets most recent :any:`Domain`/:any:`Loopout`/:any:`Extension` + to have given color. + + :param color: + color to set for :any:`Domain`/:any:`Loopout`/:any:`Extension` + :return: + self + """ + if self._strand is None: + raise ValueError('no Strand created yet; make at least one domain first') + last_domain = self._strand.domains[-1] + last_domain.color = color + return self + # remove quotes when Py3.6 support dropped def with_name(self, name: str) -> 'StrandBuilder[StrandLabel, DomainLabel]': """ @@ -3253,17 +3321,12 @@ def from_json(json_map: dict) -> 'Strand': # remove quotes when Py3.6 support d dna_sequence = optional_field(None, json_map, dna_sequence_key, legacy_keys=legacy_dna_sequence_keys) - color_str = json_map.get(color_key) + color_json = json_map.get(color_key) - if color_str is None: + if color_json is None: color = default_scaffold_color if is_scaffold else default_strand_color else: - if isinstance(color_str, int): - def decimal_int_to_hex(d: int) -> str: - return "#" + "{0:#08x}".format(d, 8)[2:] # type: ignore - - color_str = decimal_int_to_hex(color_str) - color = Color(hex_string=color_str) + color = Color.from_json(color_json) label = json_map.get(strand_label_key) diff --git a/tests/scadnano_tests.py b/tests/scadnano_tests.py index 7f4973da..0f5ecf12 100644 --- a/tests/scadnano_tests.py +++ b/tests/scadnano_tests.py @@ -67,32 +67,36 @@ def test_strand__loopouts_with_labels(self) -> None: self.assertEqual(1, len(design.strands)) self.assertEqual(expected_strand, design.strands[0]) - def test_strand__loopouts_with_labels_to_json(self) -> None: + def test_strand__loopouts_with_labels_and_colors_to_json(self) -> None: design = self.design_6helix sb = design.draw_strand(0, 0) sb.to(10) sb.loopout(1, 8) + sb.with_domain_color(sc.Color(10, 10, 10)) sb.with_domain_label('loop0') sb.to(5) sb.with_domain_label('dom1') + sb.with_domain_color(sc.Color(20, 20, 20)) sb.cross(2) sb.to(10) sb.with_domain_label('dom2') sb.loopout(3, 12) sb.with_domain_label('loop1') sb.to(5) + sb.with_color(sc.Color(30, 30, 30)) design_json_map = design.to_json_serializable(suppress_indent=False) design_from_json = sc.Design.from_scadnano_json_map(design_json_map) expected_strand = sc.Strand([ sc.Domain(0, True, 0, 10), - sc.Loopout(8, label='loop0'), - sc.Domain(1, False, 5, 10, label='dom1'), + sc.Loopout(8, label='loop0', color=sc.Color(10, 10, 10)), + sc.Domain(1, False, 5, 10, label='dom1', color=sc.Color(20, 20, 20)), sc.Domain(2, True, 5, 10, label='dom2'), sc.Loopout(12, label='loop1'), sc.Domain(3, False, 5, 10), - ]) + ], color=sc.Color(30, 30, 30)) self.assertEqual(1, len(design_from_json.strands)) self.assertEqual(expected_strand, design_from_json.strands[0]) + self.assertEqual(expected_strand.color, sc.Color(30, 30, 30)) def test_strand__3p_extension(self) -> None: design = self.design_6helix @@ -100,10 +104,11 @@ def test_strand__3p_extension(self) -> None: sb.to(10) sb.extension_3p(5) + sb.with_domain_color(sc.Color(10, 10, 10)) expected_strand: sc.Strand = sc.Strand([ sc.Domain(0, True, 0, 10), - sc.Extension(num_bases=5), + sc.Extension(num_bases=5, color=sc.Color(10, 10, 10)), ]) self.assertEqual(1, len(design.strands)) self.assertEqual(expected_strand, design.strands[0]) From fea3b32c1f3e365e41f368957d06c70df8b3a89a Mon Sep 17 00:00:00 2001 From: David Doty Date: Sun, 4 Dec 2022 08:40:11 -0800 Subject: [PATCH 14/14] added example with domain colors --- examples/domain_colors.py | 32 ++++++++++++++++++++++++ examples/output_designs/domain_colors.sc | 24 ++++++++++++++++++ 2 files changed, 56 insertions(+) create mode 100644 examples/domain_colors.py create mode 100644 examples/output_designs/domain_colors.sc diff --git a/examples/domain_colors.py b/examples/domain_colors.py new file mode 100644 index 00000000..1e8ed5e6 --- /dev/null +++ b/examples/domain_colors.py @@ -0,0 +1,32 @@ +import scadnano as sc + + +def create_design() -> sc.Design: + helices = [sc.Helix(max_offset=100) for _ in range(4)] + design = sc.Design(helices=helices, grid=sc.square) + + red = sc.Color(255, 0, 0) + dark_red = sc.Color(150, 0, 0) + green = sc.Color(0, 255, 0) + dark_green = sc.Color(0, 150, 0) + blue = sc.Color(0, 0, 255) + dark_blue = sc.Color(0, 0, 150) + black = sc.Color(0, 0, 0) + + design.draw_strand(0, 0) \ + .extension_5p(num_bases=5).with_domain_color(red) \ + .move(8).with_domain_color(green) \ + .loopout(1, 5).with_domain_color(dark_blue) \ + .move(-8).with_domain_color(dark_red) \ + .cross(2) \ + .move(8).with_domain_color(dark_green) \ + .cross(3) \ + .move(-8) \ + .extension_3p(num_bases=5).with_domain_color(black) \ + .with_color(blue) + + return design + +if __name__ == '__main__': + d = create_design() + d.write_scadnano_file(directory='output_designs') diff --git a/examples/output_designs/domain_colors.sc b/examples/output_designs/domain_colors.sc new file mode 100644 index 00000000..bb4cff75 --- /dev/null +++ b/examples/output_designs/domain_colors.sc @@ -0,0 +1,24 @@ +{ + "version": "0.17.7", + "grid": "square", + "helices": [ + {"max_offset": 100, "grid_position": [0, 0]}, + {"max_offset": 100, "grid_position": [0, 1]}, + {"max_offset": 100, "grid_position": [0, 2]}, + {"max_offset": 100, "grid_position": [0, 3]} + ], + "strands": [ + { + "color": "#0000ff", + "domains": [ + {"extension_num_bases": 5, "color": "#ff0000"}, + {"helix": 0, "forward": true, "start": 0, "end": 8, "color": "#00ff00"}, + {"loopout": 5, "color": "#000096"}, + {"helix": 1, "forward": false, "start": 0, "end": 8, "color": "#960000"}, + {"helix": 2, "forward": true, "start": 0, "end": 8, "color": "#009600"}, + {"helix": 3, "forward": false, "start": 0, "end": 8}, + {"extension_num_bases": 5, "color": "#000000"} + ] + } + ] +} \ No newline at end of file