We used all-atom molecular dynamics simulation to investigate the elastic properties of double-stranded DNA (dsDNA). We focused on the influences of temperature on the stretch, bend, and twist elasticities, as well as the twist–stretch coupling, of the dsDNA over a wide range of temperature. The results showed that the bending and twist persistence lengths, together with the stretch and twist moduli, decrease linearly with temperature. However, the twist–stretch coupling behaves in a positive correction and enhances as the temperature increases. The potential mechanisms of how temperature affects dsDNA elasticity and coupling were investigated by using the trajectories from atomistic simulation, in which thermal fluctuations in structural parameters were analyzed in detail. We analyzed the simulation results by comparing them with previous simulation and experimental data, which are in good agreement. The prediction about the temperature dependence of dsDNA elastic properties provides a deeper understanding of DNA elasticities in biological environments and potentially helps in the further development of DNA nanotechnology.

## I. INTRODUCTION

DNA carries important genetic information about the organism and has an indispensable role in biological functions. In general, DNA exhibits various mechanical elasticities in biological processes, which depend on the corresponding biological environments.^{1–7} DNA can also be used as biosensors, nanoparticles, and molecular handles owing to the elasticity, and this utilization has contributed considerably to the development of nanotechnology.^{8–11} Expectedly, the effect of biological environments on the elastic properties of DNA has attracted growing interest in recent years.^{12}

Many studies have contributed to the elasticities of double-stranded DNA (dsDNA) in various biological environments, including the salt concentrations,^{13–20} ionic environments,^{4,18,21–28} and pH values.^{27,29–32} For example, Garai *et al.* used all-atom molecular dynamics (MD) simulations to study the elastic properties of short dsDNA at different salt concentrations.^{17} The authors observed that the bending persistence length gradually decreases while the stretch modulus increases with the increase in salt concentration. The simulation evidence showed that the electrostatic repulsion on the dsDNA backbone is enhanced by greater salt concentration, which makes the dsDNA easier to bend but less likely to stretch within the moderate slat concentration ranges.^{15} The ionic environment also affects the DNA elastic properties obviously due to the polyanionic nature of DNA. For instance, Guilbaud *et al.* observed the variation in the bending persistence length of dsDNA in the presence of the monovalent and divalent metal ions, where the bending persistence length decreases with the increase in ion concentration due to the electrostatic effects.^{25} Recently, the magnetic tweezer (MT) experiment and all-atom MD simulation reported that the high-valent cations decrease the bending persistence length and stretch modulus of dsDNA, which is in contrast to those of dsRNA.^{26} For pH dependence of dsDNA elasticities, for example, Zhang *et al.* analyzed the bending persistence length of dsDNA with several pH values at a MT experiment by applying a force to the dsDNA and found that the bending persistence length increases with the buffer pH; therefore, dsDNA stability and rigidity also increase.^{31} Another experiment exhibited that the proton-driven dsDNA spring is extremely sensitive to high pH, and it suggested an important implications for the implementation of various biological functions and the recognition of nucleic acid molecules.^{33} We note that the above MD simulations and MT experiments were performed at room temperature, which ignore the effects of high and low temperatures on dsDNA elasticities.

Temperature is another crucial factor that affects dsDNA elasticities and extensively exists in biological processes and technical applications. A considerable number of experimental and simulation works have devoted to the temperature dependence of dsDNA elasticities in the last two decades.^{20,34–43} For example, an earlier MT experiment suggested that the bending persistence length decreases rapidly with increasing temperature where long dsDNA chains were used.^{35} Geggier *et al.* observed that the bending persistence length of dsDNA with 2686 base pairs (bp) becomes smaller and more flexible with the increase in temperature by measuring the j-factors and equilibrium distribution of topoisomers.^{37} Driessen *et al.* carried out a similar experiment and observed similar results that the bending persistence lengths decrease obviously with increasing temperature.^{42} Brunet *et al.* performed an experiment on 2000 bp dsDNA and observed that the bending persistence length only decreases within 4 nm when the temperatures increase from 23 to 52 K ^{◦}C, which suggested that the temperature dependence of bending persistence length is weaker.^{40} Moreover, MD simulation results provide a deeper understanding of the experimental results about the dsDNA elasticities, indicating that the entropy plays a remarkable role in dsDNA elasticity inside the base pair or even between successive base pairs.^{38} Monte Carlo (MC) and MD simulations based on various models indicated that the bending persistence length decreases with increasing temperatures.^{34,41,43} Kriegel *et al.* combined the MT measures with all-atom MD to investigate the temperature dependence of dsDNA helical twist angles, and good agreements were obtained.^{39} Recently, Zhang *et al.* performed a MT and all-atom MD combination study on the twist–diameter coupling for dsDNA, and the results showed that the twist–diameter coupling is a common driving force for temperature-induced dsDNA twist changes.^{20} Although a series of studies have focused on the temperature dependence of dsDNA bend elasticity, the temperature dependence of dsDNA twist elasticities, such as the twist persistence length and twist modulus, especially of the twist–stretch coupling, is still needed to explore.

Of course, there are other factors that affect the elasticity of dsDNA in biological environments. These factors include the dsDNA composition, such as its base composition and sequence order,^{2,10,44–48} and the dsDNA chain lengths.^{49–55} In the current work, we used a 35 bp B-type dsDNA as an example and adjusted the system temperature value over a wide range to examine the effects of temperatures on dsDNA elasticity. We use all-atom MD simulations to investigate the effects of temperature on dsDNA elasticity by analyzing elastic parameters, including the bending and twist persistence lengths, stretch and twist modulus, and the twist–stretch coupling. In Sec. II, the all-atom MD method and calculations of dsDNA elasticities are described. The results on the stretch, bend, and twist elasticities and twist–stretch coupling of dsDNA at various temperature values are then discussed in Sec. III. A summary is presented in Sec. IV.

## II. METHODS

### A. All-atom MD simulation

We selected an initial B-type dsDNA with 35 bps (5′- CGACT CTACG GAAGG GCATC TCTCG GACTA CGCGC -3′), which has been reported in a previous experiment,^{56} as the example to examine its elasticity. The melting temperature of this sequence is 79.9 ± 0.6 ^{◦}C. The model dsDNA sequence is shown in Fig. 1(a). The UCSF Chimera 1.15 software application was used to build the DNA structure file.^{57} We attempted to adjust the simulation box size by many times and then chose a box with a size of 10 × 10 × 15 nm^{3} in all simulations, where the periodic boundary conditions were used in three directions to keep the DNA as free from boundary effects. The previous experiments and simulations indicated that the electrostatic interaction between two nucleic acids decreases obviously as the distances between the two nucleic acids increases.^{58–60} These findings suggest that the simulation box size used in the current simulation is large enough to treat the nucleic acid as an isolated one even under the periodic boundary conditions. Then, the 100 mM NaCl^{61} and TIP3P water molecule models^{62} were added at random initial positions in the box to simulate the physiological environment. Given that the nucleic acid carries electronegativity, we added 68 Na^{+} with the ion model by Joung and Cheatham^{63} to ensure that the system as a whole is electrically neutral.

We pre-equilibrated the systems to ensure that the MD simulation was carried out in an isothermal and isobaric environment. The initial DNA structure is bound to change when added into the solution. Thus, we started by energy-minimizing the system so that the atoms are not too close together and to ensure the stability of the DNA structure. Then, a restrictive force needs to be applied to the DNA prior to NVT simulation so that the random movements of atoms in the DNA chain are limited. The system was then warmed up to a desired temperature *T* by the V-rescale scheme, where the temperatures are coupled by a velocity rescaling method with a random term.^{64,65} The application of a limiting force prevented damage to the DNA structure during the warming process. After the system is warmed up to temperature *T*, the system pressure was adjusted to 1 atm by NPT simulation. When these pre-equilibrium processes were finished, the restrictive force was removed, and then, the MD simulations were performed for 500 ns without the restrictive force. The processes are the standard simulation procedure using Gromacs 4.6 software with the Amber bsc1 force field.^{66–70} In the current simulations, we adjusted the temperature *T* from 280 to 320 K, with a step of 10 K. Our MD simulations with various temperatures are similar to those in the previous all-atom MD simulations, where the Amber force fields were used.^{20,38,39,43} We calculated the root mean square deviation (RMSD) for the systems to examine the validity of simulations. Two examples are shown in Fig. 1(b), where the temperature values were selected to be 290 and 310 K. The RMSD represents the degree of molecular structure change of DNA, which follows:

where *δ*_{i} is the displacement between the *i*-th atom at moment *t* and moment 0 with a step of time Δ*t* = 2 ps. We took the average value per 2 ns and plotted the black line in the middle, as shown in Fig. 1(b). The results show that the black lines reach equilibrium after 100 ns. Thus, the trajectory data from 100 to 500 ns were used in the subsequent statistical analyses in all simulations. In all data analyses, 3 bps were removed from each of the first and last strands of dsDNA during the statistical processes,^{53} and only the central 29 bps were selected to analyze the global elastic properties of the dsDNA.

### B. Elastic parameters

In the current simulation, we concentrated on the stretch, bend, and twist elasticities of dsDNA, as well as their coupling, which involve the elastic parameters of stretch and twist moduli and bending and twist persistence lengths. According to the elasticity theory, the stretch and twist moduli can be extracted from the twist–stretch elastic matrix **K**, an inverse of 2 × 2 covariance matrix.^{71} The twist–stretch elastic matrix **K** is MD-associated, which can be determined by^{12,51,61,71,72}

where *k*_{B} is the Boltzmann constant, *T* is the temperature, and *L*_{0} is an average contour length *L*. Here, Δ*L* represents the variance of the contour length *L*, and ΔΦ is the variance of the cumulative H-twist angle Φ as described in the previous work.^{61}^{,}$\Delta L2$ and $\Delta \Phi 2$ denote the mean square errors for *L* and Φ that originated from the thermal fluctuations with Gaussian distribution functions.^{54,71,73} Diagonal elements *K*_{SS} and *K*_{TT} correspond to the stretch and twist moduli, respectively, whereas elements *K*_{ST} is related to the twist–stretch coupling coefficient; *K*_{SS}, *K*_{TT}, and *K*_{ST} are also denoted as *S*, *C*, and *G* in the experiments, respectively.^{73–75} We noted that the stretch modulus, twist modulus, and the twist–stretch coupling are associated with the global contour length and twist angle.

The persistence length provides an alternative method to describe the dsDNA elasticity, as it is also experimentally accessible. According to the elasticity theory, the elastic modulus and persistence length have a simple relationship. For example, the twist elasticity can also be characterized by the twist persistence length *l*_{T} related to the twist elastic energy,^{76} and the twist persistence length *l*_{T} can be simply written as^{73,74}

where the *K*_{TT} is the twist modulus from Eq. (2).

The bending persistence length *l*_{B} was extensively used to characterize the bend elasticity of dsDNA. Based on the wormlike chain (WLC) approximation, the bending energy is quadratic to the bending angle, and then, the relationship between the bending angle and bending persistence length can be deduced.^{22} Considering the normalization constant, the bending persistence length *l*_{B} of dsDNA are determined by^{22,26,47,50,77,78}

Here, *p*(*θ*, *l*) is a probability distribution, where *θ* is the bending angle formed by a dsDNA spanning 10 bps^{44} and *l* is the constant representing the average 10 bp profile length. Equation (4) shows that −ln(*p*(*θ*, *l*)/sin *θ*) varies as a quadratic function of the curve, which can be used to fit MD simulation data to obtain *l*_{B}. We notice that Eq. (4) is accurate in the case of long chains where the contour length is much larger than the bending persistence length. However, the all-atom MD simulations showed that the short dsDNAs with excluding about 6 bps at each helix end can also be described by the WLC model with a bending persistence length of about 50 nm, which have similar flexibilities with those long dsDNA chains.^{53} In the current simulations, we analyzed the trajectory data by excluding 6 bps from the helix ends to estimate the bending persistence length and other elastic parameters, following the previous all-atom MD simulations used short dsDNA chains.^{26,28,47}

We notice that different models have been developed to describe the elastic parameters of dsDNA at a short length scale in recent years. Along with the path of WLC, the twist WLC model^{79–83} was recently extended to the non-local twist WLC model by considering the couplings between distal sites. The non-local twist WLC model is able to account for the length-scale elasticities of nucleic acids, and the results showed that the bending and twist behaviors are soft at short length scales while stiff at long length scales.^{83,84} Alternatively, the DNA molecules can be described by the mesoscopic model to investigate three-dimensional structures at the base pair level.^{54,55,85–89} In the mesoscopic model, the Hamiltonian contains both the twist and bending angles between the adjacent bps, which are able to handle the short length-scale dsDNA chains and obtain the reasonable values of bending persistence lengths for short dsDNA chains.^{54,55}

### C. Structural parameters

The elastic properties of dsDNA can be analyzed by the structural parameters because the accumulation of microstructural excursions leads to variations in the macroscopic elasticity of dsDNA. The structural parameters include a series of translation and rotation parameters to describe the chain conformations.^{12,90–93} Here, we only describe the corresponding structural parameters used in the current simulations. An intrinsic coordinate frame was constructed in the base pair where the xy plane is attached to the base pair and the z axis is along its normal direction to describe the translation and rotation parameters. Here, the z axis points to the 5′ end of the nucleic acid chain, whereas the x axis points to the major groove.^{91} The essential translation parameters are shift (*D*_{x}), slide (*D*_{y}), and rise (*D*_{z}), which denote the displacements between two adjacent base pairs along the x, y, and z directions, respectively.^{90} The translation parameters can be cumulative, exhibiting more obvious effects on the length scale, and presented in the grooved model to conveniently analyze the DNA conformations under various biological environments.^{23,26,28,94–98} Three rotation parameters, namely, roll *ρ*, tilt *τ*, and twist *ω*, were defined to describe the rotation properties of dsDNA belonging to the angles constructed between two adjacent base pairs along the x, y, and z directions, respectively.^{90–92} We focus on the helical rise *h*_{z} and helical twist *γ*, which are projection parameters on the helical axis of the nucleic acid. Helical rise *h*_{z} and helical twist *γ* are the variables that twist and stretch the dsDNA, which can be manipulated and measured in single-molecule experiments.^{93,99} During the analysis, these structural parameters of the dsDNA were obtained using the Curves+ program.^{97}

## III. RESULTS AND DISCUSSION

In the current work, we set the temperature *T* as 280, 290, 300, 310, and 320 K to investigate the temperature dependence of the elastic properties for dsDNA. In Subsection III A, we concentrated on temperature dependence of the stretch elasticity, which is described by the stretch modulus *K*_{SS} (Figs. 2 and 3). The temperature dependence of bend elasticity was characterized by the bending persistence length *l*_{B} (Figs. 4 and 5) in Subsection III B. The temperature dependence of twist elasticity was characterized by twist persistence length *l*_{T} and twist modulus *K*_{TT} (Figs. 6 and 7) in Subsection III C. Finally, we investigated the temperature dependence of twist–stretch coupling (Figs. 8 and 9) in Subsection III D.

### A. Temperature dependence of stretch elasticity

We considered the temperature dependence of stretch modulus *K*_{SS}, as shown in Fig. 2. Two typical probability distributions of contour lengths *p*(*L*) are shown in Fig. 2(a), indicating that the dsDNA has a distinct *K*_{SS} at various *T*. Our observations about *p*(*L*) agree with the fact that contour length *L* fluctuates with the Gaussian distribution in the WLC model.^{50,78} By inversing the covariance matrix according to Eq. (2), we can arrive at

where *ρ*_{LΦ} is the Pearson coefficient between contour length *L* and cumulative H-twist angle Φ, and *σ*_{L} is the standard error of contour length *L*. According to Eq. (5), we obtained *K*_{SS}, as shown in Fig. 2(b), where the available data were also inserted for comparison. In particular, we obtained *K*_{SS} = 1398.5 ± 25.2 pN at *T* = 300 K, which enabled us to compare with the available data. Bao *et al.* observed that *K*_{SS} = 1441 pN at *T* = 298 K for 40 bp dsDNA with 1M NaCl,^{61} whereas Marin-Gonzalez observed *K*_{SS} = 1280 ± 70 pN at *T* = 300 K for 16 bp dsDNA with 150 mM NaCl.^{100} Our observations are in good agreement with the results from these two all-atom MD simulations. The current results about *K*_{SS} also agree with the available data from an optical tweezer experiment, where the stretch modulus have about 1256 ± 217 pN at room temperature and various salt concentrations.^{13} However, the available experimental data^{99,101,102} showed that the dsDNA has a range about 1000–1760 pN, which is also inserted in Fig. 2(b) (more data are provided in Table S1 of supplementary material). The discrepancy between the current results and previous results are due to the various dsDNA chain length and salt solution concentrations used in the various systems. Then, we predicted the temperature dependence of *K*_{SS} as shown in Fig. 2(b). The stretch modulus *K*_{SS} decreased from 1537.0 ± 25.4 pN to 1169.1 ± 23.8 pN as *T* increased from 280 to 320 K. The linear fitting showed a linear relationship between *K*_{SS} and *T* with a slope of −8.89 pN/K. Here, our data clearly demonstrated that the stretch modulus *K*_{SS} strongly depends on temperature *T* with a linear relationship when the other parameters were constant.

We tried to explain the temperature dependence of *K*_{SS} by analyzing the structural parameters as shown in Fig. 3. According to Eq. (5), the temperature dependence of the stretch modulus is related to several parts: the Pearson correlation coefficient *ρ*_{LΦ}, the average contour length *L*_{0} and its standard error *σ*_{L}, and the temperature *T*. The Pearson correlation coefficient *ρ*_{LΦ} was ∼0.11–0.18 (see Pearson coefficient in Table S2 of supplementary material) as *T* varies from 280 to 320 K; therefore, the factor of $1\u2212\rho L\Phi 2$ can be ignored. Then, we considered the average contour length *L*_{0} at various *T*, as shown in Fig. 3(a). The average contour length *L*_{0} has a very slight upward trend in the slope of *k*_{l} = 0.002 nm/K. Indeed, the average contour length increased by Δ*L*_{0} = 0.083 nm when Δ*T* = 40 K. However, within the errors, the average contour length can be regarded as a constant value of *L*_{0} = 9.364 nm within the tested temperature range. A recent MD simulation predicted that the temperature dependence of dsDNA length is very weak because the length increase with temperature is absorbed by expanding the helix radius.^{43} This simulation confirms our observation on the average contour length. The constant *L*_{0} suggests that contour lengths are not attributed to the temperature-induced variations in the stretch modulus *K*_{SS}. We plotted the $1/\Delta hz2$ as a function of *T* as shown in Fig. 3(b), where $\Delta hz2$ denotes the mean square error for Helical rise *h*_{z} and $\sigma L2=n2\Delta hz2$. Here, *n* = 29 is the cumulative times. However, $1/\Delta hz2$ decreases with *T* in a linear relationship with a slope of *k*_{h} = −5.37 Å^{−2}/K. This result suggests that the temperature dependence of stretch modulus *K*_{SS} originates from the thermal fluctuation for the contour length.

### B. Temperature dependence of bend elasticity

We considered the temperature dependence of bending persistence length *l*_{B} for dsDNA. We fitted the simulation data by Eq. (4) to obtain the bending persistence length *l*_{B} as shown in Fig. 4. Two typical examples are shown in Fig. 4(a), where the −ln(*p*(*θ*, *l*)/sin *θ*) as functions of bending angle *θ* are quadratic curves at *T* = 290 and 310 K, respectively, and the bending angle *θ* are formed by ten consecutive base pairs in 29 base segments of the dsDNA center. We obtained the bending persistence length *l*_{B} = 57.40 ± 0.16 nm at *T* = 290 K and *l*_{B} = 52.28 ± 0.14 nm at *T* = 310 K. Here, our results about the bending persistence lengths agree with those obtained in the previous all-atom MD simulations where the short dsDNAs were used.^{22,26} In the experiments, the dsDNAs have the bending persistence lengths of about 45–50 nm, in which the long chains are usually employed.^{103–105} Actually, the dsDNA cyclization experiments have indicated that the short dsDNAs have high bending flexibilities,^{106,107} which are different from the WLC predictions.^{36,88} The high flexibilities of short dsDNA chains were also observed in the small-angle x-ray scattering experiments^{56,108} and the atomic force microscopy experiment.^{49} The deviation between the WLC prediction and experiments on short chains may result from the model used and the different conditions in the experiments and simulations. The fitting results also hint that the dsDNA chain becomes more flexible in the solution with higher *T*. We investigated the effects of *T* on bending persistence length *l*_{B} over a range of 280–320 K to observe this dependence more carefully. In particular, we plotted the bending persistence length *l*_{B} as a function of *T* with a step of 10 K as shown in Fig. 4(b). Then, we observed a linear relationship between the bending persistence length *l*_{B} and temperature *T* with the slope of −0.29 nm/K. The existing experimental data (see more data in Table S1 of supplementary material) are also inserted in Fig. 4(b) for convenient comparison, which confirmed our observations that the bending persistence length *l*_{B} decreases as *T* increases.^{37,40,42} For examples, Brunet *et al.* observed a linear dependence of the bending persistence length *l*_{B} on temperatures *T*, where Δ*l*_{B} ≈ 4 nm was obtained at Δ*T* = 29 K.^{40} A similar dependence was also observed in the experimental results from Geggier *et al.*, but the bending persistence lengths *l*_{B} are smaller at the same temperature.^{37} However, this dependence is much weaker than those observed in another experimental results from Driessen *et al.*^{42} Brunet *et al.* explained that the bias induced by detector time-averaging blurring were neglected in those observations from Driessen *et al.*, which led to the larger Δ*l*_{B} (≈22 nm) with Δ*T* = 29 K. Our observations on the temperature dependences of bending persistence lengths *l*_{B} confirmed those experimental results from Brunet *et al.*^{40} and Geggier *et al.*^{37}

We used the structural parameters involving the bend elasticity to understand the microscope mechanism about the temperature dependences of *l*_{B} shown in Fig. 4. Two typical examples for the probability distribution of the bending angle *θ* at *T* = 290 and 310 K are shown in Fig. 5(a) according to Eq. (4). Clearly, this distribution is not a Gaussian distribution because the bending angle is not an independent variable, which is related to the tilt and roll angles.^{80,109} The results show that the maximum in the corresponding curve gradually shifts to the left, that is the most probable distributions for the bending angles *θ*_{m} are equal to 12.0 and 13.5*π*/180 at *T* = 290 and 310 K, respectively, as shown in Fig. 5(a). The smaller bending angle reflects the stronger dsDNA rigidity, which can help us make initial estimations on the chain rigidity. We obtained the average bending angles ⟨*θ*⟩ over the range of *T* = 280–320 K as shown in Fig. 5(b). The results indicate that the average bending angle ⟨*θ*⟩ has a linear relationship with temperature *T* in the slope of *k*_{θ} = 0.042*π*/180 K. The weak upward trend illustrates that the dsDNA chain becomes more flexible as the temperature *T* increases because the average bending angle ⟨*θ*⟩ increases. The structural parameters showed that the average roll angle *ρ* increases, but the average tilt angle *τ* decreases as temperature *T* increases (see more data about the structural parameters in Table S2 of supplementary material), reflecting that the increase in local bending is directed toward the major groove.^{43} Then, we plotted the ratios between the thermal fluctuation for the tilt and roll angles $\Delta \rho 2/\Delta \tau 2$ as functions of temperature *T* as shown in Fig. 5(c). Here, $\Delta \rho 2$ and $\Delta \tau 2$ denote the mean square errors for the roll angle *ρ* and tilt angle *τ*, respectively. The ratio $\Delta \rho 2/\Delta \tau 2$ is about 3.8–4.0, which indicates the thermal fluctuations are greater in roll angles than those in tilt angle. This finding suggests that the bend elasticity is mainly from the roll angles for dsDNA.^{51} The $\Delta \rho 2/\Delta \tau 2$ slightly decreases with *T* in the slope of *k*_{ρτ} = −0.005/K, showing that the tilt angles have more contribution to the bending elasticity as *T* increases.

### C. Temperature dependence of twist elasticity

Two parameters, namely, the twist modulus *K*_{TT} and twist persistence length *l*_{T}, describe the twist elasticity of dsDNA, which can be measured in experiments. By calculating the matrix **K** in Eq. (2), we can arrive at

where *σ*_{Φ} is the standard error of the cumulative H-twist angle Φ. We investigated the twist modulus *K*_{TT}, according to Eq. (6) as shown in Fig. 6. We plotted the distributions of H-twist *γ* and showed two examples with *T* = 290 and 310 K in Fig. 6(a). Our simulations on *p*(*γ*) agree with the fact that the structural parameters fluctuated with Gaussian distributions in the WLC model.^{50,78} However, the difference between these two *p*(*γ*) values indicates the different twist moduli. We plotted the twist modulus *K*_{TT} with various temperature values according to Eq. (6), as shown in Fig. 6(b). We observed that the twist modulus decreases from *K*_{TT} = 455.2 ± 7.5 pN⋅nm^{2} at *T* = 280 K to *K*_{TT} = 387.0 ± 7.9 pN⋅nm^{2} at *T* = 320 K with a linear relationship of slope of −1.68 pN⋅nm^{2}/K. In particular, the present observation provides that *K*_{TT} = 402.2 ± 7.2 pN⋅nm^{2} at *T* = 300 K. Our observations are in agreement with the previous MT experiments of twist persistence lengths (i.e., *C* = 410 ± 30 pN⋅nm^{2} and *C* = 436 ± 17 pN⋅nm^{2}),^{110,111} where long-chain dsDNA were immersed in 100 mM NaCl solutions at room temperature (*T* = 296 ± 2 K). Recent coarse-grained MD studies reported that the dsDNA has *K*_{TT} = 386 ± 3 pN⋅nm^{2}^{112} and *K*_{TT} = 399 ± 1 pN⋅nm^{2}^{113} in 100 mM ion solutions at *T* = 300 K, respectively. These data were also inserted in Fig. 6(b), showing a good agreement with the current MD data. Previous MT experiments observed that the twist modulus are about 386–448 pN⋅nm^{2} at room temperatures,^{99,114} which also confirmed our MD data (see more data in Table S1 of supplementary material).

The twist persistence length *l*_{T} is another experimentally accessible parameter measured in MT experiments.^{80,110,115,116} Actually, a series of theoretical works have contributed to the twist persistence length based on the twist WLC model using rigid base pairs,^{80–84,109,117,118} as well as the simulations based on the all-atom and coarse-grained models.^{51,117} Here, we reported the twist persistence length of dsDNA by all-atom MD depending on temperatures. We considered the temperature dependence of twist persistence length *l*_{T} as shown in Fig. 7. We can statistically obtain the twist persistence lengths of dsDNA at various *T* values according to Eqs. (3) and (6) as shown in Fig. 7(a), where the data about *T* = 300 K were also inserted from the previous studies for convenient comparison. In the current MD simulations, we obtained *l*_{T} = 95.6 ± 1.7 nm at *T* = 300 K and 100 mM NaCl. Recently, Skoruppa *et al.* carried out an all-atom MD simulation and reported that *l*_{T} = 86 nm for 32 bp dsDNA with *T* = 298 K and 150 mM NaCl, which quantitatively confirmed our MD results.^{84} Caraglio *et al.* used MC method to predict that dsDNA has *l*_{T} = 118 nm.^{82} The MT experiment observed that the twist persistence length is *l*_{T} = 109 ± 4 nm for long chain dsDNA.^{99} These available data suggested that the twist persistence lengths *l*_{T} is ∼100 nm (see more available data in Table S1 of supplementary material). Our observations on the twist persistence lengths *l*_{T} are in good agreement with these data from the MC, MD, and MT methods, as well as the more simulation and experimental results. Skoruppa *et al.* also predicted that the dsDNA has *l*_{T} = 105 nm at *T* = 295 K,^{117} which suggested that exploring the temperature dependence for the twist persistence length is desired. Our results predict that the twist persistence length has an obvious decreasing trend with increasing *T* values over the range of *T* = 280–320 K in a slope of −0.76 nm/K as shown in Fig. 7(a). The theoretical predictions and numerical calculations based on the twist WLC model suggested that the twist persistence lengths have a decreasing trend with temperature within the ranges of *l*_{T} = 40–100 nm.^{80,83,109} Our observations are about *l*_{T} = 85–116 nm, which are in approximate agreement with the previous studies about the temperature-dependent twist persistence lengths.^{80,109} The experiments and simulations showed that the elasticity has length-scale dependence for short dsDNA.^{51,108,119} This is a possible reason for the difference between the results from the current simulations and those from the previous experiments. For the long dsDNA chain, conformations, such as kinks, should be considered for the effects on persistence length, which leads to the possible modification of the theoretical prediction models and all-atom simulations.^{120} Here, the current simulations indicate that twist persistence length has an approximately linear decreasing with temperature over a wide range.

Similar to stretch elasticity, we only considered the contribution of the standard error *σ*_{Φ} to the twist elasticity that described by the twist modulus and twist persistence length according to Eqs. (3) and (6). We investigated temperature dependence of $1/\Delta \gamma 2$ as shown in Fig. 7(b), where $\Delta \gamma 2$ denotes the mean square error for H-twist *γ* and $\sigma \Phi 2=n2\Delta \gamma 2$ with the cumulative times *n* = 29. The results show that $1/\Delta \gamma 2$ decreases with *T* in a linear relationship of *k*_{γ} = −0.02 × 180^{2}/*π*^{2} K. Our simulations also show that the average twist and helical twist angles decrease very slightly with temperature (see more data about the structural parameters in Table S2 of supplementary material). However, the temperature dependence of twist modulus originates from the mean square deviations of twist angles, which is caused by the thermal fluctuation. Naturally, the higher temperature, the stronger thermal fluctuation, leading to the decrease in the twist modulus and twist persistence length. However, more obvious decrease occurred in the twist persistence length than those in twist modulus due to the factor of 1/*k*_{B}*T*.

### D. Temperature dependence of twist–stretch coupling

The coupling between elastic parameters plays an important role in the elastic behavior of dsDNA. In the full parameter space, a series of coupling combinations is present between the length and angle parameters. Here, we only concentrated on the stretch couplings with the twist angle as shown in Figs. 8 and 9. In the current work, we provide two calculational methods to demonstrate the twist–stretch coupling to conveniently compare with the experimental results.

We demonstrate the temperature dependence of twist–stretch coupling in Fig. 8, where the data were extracted from the all-atom MD. The H-rise and H-twist corrections are analyzed, and an example at *T* = 300 K is shown in Fig. 8(a), in which the data between the H-rise *h*_{z} and H-twist *γ* have positive corrections with a positive slope of 0.015 × 180/*π* Å (more examples are provided in Fig. S1 of supplementary material). The positive correction indicates that the twist angle increases when the dsDNA stretches. In order to demonstrate this phenomenon clearly, we used the parameter *dL*/*dN* to characterize twist–stretch coupling, where *dL* is the change in contour length and *dN* denotes the change of helical turn.^{28,111,121} We obtained the twist–stretch coupling *dL*/*dN* = 0.53 ± 0.01 nm/turn at *T* = 300 K. Our calculations are in good agreement with several experimental measurements where the twist–stretch couplings are 0.5 ± 0.1 nm/turn,^{111} 0.44 ± 0.1 nm/turn,^{99} and 0.42 ± 0.2 nm/turn.^{122} Our results also are in agreement with the all-atom MD results with *dL*/*dN* = 0.59 ± 0.02 nm/turn^{58} and *dL*/*dN* ∼ 0.61 nm/turn,^{61} where the 16 and 30 bp dsDNA were used in 1M NaCl solution, respectively. For convenience, we also listed these comparisons between the current all-atom MD simulations and the previous MT experiments or all-atom MD results in Fig. 8(b) (see more data in Table S1 of supplementary material). This difference is reasonable because the different DNA length and concentration of ion solution can cause the changes in twist–stretch couplings.^{28} However, all these MD simulation and MT experiments were performed at a constant temperature of *T* = 300 K. We calculated the twist–stretch coupling *dL*/*dN* at various temperatures *T* as shown in Fig. 8(c) to explore the temperature dependence of twist–stretch coupling. In particular, the twist–stretch coupling has *dL*/*dN* = 0.40 ± 0.01 nm/turn at *T* = 280 K and has *dL*/*dN* = 0.59 ± 0.01 nm/turn at *T* = 320 K. Expectedly, the twist–stretch coupling increases with increasing temperature *T*. The fitting result showed the twist–stretch coupling increases linearly with a slope of *k*_{ln} = 0.004 nm/turn K, indicating that the DNA overwinds obviously with increasing temperature in an approximately linear way.

The elasticity parameters can be calculated by the perturbed from the equilibrium related to the randomly thermal fluctuations in a polymer system.^{71} The twist–stretch coupling can be expressed by *K*_{ST} according to Eq. (2). By inversing the covariance matrix, we can arrive at

where *ρ*_{LΦ} is Pearson coefficient between counter length *L* and cumulative H-twist angle Φ, *σ*_{L}, and *σ*_{Φ} are their standard errors, respectively. According to Eq. (7), we then obtained the twist–stretch coupling *K*_{ST} by analyzing the standard errors of H-rise and H-twist parameters and their Pearson correlation coefficient at various temperatures as shown in Fig. 9(a). In particular, we obtained the twist–stretch coupling *K*_{ST} = −116.7 ± 1.8 pN ⋅ nm at *T* = 300 K. This result is in agreement with the MT experiment, where the twist–stretch coupling is *g* = −90 ± 20 pN ⋅ nm at *T* = 296 ± 2 K^{111} and in good agreement with the coarse-grained MD results of *K*_{ST} ∼ −120 pN ⋅ nm.^{112,113} Actually, we can express the relationship between two coupling coefficients according to^{111}

For example, we also obtained the twist–stretch coupling *dL*/*dN* = 0.53 ± 0.01 nm/turn at *T* = 300 K according to Eq. (8), which quantitatively agrees with the data in Fig. 8(c). According to Eq. (8), we plotted the twist–stretch coupling as a function of temperature in Fig. 9(b). Clearly, the theoretical results about the twist–stretch coupling also exhibit the similar trends with the temperature, which have been demonstrated in Fig. 8(c).

To demonstrate the origin of the increasing trend of twist–stretch coupling with temperature, we investigate the correlation between the random fluctuations about the stretch and twist according to Eq. (7), as shown in Figs. 9(c) and 9(d), where both *ρ*_{LΦ} and 1/*σ*_{L}*σ*_{Φ} are shown at various temperatures. The results show that the *ρ*_{LΦ} increases obviously with temperature *T* in a linear manner, namely, *ρ*_{LΦ} increased from 0.11 to 0.18, whereas 1/*σ*_{L}*σ*_{Φ} decreased obviously from 22.9 to 15.7, as *T* varies from 280 to 320 K. Since *ρ*_{LΦ} is between 0.11 and 0.18, the factor $1/(1\u2212\rho L\Phi 2)$ can be ignored, the variation in *K*_{ST} comes from three terms: *k*_{B}*TL*_{0}, 1/*σ*_{L}*σ*_{Φ}, and *ρ*_{LΦ}. The first term directly exhibits a linear relationship between *K*_{ST} and *T* because the *L*_{0} can be treated as a constant. The last two terms, 1/*σ*_{L}*σ*_{Φ} and *ρ*_{LΦ}, are linear functions of *T*; the former is a negative correction, but the latter is a positive one. This result suggests that the increase in twist–stretch coupling *K*_{ST} is attributed to the *ρ*_{LΦ}, but *σ*_{L}*σ*_{Φ} reduces the twist–stretch coupling as the *T* increases. From a viewpoint of statistical physics, the higher *T* leads to the stronger random thermal motions of molecules; thus, the thermal fluctuations of contour length *L* and cumulative H-twist angle Φ, i.e., *σ*_{L} and *σ*_{Φ}, increase with *T*. This naturally results in the decreasing of *K*_{ST}, according to Eq. (7). However, the Pearson’s correction coefficient *ρ*_{LΦ}, which presents the corrections between the random thermal fluctuations of *L* and Φ, increases with *T*. That is to say, the *ρ*_{LΦ} and *σ*_{L}*σ*_{Φ} play contrasting roles in the *T*-dependent *K*_{ST}.

## IV. SUMMARY

In this work, we investigated the effects of temperature on the elasticity of dsDNA using an all-atom MD simulation. We selected a short dsDNA with 35 bp where only the center 29 bp were used to analyze the elastic properties of the dsDNA. We utilized the bending persistence length *l*_{B}, twist persistence length *l*_{T}, twist modulus *K*_{TT}, and stretch modulus *K*_{SS} to characterize the dsDNA bending, twist, and stretch elasticities. We used two types of twist–stretch coupling parameters, *dL*/*dN* and *K*_{ST}, to describe the twist–stretch couplings for dsDNA. We concentrated on the temperature dependence on these elastic parameters and their couplings, compared them with the available data for dsDNA.

The results showed that the stretch modulus *K*_{SS} decreases in a linear manner, with a slope of −8.89 pN/K. Our observations on *K*_{SS} were compared with the existing MT and MD data at room temperature, and the results are in good agreement. The decrease in *K*_{SS} originates from the thermal fluctuation in contour length described by $\Delta hz2$, depending on the system temperature. The MD results showed that the bending persistence length *l*_{B} decreases with temperature with a slope of −0.29 nm/K, which are in good agreement with the available data from the MT experiments. This is because that higher temperature leads to a bigger average bending angle $\theta $ in a linear relationship and subsequently leads to a smaller bending persistence length. The structural parameters showed that the bending flexibility of dsDNA mainly comes from the thermal fluctuation on roll angle $\Delta \rho 2$, but its proportion decreases linearly as *T* increases. We used the twist persistence length *l*_{T} and the twist modulus *K*_{TT} to analyze the twist elasticity for dsDNA chains by comparing with the available data at room temperature, and the results indicated that the dsDNA becomes flexible in twist elasticity in a linear manner as *T* increases because of the stronger thermal fluctuation on $\Delta \gamma 2$.

We analyzed the twist–stretch coupling for dsDNA by two methods. One is the direct analysis of *dL*/*dN* from the corrections between the corresponding structural parameters and the other is to calculate the *K*_{ST} from modulus matrix **K**. Both methods were compared with the available data at room temperatures, and good agreements were obtained. Both these two methods predicted that the twist–stretch coupling become stronger as *T* increases. The method based on the modulus matrix **K** suggested that the correlations between thermal fluctuations of the contour length and twist angles contribute to the increase in twist–stretch coupling *K*_{ST} as the *T* increases, but the *σ*_{L}*σ*_{Φ} plays a contrasting role in twist–stretch coupling *K*_{ST}. Indeed, the *T*-dependent *K*_{ST}, or even other elastic parameters, are complex and the theoretical studies about the *T*-dependent *K*_{ST} elastic parameters are still needed. These findings provide a deeper understanding for the elastic nature of DNA and the regulation of DNA elasticity under various temperature environments and the development of DNA nanotechnology.

## SUPPLEMENTARY MATERIAL

The supplementary material provided herein includes the correlation between the data about H-rise *h*_{z} and H-twist *γ* at *T* = 280, 290, 310 and 320 K. The information provided in the supplementary material includes dsDNA elasticity parameters in the previous studies and the structural parameters and Pearson correlation coefficient *ρ*_{LΦ} for dsDNA at different temperatures in present work.

## ACKNOWLEDGMENTS

We thank for the financial supports from the Program of National Natural Science Foundation of China (Grant Nos. 21973070 and 22273067).

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

### Author Contributions

**Yahong Zhang**: Data curation (equal); Formal analysis (equal); Investigation (equal); Methodology (lead); Visualization (equal); Writing – original draft (equal). **Linli He**: Conceptualization (equal); Formal analysis (supporting); Funding acquisition (equal). **Shiben Li**: Data curation (equal); Formal analysis (equal); Funding acquisition (equal); Resources (equal); Supervision (lead); Writing – original draft (equal); Writing – review & editing (lead).

## DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding author upon reasonable request.

## REFERENCES

^{2+}ions: A Raman study

^{+}/Mg

^{2+}solution