Energy Landscapes of Deoxyxylo- and Xylo-Nucleic Acid Octamers.

Artificial analogues of the natural nucleic acids have attracted interest as a diverse class of information storage molecules capable of self-replication. In this study, we use the computational potential energy landscape framework to investigate the structural and dynamical properties of xylo- and deoxyxylo-nucleic acids (XyNA and dXyNA), which are derived from their respective RNA and DNA analogues by inversion of a single chiral center in the sugar moiety of the nucleotides. For an octameric XyNA sequence and the analogue dXyNA, we observe facile conformational transitions between a left-handed helix, which is the free energy global minimum, and a ladder-type structure with approximately zero helicity. The competing ensembles are better separated in the dXyNA, making it a more suitable candidate for a molecular switch, whereas the XyNA exhibits additional flexibility. Both energy landscapes exhibit greater frustration than we observe in RNA or DNA, in agreement with the higher degree of optimization expected from the principle of minimal frustration in evolved biomolecules.

Introduction

Xeno-nucleic acids (XNAs) are a diverse family of nucleic acid structures derived from DNA or RNA by chemical modification of the sugar moiety of nucleotides.[1] Evolution and heredity have been characterized for these molecules, with associated specific folding and binding properties.[2] XNAs have emerged with important medical applications,[3] for example, as aptamers,[4−6] synthetic ribozymes,[7] artificial small interfering RNAs, and antisense oligonucleotides for the targeting of microRNAs.[8,9] The resistance of XNAs to endonucleases, a result of the inability of natural enzymes to recognize the modified nucleic acid structures, is a particularly valuable advantage of XNAs in replacing natural nucleic acids for therapeutic purposes.[10] Recently, it was demonstrated that XNA aptamers are capable of recognizing small molecules, further advancing the potential therapeutic use of such artificial nucleic acids.[11]

XNAs are also of current interest in the emergent fields of synthetic biology,[12,13] which demands the development of chemical information storage systems capable of self-replication in vitro and in vivo for artificial life and biological computation, and in systems chemistry,[14] which requires molecular switches for the control of operations in complex chemical networks.

Other potential applications of XNAs include their use as self-assembling nanomaterials, broadening the possible design scope in DNA nanotechnology, chemical sensors, and catalysts.[15,16] The study of XNAs is also motivated by the fundamental question of the origins of life, where it is important to understand the factors that led to evolution selecting ribofuranosyl nucleic acids as the genetic biopolymer for the basis of life, and where it remains unknown if an alternative nucleic acid was utilized in hypothetical organisms preceding those based on RNA.[7] Despite their importance, there are few reported structures of XNAs, and relatively little is known concerning the structural and dynamical properties of XNAs in atomistic detail.[17] Furthermore, computational studies have been thus far largely limited to molecular dynamics (MD) simulations, without the application of enhanced sampling schemes, which are required to overcome broken ergodicity.[18]

The present study focuses on nucleic acids based on xylose (XyNA) and deoxyxylose (dXyNA), referred to collectively as XyNAs, which represent some of the simplest possible perturbations to the chemical structure of natural nucleic acids. Xylose is derived from ribose by the inversion of configuration at the C3′ atom of the sugar moiety (Figure ), and likewise deoxyxylose is derived from deoxyribose. Xylose is a thermodynamic product of the formose reaction,[19] the most probable prebiotic route of sugar synthesis,[20] and so XyNAs represent arguably the most credible possibility of a genetic biopolymer adopted by prebiotic organisms that are speculative precursors to RNA-based life forms.

Figure 1

Comparison of a ribonucleotide monomer unit as it appears in canonical A-RNA (left) and a xylonucleotide monomer unit as it appears in the ladder-type conformation of the XyNA1 sequence observed by NMR (right).[36] Xylose is derived from ribose by an inversion of configuration at the C3′ chiral center, indicated by an asterisk. Note that the sugar moieties in both units exist in the C3′-endo conformation, so that the O3′ atom of XyNA units is axial, and that the glycosidic torsion angle differs between the two units.

While dXyNA-homoduplexes exhibit similar thermal stability to their DNA equivalents,[21] thermal denaturation studies have demonstrated that for dXyNA:DNA hybrid duplexes the thermodynamic stability is markedly lowered.[22] In this respect, dXyNA exhibits complementary properties to many other XNAs, which can form a stable duplex through hybridization with DNA or RNA, with strong discrimination against mismatches.[17] While this behavior precludes the use of XyNAs as an aptamer, and for other applications requiring sequence-specific binding, it is an ideal property if XyNAs are to be utilized alongside and independent of natural nucleic acids as an orthogonal information system.

Circular dichroism (CD) studies have shown that XyNA and dXyNA oligomers may adopt a left-handed helical duplex structure, or a structure with an apparent lack of helicity, and that the observed structure depends on the base sequence and sequence length, as well as external factors, including temperature and salt concentration.[21,23,24]

MD simulations of XyNA and dXyNA duplexes of length 8, 13, and 29 base pairs (bp) have revealed the existence of an intermediate between right- and left-handed helical conformations.[25,26] The observed time scale for this transition is of the order of tens of nanoseconds. For XyNA duplexes, it was observed that the left-handed helical structure is not stable, but rather that the terminal regions of the duplex undergo oscillatory movements from coiled to uncoiled states that act to repeatedly screw and unscrew the helix, suggesting structural competition between left-handed helical and “linear” ladder-type structures. This competition between ensembles makes XyNAs a potential molecular switch, and it is also a property unique to XyNAs among the known XNAs.[17] Other recent MD studies of XyNA duplexes in the presence of a carbon nanotube have demonstrated fast spontaneous unzipping as a consequence of the strained backbone, highlighting the potential of XyNAs in gene delivery for therapeutic purposes.[27,28]

Here, the computational potential energy landscape framework[29,30] is employed to investigate the structural and dynamical properties of XyNA and dXyNA duplexes. Discrete path sampling (DPS)[31,32] uses geometry optimization to locate transition states and the local minima they connect, enabling exploration of the energy landscape independent of free energy barrier heights. From the resulting databases, we can estimate and compare the relative thermodynamic stabilities of the three known major conformations, namely, left-handed helical, right-handed helical, and ladder-type structures. The free energy barriers partitioning these major conformational ensembles determine the applicability of XyNA and dXyNA duplexes as molecular switches, which requires two competing funnels to be separated by a barrier that is surmountable at ambient temperatures. Visualization of the free energy landscapes will clearly elucidate structural differences between XyNA and dXyNA duplexes, and between them and their naturally evolved counterparts. This computational framework has been successfully applied to a variety of biomolecular systems before,[30,33] including the B- to Z-DNA transition[34] and the formation of DNA mini-dumbbells.[35]

Methods

Force Field and Initial Structures

The xylo- and deoxyxylo-nucleotide monomer units were constructed with the furanose moiety in the C3′-endo (cf. canonical A-RNA) and C2′-endo conformations, respectively, using the LEAP program of AMBER.[37] Force field parameters were taken from the parm99 force field,[38] incorporating the bsc0 correction[39] for α and γ backbone torsion angles. This choice of parameters was used in previous MD studies of XyNAs[25,26] and has been employed in a number of studies of other artificial nucleic acid systems.[40−42] The bsc0 correction is acknowledged to yield a general improvement in the description of the behavior of both DNA and RNA, including noncanonical structures, and we therefore chose to use it in the present work. As it was previously shown that different parameterizations for the glycosidic torsion angle χ have little effect on the MD trajectories of XyNA sequences,[26] no such reparameterization was used here. Partial charges were obtained by the two-stage RESP fitting procedure[43] at the HF/6-31G* level of theory using ANTECHAMBER.[44] The potential function was correctly symmetrized.[45,46] To represent an aqueous solution environment, we used a generalized Born implicit solvent model with a surface area term[47,48] and an effective monovalent salt concentration of 0.1 M maintained using the Debye–Hückel approximation.[49]

The complete right-handed duplexes of the octameric sequence (5′–3′) [xG–xU–xG–xU–xA–xC–xA–xC–T] (XyNA1) and its deoxyxylose analogue (dXyNA1) were constructed based on a template of canonical B-DNA produced with NAB.[50] Ladder-type structures were obtained from the NMR solution structure of the XyNA oligomer (PDB: 2N4J).[36] From this conformation, a left-handed helical structure was obtained from a short-timescale explicit solvent MD simulation.

Exploration of the Energy Landscapes

The three distinct initial structures were used to initiate basin-hopping (BH) global optimization[51−53] runs to obtain the lowest-energy structures of the right-handed helical, left-handed helical, and ladder-type conformations. We employed group rotation moves[54,55] and simple atomic displacements, discarding all structures with inverted chirality at any center.

The potential energy landscapes were explored using discrete path sampling (DPS)[31,32] to create kinetic transition networks (KTNs),[56,57] using the low-energy structures found in the BH runs as starting points. Transition states were located with the doubly nudged[58,59] elastic band[60,61] algorithm and converged with hybrid eigenvector-following.[62−64] Local minima were then characterized by approximate steepest-descent paths, using a modified version of the L-BFGS algorithm.[65,66] More details of these methods can be found in various reviews.[29,30,33]

After the initial pathways were obtained, we continued the sampling to improve the connectivity of the landscape,[67] remove artificial kinetic traps and high energy barriers,[68] and shorten path lengths.[69] Free energy landscapes at 298 K were calculated using the harmonic superposition approximation.[70] A self-consistent recursive regrouping scheme based on a specified free energy barrier threshold was used to lump minima into free energy macrostates.[71]

Analysis of the Energy Landscapes

We visualize the free energy landscapes using disconnectivity graphs,[72,73] where each leaf corresponds to a free energy group and is colored according to the value of an appropriate order parameter for a representative potential energy minimum of that group. The chosen order parameter is the helical handedness H,[74] in the form originally proposed for describing the B to Z transition in DNA duplexes.[75] The magnitude of H is a measure of the number of turns associated with a helix. Values H > 0 and H < 0 correspond to right- and left-handed double-helical turns, respectively, while values H ≈ 0 indicate approximately zero helicity.

Processing of the data in the KTN was achieved with the CPPTRAJ module[76] of AMBER. The CURVES+ program[77] was used to extract bp-axis, inter-bp, and intra-bp geometric parameters.[78] All analyses exclude the terminal base pairs to reduce the influence of terminal effects. Molecular graphics images were produced using VMD[79] and UCSF Chimera.[80]

Results

Free Energy Landscapes

The free energy disconnectivity graphs at 298 K for XyNA1 and dXyNA1 are shown in Figure a,b, respectively. For both duplexes, the right-handed helical structures correspond to regions of the landscape with high free energy, so that the occupation probability for right-handed helices is very low. The inversion of chirality in the nucleotide units, with respect to natural nucleic acids, seeds an inversion of the preferred helical sense to favor a stable left-handed helix, and XyNAs can be thought of as effective two-state systems, with a dynamic equilibrium between left-handed helical and ladder-type structures.

Figure 2

Free energy landscape for the duplex octamers for (a) XyNA and (b) dXyNA at a temperature of 298 K, calculated using a regrouping threshold of 4 kcal mol–1 and an energy increment of 2 kcal mol–1. The branches are colored according to the helical handedness (H) of a single potential energy minimum representative of each free energy group. Some important representative structures from the different conformational ensembles are shown. The four different nucleotides are colored as follows: XA and dXA in yellow, XC and dXC in orange, XG and dXG in magenta, and XU and dXT in green. The structural representations were created with UCSF Chimera.[80]

There are striking differences between the free energy landscapes for XyNA and dXyNA duplexes. In particular, the ideal left-handed helical state of the dXyNA1 duplex is a relatively well-defined global free energy minimum on the landscape, around 5 kcal mol–1 more stable than the next lowest minimum, which is the most stable ladder-type structure. The free energy barrier for the conversion of the left-handed helical to the ladder-type structure is around 20 kcal mol–1, and for the reverse transition it is around 15 kcal mol–1. In contrast, the free energy barriers separating the two corresponding states of the XyNA1 duplex are around 10 kcal mol–1 in both directions, and the free energy difference is less than 2 kcal mol–1.

Another immediately apparent difference between the free energy landscapes of the two systems is in the distribution of values for the helical handedness order parameter. For the XyNA duplex, a much broader range of values for H is observed, from H ≈ 0 to −3.0 in the low-energy region of the landscape. In contrast, the low free energy region of the landscape for the dXyNA duplex is dominated by left-handed helical structures, with handedness in a narrower range around H ≈ −3.0. Thus, XyNA duplexes are more flexible than their dXyNA analogues.

The ideal ladder-type structure of the XyNA1 duplex compares favorably with the NMR solution structure.[36] As noted by Maiti et al.,[24] the ladder-type structures of the XyNA and dXyNA duplexes are stabilized by the interstrand stacking of adjacent bases, which arises due to the strong inclination of bases with respect to the helical axis (Figure ). However, the ladder-type structure is predicted to have a marginal left-handed helicity, as opposed to the marginal right-handed helicity in the observed structure. The overstabilization of left-handed helical states of XyNA within the force field parameterization implemented in this work is also evidenced by the fact that the left-handed helical structure is erroneously predicted to be the global free energy minimum, although it is not far below the idealized ladder-type structure. The force field is apparently more accurate in reproducing the behavior of dXyNA duplexes, where there appears to be no such bias, and the ladder-type structures are correctly predicted to have a marginal right-handed helicity.

Figure 3

Dinucleotide step in a structure that is part of the lowest free energy group in the conformational ensemble of ladder-type structures for XyNA1. The favorable interstrand stacking of adjacent bases is apparent.

The landscapes of both XyNA1 and dXyNA1 duplexes are significantly frustrated,[81,82] with a number of low-energy minima separated by high barriers. This result contrasts with the free energy landscape reported in a study of the analogous B to Z transition in CG-rich DNA sequences, where there is strong funneling to the native B-DNA state. Frustration in the free energy landscapes of XyNAs is probably attributable to the geometrical frustration of base pairs that prevents the formation of a properly extended linear duplex. Instead, it necessitates that the ladder-type conformation of XyNA duplexes partly bends inwards on itself to maintain an optimal Watson–Crick base pairing of all base pairs, which drives unwinding and rewinding of the left-handed helix. This geometric frustration may be relieved by the adoption of noncanonical base pairings at one or both of the duplex termini, which then allows the linear structures to form a properly extended conformation.

Evolved biomolecules often appear to follow the principle of minimal frustration.[81−83] The evolutionary process favors energy landscapes with a distinct bias toward a well-defined native state and promotes the elimination of kinetic traps on the pathways to this state. This principle, originally formulated for proteins, may also apply to nucleic acids[84,85] and can be extended to multistate systems.[86] The fact that the energy landscapes of XyNAs are significantly more frustrated than the landscapes of DNA and RNA represents an evolutionary argument for the adoption of DNA and RNA over XyNAs.

Pathways for Helix Transformation

Free energy pathways for helical inversion are shown in Figure a,b for the XyNA1 and dXyNA1 sequences, respectively, from a representative right-handed helical structure to a left-handed helix representing the global free energy minimum. For both duplexes, the pathways are overall downhill in energy, and the mechanism features only low barriers, leading to fast kinetics. Furthermore, the interconversion of ladder-type and left-handed helical states occurs by helix winding (unwinding) in the left (right)-handed directions, respectively, propagated inwards from one terminus, rather than from both termini simultaneously.

Figure 4

Free energy pathways for the right- to left-handed helical transitions in XyNA and dXyNA duplexes that make the single largest contribution to the steady-state rate constants. Some representative structures that are key intermediates or transition states are included. (a) Fastest pathway energy profile for the helical inversion transition of the XyNA1 duplex. The steps correspond to stationary points, i.e. minimum-transition state–minimum-transition state-···-minimum. (b) Fastest pathway energy profile for the helical inversion transition of the dXyNA1 duplex.

For the XyNA1 duplex, the helical inversion proceeds via a low-energy ladder-type structure, which subsequently evolves to the left-handed helical state via a transition state ensemble of “kinked” structures. The corresponding pathway for the dXyNA1 duplex is significantly different. In particular, structures with approximately zero helicity represent a high-energy transient state in the early stages of the pathway, and so the mechanism is not mediated by a ladder-type intermediate state, in contrast to the XyNA1 duplex. This result again reflects the greater bias toward the left-handed helical over ladder-type structures in dXyNA compared to that in XyNA duplexes. The transition then continues to progress smoothly with respect to handedness, that is, via a more regular helical structure with a partial left-handed helicity.

Observed Backbone Torsional Angles

The δ dihedral, with a characteristic value of around −40°, which effectively defines XyNA and dXyNA duplexes with respect to their natural analogues, remains relatively stable throughout the pathways. The α, β, and γ backbone dihedrals for the right-handed helical states of XyNA1 and dXyNA1 adopt values similar to those observed in canonical A-RNA and B-DNA and undergo sign inversion in the course of the transition to left-handed helical states, seeded by the sign inversion of the δ dihedral with respect to the natural nucleic acid analogues. For the α and β dihedrals, the values in the left- and right-handed helical states are of approximately equal magnitude but opposite sign, in both XyNA and dXyNA. The behavior of the ϵ and ζ dihedrals exhibits less variance along the pathways, although these angles likewise adopt values of opposite sign to the corresponding angles in canonical A-RNA and B-DNA. The glycosidic torsion angle χ takes one of two predominant values in XyNAs, around −160° (anti) or −80° (syn), and plays an important role in driving the transition.

The evolution of inter-bp, intra-bp, and bp-axis geometrical parameters[78] along the fastest potential energy pathway for the helical transition of XyNA1 exhibits large-scale changes in the bp-axis inclination angle and in the communicative parameters of the inter-bp twist angle, roll angle, and slide distance. Values of the helical rise and helical twist are also diagnostic of each of the three major conformations. Plots for the described properties are provided in the Supporting Information (Figures S1 and S2).

Conclusions

For XyNA and dXyNA duplexes, we observe an equilibrium between the left-handed helical and ladder-type structures. For dXyNA compared to XyNA, the global free energy minimum, a left-handed helical structure, is stabilized with respect to the ladder-type structures. One possible explanation for this effect is the increased solvation of the C2′ hydroxyl group, present only in XyNAs, in the ladder type compared to the helical state. The separation between the two ensembles is better defined for dXyNA than for XyNA, both with respect to the magnitude of the free energy barrier separating the basins and with respect to the helical handedness order parameter. Hence, the XyNA duplexes are more flexible than dXyNA duplexes, and the latter appear to be more promising candidates for use as a molecular switch, or as chemical information storage molecules capable of self-replication, where facile unwinding of a helical structure is undesirable.

The free energy landscapes of both XyNA and dXyNA duplexes are significantly frustrated, highlighting an important factor that may have led evolution to select ribofuranosyl nucleic acids, and not the xylose-based analogues, as the genetic basis for life. The origin of the structural competition evident in these systems is the geometrical frustration that prevents Watson–Crick base pairing without inducing strain in the nucleic acid backbone that must be relieved by bending, or else by the adoption of noncanonical base-pairing modes at one or both of the duplex termini.

Free energy pathways from a disfavored right-handed helical state to a left-handed helical state, which is the global free energy minimum, also differ significantly between XyNA and dXyNA duplexes. For dXyNA, extended linear structures represent an early and high free energy transition state along the pathway, which then proceeds via helical structures with a smooth change in the handedness order parameter. For XyNA, the transition is mediated by low-energy ladder-type structures, which transform to left-handed helical structures via a transition state ensemble of kinked structures. Left-handed helix winding and unwinding transitions of XyNA and dXyNA duplexes are driven by the highly flexible terminal base pairs. The inversion of the δ dihedral angle in XyNAs with respect to their natural nucleic acid analogues seeds a direct inversion, not only in the overall helical sense but also in the backbone dihedral angles. The glycosidic torsion angle also undergoes large-scale changes in the course of the transition, as do certain key geometric parameters, most notably the bp-axis inclination angle.

Future work could address the design and application of XyNA and dXyNA duplexes for molecular devices, for example, by investigating the sequence and length dependence of the propensity for helicity and how the equilibrium responds to environmental conditions.