Is Protein Single Molecule Dynamics under Force Described by Two or More States?

<p>The use of computational simulation to study the dynamics and interactions of macromolecules has become an important tool in the field of biochemistry. A common method to perform these simulations is to use all-atom explicit-solvent molecular dynamics (MD). However, due to the limitations in computational power currently available, this method is not practical for simulating large-scale biomolecular systems on long timescales. An alternative is to perform implicit-solvent Brownian dynamics (BD) simulations using a coarse grained (CG) model that allows for increased computational efficiency. However, if simulations using the CG model are not realistic, then the gain in computational efficiency from using a CG model is not worthwhile.</p> <p>This thesis describes the derivation of a set of bonded and nonbonded CG potential functions for use in implicit-solvent BD simulations of proteins derived from all-atom explicit-solvent MD simulations of amino acids. To determine which force field and water model to use in the MD simulations, Chapter II describes 1 Μs all-atom explicit-solvent MD simulations of glycine, asparagine, phenylalanine, and valine solutions at 50, 100, 200 and 300 mg/ml concentrations performed using eight different force field and water model combinations. To evaluate the accuracy of the force fields at high solute concentrations, the density, viscosity, and dielectric increments of the four amino acids were calculated from the simulations and compared to experimental results. Additionally, the change in the strength of hydrophobic and electrostatic interactions with increasing solute concentration was calculated for each force field and water model combination. As a result of this study, the Amber ff99SB-ILDN force field and TIP4P-Ew explicit-solvent water model were chosen for all subsequent MD simulations. Chapter III describes the derivation of CG bonded potential functions from 1 Μs all-atom explicit-solvent MD simulations of each of the twenty amino acids, including a separate simulation for protonated histidine. The angle and dihedral probability distributions sampled during the MD simulations were used to optimize the bonded potential functions using the iterative Boltzmann inversion (IBI) method. Chapter IV describes the derivation of CG nonbonded potential functions from 1 Μs all-atom explicit-solvent MD simulations of every possible pairing of the amino acids (231 different systems). The radial distribution functions calculated from these MD simulations were used to optimize a set of nonbonded CG potential functions using the IBI method. The optimized set of bonded and nonbonded potential functions, which is termed COFFDROP (COarse-grained Force Field for Dynamic Representation Of Proteins), quantitatively reproduced all of the calculated MD distributions. To determine if COFFDROP would be useful for simulations of bimolecular systems, Chapter V describes the testing of the transferability of the force field. First, COFFDROP was used to simulate concentrated amino acid solutions. The clustering of the solutes in these simulations was directly compared with results from corresponding all-atom explicit-solvent MD simulations and found to be in excellent agreement. Next, BD simulations of 9.2 mM solutions of the small protein villin headpiece were performed. The proteins aggregated during these simulations, which is in agreement with results from MD simulation but in disagreement with experiment. After scaling the strength of COFFDROP's nonbonded potential functions by a factor of 0.8 and rerunning the BD simulations, the amount of aggregation was comparable to experimental observations. Based on these results, COFFDROP is likely to be applicable in CG BD simulations of large, highly concentrated, biomolecular systems.</p>

Differential scanning calorimetry was used to measure the temperature dependence of the absolute heat capacity of the 35-residue subdomain of the villin headpiece, a protein that folds in 5 mus and is therefore assumed to have a small free-energy barrier separating folded and unfolded states. To obtain an estimate of the barrier height from the calorimetric data, two models, a variable-barrier model and an Ising-like model, were used to fit the heat capacity in excess of the folded state over the temperature range 15-125 degrees C. The variable-barrier model is based on an empirical mathematical form for the density of states, with four adjustable parameters and the enthalpy (H) as a reaction coordinate. The Ising-like model is based on the inter-residue contact map of the X-ray structure with exact enumeration of approximately 10(5) possible conformations, with two adjustable parameters in the partition function, and either the fraction of native contacts (Q) or the number of ordered residues (P) as reaction coordinates. The variable-barrier model provides an excellent fit to the data and yields a barrier height at the folding temperature ranging from 0.4 to 1.1 kcal mol(-1), while the Ising-like model provides a less good fit and yields barrier heights of 2.3 +/- 0.1 kcal mol(-1) and 2.1 +/- 0.1 kcal mol(-1) for the Q and P reaction coordinates, respectively. In both models, the barrier to folding increases with increasing temperature. Assuming a sufficiently large activation energy for diffusion on the free-energy surfaces, both models are consistent with the observation of a temperature-independent folding rate in previously published laser temperature-jump experiments. Analysis of this kinetic data, using an approximate form for the pre-exponential factor of Kramers theory and the 70 ns relaxation time for the fast phase that precedes the unfolding/refolding relaxation to determine the diffusion coefficient, results in a barrier height of 1.6 +/- 0.3 kcal mol-1 for an unspecified reaction coordinate. Although no independent test of the validity of the H, Q, or P reaction coordinates is given, the barrier-height estimates obtained with the three reaction coordinates are in quite good agreement with the value derived from a Kramers analysis of the kinetics that makes no assumptions about the reaction coordinate. However, the higher estimates obtained using Q or P appear more consistent with the finding of barrier-crossing kinetics of a villin mutant that folds in 700 ns, corresponding to a 1.3 kcal mol-1 reduction in the folding barrier relative to wild-type. All of the results suggest that the free-energy barrier to folding is sufficiently low that it should be possible to engineer this protein or find solution conditions that would eliminate the barrier to create the "downhill" folding scenario of Wolynes and Onuchic.

Adiabatic free energy dynamics (AFED) was introduced by Rosso et al. [J. Chem. Phys. 2002, 116, 4389] for computing free energy profiles quickly and accurately using a dynamical adiabatic separation between a set of collective variables or reaction coordinates and the remaining degrees of freedom of a system. This approach has been shown to lead to a significant gain in efficiency versus traditional methods such as umbrella sampling, thermodynamic integration, and free energy perturbation for generating one-dimensional free energy profiles. More importantly, AFED is able to generate multidimensional free energy surfaces efficiently via full sweeps of the surface that rapidly map out the locations of the free energy minima. The most significant drawback to the AFED approach is the need to transform the coordinates into a generalized coordinate system that explicitly contains the collective variables of interest. Recently, Maragliano and Vanden-Eijnden built upon the AFED approach by introducing a set of extended phase-space variables, to which the adiabatic decoupling and high temperature are applied [Chem. Phys. Lett. 2006, 426, 168]. In this scheme, which the authors termed "temperature accelerated molecular dynamics" or TAMD, the need for explicit coordinate transformations is circumvented. The ability of AFED and TAMD to generate free energy surfaces efficiently depends on the thermostatting mechanism employed, since both approaches are inherently nonequilibrium due to the adiabatic decoupling. Indeed, Maragliano and Vanden-Eijnden did not report any direct generation of free energy surfaces within the overdamped Langevin dynamics employed by these authors. Here, we show that by formulating TAMD in a manner that is closer to the original AFED approach, including the generalized Gaussian moment thermostat (GGMT) and multiple time-scale integration, multidimensional free energy surfaces for complex systems can be generated directly from the probability distribution function of the extended phase-space variables. The new TAMD formulation, which we term driven AFED or d-AFED, is applied to compare the conformational preferences of small peptides both in gas phase and in solution for three force fields. The results show that d-AFED/TAMD accurately and efficiently generates free energy surfaces in two collective variables useful for characterizing the conformations, namely, the radius of gyration, R(G), and number of hydrogen bonds, N(H).

The long time self-diffusion of highly charged colloidal particles consisting of polybutylacrylate was studied employing pulsed field gradient NMR (PFGNMR) and Brownian dynamics (BD) simulations. The potential parameters for the Yukawa interaction potential used in the BD simulations were derived from the static structure factor measured by means of light scattering. The long time self-diffusion coefficients resulting from NMR and BD were compared for different ionic strengths, which cause a screening of the interaction potential due to stray ions. Increasing the ionic strength of the colloidal suspensions, we observed a phase transition from a crystalline to a liquidlike structure. Whereas at high salt concentrations both long time diffusion coefficients agree within a few percent, at intermediate ionic strengths the values measured via PFGNMR are by a factor of 2 higher than those predicted by the BD simulation. This deviation is explained by hydrodynamic interactions which are expected to enhance long time self-diffusion and which are not taken into account in our BD algorithm.

Brownian Dynamics Simulations of Aldolase Binding Glyceraldehyde 3-Phosphate Dehydrogenase and the Possibility of Substrate Channeling

Multiscale Modeling of Microtubules and Actin Filaments

We carried out ab initio molecular dynamic simulations in order to determine the free energy surfaces of two selected reactions including solvents, namely a rearrangement of a ruthenium oxoester in water and a carbon dioxide addition to a palladium complex in carbon dioxide. For the latter reaction we also investigated the gas phase reaction in order to take solvent effects into account. We used two techniques to reconstruct the free energy surfaces: thermodynamic integration and metadynamics. Furthermore, we gave a reasonable error estimation of the computed free energy surface. We calculated a reaction barrier of ΔF = 59.5 ± 8.5 kJ mol−1 for the rearrangement of a ruthenium oxoester in water from thermodynamic integration. For the carbon dioxide addition to the palladium complex in carbon dioxide we found a ΔF = 44.9 ± 3.3 kJ mol−1 from metadynamics simulations with one collective variable. The investigation of the same reactions in the gas phase resulted in ΔF = 24.9 ± 6.7 kJ mol−1 from thermodynamic integration, in ΔF = 26.7 ± 2.3 kJ mol−1 from metadynamics simulations with one collective variable, and in ΔF = 27.1 ± 5.9 kJ mol−1 from metadynamics simulations with two collective variables.

In this work, we demonstrate that Linear Discriminant Analysis (LDA) applied to atomic positions in two different states of a biomolecule produces a good reaction coordinate between those two states. Atomic coordinates of a macromolecule are a direct representation of a macromolecular configuration, and yet, they are not used in enhanced sampling studies due to a lack of rotational and translational invariance. We resolve this issue using the technique of our prior work, whereby a molecular configuration is considered a member of an equivalence class in size-and-shape space, which is the set of all configurations that can be translated and rotated to a single point within a reference multivariate Gaussian distribution characterizing a single molecular state. The reaction coordinates produced by LDA applied to positions are shown to be good reaction coordinates both in terms of characterizing the transition between two states of a system within a long molecular dynamics (MD) simulation and also ones that allow us to readily produce free energy estimates along that reaction coordinate using enhanced sampling MD techniques.

All-atom Simulations of Structures and Energetics of c-di-GMP-bound and free PleD

The effect of hydrodynamic interactions (HI) on the long-time self-diffusion in quasi-two-dimensional fluids of paramagnetic colloidal particles is investigated using a combination of experiments and Brownian dynamics (BD) simulations. In the BD simulations, the direct interactions (DI) between the particles consist of a short-ranged repulsive part and a long-ranged part that is proportional to 1/r^{3}, with r the interparticle distance. By studying the equation of state, the simulations allow for the identification of the regime where the properties of the fluid are fully controlled by the long-ranged interactions, and the thermodynamic state solely depends on the dimensionless interaction strength Γ. In this regime, the radial distribution functions from the simulations are in quantitative agreement with those from the experiments for different fluid area fractions. This agreement confirms that the DI in the experiments and simulations are identical, which thus allows us to isolate the role of HI, as these are not taken into account in the BD simulations. Experiment and simulation fall onto a master curve with respect to the Γ dependence of D_{L}^{★}=D_{L}/(D_{0}Γ^{1/2}), with D_{0} the self-diffusion coefficient at infinite dilution and D_{L} the long-time self-diffusion coefficient. Our results thus show that, although HI affect the short-time self-diffusion, for a quasi-two-dimensional system with 1/r^{3} long-ranged DI, the reduced quantity D_{L}^{★} is effectively not affected by HI. Interestingly, this is in agreement with prior work on quasi-two-dimensional colloidal hard spheres.

We use a combination of Brownian dynamics (BD) simulation results and deep learning (DL) strategies for the rapid identification of large structural changes caused by missense mutations in intrinsically disordered proteins (IDPs). We used ∼6500 IDP sequences from MobiDB database of length 20-300 to obtain gyration radii from BD simulation on a coarse-grained single-bead amino acid model (HPS2 model) used by us and others [Dignon, G. L. PLoS Comput. Biol. 2018, 14, e1005941,Tesei, G. Proc. Natl. Acad. Sci. U.S.A. 2021, 118, e2111696118,Seth, S. J. Chem. Phys. 2024, 160, 014902] to generate the training sets for the DL algorithm. Using the gyration radii ⟨Rg⟩ of the simulated IDPs as the training set, we develop a multilayer perceptron neural net (NN) architecture that predicts the gyration radii of 33 IDPs previously studied by using BD simulation with 97% accuracy from the sequence and the corresponding parameters from the HPS model. We now utilize this NN to predict gyration radii of every permutation of missense mutations in IDPs. Our approach successfully identifies mutation-prone regions that induce significant alterations in the radius of gyration when compared to the wild-type IDP sequence. We further validate the prediction by running BD simulations on the subset of identified mutants. The neural network yields a (104-106)-fold faster computation in the search space for potentially harmful mutations. Our findings have substantial implications for rapid identification and understanding of diseases related to missense mutations in IDPs and for the development of potential therapeutic interventions. The method can be extended to accurate predictions of other mutation effects in disordered proteins.

We have elucidated the feasibility of trapping magnetic rod-like particles that flow in a cylindrical pipe in terms of a non-uniform magnetic field. We have performed Brownian dynamics (BD) simulations in order to clarify the influence that various factors have on the attachment or trapping characteristics of magnetic rod-like particles to the poles. The employment of a local linearisation approximation enables us to apply the BD method that was fully-developed for a simple shear flow to the present Hagen-Poiseuille flow. In the case of a sufficiently small magnetic interaction between the particles, numerous rod-like particles may attach to the magnetic poles if the influence of a non-uniform magnetic field is sufficiently more dominant than that of the flow field, otherwise, the particles flow away along the centreline of the cylindrical pipe. If the magnetic interactions are sufficiently more dominant than both the flow field and the non-uniform magnetic field, the particles aggregate to form densely-packed clusters that flow away along the centreline area of the pipe without being trapped at the poles. In the situation of insufficient magnetic interactions, thin chain-like clusters may be formed from several anchored particles trapped to the poles and incline in the downstream direction.

Absolute Protein-Protein Association Rate Constants from Flexible, Coarse-Grained Brownian Dynamics Simulations: The Role of Intermolecular Hydrodynamic Interactions in Barnase-Barstar Association

Calculating the kinetic rates for rare transitions is an important objective for molecular simulations. Here I prove equalities that relate the transition rates to the equilibrium free energy and the statistics of the transition paths. In particular, the durations of the transition paths within given intervals of the reaction coordinate provide the kinetic pre-factor in the rate formula. Based on the available free energy, the transition rates can further be rigorously calculated by initiating forward and backward simulations and evaluating the duration of each transition path. Validation on a model system confirms that the approach correctly predicts the transition rates from the simulations and demonstrates that whereas the relations here are general and valid for any chosen reaction coordinate, a good reaction coordinate will enable a more efficient sampling of the transition paths and thus a more reliable rate calculation.

Understanding how protein folds into a functional and structural configuration is arguably one of themost important and challenging problems in computational biology. Currently, the protein foldingmechanism is often characterized by calculating the free energy landscape versus the reaction coordi-nates such as the fraction of native contacts, the radius of gyration, the principal components and soon. In this paper, we present a combinatorial algorithmic approach towards understanding the globalstate changes of the configurations. The approach is based on cluster computation, each cluster beingdefined by a pattern of a combination of various reaction coordinates. We present an algorithm oftime complexity