Efficient thermal rate constant calculation for rare event systems

A skewed-momenta method to efficiently generate conformational-transition trajectories

We present a novel computational method, the skewed-momenta method (Skew'M), which applies a bias to the Maxwell distribution of initial momenta used to generate ensembles of trajectories. As a result, conformational transitions are accentuated and kinetic properties are calculated more effectively. The connection to the related puddle jumping method is discussed. A reweighting scheme permits the exact calculation of kinetic properties. Applications are presented for the rapid calculation of rate constants for molecular isomerization, and for the efficient reconstruction of free-energy profiles using a straightforward modification of the Jarzynski identity.

Efficiency analysis of reaction rate calculation methods using analytical models I: The two-dimensional sharp barrier

We analyze the efficiency of different methods for the calculation of reaction rates in the case of a simple two-dimensional analytical benchmark system. Two classes of methods are considered: the first is based on the free energy calculation along a reaction coordinate and the calculation of the transmission coefficient, the second on the sampling of dynamical pathways. We give scaling rules for how this efficiency depends on barrier height and width, and we hand out simple optimization rules for the method-specific parameters. We show that the path sampling methods, using the transition interface sampling technique, become exceedingly more efficient than the others when the reaction coordinate is not the optimal one.

Improved importance sampling distribution for rate constant calculation

An efficient method to compute the thermal rate constant for rare events within the correlation function C(t) approach is presented. This method, which is based on a modification of the sampling function used to evaluate the dynamical correlation function C(t), can be applied to high-dimensional systems having a rough energy landscape without previous knowledge on the transition states location. In this work, the sampling of a Boltzmann-like distribution for the linear momenta with a lower inverse temperature (beta(*)=1/kT(*)) than the correct one (beta=1/kT) is proposed as a way to increase the number of reactive trajectories. The mismatch between the beta(*) and beta distributions is then corrected by a reweighting procedure which allows one to obtain the exact correlation function C(t). The efficiency of this scheme in computing the rate of a particle jumping across the barrier of a simple 2D double well potential is improved by a factor 4-25 depending on the relative value of the original beta and modified beta(*) temperatures. When merged with the "puddle potential" method [S. A. Corcelli, J. A. Rohman, and J. C. Tully, J. Chem. Phys., 118, 1085 (2003)], the new importance sampling function improves the efficiency of the correlation function approach by a factor 16-800 with respect to the unbiased sampling. To test the method in a more challenging case, the previous model system was extended by adding six harmonically restrained particles, each one interacting with the diffusing particle. This model introduces both the possibility of energy exchange and a rougher energy landscape. The new sampling function alone is found to produce an improvement in efficiency of, at least, an order of magnitude when compared with the unbiased case; when merged with the puddle potential method, a 400-fold saving in computer time is found.

Tools for channels: moving towards molecular calculations of gating and permeation in ion channel biophysics

Recent X-ray structures of voltage gated potassium channels provide an exciting opportunity to connect molecular structures with measured biological function. Two of the most important connections for these channels are: first, to the molecular basis behind selectivity and the associated free energy profile underlying ionic current flow and, second, to a true molecular understanding of the large-scale conformational transitions that underlie voltage dependent gating. But, existing computational tools need to be further developed to reach these goals. In this contribution to the symposia on sampling methods we outline our dynamic importance sampling method for sampling large-scale conformational transitions as well as our studies with non-equilibrium work events and equilibrium overlap sampling (OS) methods for sampling events related to the calculation of relative free energies.

From transition paths to transition states and rate coefficients

Transition states are defined as points in configuration space with the highest probability that trajectories passing through them are reactive (i.e., form transition paths between reactants and products). In the high-friction (diffusive) limit of Langevin dynamics, the resulting ensemble of transition states is shown to coincide with the separatrix formed by points of equal commitment (or splitting) probabilities for reaching the product and reactant regions. Transition states according to the new criterion can be identified directly from equilibrium trajectories, or indirectly by calculating probability densities in the equilibrium and transition-path ensembles using umbrella and transition-path sampling, respectively. An algorithm is proposed to calculate rate coefficients from the transition-path and equilibrium ensembles by estimating the frequency of transitions between reactants and products.