All-Atom MD Simulation of DNA Condensation ... - ACS Publications

Subscriber access provided by UNIV OF NEWCASTLE

Article

All-Atom MD Simulation of DNA Condensation Using AbInitio Derived Force Field Parameters of Cobalt(III)-Hexammine Tiedong Sun, Alexander Mirzoev, Nikolay Korolev, Alexander P. Lyubartsev, and Lars Nordenskiöld J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.7b03793 • Publication Date (Web): 26 Jul 2017 Downloaded from http://pubs.acs.org on July 27, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry B is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

1

All-atom MD Simulation of DNA Condensation Using Ab-initio Derived Force Field Parameters of Cobalt(III)-Hexammine Tiedong Sun,† Alexander Mirzoev,† Nikolay Korolev,† Alexander P. Lyubartsev,*,‡ and Lars Nordenskiöld*,† †School of Biological Sciences, Nanyang Technological University, Singapore 637551 ‡Department of Materials and Environmental Chemistry, Stockholm University, 10691 Stockholm, Sweden

*E-mail: [email protected] ; [email protected]

ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 36

2

ABSTRACT: It is well established that the presence of the trivalent cobalt(III)-hexammine cation (CoHex3+) at submillimolar concentrations leads to bundling (condensation) of double-stranded DNA molecules, which is caused by DNA-DNA attraction induced by the multivalent counterions. However, the detailed mechanism of this process is still not fully understood. Furthermore, in allatom molecular dynamics (MD) simulations spontaneous aggregation of several DNA oligonucleotides in the presence of CoHex3+ have previously not been demonstrated. In order to obtain a rigorous description of CoHex3+-nucleic acids interactions and CoHex3+-induced DNA condensation to be used in MD simulations, we have derived optimized force field parameters of the CoHex3+ ion. They were obtained from Car-Parrinello molecular dynamics simulation of a single CoHex3+ ion in the presence of 125 water molecules. The new set of force field parameters reproduces the experimentally known transition of DNA from B- to A-form, and qualitatively describes changes of DNA and RNA persistence lengths. We then carried out a two microsecond long atomistic simulation of four DNA oligomers each consisting of 36 base pairs in the presence of CoHex3+. We demonstrate that in this system, DNA molecules display attractive interactions and aggregate into bundle-like structures. This behavior depends critically on the details of the CoHex3+ interaction with DNA. A control simulation with a similar setup, but in the presence of Mg2+ does not induce DNA-DNA attraction, which is also in agreement with experiment.


Page 3 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


3

1. INTRODUCTION The physical-chemical properties of solutions of DNA, which is a strongly charged polyelectrolyte, are determined largely by long range electrostatic forces. The surrounding ion environment determines the effective interactions between DNA molecules as well as between DNA and other important components of biological matter such as histone proteins, forming nucleosomes and chromatin. It is well established that in the presence of multivalent ions, the effective interaction between DNA molecules becomes attractive in some cases leading to aggregation of DNA into compact and ordered structures.1-9 Cobalt(III)-hexammine, [Co(NH3)6]3+, (abbreviated CoHex3+) is a widely used condensing agent in experimental in vitro studies of nucleic acids and other polyelectrolytes1-4,6-9 (and references cited in 6,9). Attention has been paid to formation of toroidal and rod-like DNA structures where long DNA double helices were packed as ordered parallel bundles. Molecular simulations, both on atomistic and coarsegrained (continuum solvent) levels can provide important complementary information to the experimental studies by elucidating the specific physical mechanisms of the observed phenomena and ensuring correct interpretation of the experimental data. Previously, short MD simulations have shown that two 10 bp DNA oligonucleotides exhibit attraction in the presence of CoHex3+ or the polyamines spemidine3+ and spermine4+.10 Additionally, large scale MD simulations of 64 DNA molecules constrained in parallel ordering by periodic boundary conditions have shown a DNA-DNA attraction in the presence of spermine4+.11 However, spontaneous bundling of several free DNA molecules during long (microsecond) MD simulation sampling times has not been demonstrated. The value of molecular simulation results depends crucially on the quality of the force field used. For atomistic simulations of DNA, there exist two families of force fields, AMBER12-14 and CHARMM.15-17 The ability of these force fields to reproduce DNA structure as well as interactions with other molecules, organic ligands, monovalent ions has been validated in numerous studies13,18 (and references cited therein). For interactions with multivalent ions the situation is less investigated, both due to difficulties in direct comparison of experimental and modeling data and due to strong interactions of multivalent ions with water and DNA causing slow sampling necessitating long equilibration times. In the case of CoHex3+ ions, several in silico studies have illustrated the capability of molecular dynamics simulation to probe CoHex3+ interaction with DNA.10,19-22 These studies were



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 36

4

based on the AMBER force field and they used the CoHex3+ parameters described originally in reference 19. However, these parameters once used in simulation of unrestrained DNA in the presence of CoHex3+ do not provide good agreement with WAXS structure factors describing the DNA structural changes upon addition of these ions as it was recently shown in ref.22 In the present study, we found that the CoHex3+ parameters by Cheatham and Kollman19 does not maintain a regular octahedral geometry of the complex while interacting with TIP3P water23 (see below). A later work by Zhang and Mu24 modified the previous model to study the interaction between CoHex3+ and proteins. We have tested this model in simulations of four DNA oligomers in the presence of CoHex3+ with the purpose to investigate the condensing effect of these ions, which experimentally1,5,9 are known to induce attraction between DNA molecules leading to condensation. Simulations within this CoHex3+ model24 and DNA described by the CHARMM36 force field showed only a weak attraction of DNA but did not demonstrate the ability of DNA oligomers to aggregate, contradicting the well-established experimental results showing high condensing potential of CoHex3+.2,5,9 See more information in the Results section on comparison among these models. Based on these observations we decided to develop a new set of interaction parameters for CoHex3+ that can be used in combination with either CHARMM or AMBER force fields. The new parameters were derived to reproduce CoHex3+ solution structure obtained from ab-initio molecular dynamics simulations of a CoHex3+ in water solution, and then used in a series of allatom MD simulations of DNA oligomers in presence of CoHex3+ with the ultimate goal of demonstrating the DNA bundling induced by CoHex3+.

2. METHODS 2.1. CoHex3+ model. As reported in several X-ray crystallography studies (see e.g.25), CoHex3+ has a well-defined octahedron structure as shown in Figure 1. We therefore describe CoHex3+ as a molecular ion where the cobalt atom is explicitly bonded to the nitrogen atoms of the six amminegroups, which makes it similar to the original CoHex3+ model used previosuly.19 The bonds are described by harmonic potentials in which the average bond lengths are taken from results of abinitio simulations. Furthermore, explicit harmonic N-Co-N angle potentials are applied to maintain an octahedral structure. The respective force constants for both bonds and angles are inherited from the original model.19 No torsion angles are included in the present model. All intramolecular


Page 5 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


5

non-bonded interactions are also omitted. The non-bonded interactions were parameterized from the ab-initio simulation which is described in detail below. 2.2. Car-Parrinello Molecular Dynamics simulation of cobalt(III)-hexammine in water. The reference system for the ab-initio simulations consist of one CoHex3+ and 125 water molecules placed in a 3.80 nm3 cubic periodic box (Figure 1). The system was initially equilibrated by running

standard

MD

simulation

using

AMBER-based26 CoHex3+ and TIP3P23 water models. The pressure and temperature of the system were set to 1 atm and 300 K, respectively and regulated by Berendsen weak coupling.27 The equilibrated volume and atom coordinates of the last snapshot were then employed as a starting configuration for subsequent ab-initio simulation. The ab-initio simulation was based on Car-Parrinello molecular dynamics (CPMD)28 Figure 1. CoHex3+ in solution: a system simulated by the CPMD method.

and performed using the CPMD software.29-31 The

spin-polarized

semi-local

BLYP-

functional32-33 and the plane-wave basis-set with 90 Ry cutoff were employed. The core potentials of N, C and H were described by the pseudopotentials of Troullier and Martins,34 while Goedecker pseudopotential35-37 was used to describe the core potential of cobalt. A time-step of 4 a.u. (0.1 fs) for the numerical integration of the equation of motion was employed. The temperature was kept constant at 300 K by velocity rescaling with a tolerance of 50 K. In total, 7.3 ps trajectory was produced, of which the last 5 ps was used to calculate radial distribution function between cobalt (Co) and oxygen atoms of surrounding water (Ow). In constructing the CoHex3+ model, we inherited most of the bonded parameters of Cheatham and Kollman’s CoHex3+model, with the exception of the Co-N bond length. In our CPMD simulation (see Results and Discussion section), the average bond length of Co-N was found to be 0.1950 nm, and this value was used as an equilibrium bond length in the bonded potential between Co and N atoms. From the resulting CPMD simulation, we obtained the new non-bonded CoHex3+ parameters by fitting to an all-atom MD simulation of the same system (see below). We note that the adopted Co-N bond length might



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 36

6

differ somewhat from its experimental value due to the accuracy of CPMD and the pseudopotentials chosen. The reason is that accurate non-bonded interactions are more important than accurate bond length to build a CoHex3+ model that is capable of simulating DNA condensation. Since the only target in our fitting process is the Co-Ow RDF from the CPMD simulation, using a Co-N bond length different from that of the CPMD simulation is not advisable. For instance, if we use a Co-N bond length shorter than 0.1950 nm, we need weaker non-bonded interactions between CoHex3+ and water to have a good fit of the Co-Ow RDF. Hence, we choose to use the bond length from the CPMD simulation, rather than using the experimentally derived geometry. 2.3. Molecular Dynamics simulation with new CoHex3+ interaction parameters. In order to validate the newly developed parameters for CoHex3+ and demonstrate DNA condensation induced by this cation in all-atom MD simulations, we have performed molecular dynamics studies that are described in more detail below. They were all run using the GROMACS software package38 and followed a similar protocol, in which two stages of equilibration were conducted before the production MD run. First, the system went through energy minimization with steepest decent algorithm with force tolerance of 1000.0 kJ·mol-1·nm-1. The energy minimized system, subsequently, was subject to 2 ns constant volume equilibration, during which the target temperature of 298 K was reached. Then pressure coupling was turned on to reach 1.013 bar pressure. Both thermostat and barostat were based on the Berendsen weak coupling algorithm,27 which allows to effectively reach the target temperature and pressure. Position restraints on heavy atoms of DNA and CoHex3+ were applied during these equilibration phases. Finally, production MD runs were performed with constant temperature and constant pressure. The Particle Mesh Ewald method39 (PME) was always used to treat long range electrostatic interactions with 1.0 nm real space cutoff. Velocity rescale40 and Parrinello-Rahman41-42 algorithms were used to regulate temperature and pressure in all production runs. In order to increase the integration time step to 2 fs we had constrained all bonds in both solute and water using P-LINCS43 and SETTLE44 algorithms respectively. The neighbor list is updated every 20 steps under Verlet cutoff scheme. Improved ion parameters by Yoo and Aksimentev45 are used throughout all simulations.


Page 7 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


7

One known deficiency of the adopted CHARMM3617 force field is that DNA can exhibit base-pair opening as shown by Pan et al.46 and Galinido-Murillo et al.47 This causes difficulties for simulating short nucleic acid helices, as the opening can propagate from the end of helices, resulting in their partial melting. Therefore, as in our previous work48, hydrogen bonds in the first and the last base pairs were constrained to avoid unwinding of DNA. 2.3.1 CoHex3+ induced DNA transition from B-form to A-form. As a first test we set up a reduced system consisting of a 10 bp DNA duplex with the sequence 5'-ACCCGCGGGT-3' and eight CoHex3+ to demonstrate that CoHex3+ induced the experimentally known DNA B- to A-form transition.49 Two simulations were conducted starting from A-form and B-form double helices respectively. These starting configurations were generated by "nucleic acid builder" (NAB)50 included in the AmberTools package. Both systems were set up in a periodic cubic box with volume of 6.26.26.2 nm3 filled with 8000 TIP3P51 water molecules. Six chloride ions were added to neutralize the solutes. Production runs lasted for 200 ns for each simulation. 2.3.2 Double helix DNA and RNA persistence length. To test the model performance in terms of modeling the influence of CoHex3+ on nucleic acids mechanical properties, we run MD simulation and calculated the persistence length of double helical DNA and RNA in the presence of CoHex3+. The starting configurations of dsDNA17 and dsRNA52 were also generated by NAB.50 Specifically,

the

double

helical

DNA,

with

the

sequence

5'-

CATCTGGGCTATAAAAGGGCGTCG-3', was modeled in ideal B-form, while the 25 bp homopolymeric poly(rA)poly(rU) RNA was modeled in ideal A-form. The simulation was performed with a truncated octahedron box. In simulations with CoHex3+, 16 CoHex3+ ions were added. The system was solvated by about 31,000 TIP3P water molecules. Besides neutralizing ions, extra Na+ and Cl ions were added to achieve 150 mM concentration. Both simulations were run for 100 ns. After finishing the simulations, the Curves+ software53 was used to generate axis representation of the nucleic acid in each trajectory frame. The persistence length Lp was calculated using eq 1 as described by Drozdetski et al.21

 






(1)


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 36

8

where  is the bending angle between two adjacent DNA segments, l is the average contour length of DNA segments, which are 10 bp in our case. The probability distribution P() is estimated by the histogram of  with 0.005 radian bin size. 2.3.3 DNA condensation in the presence of CoHex3+. Finally, CoHex3+ induced DNA condensation was simulated in a system consisting of four 36-bp DNA double helices with added CoHex3+ (140 ions) and a mixture of 95 Na+ and 140 K+ cations with added Cl anions to neutralize the system (salt concentration was 90 mM; 147,000 TIP3P water molecules). Double helical DNA 36-mers were modelled by CHARMM36 and AMBER Parmbsc0 in two independent simulations. They had the following sequences (sequences of the corresponding complementary strands are not shown): 5'-ATTAATGGAACGTAGCATATTCTTCAAGTTGTCACG-3' 5'-CAAAACCTGATGCACACTGTAACATGAGATCCCGCG-3' 5'-TCGGCTTATAGAGGGCCAGCTCGTATCGACGGACCG-3' 5'-GCTAGTACCCCACCAATTTAGGCGAAAGGAGTCTGC-3' The GC-content (50%) was chosen to correspond to naturally common content of DNA used in many in vitro experiments of e.g. DNA persistence length dependence on salt. The specific base pair sequence of the four DNA fragments were designed to cover all (except one) possible sequence combinations of the four-base-pair steps, which corresponds to two beads in our previously published CG model of DNA.54 This choice of DNA sequence makes it possible to use the collected all-atom MD trajectories to extract CG sequence dependent ion-DNA effective potentials for future investigation of such properties as sequence dependence of DNA helical parameters, hydration and ion binding or sequence specific DNA condensation behavior. The general features of the present result on CoHex3+ induced DNA condensation are not expected to depend on the exact choice of DNA sequences or on GC content as the DNA condensation is mainly governed by the long range electrostatic interactions between the trivalent cobalt ions and DNA phosphates. A periodic cubic box of 16.416.416.4 nm3 was used; a 2 s production MD run was performed. One simulation with each CoHex3+ ion substituted for Mg2+ was performed as a control system that experimentally is not expected to exhibit DNA condensation. The same DNA


Page 9 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


9

molecules and same number of water and Cl- ions were used. To keep charge neutrality, 70 more K+ and 70 more Na+ ions comparing to CoHex3+ simulation were added. Simulations for a corresponding system but using Zhang & Mu’s CoHex3+ parameters24 are described in the Supporting Information (see also below in the RESULTS section).

3. RESULTS AND DISCUSSION 3.1. Fitting of the Co-Ow radial distribution function and new parameters for CoHex3+. The The Co-Ow radial distribution function (RDF) obtained from the CPMD simulation is shown in Figure 2. We then optimized the values of the CoHex3+ non-bonded interaction parameters, which allows us to well reproduce the solution structure of CoHex3+, i.e. the CPMD cobalt-water RDF. The optimization was done in an iterative manner. First, the non-bonded parameter of CoHex3+ was modified. Then a 4 ns MD simulation of the same simulation box was run with the modified parameters under the same temperature and volume. The Co-Ow RDF from this all-atom simulation served as the guidance to improve the non-bonded parameters in the next iteration. Primary attention was paid to the first peak of the Co-Ow RDF curve since it is a signature of the direct interaction between CoHex3+ and water. Specifically, we aimed to fit the position of the first peak and the integral over the first peak, which reflect the average distance and coordination number of Co-Ow pairs. Due to the symmetric nature of our CoHex3+ model, the position and height of the RDF peaks are almost monotonically coupled to the Lennard-Jones (LJ) parameters of the three atom types composing CoHex3+. Hence, it was both straightforward and also effective to hand tune the parameters to achieve satisfactory result. In addition to the LJ-parameters, partial charges of CoHex3+ atoms were also tuned. We started with the charge values used in the Figure 2. Comparison of CoHex3+-Ow RDF calculated from CPMD and from the MD using 1) Cheatham & Kollman’s parameter19 (purple); 2) Zhang & Mu’s parameter24 (green) and 3) our "best fit" parameters (red).

CoHex3+ parameters by Cheatham and Kollman.19 Later in the fitting process they were altered gradually as it was necessary to



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 36

10

maintain the correct shape of the first peak of the RDF curve. The final values of CoHex3+ parameters are shown in Table 1 and Table 2. The reference and fitted CoHex3+-Ow RDFs are shown in Figure 2. Using the CoHex3+ parameters by Cheatham and Kollman19 in the all-atom MD simulation results in a RDF curve with a sharp peak at r = 0.13 nm. The small radius of cobalt and nitrogen atoms in this model allows one water molecule to be inserted among the amine groups of CoHex3+. This result is clearly artificial and lends further support for the necessity of improving these original parameters for CoHex3+. On the other hand, the parameters by Zhang and Mu24 give a too wide first peak in the Co-Ow RDF, which is resulting from a weaker interaction with water. Table 1: CoHex3+ bonded parameters in this study Bond

b0 (nm)

kb (kJ·mol-1·nm-2)

Co-N

0.1950

251040.0*

N-H

0.1037*

418400.0*

Angle

0 (degree)

k (kJ·mol-1·radian-2)

N-Co-N

90.0/180.0*

585.76*

Co-N-H

109.5*

292.88*

H-N-H

109.5*

292.88*

* same as CoHex3+ model by Cheatham and Kollman19. Table 2: CoHex3+ non-bonded interaction parameters.

 (nm)

 (kJ·mol-1)

Q (e)

Model

Co

N

H

This work

0.4122

0.3854

0.01000

Cheatham & Kollman19

0.05612

0.2047

0.06735

Zhang & Mu24

0.1412

0.3250

0.1069

This work

2.7434

0.6818

0.1925

Cheatham & Kollman19

0.04184

0.7113

0.06539

Zhang & Mu24

0.04184

0.7113

0.06569

This work

+1.368

0.439

+0.237

Cheatham&Kollman19

+1.062

0.784

+0.369

Zhang&Mu24

+1.062

0.784

+0.369


Page 11 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


11

3.2. CoHex3+ induced DNA transition from B-form to A-form. Under physiological conditions, canonical B-form DNA is the dominant conformation of double helical DNA. However it has been known for a long time (see e.g.55-56) that change of solvent conditions such as addition of ethanol or cationic species (for example polyamines) can cause the formation of other DNA structural forms, for instance A-DNA. One specific case is known that the addition of CoHex3+ induces the experimentally observed conformational transition from B-DNA to A-DNA.57 In order to test that our new parameters display the expected behavior of stabilizing A-DNA, we have followed the work of Cheatham et al.19 implementing a simulation with the newly developed CoHex3+ parameters and the CHARMM36 force field.17 Two simulations with different initial DNA conformation are performed with the protocol described in the Method section. The overall structure is compared to ideal B-DNA and A-DNA as shown in Figure 3 in terms of allatom root mead square deviation (RMSD). In the simulation started with A-DNA form (Figure 3A), the DNA remains in A-form for the whole 200 ns simulation. Although, at 38 ns and 115 ns, the conformation of DNA molecule temporarily deviated from ideal A-form, the deviation in RMSD was not larger than 0.4 nm, and it returned close to the A-form in both cases. For most of the time during the simulation, the DNA molecule exhibited a RMSD value relative to the A-form of 0.15 nm on the average. As for the second simulation, which was started from B-DNA, the transition from Bform to A-form had been observed within the first 20 ns. As shown in the inset of Figure 3B, the DNA Figure 3. All-atom root meant square deviation (RMSD) compared to ideal BDNA (black) and A-DNA (red). Two simulations which are started from A-DNA (top) and B-DNA (bottom) respectively are shown. The inset of the bottom figure shows the transition from B-DNA to A-DNA near the onset of the simulation started with BDNA.

undergoes rapid conformational change toward an intermediate conformation between B-form and Aform once the simulation had been started. This results are consistent with the AMBER-based simulation reported by Cheatham et al.19. After some conformational rearrangement the transition



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 36

12

from B-DNA to A-DNA was completed at 12 ns. DNA stayed stable in A-form for the rest of this simulation up to 200 ns. After the transition, the average RMSD to A-DNA was about 0.15 nm. To characterize A- and B- DNA conformations, quantities, such as the glycosylic torsion angle between sugar and base (χ), the base pair displacement from the helical axis (x-displacement) and the intra-strand P-P distance, are often used. However, any single criterion cannot represent a satisfactory discrimination between A- and B-DNA. The parameter Zp, defined as the displacement of the phosphorus atom from the xy-plane of the “middle frame” between neighboring base pairs,58 was proposed to function as a standalone parameter to distinguish A- and B-DNA. It was found that Zp > 0.15 nm for A-form DNA and Zp < 0.08 nm for the B-DNA, while intermediate structures give 0.08 < Zp < 0.15 nm.59 With the help of the software SCHNArP,60 we have plotted the Zp value averaged over all base pair steps as a function of simulation time (Figure 4). Simulations starting from A-DNA and from B-DNA show the value of Zp > 0.15 nm for the most of trajectories indicating that in the presence of CoHex3+ A-DNA is the dominant form. In the simulation started with A-DNA (red plot in Figure 4), the Zp value drops below 0.15 nm at ~38 ns, ~115 ns and ~183 ns, and then recovers back to Zp > 0.15 nm within 10 ns. These events correlate well with changes in RMSD values shown in the upper panel of Figure 3. In the simulation started from B-DNA (black curve in Figure 4), change of Zp values indicates that transition to A-form begins from ~12 ns and it is well correlated with the RMSD plot (Figure 3B). After the transition, DNA remains mostly in the A-form. Therefore, we can conclude that the A-form of DNA is stable in the presence of the newly developed CoHex3+ model, in agreement with experimental data. Table 3: Effect of CoHex3+ on DNA and RNA flexibility. Values of persistence length (Lp), and relative changes of persistence length (Lp/Lp)

and contour length (l/l), are given

Figure 4. Order parameter Zp, averaged over all base pair steps plotted as a function of simulation time. Values of Zp > 0.15 nm is a characteristics of an ADNA while for B-DNA Zp < 0.08 (green dashed line Zp = 0.15 nm).


Page 13 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


13

System

Lp/Lp

l/l (%)

70.7

-

-

50.6

28

14.4

80.0

-

-

116

+45

4.80

Lp (nm) (%)

DNA DNA

+

CoHex3+ RNA RNA

+

CoHex3+ 3.3. Persistence length of double helical DNA and RNA. Ionic conditions can significantly affect the mechanical properties of nucleic acids. It has been experimentally observed that the persistence length (Lp) of B-DNA is dependent on concentration of monovalent salt61-63 and can deviate strongly from the widely accepted the value Lp 50 nm observed at physiological salt conditions. Multivalent ions have an even more drastic influence on the DNA persistence length.62,64 On the other hand, the MD study by Drozdetski et al.21 reported an opposite effect of CoHex3+ on the persistence length of DNA and RNA. It turns out that due to structural differences between double stranded DNA and RNA, CoHex3+ binds inside the A-form RNA double helix whereas it cannot penetrate into the grooves of the B-form DNA. Thus, CoHex3+ makes the RNA double helix shorter and stiffer. Here, we performed calculations of the persistence length of double helical DNA11 and RNA54 using our modified model of the CoHex3+. Results of comparison of the flexibility of 25-bp-long DNA and RNA in the absence and presence of CoHex3+ are shown in Table 3 and Figure 5. The persistence length of DNA was reduced by 28%, while the RNA persistence length increased by 45%. This result is in qualitative agreement with experimental data on DNA62,64 and simulation results by Drozdetski et al.21 The value of DNA persistence length is in agreement with previous study with similar DNA length.65-66 However, in contrast to the cited computational study,21 the CoHex3+ binding to the double helices has caused a contour length shortening of both DNA and RNA.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 36

14

3.4. CoHex3+-induced DNA condensation. The main purpose of this work is to use the new model of CoHex3+ to demonstrate and investigate the ability of CoHex3+ to induce attraction between DNA

oligomers

and

condensation.

Three

independent simulations, namely CHARMM36 DNA with CoHex3+, CHARMM36 DNA with Mg2+ and AMBER Parmbsc0 DNA with CoHex3+, were conducted using similar MD protocols. DNA aggregation was observed in both simulations in the presence of CoHex3+. However, in the MD simulations using CHARMM36 force field for the DNA and the CoHex3+ parameters of Zhang and Mu24 DNA-DNA interaction did not show strong Figure 5. Harmonic fitting of bending energy used for estimation of persistence length (eq 1). DNA-related data is shown on the top plot, RNA-related data is shown on the bottom plot. Results without CoHex3+ are in black, data with CoHex3+ is in red.

attraction, which is in disagreement with the experimental data on CoHex3+-induced DNA condensation.2,5,9 See further discussion below and Figures S1 and S2 in the Supporting Information.

In the system with CHARMM36 force field for the DNA and new parameters for the CoHex3+, for the 100 ns of the simulation, DNA helices stay away from each other that is reflected by a low population of the short P-P distances in the RDF curve (Figure 6A, black curve; Figure 7A). As the simulation proceeds, CoHex3+-mediated DNA-DNA attraction is observed which is seen by increase of the RDF peaks corresponding to close DNA-DNA distances (Figure 6A, redblue-yellow). After about 500 ns (Figure 6A, blue), most DNA helices are in contact with each other forming a bundle. Following the formation of the initial condensed structure, conformational adjustments took place so that a tight DNA compaction has happened as it can be seen by comparing P-P RDF between 500 ns and 1000 ns (blue and yellow in Figure 6A), where a higher first peak is observed at 900-1000 ns (yellow curve in Figure 6A; Figure 7A). Almost identical PP RDF curves after 1000 ns (yellow, purple and green in Figure 6A) indicates that the simulations have reached equilibrium.


Page 15 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


15

Similar aggregation behavior was also observed in the MD simulations with AMBER Parmbsc0 force field. Though the trajectory is only 500 ns long, a tightly packed fiber of four DNA double helices is formed as shown in Figures 6C and 7C. The P-P RDF curves of 300-400 ns and 400-500 ns show high peaks in short distance range and are almost identical. This indicates that DNA helices are in an equilibrated, compact state. Furthermore, in a simulation with identical setup but using the CHARMM27 force field, the attraction and DNA condensation behavior was also demonstrated (data not shown). On the other hand, in the control simulation with Mg2+ replacing CoHex3+, no aggregation has been observed up to the 1000 ns. The inter-DNA phosphate-phosphate RDF shown in Figure 6B indicates that DNA double helices repel each other and are dispersed in the simulation cell. The snapshots in Figure 7B show no indication of DNA aggregation. Figure 6. Inter-DNA phosphate-phosphate radial distribution functions of simulations of DNA condensation. Top: CHARMM36 DNA + CoHex3+ simulation; Middle: DNA + Mg2+ simulation; Bottom: AMBER Parmbsc0 DNA + CoHex3+. Each curve is calculated using 100 ns trajectory window as denoted in the legend.

Close examination of the trajectories shows that the DNA oligomers are condensed in a stepwise manner (see snapshots in Figure 7 A and C). Initial DNA-DNA contacts can be observed close to the onset of simulations with CoHex3+. Soon, two DNA double helices were aggregated, before a

three-DNA bundle formed. In the simulation with CHARMM36, the three-DNA bundle can be seen from 120 ns onwards, while the three-DNA bundle appears at about 80 ns in the simulation with AMBER Parmbsc0. In the simulation with CHARMM36, the three-DNA bundle was stable for the rest of the trajectory. The fourth DNA double helix has formed occasional contacts with



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 36

16

Figure 7. Snapshots from the MD simulations examining DNA condensation. Three trajectories, namely A: CHARMM36 DNA with CoHex3+, B: CHARMM36 DNA with Mg2+ and C: AMBER Parmbsc0 DNA with CoHex3+, are presented by snapshots at selected time points. Four DNA helices are colored blue, red, yellow and green respectively. Periodic images are colored gray. The time point is indicated at the top of each snapshot.

the other three DNAs in the aggregate. At about 410 ns, the fourth DNA forms a temporary bridge with the periodic image of the bundle; however, it does not result in a long lasting configuration. Finally, at about 595 ns, a stable configuration was formed, in which the fourth DNA makes contact with the bundle and its periodic image, forming a long DNA fiber. With minor conformational variations, the fiber bundle was stable through the rest of the trajectory. In the simulation with AMBER Parmbsc0, on the other hand, the fiber was formed as early as at 100 ns. Then the three-DNA bundle went through conformational rearrangements, that included adopting a straight conformation of the DNA No 4 (colored yellow in Figure 7), sliding of the DNA No 1 (blue) to get tighter contact with the DNA No 3 (green); until a stable fiber was formed at about 300 ns. Due to the simulation setup, in which the length of four DNA double helices is larger than the dimensions of the simulation box, one might expect appearance of periodic structures for this sort of DNA- and box-size combination. However, the fact that the system aggregates to condensed


Page 17 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


17

DNA structures clearly demonstrates that this CoHex3+ model is capable of simulating cationinduced DNA attraction and condensation. It is well established by theoretical considerations and confirmed by simulations within the primitive electrolyte model since more than 20 years ago67-71 that the underlying physical reason of polyelectrolyte attraction in the presence of multivalent counterions is due to the effective attraction between like-charged polyions. This attraction is caused by dynamic correlations of counterions, which in a certain sense resembles the attractive London (van der Waals) dispersion force caused by correlation of electron clouds. The magnitude of the counterion correlations depends crucially on the ion charge and size: Higher valence and smaller size of the ions interact stronger with the polyion, which leads to stronger correlation effects and stronger attraction. The detailed atomistic simulation models modulate this picture by providing addition chemical specificity and molecular detail, while maintaining the underlying driving force of polyelectrolyte aggregation. In order to understand the difference in aggregating and non-aggregating systems and get molecular insights to the mechanism of CoHex3+-induced DNA condensation, it is useful to examine the spatial distribution function (SDF) of cations around DNA, especially near the charged phosphate groups. In Figure 8, the SDFs of CoHex3+ and Mg2+ averaged in the reference frame of the DNA phosphate group are plotted with iso-surfaces are drawn with a threshold of equal charge density. Clearly the charge-neutralizing presence of CoHex3+ is much more significant than that of Mg2+ as indicated by the volume and sizes of the orange (Figure 8A) and blue (Figure 8B) SDFs. In the MD simulations with several double helical DNAs and cationic oligopeptides mimicking fragments of the histone tails it has been shown that it is the higher occupancies of the cationic groups near the phosphates that are responsible for the

counterion-mediated

DNA-DNA

attraction.72 Figure 8. Cation spatial distribution function around phosphate group of the DNA. The volume surfaces in both plots are corresponding to the same charge density. Phosphate atom is colored brown, while oxygen atoms are colored red and other atoms are colored gray.

The above results are in contrast to the observations based on simulations using the CoHex3+ model by Zhang and Mu24 with CHARMM36 DNA. The disability of the



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 36

18

CoHex3+ in reproducing aggregation of DNA in this simulation was traced to its large size. The positions of the first peak of Co-Ow RDF curves in Figure 2 indicate that this model has the largest molecular size of all three all-atom models of CoHex3+. When interacting with DNA, as shown by the SDF presented in Figure 9A, CoHex3+ exhibits significantly less ion density around the DNA double helix compared to our newly developed model (Figure 9B). The charge density around DNA from this model is comparable to that of Mg2+, as can be seen comparing the isosufaces in Figure 9 A and C. Since CoHex3+ and DNA are oppositely charged species that interact directly, a small change in molecular size (vdW radius) can result in a large difference in the electrostatic energy. Therefore, the CoHex3+ with parameters from the work24 is not able to induce the experimentally known effect of inducing DNA aggregation. The results of our model, on the other hand, shows strong interaction of CoHex3+ with DNA, which will lead to DNA-DNA attraction in configurations where the relative position of the different oligonucleotides are correlated such that attractive interactions occurs between phosphate groups of one DNA molecule with CoHex3+ ions condensed on a neighboring DNA. The results are also in agreement with the multishell model proposed by Tolokh et al.73-74 whereby ions that are bound to the "external" shell outside the phosphate groups are key to mediation of DNA attraction. Such contacts seem to be less frequent with the old parameters, while being abundant present in the

simulations

based

on

the

new

parameters. These contacts are also less frequent in the system with Mg2+ as the counterion as shown in Figure 9.

Figure 9. Spatial distribution function of CoHex3+ or Mg2+ about DNA from simulations with (A) the Zhang & Mu’s CoHex3+ model,24 (B) our newly developed CoHex3+ model and (C) Mg2+. DNA is rendered black. Isosurfaces of SDF are plotted for the same charge density as in Figure 8.

4. CONCLUSIONS An all-atom cobalt(III)-hexammine model parameterized

from

ab-initio

(Car-

Parrinello) molecular dynamics has been developed in this work with the aim at


Page 19 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


19

simulating DNA-DNA attraction and bundling induced by this cation. The parameters of the model were tuned to reproduce the aqueous solution structure around CoHex3+ ions obtained in the CPMD simulations. Specifically, the non-bonded parameters of CoHex3+, which is the most critical aspect of this model, have been fitted to reproduce the CoHex3+-Ow RDF predicted by the CPMD simulations. MD simulations based on the new model have shown good agreement with experimental and computational results, including 1) the CoHex3+-induced B-DNA to A-DNA transition; 2) correct prediction of the CoHex3+ effect on the persistence length of both RNA and DNA; 3) ability of the model to describe CoHex3+- induced attraction of DNA double helices and formation of a condensed bundle-like structure. The model has been tested for DNA described by the CHARMM36 and AMBER Parmbsc0 force field investigating the DNA condensation behavior in the presence of CoHex3+. Since our derivation was based exclusively on ab-initio simulation of CoHex3+ ion with water, and taking into account that both CHARMM and AMBER use the TIP3P water model, we expect that this model can show comparable performance in simulations of nucleic acids within both CHARMM and AMBER force fields. This expectation is fulfilled, lastly, by one simulation with each force field that both demonstrate the experimentally observed bundling of DNA molecules in the presence of cobalt(III)-hexammine.

ASSOCIATED CONTENT Supporting Information The Supporting Information is available free of charge on the ACS Publications website at http://pubs.acs.org/. Further description of the MD simulation results using Zhang & Mu’s CoHex3+ parameters with two figures (Figure S1 and Figure S2).

AUTHOR INFORMATION Corresponding Authors *A. Lyubartsev. E-mail: [email protected] *L. Nordenskiöld. E-mail: [email protected] Notes The authors declare no competing financial interest.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 36

20

ACKNOWLEDGEMENTS This work was supported by Singapore Ministry of Education Academic Research Fund (AcRF) Tier 2 (MOE2014-T2-1-123 (ARC51/14)) and Tier 3 (MOE2012-T3-1-001) grants (to LN) and by the Swedish Research Council (to APL). We are indebted to Vishal Minhas for valuable discussions and to Prof Aatto Laaksonen for suggestions. REFERENCES (1) Wilson, R. W.; Bloomfield, V. A., Counterion-induced condensation of deoxyribonucleic acid. A light-scattering study. Biochemistry 1979, 18, 2192-2196. (2) Pelta, J.; Livolant, F.; Sikorav, J. L., DNA aggregation induced by polyamines and cobalthexamine. J. Biol. Chem. 1996, 271, 5656-5662. (3) Bloomfield, V. A., DNA condensation by multivalent cations. Biopolymers 1997, 44, 269282. (4) Deng, H.; Bloomfield, V. A., Structural effects of cobalt-amine compounds on DNA condensation. Biophys. J. 1999, 77, 1556-1561. (5) Matulis, D.; Rouzina, I.; Bloomfield, V. A., Thermodynamics of DNA binding and condensation: isothermal titration calorimetry and electrostatic mechanism. J. Mol. Biol. 2000, 296, 1053-1063. (6) Kankia, B. I.; Buckin, V. A.; Bloomfield, V. A., Hexamminecolbalt(III)-induced condensation of calf thymus DNA: circular dichroism and hydration measurements. Nucleic Acids Res. 2001, 29, 2795-2801. (7) Hud, N. V.; Vilfan, I. D., Toroidal DNA condensates: Unraveling the fine structure and the role of nucleation in determining size. Annu. Rev. Biophys. Biomol. Struct. 2005, 34, 295-318. (8) Vilfan, I. D.; Conwell, C. C.; Sarkar, T.; Hud, N. V., Time study of DNA condensate morphology: implications regarding the nucleation, growth, and equilibrium populations of toroids and rods. Biochemistry 2006, 45, 8174-8183. (9) Korolev, N.; Berezhnoy, N. V.; Eom, K. D.; Tam, J. P.; Nordenskiöld, L., A universal description for the experimental behavior of salt-(in)dependent oligocation-induced DNA condensation. Nucleic Acids Res. 2012, 40, 2808-2821.


Page 21 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


21

(10) Dai, L.; Mu, Y.; Nordenskiöld, L.; van der Maarel, J. R. C., Molecular dynamics simulation of multivalent-ion mediated DNA attraction between DNA molecules. Phys. Rev. Lett. 2008, 100, 118301. (11) Yoo, J.; Aksimentiev, A., The structure and intermolecular forces of DNA condensates. Nucleic Acids Res. 2016, 44, 2036-2046. (12) Cornell, W. D.; Cieplak, R.; Bayly, C. L.; Gould, I. R.; Merz, K. M.; Ferguson, D. M.; Spellmeyer, D. G.; Fox, T.; Caldwell, J. W.; Kollman, P. A., A second generation force field for the simulation of proteins, nucleic acids, and organic molecules. J. Am. Chem. Soc. 1995, 117, 5179-5197. (13) Cheatham, T. E.; Case, D. A., Twenty-five years of nucleic acid simulations. Biopolymers 2013, 99, 969-977. (14) Weiner, S. J.; Kollman, P. A.; Case, D. A.; Singh, U. C.; Ghio, C.; Alagona, G.; Profeta, S.; Weiner, P., A new force field for molecular mechanical simulation of nucleic acids and proteins. J. Am. Chem. Soc. 1984, 106, 765-784. (15) Foloppe, N.; MacKerell, A. D., All-atom empirical force field for nucleic acids: I. Parameter optimization based on small molecule and condensed phase macromolecular target data. J. Comp. Chem. 2000, 21, 86-104. (16) MacKerell, A. D.; Banavali, N., All-atom empirical force field for nucleic acids: II. Application to molecular dynamics simulations of DNA and RNA in solution. J. Comp. Chem. 2000, 21, 105-120. (17) Hart, K.; Foloppe, N.; Baker, C. M.; Denning, E. J.; Nilsson, L.; MacKerell, A. D., Optimization of the CHARMM additive force field for DNA: Improved treatment of the BI/BII conformational equilibrium. J. Chem. Theory Comput. 2012, 8, 348-362. (18) MacKerell, A. D.; Nilsson, L., Molecular dynamics simulations of nucleic acid–protein complexes. Curr. Opin. Struct. Biol. 2008, 18, 194-199. (19) Cheatham, T. E.; Kollman, P. A., Insight into the stabilization of A-DNA by specific ion association: spontaneous B-DNA to A-DNA transitions observed in molecular dynamics simulations of d[ACCCGCGGGT]2 in the presence of hexaamminecobalt(III). Structure 1997, 5, 1297-1311. (20) Wu, Y. Y.; Zhang, Z. L.; Zhang, J. S.; Zhu, X. L.; Tan, Z. J., Multivalent ion-mediated nucleic acid helix-helix interactions: RNA versus DNA. Nucleic Acids Res. 2015, 43, 6156-6165.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 36

22

(21) Drozdetski, A. V.; Tolokh, I. S.; Pollack, L.; Baker, N.; Onufriev, A. V., Opposing effects of multivalent ions on the flexibility of DNA and RNA. Phys. Rev. Lett. 2016, 117, 028101. (22) Pabit, S. A.; Katz, A. M.; Tolokh, I. S.; Drozdetski, A.; Baker, N.; Onufriev, A. V.; Pollack, L., Understanding nucleic acid structural changes by comparing wide-angle x-ray scattering (WAXS) experiments to molecular dynamics simulations. J. Chem. Phys. 2016, 144, 205102. (23) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L., Comparison of simple potential functions for simulating liquid water. J. Chem. Phys. 1987, 79, 926-935. (24) Zhang, T.; Mu, Y., Initial binding of ions to the interhelical loops of divalent ion transporter CorA: replica exchange molecular dynamics simulation study. PLoS One 2012, 7, e43872. (25) Brennan, R. G.; Westhof, E.; Sundaralingam, M., Structure of a Z-DNA with two different backbone chain conformations. Stabilization of the decadeoxyoligonucleotide d(CGTACGTACG) by [Co(NH3)6]3+ binding to the guanine. J. Biomol. Struct. Dyn. 1986, 3, 649-665. (26) Cheatham, T. E.; Kollman, P. A., Molecular dynamics simulations highlight the structural differences among DNA:DNA, RNA:RNA, and DNA:RNA structures. J. Am. Chem. Soc. 1997, 119, 4805-4825. (27) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A.; Haak, J. R., Molecular dynamics with coupling to an external bath. J. Chem. Phys. 1984, 81, 3684-3690. (28) Car, R.; Parrinello, M., Unified approach for Molecular Dynamics and density-functional theory. Phys. Rev. Lett. 1985, 55, 2471-2474. (29) Marx, D.; Hutter, J., Ab initio molecular dynamics. Theory and Implementation. In Modern Methods and Algorithms of Quantum Chemistry, Grotendorst, J., Ed. John von Neumann Institute for Computing: Julich, 2000; Vol. 1, pp 301-449. (30) Andreoni, W.; Curioni, A., New advances in chemistry and materials science with CPMD and parallel computing. Parallel Computing 2000, 26, 819-842. (31) Marx, D.; Hutter, J., Ab Initio Molecular Dynamics: Basic Theory and Advanced Methods. Cambridge University Press: Cambridge, 2009; p 578.


Page 23 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


23

(32) Becke, A. D., Density-functional exchange-energy approximation with correct asymptotic behavior. Phys. Rev. A 1988, 38, 3098-3100. (33) Lee, C.; Yang, W.; Parr, R. G., Development of the Colle-Salvetti correlation-energy formula into a functional of the electron density. Phys. Rev. B 1988, 37, 785-789. (34) Troullier, N.; Martins, J. L., Efficient pseudopotentials for plane-wave calculations. Phys. Rev. B 1991, 43, 1993-2006. (35) Goedecker, S.; Teter, M.; Hutter, J., Separable dual-space Gaussian pseudopotentials. Phys. Rev. B 1996, 54, 1703-1710. (36) Hartwigsen, C.; Goedecker, S.; Hutter, J., Relativistic separable dual-space Gaussian pseudopotentials from H to Rn. Phys. Rev. B 1998, 58, 3641-3662. (37) Krack, M., Pseudopotentials for H to Kr optimized for gradient-corrected exchangecorrelation functionals. Theor. Chem. Acc. 2005, 114, 145-152. (38) Abraham, M. J.; Murtola, T.; Schulz, R.; Páll, S.; Smith, J. C.; Hess, B.; Lindahl, E., GROMACS: High performance molecular simulations through multi-level parallelism from laptops to supercomputers. SoftwareX 2015, 1–2, 19-25. (39) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G., A smooth particle mesh Ewald method. J. Chem. Phys. 1995, 103, 8577-8593. (40) Bussi, G.; Donadio, D.; Parrinello, M., Canonical sampling through velocity rescaling. J. Chem. Phys. 2007, 126, 014101. (41) Parrinello, M.; Rahman, A., Crystal structure and pair potentials: A Molecular-Dynamics study. Phys. Rev. Lett. 1980, 45, 1196-1199. (42) Parrinello, M.; Rahman, A., Polymorphic transitions in single crystals: A new molecular dynamics method. J. Appl. Phys. 1981, 52, 7182-7190. (43) Hess, B., P-LINCS: A parallel linear constraint solver for Molecular Simulation. J. Chem. Theory Comput. 2008, 4, 116-122. (44) Miyamoto, S.; Kollman, P. A., Settle: An analytical version of the SHAKE and RATTLE algorithm for rigid water models. J. Comput. Chem. 1992, 13, 952-962. (45) Yoo, J.; Aksimentiev, A., Improved parametrization of Li+, Na+, K+, and Mg2+ oons for all-atom molecular dynamics simulations of nucleic acid systems. J. Phys. Chem. Lett. 2012, 3, 45-50.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 36

24

(46) Pan, Y.; MacKerell, A. D., Altered structural fluctuations in duplex RNA versus DNA: a conformational switch involving base pair opening. Nucleic Acids Res. 2003, 31, 7131-7140. (47) Galindo-Murillo, R.; Roe, D. R.; Cheatham, T. E., Convergence and reproducibility in molecular dynamics simulations of the DNA duplex d(GCACGAACGAACGAACGC). Biochim. Biophys. Acta 2015, 1850, 1041-1058. (48) Korolev, N.; Nordenskiöld, L., H4 histone tail mediated DNA-DNA interaction and effects on DNA structure, flexibility and counterion binding. A molecular dynamics study. Biopolymers 2007, 86, 409–423. (49) Robinson, H.; Wang, A. H., Neomycin, spermine and hexaamminecobalt (III) share common structural motifs in converting B- to A-DNA. Nucleic Acids Res. 1996, 24, 676-682. (50) Macke, T. J.; Case, D. A., Modeling unusual nucleic acid structures. In Molecular modeling of nucleic acids, Leontis, N. B.; SantaLucia, J., Eds. American Chemical Society: Washington, DC, 1997; Vol. 682, pp 379–393. (51) Neria, E.; Fischer, S.; Karplus, M., Simulation of activation free energies in molecular systems. J. Chem. Phys. 1996, 105, 1902-1921. (52) Denning, E. J.; Priyakumar, U. D.; Nilsson, L.; MacKerell, A. D., Impact of 2′-hydroxyl sampling on the conformational properties of RNA: Update of the CHARMM all-atom additive force field for RNA. J. Comput. Chem. 2011, 32, 1929-1943. (53) Lavery, R.; Moakher, M.; Maddocks, J. H.; Petkeviciute, D.; Zakrzewska, K., Conformational analysis of nucleic acids revisited: Curves+. Nucleic Acids Res. 2009, 37, 59175929. (54) Korolev, N.; Di, L.; Lyubartsev, A. P.; Nordenskiöld, L., A coarse-grained DNA model parameterized from atomistic simulations by inverse Monte Carlo. Polymers 2014, 6, 1655-1675. (55) Ivanov, V. I.; Minchenkova, L. E.; Minyat, E. E.; Frank-Kamenetskii, M. D.; Schyolkina, A. K., The B to A transition of DNA in solution. J. Mol. Biol. 1974, 87, 817-833. (56) Minyat, E. E.; Ivanov, V. I.; Kritzyn, A. M.; Minchenkova, L. E.; Schyolkina, A. K., Spermine and spermidine induced B to A transition of DNA in solution. J. Mol. Biol. 1978, 128, 397-409. (57) Xu, Q.; Shoemaker, R. K.; Braunlin, W. H., Induction of B-A transitions of deoxyoligonucleotides by multivalent cations in dilute aqueous solution. Biophys. J. 1993, 65, 1039-1049.


Page 25 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


25

(58) Lu, X.-J.; Shakked, Z.; Olson, W. K., A-form Conformational Motifs in Ligand-bound DNA Structures. J. Mol. Biol. 2000, 300, 819-840. (59) Kulkarni, M.; Mukherjee, A., Understanding B-DNA to A-DNA transition in the righthanded DNA helix: Perspective from a local to global transition. Prog. Biophys. Mol. Biol. 2017, in press, doi: 10.1016/j.pbiomolbio.2017.05.009. (60) Lu, X.-J.; El Hassan, M. A.; Hunter, C. A., Structure and conformation of helical nucleic acids: rebuilding program (SCHNArP)11Edited by K. Nagai. J. Mol. Biol. 1997, 273, 681-691. (61) Hagerman, P. J., Investigation of the flexibility of DNA using transient electric birefringence. Biopolymers 1981, 20, 1503–1535. (62) Baumann, C. G.; Smith, S. B.; Bloomfield, V. A.; Bustamante, C., Ionic effects on the elasticity of single DNA molecules. Proc. Natl. Acad. Sci. U.S.A. 1997, 94, 6185-6190. (63) Lu, Y.; Weers, B.; Stellwagen, N. C., DNA persistence length revisited. Biopolymers 2002, 61, 261-275. (64) Baumann, C. G.; Bloomfield, V. A.; Smith, S. B.; Bustamante, C.; Wang, M. D.; Block, S. B., Stretching of single collapsed DNA molecules. Biophys. J. 2000, 78, 1965-1978. (65) Bomble, Y. J.; Case, D. A., Multiscale modeling of nucleic acids: Insights into DNA flexibility. Biopolymers 2008, 89, 722-731. (66) Mazur, A. K., Evaluation of elastic properties of atomistic DNA models. Biophys. J. 2006, 91, 4507-4518. (67) Kjellander, L.; Marcelja, S., Correlation and image charge effects in electric doublelayers. Chem. Phys. Lett. 1984, 112, 49-53. (68) Guldbrand, L.; Nilsson, L. G.; Nordenskiöld, L., A Monte Carlo simulation study of electrostatic forces between hexagonally packed DNA double helices. J. Chem. Phys. 1986, 85, 6686-6698. (69) Lyubartsev, A. P.; Nordenskiöld, L., Monte Carlo simulation study of ion distribution and osmotic pressure in hexagonally oriented DNA. J. Phys. Chem. 1995, 99, 10373-10382. (70) Gronbech-Jensen, N.; Mashl, R. J.; Bruinsma, R. F.; Gelbart, W. M., Counterion-induced attraction between rigid polyelectrolytes. Phys. Rev. Lett. 1997, 78, 2477-2480. (71) Ray, J.; Manning, G. S., Effect of counterion valence and polymer charge density on the pair potential of two polyions. Macromolecules 1997, 30, 5739-5744.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 36

26

(72) Korolev, N.; Yu, H.; Lyubartsev, A. P.; Nordenskiöld, L., Molecular dynamics simulations demonstrate the regulation of DNA-DNA attraction by H4 histone tail acetylations and mutations. Biopolymers 2014, 101, 1051-1064. (73) Tolokh, I. S.; Drozdetski, A. V.; Pollack, L.; Baker, N. A.; Onufriev, A. V., Multi-shell model of ion-induced nucleic acid condensation. J. Chem. Phys. 2016, 144, 155101. (74) Tolokh, I. S.; Pabit, S. A.; Katz, A. M.; Chen, Y.; Drozdetski, A.; Baker, N.; Pollack, L.; Onufriev, A. V., Why double-stranded RNA resists condensation. Nucleic Acids Res 2014, 42, 10823-10831.


Page 27 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32



The Journal of Physical Chemistry Page 28 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16


Page 29 TheofJournal 36 of Physical Chemistry

1 2 3 4 5 6 7 8 9 10



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24



1 2 3 4 5 6 7 8 9 10



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26


Page 34 of 36


1 2 3 4 5



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20


All-Atom MD Simulation of DNA Condensation ... - ACS Publications

Recommend Documents