Markov State Model of Ion Assembling Process - The Journal of

Article pubs.acs.org/JPCA

Markov State Model of Ion Assembling Process Roman Shevchuk* Department for Mathematics and Computer Science, Freie Universität Berlin, Arnimallee 6, 14195 Berlin, Germany Institute for Microbiology and Genetics, Georg-August-University Göttingen, Justus-von-Liebig-Weg 11, 37077 Göttingen, Germany ABSTRACT: We study the process of ion assembling in aqueous solution by means of molecular dynamics. In this article, we present a method to study many-particle assembly using the Markov state model formalism. We observed that at NaCl concentration higher than 1.49 mol/kg, the system tends to form a big ionic cluster composed of roughly 70−90% of the total number of ions. Using Markov state models, we estimated the average time needed for the system to make a transition from discorded state to a state with big ionic cluster. Our results suggest that the characteristic time to form an ionic cluster is a negative exponential function of the salt concentration. Moreover, we defined and analyzed three different kinetic states of a single ion particle. These states correspond to a different particle location during nucleation process.

I. INTRODUCTION Aqueous ion solutions are natural environment of many biological systems, and therefore they affect many biological processes.1,2 Solvated salts affect the strength of hydrogen bond network of waters as well as the strength of interatomic interactions in the solutions.3,4 If salt ions are added into the water, there is a certain limit at which further salt addition will not increase the amount of dissolved salt. This additional amount of ions will aggregate into a solid phase. However, not much is known about kinetical properties of ions at such supersaturated concentrations. It was found that the increase of the salt concentration decreases the mobility of both water and ions.5,6 Several studies were made to estimate the effect of ions on water structure.7−9 Recently, Soniat et al. employed a polarizable force field to study ionic clusters. For the case of NaCl they found that the amount of ion pairing is not larger than random. Another interesting result is the fact that ionic clusters tend to be positively charged.10 Several methods for calculating the solubility point from molecular perspective were proposed.11−14 Experiments suggest that in ionic aqueous solutions a significant amount of ions are present in clusters.15,16 Several models were used to describe nucleation growth of salts and other materials.16 Hassan observed ion clusters of around 10 ions in undersaturated solution, however without any indications of further aggregation.17,18 He suggested that the size and morphology of initial nucleus and its hydration context are important factors for the nucleation process. Direct nucleation of ions in molecular dynamics simulation was observed by Zahn.19 His observations suggest that aggregation starts with the single sodium ion, which lacks water molecules in its solvation shell. However, his study did not consider the time scales of the aggregation process. Recently, the formation of ionic clusters gradually increasing salt concentration was studied by Hassan.20 Up to a certain concentration he observed © 2016 American Chemical Society

relatively small unstable clusters, whereas a further increase of concentration led to formation of a nanocrystal. Lanaro et al. studied a dissolution of ionic crystals with different starting morphology. They found that dissolution of ion clusters occurs in three distinct stages. Interestingly, they found that the shape of the microscopic crystal determines the dissolution rate.21 Markov state models (MSM) is a powerful tool to study complex molecular processes. It was used to study proteins,22−25 liquid dynamics.26−28 Recently, it was employed to study assembling of nanoparticles.29 The essential idea of this method is to construct a kinetic model from the molecular dynamics trajectories. From the obtained model one can extract all important information about the system, such as the population of metastable states and the characteristic times of the transitions.30 Here we employ MSM to study the effect of ion concentration to the ion assembling process.

II. METHODS II.A. Molecular Dynamics Setup. Molecular dynamics simulations were run with the program GROMACS31 with an integration time step of 2 fs. The simulation box consisted of 1024 molecules in the NPT ensemble with pressure of 1 atm and temperature 300 K. To mimic an infinite system, we used cubic simulation box where periodic boundary conditions were applied in all three directions. Different salt concentrations were obtained by replacing the pair of water molecules with ion pair. Here we studied the molalities from 1.37 to 1.85 mol/kg with the step of 0.12 mol/kg, where the highest concentration corresponds to 32 NaCl pairs. The Berendsen barostat,32 velocity rescale thermostat,33 and PME34 were used for pressure coupling, temperature coupling, and long-range Received: November 29, 2015 Revised: April 21, 2016 Published: April 21, 2016 2783

DOI: 10.1021/acs.jpca.5b11650 J. Phys. Chem. A 2016, 120, 2783−2788

Article

The Journal of Physical Chemistry A

time scales (ITS) analysis. By definition, the ith implied time scale of the transition matrix is given by τ ITSi = − ln λi (3)

electrostatics, respectively. The cutoff of 1 nm was applied to both direct-space electrostatic and Lennard-Jones interactions. The TIP4P/2005 model was used for water molecules.35 To model interatomic interactions, we used the AMBER03 force field36 applying Lorentz−Berthelot rules for the ion interactions.37−39 To test possible system size effects, we simulated a system where the numbers of both waters and ions were increased by a factor of 10. In this system we also observed a qualitatively similar ion aggregation process. Though it is known that the AMBER force field tends to overestimate ion pairing,40 our aim here is to show the possibility of studying the aggregation process with the usage of MSM. To build MSM, which accurately describes the kinetics of the process, one needs good statistics and using the current force field allows us to get nucleation molecular dynamics trajectories relatively fast in terms of the computational speed. Therefore, it is important to note that because the used force field overestimates ion pairing, the nucleation may not be observed in experiments at concentrations studied here. The analysis was done for each concentration value over 9 trajectories of 10 000 snapshots obtained from a 200 ns long run. II.B. Order Parameters. To build the Markov model, which can separate different metastable states of the system one requires a good order parameter. To measure the local order around a single ion, we used the order parameter proposed by Steinhardt et al.41 First, we calculate the complex vector qlm for each ion i as qlm(i) =

1 Nb(i)

Nb(i)

⎛ rij ⎞ ⎟⎟ ⎝ |rij| ⎠

∑ Ylm⎜⎜ j=1

where λi is the ith eigenvalue of the transition matrix and τ is the lag time. By increasing lag time τ, one will reach the point where all implied time scales converge to the certain value tc. The converged value corresponds to the time of transition from one state to another. We validated that in the studied system the slowest ITS corresponds to the transition between the disordered and ordered states of the system by checking the values of corresponding eigenvectors of the transition matrix along chosen order parameter.43 To check the obtained characteristic time of transition between the disordered and ordered state of the system, we calculated the distribution of the mean first passage time between two states. Arrival times were obtained by running 10 000 random walks over the transition matrix. As a starting configuration we chose the state with the lowest value of an order parameter, and as a target configuration we picked the states in which the biggest cluster consists of at least 80% of the total number of ions. Because single ion transitions between metastable states occur at a much faster time scale than aggregation of the cluster, we used the Markov Clustering Algorithm (MCL) in combination with Kolmogorov−Smirnov test to build a Markov model from the single ion transition network.44−46 The MCL algorithm is based on the properties of the random walker on transition matrix and was successfully used in numerous applications.47−49 The essential part of this analysis is that it takes into account not only the instant value of an order parameter but also its history. To assign a specific frame at time t to a microstate, we consider not only its value Xt but also the distribution in a time window Xt−τ ... Xt+τ. After obtaining distributions of order parameter for two different frames, we applied the Kolmogorov−Smirnov test.45 If the difference of two distributions is smaller than a certain tolerance value, we assign the values of both frames to the microstate Mi. Thus, we obtained the transition matrix with the time-dependent microstates and grouped them into metastable states (groups of microstates) with the help of MCL algorithm.44 For our analysis of the distributions we used the time window of 100 frames, which corresponds to the time of 2 ns whereas for the MCL algorithm we set the granularity value to 1.2. Because this analysis requires equilibrium trajectories, we skipped the part of the trajectory until the point where the size of the biggest cluster converged. The detailed description of this method can be found elsewhere.46 Markov state model analysis was done with the EMMA and Aqualab software packages.46,50

(1)

where Nb is the number of neighboring ions, Ylm are the spherical harmonics where l is a free integer parameter and m is the integer that runs from −l to +l. For simplicity we always put the number of nearest neighbors equal to six to avoid an uncertainty in the cases where the ion is floating in bulk without any neighbors inside its solvation shell. Thus, one can obtain the value of Ql for single particle as Ql =

1 Nb(i)

Nb(i)

l

∑ ∑ j = 1 m =−l

* qlmqlm

(2)

where * denotes complex conjugate. II.C. Markov State Models. We constructed transition matrix in the following way. First, we discretized the space of the order parameter to the finite number of microstates rounding the value of an order parameter to the integer value. We increase the count matrix element Cij if the system was in state i at time t and in state j at time t + τ where τ is the lag time of the Markov model. Subsequently, the transition matrix T was obtained by normalizing every row of the count matrix to one. One has to check that the transition matrix holds Markovity. The easy way to do this is to calculate several different transition matrices decreasing lag time τ and to estimate the slowest transition rate of the system. Once the slowest transition rate converges, the transition matrix at that particular lag time is Markovian. This transition matrix contains all the kinetical information on the system. The detailed description of the transition matrix construction and its properties were presented elsewhere.42 To obtain the characteristic transition times from the disordered to the ordered state of the system, we used implied

III. RESULTS III.A. Aggregation Process. We started the simulations at all concentrations from the state where all the ions were placed at random positions and no big ion cluster was present. During the first few nanoseconds of the simulation, the ion cluster that consists typically of around 30% of the total number of ions is formed. After the formation, this cluster did not dissociate and kept its relative size with some single ion binding and unbinding events. Then, after several tens of nanoseconds, a spontaneous aggregation occurred and the size of the cluster increased up to around 90% of the total number of ions. The 2784


Article


comparable with the total number of ions, its value of Q6 increases, which shows that the cluster is not only increased in size but also increased in local order (Figure 2, central and right panels). III.B. Implied Time Scales of Different Order Parameters. A recent study by Beckham et al. suggested that, to study the nucleation processes, the combination of size of the cluster and its local order parameter Q6 works better than other order parameters.51 An ion cluster (Ncl) here is defined as an structure of ions where every ion has a connection with another ions via a path of bonds. Two ions were considered as bonded if the distance between them is smaller than 0.4 nm, which is the value of first minimum of the radial distribution function. The same process can be described with different sets of order parameters, and one can obtain different transition times from the Markov model using different order parameters. Due to variational principle, the order parameter that gives as an output the slowest implied time scale has to be picked.42 To check that this particular order parameter produces a good description of the aggregation process, we compared the converged implied time scale of transition matrix constructed by Ncl·Q6cl with other possible parameters. In particular, we tested Ncl, Q̅ 6, Q4cl, Q̅ 4, Q4cl·Ncl and found that the Markov model produced the slowest time scale if Ncl·Q6cl is used as an order parameter, whereas in the cases of other order parameters the implied time scale is 10−20% smaller (Table 1). To check

fully aggregated cluster is shown in Figure 1. The time for this aggregation varied a lot among different concentrations and

Figure 1. Ion cluster at the end of the aggregation process. Sodium and chloride ions are shown as blue and light green spheres, respectively. Water molecules are not shown for clarity.

Table 1. Values of ITS for Different Order Parameters Calculated from the Simulation at m = 1.74 mol/kg

runs. Moreover, in some of the trajectories we did not observe any cluster growing after it reached the size of around 30%. For example, from nine trajectories of 200 ns at m = 1.74 mol/kg, the fastest growing occurs after only 20 ns whereas the slowest after around 150 ns. The typical trajectory of the ion cluster size is shown on the left panel of Figure 2. Likewise, at the lowest studied molality, m = 1.37 mol/kg, only four of the nine simulations led to the growing of the stable cluster with size of 90% of the system size, whereas in other simulations such a cluster was formed for a few nanoseconds and then dissociated back to the smaller size. We monitored the value of the Q6 order parameter and did not find any significant differences in its value when the cluster size is below 50% of the total number of ions. From this observation we concluded that although the cluster has a relatively big size it does not have any crystal order inside it. On the contrary, at the point the cluster size is

order parameter

ITS, ns

Ncl Q̅ 6 Q6cl Ncl·Q6cl Q̅ 4 Q4cl Ncl·Q4cl

31.26 29.12 31.98 33.96 31.68 30.52 32.56

that the chosen order parameter captures the transition, we calculated the second eigenvector of the transition matrix used in MSM.43 In Figure 3 the typical distribution of eigenvector values for Q6·Ncl is shown. It clearly demonstrates the transition from disordered state to a crystal. III.C. MSM of the Aggregation Process. The only direct observation we can make strictly from the results of the

Figure 2. Time series of the biggest cluster size (left panel), Q6 order parameter (central panel), and their relative counts (right panel). 2785


Article


Using transition matrices calculated at different lag times τ and eq 3, we obtained the converged implied time scale for all molalities except the lowest studied. The results differ from the direct trajectory observation times and show exponential dependence of molality; e.g., the lower the molality, the more time is needed to form a big cluster. To verify this estimation, we calculated average time to form big cluster for the system via mean-first-passage-time analysis (see Methods for details). The obtained result is in a good agreement with the implied time scale values, moreover the point at m = 1.37 mol/kg has the same trend, although it was unreachable via implied time scale analysis due to relatively short length of our trajectory. D. MSM of a Single Ion. After studying the behavior of the whole system, we tried to describe the kinetics of the single ion particles. First, we calculated the lifetime of ion contact. We found that an increase of the molality does not increase the lifetime of ionic contact in a direct way. The lifetime of a contact between two ions from the bulk is the same at all studied molalities. On the contrary, if a single ion binds to the big cluster, it stays there for a much longer time than in a case of a single pair association (Figure 5). These results suggest

Figure 3. Typical profile of the second eigenvector for order parameter Q6·Ncl that was used for MSM construction.

simulation is that the system goes through the transition. However, our simulation data are not enough to directly obtain kinetic rates for this process from raw molecular dynamics trajectories. We tried to estimate the average time going from disordered state to the big cluster (Figure 4, blue line), but we did not find any strict dependency of the aggregation time on the molality.

Figure 5. Lifetime of ion contact distribution. The distribution from the simulation fragments where big ionic cluster is present or absent are shown in magenta and blue, respectively. Black lines show exponential fit.

that there are at least two different types of ionic contact. To check the latter, we constructed a MSM for a single ion particle. To detect the metastable states of a single ion, we used the same order parameter as in the description of the whole system but calculated it for every ion independently. After calculating an order parameter, we constructed the transition matrix using every particle trajectory. Because forming and breaking atomic contacts means fast recrossings over the free-energy barrier for detecting metastable states we used Kolmogorov−Smirnov algorithm46 (see Methods for details). Applying this algorithm, we found three different metastable states for a single particle (Figure 6, left panel). The first state corresponds to the situation where an ion is in a bulk or attached to a small disordered cluster. The second and third metastable states represent an ion being involved in the formation of the biggest ionic cluster of the system. They differ only in the value of the local order of the ions first solvation shell. We found that the population of the state with the higher value of order parameter is increasing with increase of the salt concentration. This can be

Figure 4. Characteristic time to go from disordered to ordered state for the system. The blue curve represents the values obtained directly from MD simulations, and red and black curves show the time obtained from Markov model via mean-first-passage-time and implied time scale convergence. For m = 1.37 mol/kg we did not obtain the converged implied time scale due to relatively short length of our trajectory.

To tackle this aim, we employ the Markov state model. As an order parameter for the aggregation process, we used Q6·Ncl. The value of the order parameter was calculated at every snapshot of the trajectory and then was uniformly split into bins with size of 0.1. Following this procedure we obtain approximately 250 nonzero bins. Every bin was accounted as a microstate for the Markov state model. If two microstates occurred sequentially in the trajectory with lag time τ, the corresponding element of the count matrix was increased by one. Thus, a transition matrix was obtained from the count matrix by normalization its every column. 2786


Article


Figure 6. (Left panel) separated distribution of three kinetically different metastable states at m = 1.74 mol/kg. (right panel) population of three metastable states versus molality.

semiordered state, an ion is the part of the big ionic cluster; however, its neighbors still lack perfect crystal order; (iii) fully ordered state, an ion is at the center of ionic cluster and its surrounding follows crystal pattern.

explained by the fact that the increase of the concentration allows the ion to attach to the bigger clusters, and as a consequence, its value of an order parameter is also increasing. Another interesting result is that at the two highest studied molalities, m = 1.74 mol/kg and m = 1.85 mol/kg, the population of a fully ordered state significantly increases with respect to other molalities. This means that some ions in the cluster have a significantly higher value of an order parameter. We found that these ions are typically located in the center of the cluster. Overall, our observations suggest that a single ion makes two kinetic steps during formation of the cluster. First, it connects to the cluster, without having a maximum possible value of an order parameter. Second, the ion appears in the core of the cluster, being the part of perfect crystal grid.

■

AUTHOR INFORMATION

Corresponding Author

*E-mails: [email protected], [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS The author thanks the European Research Council (ERC) starting grant pcCell for funding. Frank Noé is acknowledged for helpful discussions.

IV. CONCLUSIONS In summary, we carried out a series of classical molecular dynamics simulation of sodium chloride in aqueous solution at high salt concentrations. Most of our simulations led to the aggregation of an ion cluster of size of around 90% of the total number of ions. Various order parameters were tested by means of MSM. As a result of our assessment, we found the combination of order parameters which yields the slowest implied time scale and describes aggregation kinetics better that the other order parameters. Using Ncl·Q6 order, we constructed a Markov state model from which we extracted the characteristic time for an aggregation. The time needed to form an ordered ionic cluster at different ion concentrations shows a negative exponential dependence a on salt concentration. The exact value of aggregation times calculated here may be different from the experimental ones due to the fact that the force field we used here overestimates ion pairing. However, our aim here is to show that these times can be calculated with the usage of MSM. Moreover, results obtained with MSM are in better agreement with mean-first-passage times than bruteforce estimation. By studying the properties of a single ion inside the ionic cluster, we found that it forms much more stable contacts than a single ion pair. Also, we built MSM of a single ion, which suggests that three metastable states of a single ion are involved in the assembly. (i) bulk state, an ion is surrounded by water molecules or form small clusters with size of few particles; (ii)

■

REFERENCES

(1) Cacace, M. G.; Landau, E. M.; Ramsden, J. J. The Hofmeister Series: Salt and Solvent Effects on Interfacial Phenomena. Q. Rev. Biophys. 1997, 30 (03), 241−277. (2) Ru, M. T.; Hirokane, S. Y.; Lo, A. S.; Dordick, J. S.; Reimer, J. A.; Clark, D. S. On the Salt-Induced Activation of Lyophilized Enzymes in Organic Solvents: Effect of Salt Kosmotropicity on Enzyme Activity. J. Am. Chem. Soc. 2000, 122 (8), 1565−1571. (3) Cappa, C. D.; Smith, J. D.; Wilson, K. R.; Messer, B. M.; Gilles, M. K.; Cohen, R. C.; Saykally, R. J. Effects of Alkali Metal Halide Salts on the Hydrogen Bond Network of Liquid Water. J. Phys. Chem. B 2005, 109 (15), 7046−7052. (4) Thomas, A. S.; Elcock, A. H. Molecular Dynamics Simulations of Hydrophobic Associations in Aqueous Salt Solutions Indicate a Connection Between Water Hydrogen Bonding and the Hofmeister Effect. J. Am. Chem. Soc. 2007, 129 (48), 14887−14898. (5) Chowdhuri, S.; Chandra, A. Molecular Dynamics Simulations of Aqueous NaCl and KCl Solutions: Effects of Ion concentration on the Single-Particle, Pair, and Collective Dynamical Properties of Ions and Water Molecules. J. Chem. Phys. 2001, 115 (8), 3732−3741. (6) Lyubartsev, A. P.; Laaksonen, A. Concentration Effects in Aqueous NaCl Solutions. A Molecular Dynamics Simulation. J. Phys. Chem. 1996, 100 (40), 16410−16418. (7) Uchida, H.; Matsuoka, M. Molecular Dynamics Simulation of Solution Structure and Dynamics of Aqueous Sodium Chloride Solutions from Dilute to Supersaturated Concentration. Fluid Phase Equilib. 2004, 219 (1), 49−54. 2787


Article

The Journal of Physical Chemistry A (8) Galamba, N. Mapping Structural Perturbations of Water in Ionic Solutions. J. Phys. Chem. B 2012, 116 (17), 5242−5250. (9) Galamba, N. On the Effects of Temperature, Pressure, and Dissolved Salts on the Hydrogen-Bond Network of Water. J. Phys. Chem. B 2013, 117 (2), 589−601. (10) Soniat, M.; Pool, G.; Franklin, L.; Rick, S. W. Ion Association in Aqueous Solution. Fluid Phase Equilib. 2016, 407, 31. (11) Ferrario, M.; Ciccotti, G.; Spohr, E.; Cartailler, T.; Turq, P. Solubility of KF in Water by Molecular Dynamics Using the Kirkwood Integration Method. J. Chem. Phys. 2002, 117 (10), 4947−4953. (12) Lísal, M.; Smith, W. R.; Kolafa, J. Molecular Simulations of Aqueous Electrolyte Solubility: 1. The Expanded-Ensemble Osmotic Molecular Dynamics Method for the Solution Phase. J. Phys. Chem. B 2005, 109 (26), 12956−12965. (13) Sanz, E.; Vega, C. Solubility of KF and NaCl in Water by Molecular Simulation. J. Chem. Phys. 2007, 126 (1), 014507. (14) Aragones, J. L.; Sanz, E.; Vega, C. Solubility of NaCl in Water by Molecular Simulation Revisited. J. Chem. Phys. 2012, 136 (24), 244508. (15) Bian, H.; Wen, X.; Li, J.; Chen, H.; Han, S.; Sun, X.; Song, J.; Zhuang, W.; Zheng, J. Ion Clustering in Aqueous Solutions Probed with Vibrational Energy Transfer. Proc. Natl. Acad. Sci. U. S. A. 2011, 108 (12), 4737−4742. (16) Wang, F.; Richards, V. N.; Shields, S. P.; Buhro, W. E. Kinetics and Mechanisms of Aggregative Nanocrystal Growth. Chem. Mater. 2014, 26 (1), 5−21. (17) Hassan, S. A. Computer Simulation of Ion Cluster Speciation in Concentrated Aqueous Solutions at Ambient Conditions. J. Phys. Chem. B 2008, 112 (34), 10573−10584. (18) Hassan, S. A. Morphology of Ion Clusters in Aqueous Electrolytes. Phys. Rev. E 2008, 77 (3), 031501. (19) Zahn, D. Atomistic Mechanism of NaCl Nucleation from an Aqueous Solution. Phys. Rev. Lett. 2004, 92 (4), 040801. (20) Hassan, S. A. Microscopic Mechanism of Nanocrystal Formation from Solution by Cluster Aggregation and Coalescence. J. Chem. Phys. 2011, 134 (11), 114508. (21) Lanaro, G.; Patey, G. N. Molecular Dynamics Simulation of NaCl Dissolution. J. Phys. Chem. B 2015, 119 (11), 4275−4283. (22) Rao, F.; Caflisch, A. The Protein Folding Network. J. Mol. Biol. 2004, 342 (1), 299−306. (23) Noé, F.; Horenko, I.; Schütte, C.; Smith, J. C. Hierarchical Analysis of Conformational Dynamics in Biomolecules: Transition Networks of Metastable States. J. Chem. Phys. 2007, 126 (15), 155102. (24) Bowman, G. R.; Pande, V. S. Protein Folded States are Kinetic Hubs. Proc. Natl. Acad. Sci. U. S. A. 2010, 107 (24), 10890−10895. (25) Buch, I.; Giorgino, T.; De Fabritiis, G. Complete Reconstruction of an Enzyme-Inhibitor Binding Process by Molecular Dynamics Simulations. Proc. Natl. Acad. Sci. U. S. A. 2011, 108 (25), 10184− 10189. (26) Rao, F.; Garrett-Roe, S.; Hamm, P. Structural Inhomogeneity of Water by Complex Network Analysis. J. Phys. Chem. B 2010, 114 (47), 15598−604. (27) Prada-Gracia, D.; Shevchuk, R.; Hamm, P.; Rao, F. Towards a Microscopic Description of the Free-Energy Landscape of Water. J. Chem. Phys. 2012, 137 (14), 144504. (28) Shevchuk, R.; Agmon, N.; Rao, F. Network Analysis of Proton Transfer in Liquid Water. J. Chem. Phys. 2014, 140 (24), 244502. (29) Perkett, M. R.; Hagan, M. F. Using Markov State Models to Study Self-Assembly. J. Chem. Phys. 2014, 140 (21), 214101. (30) Bowman, G. R.; Pande, V. S.; Noé, F. An Introduction to Markov State Models and Their Application to Long Timescale Molecular Simulation; Springer Science & Business Media: Berlin, 2013; Vol. 797. (31) Berendsen, H. J. C.; van der Spoel, D.; van Drunen, R. GROMACS: A Message-Passing Parallel Molecular Dynamics Implementation. Comput. Phys. Commun. 1995, 91 (1), 43−56. (32) Berendsen, H. J. C.; Postma, J.; van Gunsteren, W. F.; DiNola, A.; Haak, J. Molecular Dynamics with Coupling to an External Bath. J. Chem. Phys. 1984, 81 (8), 3684−3690.

(33) Bussi, G.; Donadio, D.; Parrinello, M. Canonical Sampling Through Velocity Rescaling. J. Chem. Phys. 2007, 126 (1), 014101. (34) Darden, T.; York, D.; Pedersen, L. Particle Mesh Ewald: An N log (N) Method for Ewald Sums in Large Systems. J. Chem. Phys. 1993, 98 (12), 10089−10092. (35) Abascal, J. L. F.; Vega, C. A General Purpose Model for the Condensed Phases of Water: TIP4P/2005. J. Chem. Phys. 2005, 123 (23), 234505. (36) Duan, Y.; Wu, C.; Chowdhury, S.; Lee, M. C.; Xiong, G.; Zhang, W.; Yang, R.; Cieplak, P.; Luo, R.; Lee, T.; et al. A Point-Charge Force Field for Molecular Mechanics Simulations of Proteins Based on Condensed-Phase Quantum Mechanical Calculations. J. Comput. Chem. 2003, 24 (16), 1999−2012. (37) Lorentz, H. Ü ber die Anwendung des Satzes vom Virial in der kinetischen Theorie der Gase. Ann. Phys. 1881, 248 (1), 127−136. (38) Berthelot, D. Sur le Mélange des Gaz. CR Chim. 1898, 126, 1703−1706. (39) Allen, M. P.; Tildesley, D. J.; et al. Computer Simulation of Liquids; Oxford Science Publishers: Oxford, U.K., 1987. (40) Auffinger, P.; Cheatham, T. E.; Vaiana, A. C. Spontaneous Formation of KCl Aggregates in Biomolecular Simulations: a force field issue? J. Chem. Theory Comput. 2007, 3 (5), 1851−1859. (41) Steinhardt, P. J.; Nelson, D. R.; Ronchetti, M. BondOrientational Order in Liquids and Glasses. Phys. Rev. B: Condens. Matter Mater. Phys. 1983, 28 (2), 784. (42) Prinz, J.-H.; Wu, H.; Sarich, M.; Keller, B.; Senne, M.; Held, M.; Chodera, J. D.; Schütte, C.; Noé, F. Markov Models of Molecular Kinetics: Generation and Validation. J. Chem. Phys. 2011, 134 (17), 174105. (43) Noé, F.; Fischer, S. Transition Networks for Modeling the Kinetics of Conformational Change in Macromolecules. Curr. Opin. Struct. Biol. 2008, 18 (2), 154−162. (44) Van Dongen, S. Graph Clustering by Flow Simulation. Ph.D Thesis, University of Utrecht, 2000. (45) Smirnov, N. V. On the Estimation of the Discrepancy Between Empirical Curves of Distribution for Two Independent Samples. Bull. Math. Univ. Moscou 2 (2), 1939. (46) Berezovska, G.; Prada-Gracia, D.; Mostarda, S.; Rao, F. Accounting for the Kinetics in Order Parameter Analysis: Lessons from Theoretical Models and a Disordered Peptide. J. Chem. Phys. 2012, 137 (19), 194101. (47) Gfeller, D.; Chappelier, J.-C.; De Los Rios, P. Finding Instabilities in the Community Structure of Complex Networks. Phys. Rev. E 2005, 72 (5), 056135. (48) Gfeller, D.; De Los Rios, P.; Caflisch, A.; Rao, F. Complex Network Analysis of Free-Energy Landscapes. Proc. Natl. Acad. Sci. U. S. A. 2007, 104 (6), 1817−1822. (49) Vlasblom, J.; Wodak, S. J. Markov Clustering Versus Affinity Propagation for the Partitioning of Protein Interaction Graphs. BMC Bioinf. 2009, 10 (1), 99. (50) Senne, M.; Trendelkamp-Schroer, B.; Mey, A. S. J. S.; Schütte, C.; Noé, F. EMMA: a Software Package for Markov Model Building and Analysis. J. Chem. Theory Comput. 2012, 8 (7), 2223−2238. (51) Beckham, G. T.; Peters, B. Optimizing Nucleus Size Metrics for Liquid-Solid Nucleation from Transition Paths of Near-Nanosecond Duration. J. Phys. Chem. Lett. 2011, 2 (10), 1133−1138.

2788


Markov State Model of Ion Assembling Process - The Journal of

Recommend Documents