The Chirality of DNA: Elasticity Cross-Terms at Base-Pair Level

May 24, 2010 - A systematic analysis of B-DNA elasticity cross-terms was performed using molecular .... were subject to 150 ns of productive MD simula...
0 downloads 0 Views 2MB Size
Download PDF
8022

J. Phys. Chem. B 2010, 114, 8022–8031

The Chirality of DNA: Elasticity Cross-Terms at Base-Pair Level Including A-Tracts and the Influence of Ionic Strength Agnes Noy* and Ramin Golestanian* Department of Physics and Astronomy, UniVersity of Sheffield, Sheffield S3 7RH, U.K. ReceiVed: May 6, 2010

A systematic analysis of B-DNA elasticity cross-terms was performed using molecular dynamics simulations of three different duplexes designed to contain all dinucleotide steps including a 6-mer A-tract. The influence of ionic strength was also evaluated by several trajectories of each molecule with different NaCl concentrations at physiological rank. Simulations show DNA flexibility is independent of salt, in agreement with the Odijk-Skolnick-Fixman model. In addition, our results demonstrate DNA asymmetry at this scale is more complex than predicted by long-scale DNA models, with the cross-terms relating twist, slide, roll, and twist on the one hand and tilt and shift on the other being most essential. We find the rest of the coupling terms can be generally discarded due to their lack of correlation, with the exception of purine-purine’s rise-tilt and shift-tilt. More specifically, A-tracts do not present any specific features in terms of their flexibility and chirality properties within those expected for purine-purine steps. Finally, some hints about coupling mechanisms are provided; we suggest cup deformation is mostly responsible for the positive twist-rise correlation at step level, whereas roll-rise and tilt-rise correlations can be understood via changes in stagger. Introduction Although our current knowledge of genome sequence is extensive, the regulation of gene expression is not only coded by the sequence but also by the DNA structure and flexibility. Thus, the deformability of DNA is essential during many life processes inside cells including protein-DNA interactions1,2 and loop-DNA making contact with nonadjacent binding sites of related transcription factors.3 More specifically, several recent experimental studies have indicated that sequence-dependent mechanical properties are crucial in nucleosome positioning.4 However, how these specific interactions are achieved and the way different sequences can construct a global curvature are still poorly understood. A well-known example of such an effect is the curvature caused by the so-called A-tract sequence, a succession of 4-6 consecutive As. Since the 1980s, when gel electrophoresis experiments observed an anomalous retardation behavior in a series of helix-turn separated A-tracts,1 several theories have been developed to explain the origin of this deformation. The main question have been whether the origin of the curvature is in the A-tract itself (wedge model),5 on the junction between the A-tract and the flanking sequences (junction model),6 or in the pyrimidine-purine (YpR) steps, which were seen as the primary source of DNA deformation in X-ray and NMR experiments (non-A-tract model).7 Despite several studies that have tried to combine data from both kinds of experiments8 and perform atomistic molecular dynamics simulations (MD),9 no clear explanation of the effect has yet been found. DNA deformability and curvature at the base-pair (bp) level have been studied by analyzing an ensemble of X-ray structures.2 In a seminal work, Olson et al.10 calculated the elastic constants for 10 different bp steps associated with the geometric parameters that relate successive bp planes: three rotations (tilt, roll, * To whom correspondence should be addressed. E-mail: a.noy@ sheffield.ac.uk; [email protected]. Telephone: +44 (0)114 2224273. Fax: +44 (0)114 2223555.

twist) and three translations (shift, slide, rise) using crystallographic data. Unfortunately, lack of structures in some sequences in the databases makes the DNA experimental conformational map incomplete. Since then, several MD studies have dealt with the calculation of elastic constants in different circumstances including the dependence on main DNA force fields11 or the trajectory length.12 Most notably, the ABC consortium have been working toward characterizing all the different 4-bp sequences in a systematic way.13 The role of the electrostatic energy is another important factor in determining the DNA flexibility due to its polyelectrolyte nature. One might expect that the screening of the phosphate backbone negative charges will enhance the flexibility and the extent of the conformational changes. However, nonlinear screening of the DNA negative charges by monovalent cations present in the solution, as described by the counterion condensation (CC) theory, determines the number of bound counterions per phosphate to be 0.76 (almost independently of the ionic strength).14 As a consequence, the electrostatic contribution to the persistence length predicted by Odijk-Skolnick-Fixman (OSF) model when the CC theory is included15,16 is rather small in the physiological ionic range (0.1-0.5 M) as was confirmed by single-molecule experiments.17 On the other hand, several studies have shown a certain dependence between the widths of the grooves and their Na+ occupancy18 whereas other authors did not find any solid relation.19,20 Because none of these studies have related ion concentration with DNA mechanical properties at this short scale, in this article we evaluated the elastic constants of three different dodecamers in seven different NaCl concentrations (from 0.22 to 1.2 M) to see whether they change with ionic strength. Obtaining the elastic constants allows us to analyze indirect readout of DNA by proteins. Moreover, the bp elastic constants could be used to construct a coarse-grained mesoscopic model of DNA that can be applied in specific DNA fragments with biological interest. MD-derived stiffness analysis has already

10.1021/jp104133j  2010 American Chemical Society Published on Web 05/24/2010

DNA Elasticity Cross-Terms at Base-Pair Level given promising results in tracing promoter location21 or understanding Lac operon regulation.22 However, the number of parameters that should be introduced in a mesoscopic model is rather large if we use the full (6 × 6) elastic matrix for all the different steps (10), or even more if the 136 different 4-bp sequences were to be considered. This has raised concern since it causes a sizable complexity in multiscale simulations. Using general geometrical arguments, Marco and Siggia suggested that twist-roll and twist-stretch are the most significant (off-diagonal) elastic coupling constants.23 They predicted a positive value for the latter: DNA should unwind under tension; nevertheless recent single-molecule experiments in a low-force regime,24,25 all-atom simulation26 as well as X-ray analysis27 suggest a relation in the other direction: surprisingly, DNA overwind under tension. On the other hand, preliminary studies of X-ray structures2,10 only showed twist-roll, roll-slide, and twist-slide as measurable couplings. These results suggest that, whereas it could be the case that not all of the elastic matrix elements are necessary to describe DNA elasticity, the complexity of the structural properties of relatively short DNA fragments necessitates a more elaborate description in terms of an effective elastic model. In this article, we aim to reexamine the issue of DNA flexibility with a systematic analysis of the cross-terms in the elastic matrix and attempt to help draw a more rigorous picture of DNA deformability. Our main purpose is to determine whether each cross-term is important or can be eliminated due to lack of systematic correlation between the two relevant parameters. Moreover, we intend to understand the mechanism by means of which a coupling is achieved as well as the resulting sign of the correlation. In particular, we determine whether cross-terms present a uniform behavior along all steps, which would imply a generic asymmetry of DNA at step level, or acquire different signs depending of the type of the base step, which would cause sequence-dependent variation in mechanical properties at this scale. To this end we recalculate the elastic matrix at the base step level using a reliable ensemble determined through 150 ns MD trajectories for three different dodecamers containing several copies of each nonredundant dinucleotide step in seven different NaCl concentrations (covering in toto more than 3 µs). We find that the elastic constants as well as cross-terms are independent of the salt concentration in the considered range of ionic strengths, which means that they can be considered as independent measures. This gives us the possibility to perform statistical analysis of the stiffness matrix and thus obtain a more robust picture of DNA flexibility. Finally, we calculate, for the first time according to our knowledge, the elastic constants and cross-terms of an A-tract to look for special features that could shed some light on the origin of the A-tract curvature. Methods Molecular Dynamics. Three different DNA dodecamers with at least one copy of the 10 unique dinucleotide steps and a 6-mer A-tract were used to calculate the elastic matrix: Dickerson’s dodecamer (DD,28 (CGCGAATTCGCG)2), (CATAGGCCTATG)2 and (CACAAAAAACAC)2. Duplexes were built from Arnott B-DNA models29 using the NUCGEN module of AMBER830 and with the force field AMBER parm99 and parmbsc0 corrections.31–33 All structures were solvated in ∼60 Å truncated octahedron boxes (∼5300 water molecules), 22 Na+ counterions to balance the DNA charge and different number of additional Na+/Cl- ion pairs (0, 10, 20, 30, 50, 70, 100) corresponding to the following Na+ concentrations respectively:

J. Phys. Chem. B, Vol. 114, No. 23, 2010 8023 0.22, 0.33, 0.42, 0.54, 0.75, 0.98, 1.22 M. To avoid initial bias, ion positions were randomized by swaps with water molecules such that no ion was closer than 5 Å to DNA and 3.5 Å to any other ion. SPC/E34 water molecules together with Smith and Dang ions models35 were used to represent solvent molecular interactions due to their capacity to describe properly high ionic strength solutions.36 The systems were minimized, thermalized (T ) 298 K) and equilibrated using our standard equilibration protocol37 doubling the length of the individual periods followed by a final re-equilibration of 10 ns. The equilibrated structures were subject to 150 ns of productive MD simulation at constant temperature (298 K) and pressure (1 atm) using periodic boundary conditions and particle mesh Ewald.38 SHAKE39 was used for hydrogen atoms in conjunction with an integration time step of 2 fs. All MD simulations were performed using the PMEMD module of the AMBER8.0 computer program.30 Elastic Model at Base-Pair Step level. DNA trajectories of 150 ns were analyzed at the bp step level by three rotations (roll, tilt, and twist) and three translations (shift, slide, and rise) obtained using the 3DNA code.40 Then, the elastic matrix F was calculated by:10

F ) kBTC-1

(1)

where C stands for the covariance matrix in the helical space. Thus, elastic energy for the deformation of DNA bp step level in units of thermal energy was written as:26

E b ) kBT 2

[∑

Kiω20(Xi - X0)2 +

i(tra)

∑ Ki(Xi - X0)2 +

i(rot)

∑

2Kijω20(Xi - X0)(Xj - X0) +

i>j(tra, tra)

∑

2Kij(Xi - X0)(Xj - X0) +

i>j(rot, rot)

∑

]

2Kijω0(Xi - X0)(Xj - X0)

i>j(rot, rot)

(2)

where b stands for the dinucleotide step length, Ki are the elastic constants, and Kij the corresponding coupling terms. The distances and angles are interconverted by the factor ω0 ) 2π/p where p ) 3.4 nm is the standard helical pitch of B-DNA. Persistence Length. The persistence length at the step level -1 -1 + Ktilt ). Because DNA is a was evaluated via P-1 ) 1/2(Kroll polyelectrolyte, its persistence length can be decoupled into the nonelectrostatic (P0) and the electrostatic (Pe) contributions. Using a generalization of the OSF model that includes the CC theory, we can write the persistence length as:15,16

P ) P0 + Pe ) P0 +

1 4κ2l B

(3)

where 1/κ is the Debye screening length and l B is the Bjerrum length (7.14 Å in water at 25 °C). Base-Pair Parameters. Base-pair parameters (propeller twist, buckle, open, stretch, stagger, and shear), which describes the geometrical relationship between the bases forming the canonical Watson-Crick bp, were calculated using the 3DNA code40 along all DNA trajectories. We subsequently computed the sum and the difference of these parameters (Y) considering two consecutive bp (i and j ) i + 1) to achieve equivalence with step parameters:

8024

J. Phys. Chem. B, Vol. 114, No. 23, 2010

Yi+j ) Yi + Yj Yi-j ) Yi - Yj

Noy and Golestanian

(4)

Note that Yi+j becomes significant when the deformations of two consecutive bp are similar whereas Yi-j is important when the two bp move in the opposite sense (see Figure 1). Finally, the correlation between these new parameters Y′ and bp step ones X denoted as FX, Y′ were evaluated to understand the local origin of the bp step couplings:

FX,Y' )

cX,Y' σXσY'

(5)

where cX, Y′ is the covariance and σX and σY are the corresponding standard deviations. F2 was also probed as a measure of the percentage of the variance of the variables that are correlated with each other. Technical Details. Central 8-mers of the three sequences were used to avoid end effects and average over all the occurrences of a dinucleotide step was done. In the case of the A-tract, only three internal AA steps of the 6-mer A-tract were taken into account to probe the specific elastic profile. Besides, elastic constants, cross-terms, and correlations from seven different simulations of the same sequence were considered as independent measures due to the results being insensitive to the ionic strength in the range considered. Hence, average (avg) and standard deviation (sd) were calculated. Results and Discussion General Structure. RMSd and backbone analysis (Supporting Information, Figure S1 and Table S1) show different DNA trajectories sample regions of conformational space that are expected for a B-DNA. In agreement with experimental data and previous MD simulations,12,41 noncanonical North puckering is present around 5% of the time, B-II ε/ζ torsion conformation appears in 15% of snapshots approximately, and noncanonical R/γ torsion values are taken around 3%. Extra analysis were done over A-tract trajectories due to the fact that this sequence adopts a peculiar conformation. Afterward, results were compared with NMR (1FZX42) and crystallographic (1D98,43 1D8944) 6-mer A-tract structures (Figure S2, Supporting Information). Similar to NMR structures, MD ApA steps present slightly negative values of tilt and an almost negligible bend in terms of roll compared to the 5° average value of a random sequence.11 In addition, trajectories reproduce quite well the strong negative propeller twist, the decrease of buckle observed from the 5′ to the 3′ junction, and the difference in puckering between adenines and thymines observed in NMR data.42 In summary, these results suggest our DNA duplexes show a stable but flexible behavior in reasonable agreement with experimental data as well as previous simulation studies. DNA Deformation and Ionic Strength. As has been previously discussed,12 the actual time scale of the MD simulations (100-200 ns) is adequate to provide reliable estimations of the bp step elastic matrix. Although sequencedependent effects have been found to be significant beyond the nearest neighbor bp level, the dinucleotide model will be sufficient to treat the questions that we aim to address here in this study. In general terms, MD simulations reveal that elastic constants are quite stable with regards to increasing ionic concentration (see Figure 2). This finding is in general agreement with the

Figure 1. Schematic diagram of GC where subtracted-buckle or cup (A) and added-buckle (B) are significant (see Methods). (A) Images are generated with 3DNA building module using average values of step parameters and buckle from snapshots with rise and buckle higher/ lower than corresponding avg ( 2 · sd, whereas in (B) standard step parameters were used with sizable buckle.

single-molecule experiments17 and the OSF model, which evaluates the electrostatic contribution as at most 2% of the persistence length in the studied ionic conditions. Thus, MD simulations with explicit water and ions suggest the OSF model, with the CC theory included, would be valid in step level. This reflects the atomic nature of the CC theory, which states that cations are trapped around DNA-phosphates for sufficiently low ionic strength. This minimal concentration is determined by the balance between the loss of entropy and the gain in electrostatic potential energy for the cations (approximately 0.05 M for sodium).14 General Properties of Flexibility. A further inspection of the results in Figure 2 confirms the exceptional softness of YR steps in all helical parameters as has been seen in previous studies.10,26 One can note that the A-tract elastic constants have intermediate values without presenting any special feature: they are neither the most flexible as implied by the wedge model, nor are they the most rigid as the non-A-tract theory suggests. Note also that the junction steps do not seem to stand out with any particular property. If individual steps are observed in experimental structures,42 only a small roll and tilt is found in accordance with our trajectory values (Figure S2, Supporting Information). Therefore, our results suggest that it is not possible to explain the origin of the strong curvature at dinucleotide scale. This is not in contradiction with other theoretical and experimental conformations, since the 19° bend for A-tract predicted from gel electrophoresis data is obtained when the global bend of a full DNA turn is considered.42,45,46 How minor roll and tilt values build the global bend remains unclear. That is why we pay attention to the dinucleotide cross-terms, to see if there are special chiral properties that confer a characteristic anisotropy capable of explaining the strong curvature (see below). Despite the difference in the length scale between DNA fragments used in single-molecule experiments and the dinucleotides step levels in our case, rotational elastic constants are notably similar. In the case of twist, experimental values range from 80 to 120 nm,47,48 whereas our values range from 20 to 120 nm (see Figure 2). We find that this parameter is significantly more sequence dependent than the other rotational elastic couplings, which is consistent with previous reports in the tetranucleotide context.13 Considering the persistence length, an average value around 40 nm is found in our simulations, whereas the experimental dynamical persistence length has been found to be slightly more rigid (70 nm with magnetic beads49 and 80 nm using cryo-electron microscopy50). However, recent atomic force microscopy,51 FRET/SAXS measures52 as well as several cyclization experiments53,54 suggest that DNA might be

DNA Elasticity Cross-Terms at Base-Pair Level

J. Phys. Chem. B, Vol. 114, No. 23, 2010 8025

Figure 2. Force constants (in nm) of helical coordinates associated with the deformation of different bp steps increasing the ionic strength (in M). A-tract is defined by considering the three internal AA steps of the 6-mer A-tract (in yellow). Values reported here are averages over all occurrences of a dinucleotide step and over all different steps (in black).

more flexible at smaller scales, putting our values inside the relevant experimental range (see Figure 2). If we calculate the sequence-average rise elastic constant, we obtain a value of 1950 pN that is almost twice the stretch modulus found in single-molecule experiments (1100 pN;55). The restricted variation of rise was previously noted by X-ray analysis,2 and it is usually related to the strong base-stacking interactions.56 Nevertheless, rise is stiffer than movements in the other two local directions, namely, shift and slide, which have average values of 382 and 558 pN, respectively. Because bases present a certain inclination and displacement with respect to a global molecular axis, the movements perpendicular to stacking interactions would contribute toward enhancing fluctuations of the global axis. Thus, if we substitute the three local translational parameters of eq 2 by a single global axis estimated by the CURVES program,57 we obtain a stretch elastic constant of 1360 pN which might be described as a kind of average between the three local force constants and is reasonably close to the experimental value of the stretch modulus. The magnitudes of the coupling terms not only give us information about the independence or correlation between two parameters but also reflect the degree of chirality of DNA. Each cross-term effectively shows how different it is to deform the polymer along a given geometrical variable in opposite directions.23 This can be important to build an overall curvature from short-scale flexibility since a certain anisotropy would be induced on the DNA molecule. Thus, they become essential for a complete description of the DNA elastic scenario and for a construction of a notable mesoscopic model. We computed them for all of the 10 unique steps as well as the A-tract using the three different dodecamers over all ionic strength range simulations. Results are presented as follows: Figure 3 shows the coupling terms between the translation and the rotation

parameters, whereas the translation and the rotation cross-terms among themselves are presented in Figure 4. Twist, Slide, and Roll: B-A Motion. Twist-slide and roll-twist coupling terms present a uniform behavior along all different steps (Figure 3 and Figure 4). The sequence-average roll-slide cross-term might appear to suggest that these two parameters are independent. However, if attention is paid to the different steps, a positive value is observed for most of the dinucleotides except for AG, GA, and AA steps that have negative values (Figure 3). As was already noted by X-ray pioneering studies,2,10 all these values would reflect a coordinated motion similar to the B-A DNA transition. Especially, El Hassan and Calladine2 used this model to classify the deformability of different steps, suggesting that YR, GG, and GC presented a spread of structures between B-DNA and A-DNA conformations, whereas AG, GA, and AA were mostly considered rigid. It is worth mentioning that A-tracts exhibit a big negative roll-slide term, which is in stark contrast to the behavior of the other steps. In summary, we can assert that these cross-terms are crucial for explaining the chirality of DNA at the step level since they were described as significant in crystal data analysis,2,10 and they build a global motion similar to B-A transition that was identified as the first essential mode of a full turn of DNA.58 Twist-Rise: B-A or Buckle Movements? It has been recently suggested that24,25 B-A movements can provide an explanation to understand the counterintuitive relation between stretch and twist detected in recent single-molecule experiments:24,25 DNA overwinds under tension. The same sense of correlation has been found between rise and twist parameters in both allatom simulations26 and X-ray analysis.27 In accordance with the above model, simulations present negative rise-twist coupling terms for all steps (Figure 3). However, in the context of the

8026

J. Phys. Chem. B, Vol. 114, No. 23, 2010

Figure 3. Cross-terms (in nm) between translational and rotational parameters associated with the deformation of different bp steps. Values reported here are averages over all trajectories with different ionic strengths, and corresponding standard deviations are given as a error bars.

B-A transition model all the rise-slide cross-terms also need to be negative, which is contrary to our findings as Figure 4 demonstrates, in agreement with previous all-atom simulations26 and X-ray analysis.27 To find out the origin of this coupling from the perspective of the local movements, we looked for measurable correlations between slide or rise and bp parameters (see Methods). Figure 5A indicates a high positive correlation between rise and buckle (Fj ) 0.83) and a sizable anticorrelation between slide and buckle (Fj ) -0.21), providing a hint for the anticorrelation between rise and slide since no other significant values were found (Figure S3 in Supporting Information contains results for step and bp parameters in a systematic way). Specifically, the crucial bp parameter seems to be the subtraction of two consecutive bp buckles, which is usually called cup59 which reaches big values when bp’s deform in the opposite sense (see Methods and Figure 1). Nonetheless, caution is needed because Fj2 between slide and cup is 0.05 which means only 5% of the slide variation is explained by the cup parameter. Therefore, to assess whether rise and slide are anticorrelated due to cup movement, we analyzed the relevant configurations by choosing the snapshots with rise and cup higher/lower than avg ( 2 · sd, respectively,

Noy and Golestanian from all simulations. We observe significant slide differences between both sets of structures, hence confirming this mechanism with the only exception of CG which is the only step with positive cup-slide correlation (Figure 5B). Interestingly, cup also correlates very well with twist (Fj ) 0.35, Figure 5A), which would cause an alternative way to achieve positive correlation between rise and twist: when two neighboring bp buckle in different senses, widening the center of the step, rise and twist increase, whereas slide decreases (see Figure 1). If we take the same snapshots as before (rise and cup higher/lower than avg ( 2 · sd, respectively), a difference in twist is obtained between the two groups (see Figure 5C). Although there is no clear distinction for the case of RY steps, in general we observe that rise-cup is a mechanism for changing the twist. Thus, what is the DNA motion that causes positive correlation between twist and rise at this scale; the B-A transition, the buckle movement, or both? To clarify this issue, we selected the snapshots in all simulations that had rise and/or twist values higher/lower than avg ( 2 · sd. Differences of cup between those sets of snapshots with rise or twist higher/lower than avg ( 2 · sd were found to be significant and consistent along all steps (see Figure 5E). When we took into account both parameters together, only a slightly bigger difference was observed. This indicates redundancy in the chosen snapshots and thus a similar mechanism to achieve cup differences independently of the used criterion. As a consequence, all of these results confirm the participation of this mode of deformation in the positive rise-twist correlation. On the other hand, we found a sizable increase in slide for extreme values of twist, whereas we found a slide decrease for extreme values of rise, as could be predicted by correlation signs (see Figure 5D). The interesting point here is to observe that, when both parameters are considered together, differences are smaller than when only twist is used to build the threshold. This suggests that the set of snapshots with high/ low twist presents the same slide-rise relation as the full trajectories, not showing a hidden positive slide-rise correlation. In summary, we find that the cup mode is likely to be responsible for the positive rise-twist correlation at this short scale. Interestingly, both Dickerson’s59 and Jame’s60 groups found the same rise-cup correlation when analyzing individual X-ray and NMR structures respectively, whereas Sponer and Kypr61,62 confirmed it by means of base-stacking calculations. In addition, some steps in DNA molecules that bind to different proteins present low cup-twist-rise or vice versa, which shows the biological importance of this deformation (see Figure S6, Supporting Information). It is worth mentioning that a different mechanism will cause the stretch-twist correlation observed in the single-molecule experiments since we analyzed flexibility of DNA at step scale, whereas these experiments are on 1-10 kilo-bp DNA fragments. The cup deformation cannot contribute to this large-scale deformation because it cannot be propagated along neighboring steps. In any case, trajectories show twist-rise coupling is important for understanding the asymmetry of DNA at the 2 bp range. Roll-Rise. From Figure 3, we find that simulations suggest a varied picture for the coupling between roll and rise coordinates. Despite the average value indicating a negative coupling, this behavior is not uniform over all steps. Specifically, YR and RR steps present negative values, whereas RY crossterms are positive. Considering that, in general, roll average values are positive for YR and RR steps (bend to major groove) and negative for RY (bend to minor groove), we could deduce

DNA Elasticity Cross-Terms at Base-Pair Level

J. Phys. Chem. B, Vol. 114, No. 23, 2010 8027

Figure 4. Cross-terms (in nm) of translations among themselves as well as rotations among themselves associated with the deformation of different bp steps. Values reported here are averages over all trajectories with different ionic strengths, and corresponding standard deviations are given as a error bars.

that rise increases when the absolute value of roll increases. Although this finding is rather surprising as polymer bending usually means shortening, it agrees with the crystallographic data and other MD simulations.26,27 This could be understood by taking into account the fact that the elastic matrix has been calculated at the bp step scale, which implies that roll deformation could only be achieved by an increase of the mean separation between both bp to avoid steric clashes. To confirm this hypothesis, we looked for local correlations between rise and roll and bp parameters (see Figure 6). Interestingly, we found a high correlation between rise and added-stagger with negative sign for YR steps (FjYR ) -0.54), positive for RY (FjRY ) 0.34), and negligible for RR (FjRR ) -0.05), reflecting the pattern of roll--rise cross-terms (see Figure 3). This finding is supported by the significant correlation between roll and added-stagger (Fj ) -0.23, Figure 6) although we need to keep in mind that only 5% of the variation is shared between both parameters. Thus, to determine whether the roll-rise correlation is reached through stagger, we used a strategy similar to that before: we choose snapshots with roll and stagger higher/lower than the corresponding avg ( 2 · sd, and examined whether we would obtain significant differences in rise. The results in Figure 6 support this view: when roll increases, rise increases for YR steps, whereas it decreases for RY, and the differences are small for RR steps. Figure 7 and Table 1 help us understand why stagger is involved in roll-rise coupling. In the case of rise, it is rather intuitive: when R bases are near, interstrand steric clashes are important, hence, rise needs to augment (see Figure 7A,C); whereas when Y bases are near, bp can approach more closely, hence, decreasing rise (see Figure 7B,D). On the other hand, to interpret roll-stagger coupling we use Calladine’s steric clashes model that was used to rationalize why YR steps tend to bend

toward the major groove and RY to the minor groove. In YR steps, when R bases are near due to stagger (see Figure 7A), R interstrand steric clashes in the minor groove cause bending toward the major groove (see Figure 7E and Table 1). However, if stagger is in the opposite sense and Y bases are close (Figure 7B), distances between Rs augment, and C-N4/T-O4 with G-O6/ A-N6 distances become the smallest ones, which now raises the likelihood of a steric clash in the major groove instead (Figure 7E and Table 1). In the case of RY steps, major groove steric clashes between interstrand Rs are the important ones, which prevent a positive roll (Figure 7C,F and Table 1). Nonetheless, when stagger changes (Figure 7D), distances among intrastrand R-Y atoms in the minor groove become the smallest ones, and this switches the sense of roll (Figure 7F and Table 1). It is worth mentioning that Sponer and Kypr62 deduced the same roll-stagger correlation as well as similar scheme of Figure 7 only with stacking calculations, which would agree with the steric clashes explanation provided above. Tilt-Shift. The tilt-shift coupling appears to be homogeneous for all steps, as has been previously found by other authors.10 Although no clear correlation is found between the two parameters and the set of bp ones (Figure S3, Supporting Information), the mechanism can be readily understood: when a bp moves to a backbone (absolute value of shift increases), a bend toward the same backbone is achieved more easily (tilt presents same sign as shift). In the case of YR and RY steps, tilt and shift parameters remain independent from the other step parameters (slide, twist, roll, and rise) as can be seen in Figures 3 and 4 where the couplings between the two groups are negligible and the measures are spread over positive and negative values. This suggests that both strands can be considered equivalent so that a rotation or translation toward any backbone should not

8028

J. Phys. Chem. B, Vol. 114, No. 23, 2010

Noy and Golestanian

Figure 5. General characteristics of twist, rise, and slide coupling. (A) Correlation between these step parameters and cup. Slide (B) and twist (C) differences between two sets of snapshots where rise and cup are higher/lower than the corresponding avg ( 2 · sd. Slide (D) and cup (E) differences between two sets of snapshots where twist and/or rise are higher/lower than the corresponding avg ( 2 · sd. Cup refers to the subtraction of buckles from neighboring bp (see Methods and Figure 1). Values reported here are averages over all trajectories with different ionic strengths, and corresponding standard deviations are given as a error bars.

TABLE 1: Roll (deg), Rise and Stagger (Å) from structures of Figure 7 as well as Distances over Major and Minor Groove (Å) from Structures Generated by 3DNA Building Program with Roll ) 0°, Rise ) 3.3 Å, but Original Stagger to See Possible Steric Clashes When Roll Is Imposed CG (Figure 7A) roll rise stagger Major GrooVe Gua O6/Gua O6 Cyt N4/Cyt N4 Cyt N4/Gua O6 Minor GrooVe Gua N2/Gua N2 Cyt O2/Cyt O2 Gua N2/Cyt O2

23 3.9 -0.6

CG (Figure 7B) -5 2.7 0.8

GC (Figure 7C)

GC (Figure 7D)

-8 3.9 0.8

14 3.0 -0.7

3.8 4.8 3.5/3.5

4.9 3.8 3.5/3.5

2.7

4.1

3.0

4.3

3.2 5.5 3.5/3.4

4.4 4.5 3.4/3.5

influence the flexibility of DNA among all other step parameters. In the case of AT, TA, CG, and GC, the previous statement holds, and it can be easily extended to TG and GT. When constructing a mesoscopic model, the 6 × 6 elastic matrix could be considerably simplified by discarding the following eight cross-terms: tilt-slide, tilt-twist, tilt-roll, tilt-rise, shift-slide, shift-twist, shift-roll, and shift-rise. On the other hand, if we consider RR steps, we realize that both strands are quite different since now the two purines and pyrimidines are separated in opposite strands. As a consequence, an axis of symmetry is broken, leaving this kind of step more

chiral than the others, as reflected by the sizable cross-terms between tilt-rise (FjRR ) 0.38) and shift-twist (FjRR ) -0.40). The other terms (shift-slide, shift-rise, shift-roll, tilt-roll, tilt-twist, tilt-slide) could also be discarded since these couplings mostly represent correlations of less than 0.2, which means that only 4% of the variance of one quantity is explained by the other (Figures S4 and S5, Supporting Information). Whereas there is no significant correlation between shift/twist and bp parameters that help us to explain a possible mechanism for this coupling (Figure S3, Supporting Information), rise and tilt present a high correlation with subtracted-stagger (FjRR )

DNA Elasticity Cross-Terms at Base-Pair Level

Figure 6. General characteristics of roll and rise coupling. (A) Correlation between these step parameters and stagger. (B) Rise differences between two sets of snapshots where roll is higher/lower than avg ( 2 · sd and stagger lower/higher than avg ( 2 · sd. Stagger refers to the addition of values from neighboring bp (see Methods and Figure 7). Values reported here are averages over all trajectories with different ionic strengths, and corresponding standard deviations are given as a error bars.

Figure 7. Side views of rise-, roll-, and stagger-coupled motion in the case of CG as an example of YR steps (A and B) and GC as an example of RY steps (C and D). Images are generated with 3DNA building module using average values of step parameters and stagger from snapshots with roll and stagger higher/lower than corresponding avg ( 2 · sd (see Table 1). Axial view of CG (E) and GC (F) to see possible steric clashes along minor (N2 and O2 atoms) or major grooves (O6 and N4) when positive or negative roll is introduced. Stagger refers to the addition of values from neighboring bp (see Methods). C bases in green, G in gray, C atoms in black, N in blue, and O in red.

0.52 and Fj ) 0.67, respectively, Figure 8). We proceed as before by selecting snapshots with tilt and stagger higher/lower than avg ( 2 · sd correspondingly, and find significant differences in the rise of RR steps, thus confirming this model. Figure 9 shows how stagger facilitates tilt and hints at the possible role of rise: when DNA bends toward the RR strand, intrastrand purines would present more overlap than pyrimidines, making it possible for the latter ones to get closer.

J. Phys. Chem. B, Vol. 114, No. 23, 2010 8029

Figure 8. General characteristics of tilt and rise coupling. (A) Correlation between these step parameters and stagger. (B) Rise differences between two sets of snapshots where tilt and stagger are higher/lower than the corresponding avg ( 2 · sd. Stagger refers to the subtraction of values from neighboring bp (see Methods and Figure 9). Values reported here are averages over all trajectories with different ionic strengths, and corresponding standard deviations are given as a error bars.

Figure 9. Side views of rise-, tilt-, and stagger-coupled motion in the case of GA as an example of RR steps (A and B). Images are generated with 3DNA building module using average values of step parameters and stagger from snapshots with tilt and stagger higher/lower than corresponding avg ( 2 · sd. Axial view of GA (C) to see the different overlaps Rs and Ys. Stagger refers to the subtraction of values from neighboring bp (see Methods). C bases in green, G in gray, C atoms in black, N in blue, and O in red.

To sum up, cross-terms relating shift-tilt with other parameters could be eliminated from a mesoscopic model with the exception of twist-shift and tilt-rise for RR steps. Tilt-rise together with roll-rise are couplings with local nature, where stacking between two bp are involved by means of stagger parameter. Because both cross-terms depend on the nature of the step, they cause a specific pattern of deformation that might be important for discriminating different sequences in some DNA-protein interactions, especially those with any intercalation. On the other hand, tilt-shift coupling was detected as one of the main motions of a DNA turn,10 thus being involved in the general asymmetry of short fragments of DNA. A-Tracts. We find that A-tracts have stiffness properties in accordance with the RR properties: they have intermediate values for their elastic constants and corresponding averages, as well as relatively high tilt-rise and twist-shift couplings.

8030

J. Phys. Chem. B, Vol. 114, No. 23, 2010

In addition, the high propeller-twist, which was described as a unique property of this sequence, does not seem to be correlated with any step parameter (Figure S3, Supporting Information) in agreement with other works.63 Hence, our results do not reveal any specific feature that could explain the origin of the curvature anomaly from the step scale, although they do not contradict other theoretical and experimental studies since always an overall bend of a full DNA turn is observed and roll/ tilt values found are quite small.42,45,46 Specifically, cross-terms do not reveal any extraordinary chiral properties that could cause a remarkable anisotropy. We note, however, that the roll-slide coupling for A-tracts is quite different from all the other steps, although the difference may not be enough to explain the anomalous properties. Therefore, it is likely that the curvature would be caused by a combination of subtle differences, rather than by an exceptional deformability quality. Concluding Remarks Our simulations strongly suggest that DNA flexibility is stable at the physiological range of monovalent ion concentrations in agreement with the OSF model including the CC theory. In addition, bending and twisting elastic constants are reasonably close to the experimental values despite the different length scales. Rise elastic constant is very stiff because of the solid stacking interactions, whereas deformations along the other two directions (shift and slide) are easier. If we consider a global molecular axis, an intermediate stretch elastic constant close to the single-molecule value is obtained. The inclination and displacement of bases with respect to this axis invokes the participation of the three directions in stretch fluctuations. Generally, cross-terms relating twist/slide/roll/rise on one side and tilt/shift on the other are significant and uniform along all steps, which is what defines the generic asymmetry of DNA at the step level. As a consequence, the coupling terms between the two groups of helical parameters could be discarded from mesoscopic models due to the low correlation (less than 0.2) with the exception of RR’s tilt-rise and twist-shift. These couplings together with the roll-rise are the only ones that depend on the step nature, and could thus create a specific pattern of deformation. Some hints about coupling mechanisms are provided with the help of bp parameters. Simulations suggest dinucleotide twist and rise are correlated by means of cup, although this deformation mode will not participate in the long-scale stretch-twist coupling seen in single-molecule experiments due to cup not being propagated along neighboring steps. Moreover, our results show that roll-rise and tilt-rise are also couplings with local nature, where stacking between 2 bp are involved by means of the stagger parameter. Finally, neither the elastic constants and average values nor the coupling terms of A-tract exhibit any extreme value at dinucleotide scale that could readily explain their peculiar curvature properties. Elastic and chiral properties of this sequence are inside the framework of RR steps deformability. Thus, to understand the origin of their singular curvature, some calculations combining the subtle mechanical differences would be necessary. Acknowledgment. We thank Prof. J. Sponer and Prof. M. Orozco for useful discussions and comments, and we thank the Red Espanola de Supercomputacio´n for computational resources. A.N. is a EMBO long-term fellow. Supporting Information Available: Complete refs 13 and 30, Table S1 and Figures S1-S6 which contain general analysis

Noy and Golestanian of simulations, correlations of step parameters among themselves and with bp parameters, and some examples of protein-DNA complexes. This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes (1) Hagerman, P. J. Annu. ReV. Biochem. 1990, 59, 755–781. (2) Hassan, M. A. E.; Calladine, C. R. Philos. Trans. R. Soc. London, Ser. A 1997, 355, 43–100. (3) Haldorf, S. E.; Welsh, A. J.; Szczelkun, M. D. Ann. ReV. Biophys. Biomol. Struct. 2004, 33, 1–24. (4) Segal, E.; Fondufe-Mittendorf, Y.; Chen, L.; Thåstro¨m, A.; Field, Y.; Moore, I. K.; Wang, J.-P. Z.; Widom, J. Nature 2006, 442, 772–778. (5) Trifonov, E. N.; Sussman, J. L. Proc. Natl. Acad. Sci. U.S.A. 1980, 77, 3816–3820. (6) Koo, H.-S.; Crothers, D. M. Proc. Natl. Acad. Sci. U.S.A. 1988, 85, 1763–1767. (7) Maroun, R. C.; Olson, W. K. Biopolymers 1988, 27, 585–603. (8) Zhurkin, V. B.; Ulyanov, N. B.; Gorin, A. A.; Jernigan, R. L. Proc. Natl. Acad. Sci. U.S.A. 1991, 88, 7046–7050. (9) Dixit, S. B.; Pitici, F.; Beveridge, D. L. Biopolymers 2004, 75, 468–479. (10) Olson, W. K.; Gorin, A. A.; Lu, X.-J.; Hock, L. M.; Zhurkin, V. B. Proc. Natl. Acad. Sci. U.S.A. 1998, 95, 11163–11168. (11) Perez, A.; Lankas, F.; Luque, F. J.; Orozco, M. Nucleic Acids Res. 2008, 36, 2379–2374. (12) Perez, A.; Luque, F. J.; Orozco, M. J. Am. Chem. Soc. 2007, 129, 14739–14745. (13) Lavery, R.; et al. Nucleic Acids Res. 2010, 38, 299–313. (14) Manning, G. S. Q. ReV. Biophys. 1978, 11, 7796–7807. (15) Skolnick, J.; Fixman, M. Macromolecules 1977, 10, 944–948. (16) Odijk, T. Macromolecules 1979, 12, 688–693. (17) Baumann, C. G.; Smith, S. B.; Bloomfield, V. A.; Bustamante, C. Proc. Natl. Acad. Sci. U.S.A. 1997, 94, 6185–6190. (18) Hud, N. V.; Plavec, J. Biopolymers 2003, 69, 144–159. (19) Va´rnai, P.; Zakrzewska, K. Nucleic Acids Res. 2004, 32, 4269– 4280. (20) McConnell, K. J.; Beveridge, D. L. J. Mol. Biol. 2000, 304, 803– 820. (21) Goni, J. R.; Pe´rez, A.; Torrens, D.; Orozco, M. Gen. Biol. 2007, 8, 263–273. (22) Swigon, D.; Coleman, B. D.; Olson, W. K. Proc. Natl. Acad. Sci. U.S.A. 2006, 103, 9879–9884. (23) Marko, J. F.; Siggia, E. D. Macromolecules 1994, 27, 981–988. (24) Gore, J.; Bryant, Z.; No¨llmann, M.; Le, M. U.; Cozzarelli, N. R.; Bustamante, C. Nature 2003, 442, 836–840. (25) Lionnet, T.; Joubaud, S.; Lavery, R.; Bensimon, D.; Croquette, V. Phys. ReV. Lett. 2006, 96, 178102–178106. (26) Lankas, F.; Sponer, J.; Langowski, J.; Cheatham, T. E., III. Biophys. J. 2003, 85, 2872–2883. (27) Pe´rez, A.; Noy, A.; Lankas, F.; Luque, F. J.; Orozco, M. Nucleic Acids Res. 2004, 32, 6144–6151. (28) Drew, H. R.; Wing, R. M.; Takano, T.; Broka, C.; Tanaka, S.; Itakura, K.; Dickerson, R. E. Proc. Natl. Acad. Sci. U.S.A. 1981, 78, 2179– 2183. (29) Arnott, S.; Chandrasekaran, R.; Hall, I. H.; Puigjaner, L. C. Nucleic Acids Res. 1983, 11, 4141–4155. (30) Case, D. A.; et al. AMBER8; University of California: San Francisco, 2004. (31) Cornell, W. D.; Cieplak, P.; Bayly, C. I.; Gould, I. R.; Kenneth, M.; Merz, J.; Ferguson, D. M.; Spellmeyer, D. C.; Fox, T.; Caldwell, J. W.; Kollman, P. A. J. Am. Chem. Soc. 1995, 17, 5179–5197. (32) Cheatham, T. E., III; Cieplak, P.; Kollman, P. A. J. Biomol. Struct. Dyn. 1999, 16, 845–862. (33) Pe´rez, A.; Marcha´n, I.; Svozil, D.; Cheatham, T. E., III; Laughton, C. A.; Orozco, M. Biophys. J. 2007, 92, 3817–3829. (34) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. J. Phys. Chem. 1987, 91, 6269–6271. (35) Smith, D. E.; Dang, L. X. J. Chem. Phys. 1994, 100, 3757–3766. (36) Noy, A.; Soteras, I.; Luque, F. J.; Orozco, M. Phys. Chem. Chem. Phys. 2009, 11, 10596–10607. (37) Shields, G. C.; Laughton, C. A.; Orozco, M. J. Am. Chem. Soc. 1999, 119, 7463–7469. (38) Darden, T.; York, D.; Pedersen, L. J. Chem. Phys. 1993, 98, 10089– 10092. (39) Ryckaert, J.-P.; Ciccotti, G.; Berendsen, H. J. C. J. Comput. Phys. 1977, 23, 327–341. (40) Lu, X.-J.; Olson, W. Nucleic Acids Res. 2003, 31, 5108–5121. (41) Va´rnai, P.; Djuranovic, D.; Lavery, R.; Hartmann, B. Nucleic Acids Res. 2002, 30, 5398–5406.

DNA Elasticity Cross-Terms at Base-Pair Level (42) MacDonald, D.; Herbert, K.; Zhang, X.; Polgruto, T. J. Mol. Biol. 2001, 306, 1081–1098. (43) Nelson, H. C. M.; Finch, J. T.; Luisi, B. F.; Klug, A. Nature 1987, 330, 221–226. (44) DiGabriele, A. D.; Steitz, T. A. J. Mol. Biol. 1993, 231, 1024– 1039. (45) Curuksu, J.; Zacharias, M.; Lavery, R.; Zakrewska, K. Nucleic Acids Res. 2009, 37, 3766–3773. (46) Lankas, F.; Spackova, N.; Moakher, M.; Enkhbayar, P.; Sponer, J. Nucleic Acids Res. 2010. In press. (47) Strick, T. R.; Allemand, J. F.; Bensimon, D.; Bensimon, A.; Croquette, V. Science 1996, 271, 1835–1837. (48) Moroz, J. D.; Nelson, P. Proc. Natl. Acad. Sci. U.S.A. 1997, 94, 14418–14422. (49) Smith, S. B.; Finzi, L.; Bustamante, C. Science 1992, 258, 1122– 1126. (50) Bednar, J.; Furrer, P.; Katritch, V.; Stasiak, A. Z.; Dubochet, J.; Stasiak, A. J. Mol. Biol. 1995, 254, 579–594. (51) Wiggins, P. A.; Heijden, T. V. D.; Moreno-Herrero, F.; Spakowitz, A.; Phillips, R.; Widom, J.; Dekker, C.; Nelson, P. C. Nat. Nanotechnol. 2006, 1, 137–141.

J. Phys. Chem. B, Vol. 114, No. 23, 2010 8031 (52) Yuan, C.; Chen, H.; Lou, X. W.; Archer, L. A. Phys. ReV. Lett. 2008, 100, 1018102–1018106. (53) Cloutier, T. E.; Widom, J. Mol. Cell 2004, 14, 355–362. (54) Du, Q.; Smith, C.; Shiffeldrim, N.; Vologodskaia, M.; Vologodskii, A. Proc. Natl. Acad. Sci. U.S.A. 2005, 102, 5397–5402. (55) Bustamante, C.; Smith, S. B.; Liphardt, J.; Smith, D. Curr. Opin. Struct. Biol. 2000, 10, 279–285. (56) Hobza, P.; Sponer, J. Chem. ReV. 1999, 11, 3247–3276. (57) Lavery, R.; Moakher, M.; Maddocks, J. H.; Petkeviciute, D.; Zakrzewska, K. Nucleic Acids Res. 2009, 17, 5917–5929. (58) Noy, A.; Pe´rez, A.; Laughton, C. A.; Orozco, M. Nucleic Acids Res. 2007, 35, 3330–3338. (59) Yanagi, K.; Prive, G. G.; Dickerson, R. E. J. Mol. Biol. 1991, 217, 201–214. (60) Ulyanov, N. B.; James, T. L. Appl. Magn. Reson. 1994, 7, 21–42. (61) Sponer, J.; Kypr, J. J. Biomol. Struct. Dyn. 1990, 7, 1211–1220. (62) Sponer, J.; Kypr, J. J. Biomol. Struct. Dyn. 1993, 11, 27–41. (63) Sherer, E. C.; Harris, S. A.; Soliva, R.; Orozco, M.; Laughton, C. A. J. Am. Chem. Soc. 1999, 121, 5981–5991.

JP104133J