# Controlling the temperature sensitivity of DNA-mediated colloidal interactions through competing linkages

## Abstract

We propose a new strategy to improve the self-assembly properties of DNA-functionalised colloids. The problem that we address is that DNA-functionalised colloids typically crystallize in a narrow temperature window, if at all. The underlying reason is the extreme sensitivity of DNA-mediated interactions to temperature or other physical control parameters. We propose to widen the window for colloidal crystallization by exploiting the competition between DNA linkages with different nucleotide sequences, which results in a temperature-dependent switching of the dominant bond type. Following such a strategy, we can decrease the temperature dependence of DNA-mediated self assembly to make systems that can crystallize in a wider temperature window than is possible with existing systems of DNA functionalised colloids. We report Monte Carlo simulations that show that the proposed strategy can indeed work in practice for real systems and specific, designable DNA sequences. Depending on the length ratio of the different DNA constructs, we find that the bond switching is either energetically driven (equal length or ‘symmetric’ DNA) or controlled by a combinatorial entropy gain (‘asymmetric’ DNA), which results from the large number of possible binding partners for each DNA strand. We provide specific suggestions for the DNA sequences with which these effects can be achieved experimentally.

(a) Electronic mail: bm411@cam.ac.uk, m.e.leunissen@amolf.nl, df246@cam.ac.uk

## I Introduction

The high selectivity of hybridisation of DNA makes it an interesting molecule to be used as “smart glue” in the self-assembly of complex, nano-structured materials. Some of the advantages of DNA as a selective linker are that it can code for a large variety of specific linkages, it is readily available, it can be used under near ambient conditions and the binding process is reversible. Consequently, applications of DNA-mediated self-assembly range from computational biology (1) and the assembly of scaffolded “DNA-origami” structures (2) to the development of targeting strategies (see e.g. (3); (4); (5)). Much experimental work has focused on the application of DNA as a selective linker that enables the self-assembly of complex colloidal structures (6); (7); (8); (9); (10); (11); (12); (13); (14); (15); (16); (17); (18); (19); (20); (21); (22); (23); (24); (25). Several examples have been reported in the literature of (relatively simple) DNA-linked colloidal crystals, either consisting of nano colloids (see e.g. (12); (13); (14); (15); (16)) or of micron–size colloids with complementary DNA coating (see e.g. (11)). In parallel, theoretical investigations have provided insight into the factors that influence DNA-mediated self-assembly (26); (27); (28); (29); (30).

There are, however, some disadvantages associated with the use of DNA as a tool to link (nano)colloids. In particular, the very factors that lead to the specificity and reversibility of DNA linkages also cause the strength of DNA mediated interactions to depend strongly on the external conditions, such as the temperature or ionic strength (32); (8); (31); (10); (33). This sensitivity results in an abrupt onset of aggregation as the temperature is lowered - a phenomenon that can be problematic for self-assembly because it narrows the “window” of conditions within which reversible self-assembly of ordered structures is possible. If this window is missed, self-assembly either does not take place at all or results only in disordered aggregates. It is therefore important to explore possible approaches to combine the high selectivity of DNA-mediated interactions with a more gradual response to external conditions.

In the present paper, we use numerical simulations to demonstrate that ‘competing’ DNA interactions can be used to create colloidal systems with a more gradual temperature response. Specifically, we consider a binary mixture of colloids ( and ) that, unlike most systems studied so far, are functionalised by a mixture of different DNA strands ( and , respectively). These DNA sequences are chosen such that can bind to both and (and similarly to and ), but cannot bind to (we will show that such DNA sequences can be readily designed). Importantly, the – and – linkages are weaker than those between and . As a result, – can form at a higher temperature than the linkages involving or . However, whereas a given , pair can form only one – linkage, it can form two weaker linkages (– and –). This leads to a competition between the different types of linkages in which both their free energy of formation and combinatorial entropy effects,which depend on the number of different ways in which the linkages can form, play a role - as is characteristic for systems with multivalent binding (34); (35); (36); (37). Here, we investigate how the majority of the linkages can switch from one type to the other, not only as a function of the temperature and the difference in binding strength of the two linkage types, but also as a function of the surface coverage and length ratio of the DNA constructs.

The remainder of this paper is organised as follows: in Sec. II we introduce the model and the Monte Carlo algorithm. Sec. III reports our results. In Sec. III.1 we consider the symmetric model, in which the length of the competing DNA constructs is the same (). We show how the effective pair interaction depends on the hybridisation free energy of the DNA sticky ends and how the energetically driven switching from one linkage type to another broadens the association-dissociation transition of the particles; we rationalise these findings in the context of a Mean-Field model. In Sec. III.2 we consider the asymmetric model (). We show how the length ratio of the different DNA constructs can be used to enhance the bond switching through combinatorial entropy effects. Finally, in Sec. IV we indicate how our approach could be implemented experimentally.

## Ii The Model

For the sake of simplicity, we use the relatively simple but well-tested model of ref. (31) to describe the interaction between DNA-functionalised colloids.

Following ref. (31), colloids are coated with double-stranded (ds) DNA that is shorter than its persistence length of (38) (corresponding to 150 nucleotides (39)). The dsDNA is terminated by a short single-stranded (ss) DNA sequence (’sticky end’) that can hybridise to complementary ssDNA. Binding is only allowed between colloids of type and . DNA is modelled as thin rods (see Fig. 1) randomly tethered to the colloidal surface. For typical surface coverages, the steric repulsion due to self-avoidance of the DNA strands ( diameter) is negligible compared to the entropic effects involved in the binding of tethered strands (31) and we therefore treat the DNA as non-self-avoiding. Except for the DNA funtionalization, the colloids are assumed to be smooth and hard. We consider two types of DNA (, on and , on , Fig. 1) defined by different ssDNA end sequences, and (later on) different lengths , . The characteristic distance between DNA strands on a given colloid is , where is the total surface area of that colloid (say type ) and and are the total number of and strands on (similar expressions apply to ).

Here we focus on the interaction between two parallel, planar surfaces. The pair potential between spherical colloids (of radius ) can be computed from these planar surface interactions via the Derjaguin approximation (40). In Ref. (31) we have shown that this is a reliable approximation when . From geometry, it further follows that the interactions are strictly pairwise additive if , assuming that the hard cores of the particles can come into contact. For smaller particles, curvature and so-called three-body effects may cause deviations of the exact interactions from what is predicted here. Nevertheless, we expect that bond switching can also occur in these systems, albeit at different values for the relevant parameters.

In the simulations, we consider two square planes with side and periodic boundary conditions in the directions parallel to the planes. We fix (31). In real units, the width of the simulation box is . Varying and , we obtain characteristic inter-chain separations . These values span typical experimental values, e.g. (10); (11). We verified that averaging over different random realisations of tethering points does not alter our results within a 2% of tolerance.

The hybridisation of pairs of ssDNA depends on the specific nucleotide sequences of the individual strands and on the solvent properties (temperature and salt concentration). For the present paper it is important that ssDNA sequences can be designed in such a way that only the following linkages are possible: –, but also – and – (see Fig. 1a and b). In Sec. IV we show that it is possible to design nucleotide sequences that will yield this behaviour.

We denote the hybridisation free energy of two free ssDNA by , where for an – linkage and or for an – or – linkage. We assume that . The binding free energy between two surface-bound DNA strands must also include a configurational entropy term (31); (37), which accounts for the reduced freedom of motion of the strands upon binding (Fig. 1). If the dsDNA can swivel freely on the surface and the ssDNA is connected flexibly to the dsDNA (13), then the configurational space of the DNA rods is a circle ( and in Fig. 1), that is bounded by the two surfaces. In App. A we give an explicit expression for as a function of and . As can be seen from Figs. 1c and d, and are equal to if the distance between the planes is larger than and equal to if . As a result (see (31) and App. B) the hybridisation free energy of two tethered strands ( or ) is given by

(1) |

where is the number density corresponding to one molar . is a non–universal factor that depends on the details of the linkage formed between the two dsDNA strands. As explained in App. B, slightly different expressions for result depending on how the sticky ends are coarse grained. However, the resulting differences in the predicted binding strengths are tiny and therefore irrelevant for the present discussion.

In each MC move (31), we randomly select one of the (or ) functional arms and attempt to make, break or switch a linkage. First we consider the list of all the free () that can bind to , possibly including the partner to which is already bound. A linkage – is created with probability . The probability that no linkage is formed is given by , where . is computed using Eq. (1). In order to enhance the sampling of the model we also implemented a ‘linkage-swapping’ MC move that attempts to switch between a single – linkage and two weaker –+– linkages. Details are given in App. C.

## Iii Results

### iii.1 Symmetric DNA model ()

We first consider systems in which the competing DNA constructs have equal length . Fig. 2 shows , the average number of – linkages, and , which denotes the number of –+– linkages. and are plotted as a function of . There is a linear relation between the hybridisation free energy and the temperature , where the hybridisation enthalpy/entropy () are negative constants. Hence decreasing corresponds to a decrease in . For real DNA sequences, also depends on (see Sec. IV). However for the sake of simplicity in Fig. 2 we sketch our results at constant .

Upon decreasing from an initial situation where all DNAs are unbound, the stronger – linkages form first. For the specific choice of simulation parameters listed in Fig. 2, this happens when . As we lower even more, the linkages disappear in favour of the weaker bonds between – and –. This may seem counter-intuitive at first, but the reason is simple: replacing a single – linkage by two weaker ones (–+–) is favorable as long as the total gain in binding free energy upon forming these two linkages outweighs the loss in binding free energy when breaking an – linkage. The switching between the different bond types in Fig. 2 is driven by the difference in the hybridization free energies of the individual linkages and the main effect of increasing is to decrease the value of the crossing point where . As will be shown below, the combinatorial entropy of the overall system (37); (31) plays a less important role in this particular case.

The behaviour shown in Fig. 2 can be understood in a more quantitative way using a mean-field theory. We define and as the average configurational entropy cost of an – and – linkage

where and are the circles on enclosing all the possible and (tethered at position ) which could be hybridised by on –DNA coated colloids represent a multivalent binding system in which each DNA strand can form several different linkages. In the low binding regime () this combinatorial entropy contribution is accounted for by noticing that the number of possible – and – linkages are and , where () is the mean distance between () strands. It follows that the partition function and the average fraction of the two linkages ( ) is:

(3) | |||||

(4) |

Fig. 3 compares the mean-field predictions based on Eq. (4) for the single linkage model (31) (part a) and for the competing linkages model of Sec. II (part b) with the results of the MC simulations. There is qualitative, though not quantitative agreement between theory and simulations.

A better agreement is obtained when we improve our estimate of the combinatorial entropy terms. If linkages are already present, the number of ways an extra linkage can be added is smaller than if . This is, because the mean distance between un–hybridised strands increases like (and similarly for ). Instead of exactly dealing with this correction we used Eq. (4), while correcting and by the average number of linkages ( and ) which are then computed self–consistently in the following way

Dotted curves of Fig. 3 show the predictions based on the solution of Eq. LABEL:mean2. As can be seen from this figure, the agreement between theory and simulation is now satisfactory. Although producing different profiles, Fig. 3b shows that the two different estimates of the combinatorial entropy (Eqs. 4 and LABEL:mean2), place the bond switching transition (where ) at the same . Indeed, using Eq. (4) without any combinatorial prefactors (but allowing for each / strand only one and one partner) we find that at the bond switching transition decreases by . This indicates that the bond switching observed here is mainly an ‘energetic’ effect, having to do with the different hybridization free energies of the individual bonds.

The fact that the dominant bond-type between colloids switches as the temperature is decreased, gives rise to a more gradual temperature dependence of the attractive interaction between the colloids. To see this, we consider the effective pair potential , where is the separation between the two surfaces. Following Ref. (31), can be computed using thermodynamic integration starting from the case where the DNAs do not hybridise. In that case the pair potential is purely repulsive and is given per unit area by:

(6) |

If we now consider colloids coated with two different types of linkers with constant , we can generalise the arguments of Ref. (31) to show that the full potential is given by

(7) |

and are the mean spacing of the and strands (, ) and denotes the unit area inside which linkages are counted. In Eq. (7) the average is taken with . For the lower limit of the integration, we choose a value for such that . Fig. 4a shows the dependence of the depth of the attractive well of (usually located near due to the steric DNA–surface repulsion (31)) as a function of for three different values of . For sake of comparison we also show the results for a system with the same distribution of and sticky ends, but in which the – linkages are forbidden (equivalent to taking ). This choice allows us to compare between a competing and a single linkage model with equal maximum number of possible bonds. At low only linkages are present, and all the curves approach the same binding strength. Importantly, the figure illustrates that the strength of the attraction changes more gradually with temperature when competing DNA linkers are present - this is especially evident at high . Experimentally, the more gradual response to temperature means that the DNA coated colloids display a broader association-dissociation transition and have a larger range of conditions under which they can form ordered assemblies (instead of kinetically disordered aggregates).

Fig. 4b shows versus the fraction of hybridised bonds , where is the maximum number of available linkages. Compared to single linkage models (black curves), the competing linkages model acquires a reasonable attractive well at lower . It is worth remembering that the strong, high temperature – linkages are replaced with weaker – linkages as the temperature is lowered. The fact that the – linkages disappear before they become prohibitively strong means that kinetically trapped configurations should be automatically avoided, which represents a great experimental benefit. Indeed the formation of an ’irreversible’ – linkage starting from two ’dynamic’ –– ones needs the breaking of two independent linkages. Under bond–switching conditions, this requires an average time () even longer than the life–time of – (). Nevertheless, Fig. 2 does show signs of equilibration problems. Using only sequential single linkage MC moves (full lines in Fig. 2) we found that, at low and high , the system was not able to completely remove the – linkages. The problem is due to the fact that although the formation of two weaker bonds is thermodynamically favourable the system needs to first break a strong linkage. To equilibrate the system (dotted lines in Fig. 2) it was necessary to use the MC move described in App. C, which swaps between a strong linkage and two favourable weaker ones. This implies that, in order to engineer kinetically accessible experiments, it is important to design nucleotide sequences with a small , and to think about possible strategies to push the bond switching transition to higher . This will be investigated in the next section.

### iii.2 Asymmetric DNA model ()

In the previous section we have discussed symmetric systems () for which the switching from one strong bond to two weak bonds was mainly driven by the gain in binding free energy, due to the difference in hybridization free energy of the different linkage types. However, the bond switching can be further enhanced and, more importantly, kinetically enabled, by modifying the number of possible binding partners that each type of linker ‘sees’ on the opposing surface. In doing so, we take full advantage of the multivalent nature of the system, going from a mostly ‘energetically’ driven bond switching mechanism to one that is driven by combinatorial entropy effects.

Eqs. 4 show that multivalence controls the appearance of – or –+– linkages by and , where and is the area of the surface that enclose the tethering points of all the possible targets of ( and ), while is the mean distance between DNA strands. Using simple algebra we find

(8) |

where we have introduced the number of strands () for unit area . One obvious way to favour the weaker – over the stronger – bonds is to decrease the concentration of -terminated DNA () relative to that of the more weakly binding -DNA (). However Eq. (8) suggests that – linkages may be suppressed even more efficiently by decreasing the length . For instance, at , the same effect of having a relative concentration of strands equal to 1%, can be obtained by an asymmetric system with .

Fig. 5a shows and under the same conditions as those considered in Fig. 2, but using (). As compared to Fig. 2 the value for which is higher if , higher if , and higher if . We used both local and linkage swapping MC moves (App. C). Interestingly, the linkages first saturate at a fraction well below (). This means that many and sticky ends remain free and that the linkages are formed much earlier than in the symmetric case (see Fig. 2). The plateau value where saturates depends on the DNA mean distance (Eq. 8) and ranges from if to if . The value of where also has a dependence on the coverage density, when using two lengths for the DNA constructs, but this dependence is found to be small ( for the tested cases).

Fig. 5b shows the effective pair potential per unit area at different values of , keeping fixed at . For large values of , there is a local minimum in the potential at . This minimum becomes a global minimum as is lowered. Important features of Fig. 5b are the presence of repulsive tails that become shoulders at low . By using a Derjaguin approximation (40); (31) it is possible to design sensible interactions that are useful to stabilise complex assemblies.

## Iv Designing the sticky–end sequences

The analysis described in the previous sections allows us to predict the relation between the formation of –, – or – bonds and the hybridisation free energies of the individual strands. In addition, we can now understand the effect of surface coverage and chain length. We find that the dsDNA lengths ( and ) are the most relevant parameters that control the bond switching transition.

We now consider possible choices for the sticky end sequences that in experiments would result in a broadened dissociation transition of the colloids. We need to design our ssDNA sequences such that the sequence binds more strongly to than to (and similarly for ). If we limit the discussion to Watson-Crick base pairs, then different domains of hybridise with and with . In Fig. 6 we consider two architectures for the ssDNA in which the – linkage involves (a) the full length of the reactive fragments or (b) pieces of ssDNA which are not used in the weak bonds. As a proof of concept in Tab. 1, we report possible sequences for these two families (A for Fig. 6a and B for Fig. 6b).

seq. | = | = | ||||
---|---|---|---|---|---|---|

= | = | |||||

seq. | = | = | ||||

= | = | |||||

seq. | = | = | ||||

= | = | |||||

seq. | = | = | ||||

= | = |

Sequences of type have been chosen with (Fig. 6a). This straightforwardly balances the hybridisation free energy of the – and the – linkages (42). The bottleneck is the possibility of – linkages, weaker than –. This problem is avoided by using the scheme of Fig. 6b and less degenerate sequences ( and ) to further enforce the selectivity of each subdomain. In particular has been assembled using for the weak linkages two couples of sequences with nearly equal hybridisation free energy (41).

As input for our model we use the hybridisation free energies of ssDNA fragments in solution. We used the “DINAMelt’ estimates of the ssDNA binding free energies (43). For more details on the procedure, see Appendix B. Fig. 7 shows the relation between and that is predicted for and listed in Tab. 1 and different values of . Based on the MC results of the previous section, we can predict the different bonding regimes both for DNA strands of equal length (Fig. 7a) and for the asymmetric case (Fig. 7b).

Fig. 7 can be used as the starting point to design possible experiments. First, given a certain colloidal architecture, MC results allow to draw the region where none, or linkages are expected. For given ssDNA sequences, it is then possible to predict the transition temperatures computing and and overlapping them with MC results like those in Fig. 7. Sequences can be optimised (e.g. changing in Tab. 1) to avoid kinetically trapped configurations.

## V Conclusions

In this paper we have proposed and tested a strategy that could make the crystallization transition of DNA-functionalised colloids less sensitive to external conditions.

Specifically, we have considered a binary system of colloids ( and ) covered by two families of reactive sticky ends (, on and , on ). Exploiting the selectivity of DNA, we have shown that it is possible to design sequences that only allow for binding between –, – and –. We choose the hybridisation free energy of and in solutions stronger than the other two possible pairings with equal hybridisation free energy. Because or participate in all the possible linkages, the strong – linkages compete with the weaker ones (–+–). However, the number of possible weaker linkages is twice the possible number of the strong ones.

We have demonstrated that a bond switching transition is possible: upon decreasing the temperature the first linkages to appear are the strongest – ones while at lower temperatures the system can gain free energy by replacing – with –+–. For symmetric systems (in which the length of the DNA constructs is equal), the bond switching is energetically driven. This means that the hybridisation free energy of having two weaker linkages is lower than the hybridisation free energy of a single strong one. When choosing the length of the and linkers shorter than and (asymmetric model) the bond switching transition is enhanced and occurs at higher temperatures than in the symmetric case. Here the transition is driven by the combinatorial entropy gain related to the fact that an strand, for instance, can bind more than .

The main effect of the competing DNA linkages is that the resulting effective inter-colloid pair potential is less strongly temperature dependent than is observed in the conventional case where only a single type of linkage is possible. This enhances the experimental control over the self–assembly of DNA functionalized colloids. A further advantage of the proposed strategy is that the strong linkages are replaced as temperature is lowered, which could be used in in step–wise assembly schemes.

Our procedure predicts the temperature range where transitions between different kinds of linkages are to be expected in experiments. Moreover, we have shown how one can optimise the DNA sequences and colloid architectures such that the linkage switching transition remains kinetically accessible.

ACKNOWLEDGMENTS: This work was supported by the ERC (Advanced Grant agreement 227758). DF acknowledges support from a grant of the Royal Society of London (Wolfson Merit Award). The work of ML is part of the research programme of the Foundation for Fundamental Research on Matter (FOM), which is part of the Netherlands Organisation for Scientific Research (NWO). We acknowledge S. Angioletti-Uberti for a critical reading of the manuscript; M. S. Llewellyn-Jones, B. M. Mladek, and F. J. Martinez-Verachocea for useful discussions.

## Appendix A Computation of

Here, we derive expressions for the configurational space of the hybridised strands, , (Eq. 1) as a function of the distance between the tethering points of the hybridised sticky ends , the rods lengths , (), and the angle between and the vector normal to the colloid surface (see Fig. 8). If is the angle^{1}

(9) | |||

(10) | |||

(11) |

where we have defined as the planar angle of a cone (of amplitude ) which is cut by a plane tilted by an angle

If (Fig. 8b) there might be a single cut by the plane where is tethered (Eq. 12), or the second plane could also cut (Eq. LABEL:A1two)

(12) | |||

Notice that when . This happens only if , implying that for the hybridization of asymmetric linkages needs to stay inside a narrower region than when .

## Appendix B Estimate of the configurational entropy

Following (31), we first consider ssDNA fragments in solution (, and ) in thermodynamic equilibrium . In the ideal limit, the partition functions of the three species are

(14) |

, and are the contributions to the partition function due to the internal degrees of freedom (31), are the de-Broglie wave lengths, and the volume and the number of molecules, while is related to the specific potential () by which the two fragments are hybridised. Using classic pair potentials would give, for instance, in the case of a harmonic potential with spring constant , and for a square-like potential with amplitude . At equilibrium the relation links the ratio of the internal partition functions () with the number densities () which can be expressed in terms of the equilibrium constant

(15) |

We now consider the case of tethered rods. When un–hybridised the partition functions of the free fragments are and , where is the configurational space available to the sticky ends (Fig. 1). When hybridised

(16) |

where and are the coordinates of the two sticky ends. Taking the large limit (for spring potentials) or the small limit (for square–well potentials) we find

(17) |

where and are the angles between the rods and the vector linking the tethering points ( in Fig. 8), and as defined in Fig. 1. Finally using Eqs. (17) and (15) we can compute

(18) |

It is important to note that the term which appears in Eq. (18) is specific to the case in which the sticky ends are modelled as point particles, but that it has no physical meaning. For the purposes of the present work, different prefactors give tiny differences which do not affect any of the presented results (compare black symbols in Fig. 7).

## Appendix C Linkage swapping MC moves

Here, we discuss a Monte Carlo scheme which switches between a strong linkage – and two weak ones –+– (see Fig. 9). In the ––+– move, we randomly select an – linkage and two unhybridised strands (, ) from the that could be connected to and (Fig. 9). For the reverse move –+–– we use one of two similar schemes: A) We randomly choose two hybridised couples –– within the possible with and neighbours or B) we randomly choose an – (–) hybridised pair from the () set () and, subsequently, an – (–) from the () set with and that could be connected. In the two cases the acceptance rules are given by:

(19) | |||||