A Method for Constructing Large DNA Codesets

A Method for Constructing Large DNA CodesetsVinhthuy PhanAbstractThe word design problem is an important problem in DNA computing.The goal is to design a set of single-stranded DNA molecules that are structure-free and mutually non-crosshybridizing.We present an innovative methodof constructing such large sets,so called codesets,whose degree of structure-freeness and non-crosshybridization is controlled by a given ing awell-known,simplified model of hybridization affinity,known as the h-distancemodel,our method produces sets larger than those produced under similarmodels,for example,a similar(slightly less realistic)model can construct1088-mers,while our construction produces2568-mers with the same constraints interms of the strands’mutual hybridization affinity and GC-content.We furtherprovide a justification for the claimed large sizes of produced sets by estimatingtheir relative sizes to those of theoretically maximal sets.The existence of suchlarge sets of structure-free,non-crosshybridizing strands is necessary for large-scale DNA computations and approaches.Current methods that employ morerealistic estimates of DNA hybridization affinity have been unable to producesets large enough to encode more than several thousands or even hundred pa-rameters.It seems feasible to extend our methodology to include more realisticassumptions,such as making use of known thermodynamic parameters.Keywords:Data representation,codeword design,Gibbs energy,h-distance,insertion-deletion block-isomorphic distance.1IntroductionSince the discovery of the DNA double helix by Watson and Crick[46]and the sub-sequent discovery of the base-pairing rule that provides a mechanism for the cell to replicate and correct mistakes made during various biological processes,DNA mole-cules have been known as a very powerful storage medium of genetic information. Recently,this genetic material and its base-pairing characteristic have served as a basis for computation and the building of nano-structures in test tubes.Adelman [1]designed an experiment using DNA molecules to encode and solve a computa-tionally hard problem,known to be NP-complete.This seminal work effectively proved that DNA materials can be used to compute,and in fact solve very hard computational problems.Since then,researchers have proposed solutions that use DNA materials to solve and model various types of problems,including models of DNA memories[4,10,24,25,35],the design and construction of nano-structures1[15,18,26,28,31,39,44,49],the construction of sophisticated protein arrays [27,30],and other models and approaches in DNA computing[6,32,41,42,45,48]. The usage of DNA materials in bio-computing is attractive mainly because the inher-ent massively parallelism of chemical reactions that take place once single-stranded DNAs are placed in a test tube.In fact,this was the reason Adelman argued that an NP-complete problem can be solved tractably in a biocomputing framework.At the heart of such an experiment is the process that transforms a computational problem into a biological experiment.The obvious approach to encode information in DNA is to encode symbolic strings in DNA strands.Direct encoding is,how-ever,not very efficient for storage or processing of massive amounts of abiotic data because of the enormous implicit cost of DNA synthesis to produce the encoding se-quences.Indirect and more efficient methods have been proposed[19,23],assuming the existence of a large set of structure-free and non-cross-hybridizing DNA mole-cules.Cross-hybridization happens when two different strands interact and form a secondary-structure duplex,or when even a strand interacts with itself,forming a secondary structure.For the encoding to work,there must be a way to rep-resent parameters of the computational problem as single-stranded DNAs so that no cross-hybridization(hybridization among unintended molecules)is allowed.In other words,a strand should hybridize only to its designated strands.When many strands do not hybridize to designated strands,the experiment will not work as intended.Hence,the modeling of the original problem using DNA strands must be very carefully done,to discourage as few cross-hybridizations as possible.The codeword design problem[7,12,16,18,21,22]is to design a set of single-stranded DNAs in such a way that cross-hybridizations and freeness of structures can be controlled effectively.The goal is to produce sets of strands that are likely to hybridize to targeted strands while minimizing the possibility of unintended hy-bridizations that may induce erroneous outcomes.The ability to construct large sets of strands with such properties gives hope to the modeling of large problem inputs using DNA molecules without worrying much about these side effects as a result of cross-hybridization.Without such large sets,computing in a test tube is restricted to only solving small problems.Stability of DNA and RNA structures has been described using simple measures of such as the number of matching pairs,the GC-content[34]as well as more real-istic estimates of free-energy using thermodynamic models[2,8].Current methods based on latest thermodynamic parameters[37]have done an admirable job in pre-dicting the stability of DNA/RNA single strands and ing these methods to create large sets of non-crosshybridized,structure-free single-stranded DNAs is, however,a computationally expensive and often infeasible task.The complexity of the estimation of free energy used by these model has made it difficult to design a simple,systematic,well-parametrized procedure to create such large sets,other than searching through the space of single-stranded DNAs or performing a limited heuristic.Although various algorithms have been proposed for testing the quality ofcodeword sets in terms of being free of secondary structure[7,14],very few methods have been proposed to systematically produce large sets of high enough quality to guarantee good performance in test tube protocols[19].Similarly,although bounds on the cardinality of maximal sets of given parameters have been established using very simple models[29,33],no such estimates on the cardinality of maximal sets are known under more realistic measures of strand stability.Attempts to construct non-crosshybridizing sets in silico using various methods such as templates[3,16] or heuristic search[40,43]have produced relatively small sets.Other attempts to create non-crosshybridizing sets[9]in vitro are attractive in being able to construct a physical set but suffer from the challenge of knowing precisely its size and more seriously its content,which is necessary to encode a problem.In this work,we present a methodology to construct very large sets of non-crosshybridizing singled-stranded DNAs.We employ a measure of duplex stability known as the h-distance[21,22,35]that is has been shown in vitro to be a good approximation of free energy[5,19,21,22,24].We use the notion of a codeset,which is a set of single-stranded DNAs satisfying a hybridization constraint,to define the distance between two single-stranded DNAs,which relates negatively to the number of complementary matches in any rigid alignment between the two strands themselves and between their complements.This ensures that the strands in the set and their Watson-Crick complements do not cross-hybridize.This definition will be articulated more precisely in the following sections.Our methodology shares similar characteristics to such methods as the template methods[3,16],where a defined product between a base set with another set is taken to produce the desired codeset. One difference is that we use a base set whose strands are short and thus can be constructed in anyway with stringent thermodynamic properties(computationally feasible due to the presumed short strand length).From these base set,we then construct much larger codesets.As such,our methodology produces very large codesets in comparison to other methods using similar models of strand stability. For instance,a model that is less realistic than ours[16]produces a codeset of108 8-mers with Hamming distance of4and GC-content of4;our methodology can produce a codeset of2568-mers with h-distance(more stringent than Hamming distance)of4and GC-content of4.Our methodology can produce easily very large codesets with a given x%of h-distance and x%GC-content.Table1shows the sizes of codesets with various constraints in the h-distanceτand GC-content.2Hybridization of Single-stranded DNAsBio-computation relies on the base-pairing mechanism of DNA molecules.A single stranded DNA–and we will simply call a strand–is a directed sequence of Adenine (A),Thymine(T),Cytosine(C)and Guanine(G),conventional read from the5’end to3’end.Further,when two nucleotides(a.k.a bases)of the same strand or of different strands come into contact,the base-pairing rule dictates that A pairs withLengthτ=Lengthτ=Length/2GC-content=Length/2GC-content=Length/28-mer1625612-mer644,09616-mer25665,53620-mer10241,048,57624-mer409616,777,216Table1:Sizes of produced codesets with stringency,τat100%and50%of lengths, respectively and GC-content at50%of lengths.The larger the stringencyτ,the more non-crosshybridizing the strands in a codeset are.τ=100%of length is the most stringent condition.T and C pairs with G by forming hydrogen bonds.Whether two bases A and T or C and G pair together also depend on the environment’s temperature and–according to models such as that of[2,8,36,37,38]–also what the neighboring bases are. As a consequence,when two strands come to close contact,if the number of hydro-gen bonds between them formed via A-T,C-G pairs is strong enough at the given temperature,they will hybridize and form a secondary-structure duplex.A strand can even hybridize to itself to fold into a secondary structure,if it has enough com-plementary bases arranged in an appropriate manner.When two strands hybridize with sufficiently many hydrogen bonds,small variations in strand composition will not cause major changes in hybridization events.In a bio-computing experiment,the experiment designer typically aims to have a DNA strand(representing part of the original problem)hybridize only to a small set of designated strands,most often only its perfect complement.When these strands float freely in the test tube,they however can come into contact with any other strand,and possibly hybridize to any of them if the two form enough hydrogen bonds.When such a cross-hybridization occurs,the designated strand will not be able to hybridize to its target and the experiment does not occur as the experi-ment designer desires.To avoid these undesirable effects,input strands should be by design”far apart”from one another in hybridization affinity.Unfortunately,the hybridization affinity between DNA strand sets is difficult to quantify and the most precise models of approximating the Gibbs free-energy of a DNA duplex is compu-tationally expensive[11,37,38].Further,the problem offinding a maximal set of strands that are mutually far apart in hybridization affinity can be abstracted as the maximum independent set problem,which is computationally infeasible[17]when the problem size gets large.In fact,evenfinding a very large set is computationally difficult[11],even for the small oligo-nucleotide sizes useful in DNA computing.Re-cent efforts based on the thermodynamic estimate of free energy has only been able to produce sets of no larger than400strands of12-to20-mers[40,43].This greatlyrestrict the ability to model and solve large problem instances in the test tube.As a result,researchers have either used most accurate models to design very small sets of non-crosshybridizing strands or relied on simpler models to design larger sets. Our work presented here uses a good approximation to the most accurate model of free energy to produce very large sets of strands that are mutually far apart in hybridization affinity.2.1Estimating Hybridization AffinityThe Gibbs free energy between two strands is the most appropriate criterion of quality for sets of non-crosshybridizing strands or experiments in-vitro[47].And al-though hybridization reactions in vitro are governed by well established rules of local interaction between base pairings in the nearest-neighbor model and the staggered-zipper mode[36,37,38],an exhaustive search of strand large sets of strands that are non-crosshybridizing and structure-free is infeasible,even for the small size of oligo-nucleotides useful in DNA computing.Even computing large sets of non-crosshybridizing strands is too computationally expensive.Thus,we need to address this serious problem,if we wish to carry out large-scale bio-computing tasks.Several approaches have been proposed to estimate hybridization affinity by making assumptions about the type of interactions between neighboring bonds,for the range of oligo-nucleotides used in DNA-based computing.These approaches include simple Hamming distance[29,33],Hamming distance with limited overlap regions[3],h-distance[21,22,35],insertion-deletion block-isomorphic[13].An ex-tension of the nearest neighbor model proposed by[11]computes optimal alignments between DNA oligo-nucleotides using a dynamic programming algorithm.Of these, models based strictly on the Hamming distance is most efficient as the distance between any two strands can be computed in O(n)time(n being the length of the two strands);the h-distance model is less efficient(O(n2)time);and the insertion-deletion block-isomorphic model and nearest neighbor models are least efficient (O(n3)time).The simple Hamming distance model is too crude as an approximation about hybridization affinity as it exclude the possibility of two strands hybridizing in shifted alignments.The h-distance model appears to be a good enough estimation to hybridization affinity and in vitro protocols(such as PCR selection)based on it has appeared to produce decent-quality sets of mutually non-crohybridizing strands [5,19,21,22,24].At the same time,the h-distance model isflexible enough to allow construction of very large sets of mutually non-crohybridizing strands,as shown in this work.In the model using the h-distance,the Gibbs energy is approximated by counting the maximal number of basepair matches in all possible alignments of two strands. More precisely,Let x,y be two strands(written from the5 -to the3 -end.)Defineh(x,y):=min−n<k<n|k|+H(x,σk(y ))whereσk is the shift by k positions(k>0means right shift;k<0means left shift); y is the Watson-Crick complement of y obtained by reversing y and exchanging A s for T s and vice versa,and C s for G s and vice versa;and H(·,·)is the ordinary Hamming distance measuring the number of different bases in the specified shift. The h-distance considers hybridization in all possible frame-shifts and is therefore more realistically restrictive than simpler models[3,16,29,33].A distance of0 indicates perfect complementarity.A large distance indicates that even when x finds itself in the proximity of y,they contain few complementary basepairs,and are less likely to hybridize.This distance measure h can be precisely related to the maximum number of complementary base pairs in all frame shifts as follows: Lemma2.1.h(x,y)=n−max−n≤i≤n{m(x,σi(y ))}where m(x,σs(y ))is the number of agreed bases of x and y when they are aligned in frame shift i.Proof.Let s be the shift in which h(x,y)=τ=s+H(u,σs(v ))is minimum.Since s+H(u,σs(v ))+m(u,σs(v ))=n,we get h(x,y)=n−m(u,σs(v )).Further,since s+H(u,σs(v ))is minimum over all shifts s,m(u,σs(v ))is maximum over all shifts s.Example:Suppose we have x=AGC,y=T GG(and y =CCA).In the computation of h(x,y),we look at the number of differing characters in all shift of x and y .At shift-2,the distance is2+H(A,A)=2.At shift-1,the distance is 1+H(AG,CA)=3.At shift0,the distance is H(AGC,CCA)=3.At shift1,the distance is1+H(GC,CC)=2.At shift2,the distance is2+H(C,C)=2.So, h(ACG,CGG)=2.Lemma2.1shows that the h-distance is a special case of the insertion-deletion like block-isomorphic metric introduced by[13]to estimate the free energy between two strands.The notion of block-isomorphic subsequences of two strands allows alignment of the two strands with insertions and deletions of nucleotides.This model is more general than the h-distance model because it allows two DNA strands hybridizing in blocks with non-equal gaps in between.For short strands(12-to20-mers),the difference does not seem to be substantial.In the insertion-deletion block-isormorphic model,the distance between two strands x and y is defined as:ψ(x,y)=min{ x , y }−max{b(x,y )}where x =x i is the weight of the strand x,which is a sum of all weightsof its bases.This means each base can have a different weight,for instance indi-cating the number of hydrogen bonds it could form with its complementary base. b(x,y )=max{ σ(x) :σ(x)=σ(y)}is the weight of the maximum isomorphic block that are identical in both x and y .As an example,one of such blocks(in-dicated by capitalized letters)in x=aCGT ggAA and y =aatCGT aaaggAAccis CGT-AA.Each corresponding pair of subblocks in the block must have equal length;for example,the blocks(indicated by capitalized letters)in aCaaGT ggAA and aatCGT aaaggAAcc are not isomorphic.3Construction of Large DNA CodesetsThe notion of a codeset is borrowed from information theory,where codewords of a codeset have the property of being mutually far from one another so that errors that occur during the transmission of codewords through a noisy channel can be best detected.In the DNA word design problem,a codeset[24,35]is a set of single-stranded DNAs that are mutually far from one another,with respect to a distance measure that is an estimate of free energy,which in our case is the h-distance.The larger their mutual distances are,the less likely they cross-hybridize. Stronger conditions for a codeset require that not only these strands but also their complements are far from one another so that they and their complements are less likely to cross-hybridize[40,43].In this work,we define a(m,τ)-codeset using these requirements,which are stronger than other simpler models[3,16,29,33]. Our codeset with parameterτis very similar to the set of strands with parameter δof another model that employs a nearest-neighbor thermodynamic model[40] to construct(much smaller)codesets.Specifically,given a distance measure d,an (m,τ)-codeset S is a set of strands with length m,satisfying the following properties:1.d(u,v)≥τ,for all u,v∈S possibly u=v.2.d(u,v )≥τ,for all u,v∈S,u=v.Note that d(u,u )=0for all u.3.d(u ,v )≥τ,for all u,v∈S,u=vThefirst property ensures that cross-hybridization between strands within the code-set is unlikely.The second property ensures that cross-hybridization between strands within the codeset and complements of strands within the codeset(other than their own complements)is unlikely.This property maximizes the chance that the com-plement of a strand within the codeset will hybridize to it,not other strands within the codeset.The third property ensures that the complements of strands within the codeset will not cross-hybridize,again to maximize the chance that complements of strands within the codeset will hybridize to their targeted strands(in the set)and not among themselves.It is desirable to require a codeset to have an additional property of having a similar GC-content[3,29].If each strand has the same number of G’s and C’s, they are expected to have a similar melting temperature,due to the fact a GC bond with three hydrogen bonds is stronger than an AT bond with only two hydrogen bonds.Similar melting characteristics are particular important in bio-molecular experiments that rely heavily on repeated melting and annealing DNA strands,such as PCR selection.To incorporate this requirement,we define an(m,τ,w)-codesetto be the same as an(m,τ)-codeset with the additional property that each strand in the set has exactly w G’s or C’s in totality.Our methodology of constructing large codesets is a two-step process and can be summarized as follows:•First,a base(m,k)-codeset(or(m,k,w)-codeset)S is constructed.•Second,the base codeset is shuffled l times to create an(ml,kl)-codeset(or (ml,kl,wl)-codeset)of cardinality|S|l.This methodology assumes the existence of a base codeset,which can be con-structed by any method.For example,it is easy to construct an(m,1)-codeset (or(m,1,w)-codeset)under the model of the h-distance,or the insertion-deletion block-isomorphic distance for that matter.The attractiveness of this construction is that it scales linearly the stringency l and scales multiply the cardinality of the base codeset,consequently producing very large codesets.To summarize our specific results,we can construct•An(ml,l)-code of size4m−|P m|2l,where|P m|=|{x:h(x,x)=0,|x|=m}=4m2if m is even,and0if m is odd.•An(ml,l,wl)-code(each strand having GC-content wl)of size(mw)2m−|Q m|2l,where|Q m|= m/2w/22m/2if m and w are even and|Q m|=0if m or w is odd.•In general,given a(m,k)-codeset(or(m,k,w)-codeset with GC-content of w) S that is amplifiable,we can construct an(ml,kl)-code(or(ml,kl,wl)-codeset respectively)of size|S|l.Table1shows the sizes of a few sample codesets that we can construct.Its sizes compare favourably to others such as[16].3.1The Base CodesetThefirst step of constructing a codeset using our methodology is to create a base codeset.A base codeset is a codeset that already exists or can be easily constructed. This means that we can make use of existing codesets to construct larger codesets with longer strands or if we choose to construct it expensively,it should have short strands and/or a small parameterτ.Here,we describe a simple construction of an (m,1)-codeset and respectively an(m,1,w)-codeset(with GC-content w).Lemma3.1.Algorithm1constructs a set S m,which is an(m,1)-codesetProof.By construction it consists of no pair u,v such that d(u,v)=0.The distance d can be either the h-distance or the insertion-deletion distance or any distance function with which d(u,v)=0if and only if u=v .More specifically,1.d(u,v)≥1for all u,v∈S m(by construction).2.d(u,v )≥1for all u,v∈S m,u=v.This is because d(u,v )=0if and only ifv =u or v=u.Since the algorithm removes all strand u from S m for which d(u,u)=0,it is impossible to have d(u,v )=0.3.d(u ,v )≥1for all u,v∈S m,u=v .This is because d(u ,v )=0if and onlyif v =(u ) =u,which is again impossible for the same reason as above.We can further determine the cardinality of S m,|S m|,to be equal to4m−|P m|2,where P m is the set of strands that are their own Watson-Crick reverse complements, i.e.strands u’s such that d(u,v)=0;The cardinality of this set is|P m|=4m2if m is even,and0if m is odd.After removing these strands,there is exactly one Watson-Crick reverse complement v for every strand v left,such that v =v.Therefore, |S m|=4m−|P m|2.Algorithm1Input:n=m.Output:(m,1)-codesetConstruct S m,the set of all possible DNA strands of length m.for all strand v∈S m doLet v be the unique complement of v.S m=S m−{v }end forAlgorithm2Input:m,w.Output:an(m,1,w)-code(w is the GC-content)S m=Ø.for all configuration f of all nwpositions containing G or C dofor all strand v in configuration f doS m,w=S m,w∪vend forend forfor all strand v∈S m,w doConstruct v ,the unique complement of v.S m,w=S m,w−vend forConstructing an(m,1,w)-codeset follows the same ideas as in Algorithm1.Lemma3.2.Algorithm2constructs an(m,1,w)-codeset S m,w,which is an(m,1)-code with the additional property that every strands of S m,w has exactly w G’s and C’s in totality.Proof.Similar to the proof above.Using the same argument,we can show that the set S m,w constructed is an(m,1,w )-code.Further,the cardinality of S m,w is |S m |=(n w )2n −|Q m |2,where Q m ={x :d (x,x )=0,x has exactly w G’s or C’s in total.},and |Q m |= m/2w/2 2m/2if mand w are even,and 0if either m or w is odd.More specifically,the set S m,w starts out with n w 2n strands,because there are n w configurations for G’s and C’s,in each of which there are exactly two choices for each position.Then,complements of strands in S m,w are removed from the set.These strands consists of the set Q m ,all strands v ’s that are its own reverse complements.|Q m |= m/2w/22m/2because the i th and (m −i +1)th characters of these strands must be complementary,which implies m and w must be both even and further the GC-content of w implies that each position has exactly 2choices.After strands in Q m are removed from S m,w ,as before thereis exactly one complement for each strand left.Hence,|S m |=(n w )2n −|Q m |2.3.2Amplifiable Functions and the Shuffle FunctionAfter constructing a base (m,1)-codeset or using a good,existing codeset,the second step of our construction is apply the shuffle function to this base codeset l times to create a larger (ml,l )-codeset.Similarly,we can shuffle an (m,1,w )-codeset (with GC-content w )to construct an (ml,l,wl )-codeset (with GC-content wl ).The reason this construction works as claimed is based on the property of the shuffle function being amplifiable .A function f :D m ×D m ×···D m (l times)→D ml (D m is a set of strands of length m )is amplifiable with respect to the distance function d ifd (f (u 1,···,u l ),f (v 1,···,v l ))≥min 1≤i,j ≤l{l ·d (u i ,v j )}Theorem 3.3.Let S m be (m,k )-codeset (or (m,k,w )-codeset),f be an amplifiable function,andS ml ={f (u 1,···,u l ):u i ∈S m ,1≤i ≤l }Then,S ml is an (ml,kl )-codeset (or (ml,kl,wl )-codeset,respectively).Proof.for all u,v ∈S ml ,where u =f (u 1,···,u l ),v =f (v 1,···,v l ),u i ,v j ∈S m ,1≤i,j ≤l :d (u,v )≥min 1≤i,j ≤l{l ·d (u i ,v j )}≥lkOur construction is based on the theorem above and an amplifiable function called Shuffle .Let u 1,u 2,···,u l be l strands,each of length m .Specifically,let u i =u i 1u i 2···u im .Define:Sh (u 1,···,u l )=u 11u 21···u l 1u 12u 22···x l 2···u 1m u 2m ···u lmIn other words,Sh (u 1,···,u l )is a string of length lm ,whose first l characters are the first characters of u 1,···,u l concatenated together,and whose next l characters are the second characters of u 1,···,u l concatenated together,and so on.This function is analogous to shuffling l decks of cards (each having m cards)together in such a way that the cards are perfectly interleaved.For example,suppose u 1=ABC,u 2=LMN,u 3=XY Z ,thenSh (u 1,u 2,u 3)=ALXBMY CNZTheorem 3.4.Sh is amplifiable with respect to the h -distance function.Proof.According to Lemma 2.1,the h -distance function defines the distance be-tween two strands u ,v as:h (u,v )=n −max −n ≤i ≤n{m (u,σi (v ))}where m (u,σi (v ))is the number of agreed bases of u and v when they are aligned in frame shift i .Suppose we construct a set S by applying Sh to S m ,an (m,k )-code l times.Let u =Sh (u 0,···,u l −1),v =Sh (v 0,···,v l −1)∈S and u i ,v i ∈S m ,0≤i ≤l −1.Note that the Sh function is designed in such a way that when u and v are aligned in any frame shift,each u j is aligned perfectly to some v f (f is uniquely determined by the frame shift).Let i be the frame shift where h (u,v )attains its value,i.e.h (u,v )=lm −m (u,σi (v )).In this case,each u j is aligned perfectly with v j +i mod l .And we have,m (u,σi (v ))≤ l j =1max −m ≤k ≤m {m (u j ,σk (v j +i mod l ))}.Hence,h (u,v )=lm −m (u,σi (v ))≥lj =1m −max −m ≤k ≤m {m (u j ,σk (v j +i mod l ))}=lj =1h (u j ,v j +i mod l )≥l ·min 1≤i,j ≤l{h (u i ,u j )}Therefore,Sh is amplifiable. Corollary 3.5.Under the h -distance model,our construction can construct•An (ml,l )-code of size 4m −|P m |2l ,where |P m |=|{x :h (x,x )=0,|x |=m }=4m 2if m is even,and 0if m is odd.•An (ml,l,wl )-code (each strand having GC-content wl )of size(m w )2m −|Q m |2 l,where |Q m |= m/2w/2 2m/2if m and w are even and |Q m |=0if m or w is odd.。