Search

The proper balancing of information from experiment and theory is a long-standing problem in the analysis of noisy and incomplete data. Viewed as a Pareto optimization problem, improved agreement with the experimental data comes at the expense of growing inconsistencies with the theoretical reference model. Here, we propose how to set the exchange rate a priori to properly balance this trade-off. We focus on gentle ensemble refinement, where the difference between the potential energy surfaces of the reference and refined models is small on a thermal scale. By relating the variance of this energy difference to the Kullback–Leibler divergence between the respective Boltzmann distributions, one can encode prior knowledge about energy uncertainties, i.e., force-field errors, in the exchange rate. The energy uncertainty is defined in the space of observables and depends on their type and number and on the thermodynamic state. We highlight the relation of gentle refinement to free energy perturbation theory. A balanced encoding of prior knowledge increases the quality and transparency of ensemble refinement. Our findings extend to non-Boltzmann distributions, where the uncertainty in energy becomes an uncertainty in information. [ABSTRACT FROM AUTHOR]

Proton-powered c-ring rotation in mitochondrial ATP synthase is crucial to convert the transmembrane protonmotive force into torque to drive the synthesis of adenosine triphosphate (ATP). Capitalizing on recent cryo-EM structures, we aim at a structural and energetic understanding of how functional directional rotation is achieved. We performed multi-microsecond atomistic simulations to determine the free energy profiles along the c-ring rotation angle before and after the arrival of a new proton. Our results reveal that rotation proceeds by dynamic sliding of the ring over the a-subunit surface, during which interactions with conserved polar residues stabilize distinct intermediates. Ordered water chains line up for a Grotthuss-type proton transfer in one of these intermediates. After proton transfer, a high barrier prevents backward rotation and an overall drop in free energy favors forward rotation, ensuring the directionality of c-ring rotation required for the thermodynamically disfavored ATP synthesis. The essential arginine of the a-subunit stabilizes the rotated configuration through a salt bridge with the c-ring. Overall, we describe a complete mechanism for the rotation step of the ATP synthase rotor, thereby illuminating a process critical to all life at atomic resolution. [ABSTRACT FROM AUTHOR]

Single-molecule localization microscopy allows practitioners to locate and track labeled molecules in biological systems. When extracting diffusion coefficients from the resulting trajectories, it is common practice to perform a linear fit on mean-squared-displacement curves. However, this strategy is suboptimal and prone to errors. Recently, it was shown that the increments between the observed positions provide a good estimate for the diffusion coefficient, and their statistics are well-suited for likelihood-based analysis methods. Here, we revisit the problem of extracting diffusion coefficients from single-particle tracking experiments subject to static noise and dynamic motion blur using the principle of maximum likelihood. Taking advantage of an efficient real-space formulation, we extend the model to mixtures of subpopulations differing in their diffusion coefficients, which we estimate with the help of the expectation–maximization algorithm. This formulation naturally leads to a probabilistic assignment of trajectories to subpopulations. We employ the theory to analyze experimental tracking data that cannot be explained with a single diffusion coefficient. We test how well a dataset conforms to the assumptions of a diffusion model and determine the optimal number of subpopulations with the help of a quality factor of known analytical distribution. To facilitate use by practitioners, we provide a fast open-source implementation of the theory for the efficient analysis of multiple trajectories in arbitrary dimensions simultaneously. [ABSTRACT FROM AUTHOR]

Gasdermin-D (GSDMD) is the ultimate effector of pyroptosis, a form of programmed cell death associated with pathogen invasion and inflammation. After proteolytic cleavage by caspases, the GSDMD N-terminal domain (GSDMDNT) assembles on the inner leaflet of the plasma membrane and induces the formation of membrane pores. We use atomistic molecular dynamics simulations to study GSDMDNT monomers, oligomers, and rings in an asymmetric plasma membrane mimetic. We identify distinct interaction motifs of GSDMDNT with phosphatidylinositol-4,5-bisphosphate (PI(4,5)P2) and phosphatidylserine (PS) headgroups and describe their conformational dependence. Oligomers are stabilized by shared lipid binding sites between neighboring monomers acting akin to double-sided tape. We show that already small GSDMDNT oligomers support stable, water-filled, and ion-conducting membrane pores bounded by curled beta-sheets. In large-scale simulations, we resolve the process of pore formation from GSDMDNT arcs and lipid efflux from partial rings. We find that high-order GSDMDNT oligomers can crack under the line tension of 86 pN created by an open membrane edge to form the slit pores or closed GSDMDNT rings seen in atomic force microscopy experiments. Our simulations provide a detailed view of key steps in GSDMDNT-induced plasma membrane pore formation, including sublytic pores that explain nonselective ion flux during early pyroptosis. [ABSTRACT FROM AUTHOR]

The demands on the accuracy of force fields for classical molecular dynamics simulations are steadily growing as larger and more complex systems are studied over longer times. One way to meet these growing demands is to hand over the learning of force fields and their parameters to machines in a systematic (semi)automatic manner. Doing so, we can take full advantage of exascale computing, the increasing availability of experimental data, and advances in quantum mechanical computations and the calculation of experimental observables from molecular ensembles. Here, we discuss and illustrate the challenges one faces in this endeavor and explore a way forward by adapting the Bayesian inference of ensembles (BioEn) method [Hummer and Köfinger, J. Chem. Phys. (2015)] for force field parameterization. In the Bayesian inference of force fields (BioFF) method developed here, the optimization problem is regularized by a simplified prior on the force field parameters and an entropic prior acting on the ensemble. The latter compensates for the unavoidable over simplifications in the parameter prior. We determine optimal force field parameters using an iterative predictor–corrector approach, in which we run simulations, determine the reference ensemble using the weighted histogram analysis method (WHAM), and update the force field according to the BioFF posterior. We illustrate this approach for a simple polymer model, using the distance between two labeled sites as the experimental observable. By systematically resolving force field issues, instead of just reweighting a structural ensemble, the BioFF corrections extend to observables not included in ensemble reweighting. We envision future force field optimization as a formalized, systematic, and (semi)automatic machine-learning effort that incorporates a wide range of data from experiment and high-level quantum chemical calculations, and takes advantage of exascale computing resources. [ABSTRACT FROM AUTHOR]

We describe different Bayesian ensemble refinement methods, examine their interrelation, and discuss their practical application. With ensemble refinement, the properties of dynamic and partially disordered (bio)molecular structures can be characterized by integrating a wide range of experimental data, including measurements of ensemble-averaged observables. We start from a Bayesian formulation in which the posterior is a functional that ranks different configuration space distributions. By maximizing this posterior, we derive an optimal Bayesian ensemble distribution. For discrete configurations, this optimal distribution is identical to that obtained by the maximum entropy "ensemble refinement of SAXS" (EROS) formulation. Bayesian replica ensemble refinement enhances the sampling of relevant configurations by imposing restraints on averages of observables in coupled replica molecular dynamics simulations. We show that the strength of the restraints should scale linearly with the number of replicas to ensure convergence to the optimal Bayesian result in the limit of infinitely many replicas. In the "Bayesian inference of ensembles" method, we combine the replica and EROS approaches to accelerate the convergence. An adaptive algorithm can be used to sample directly from the optimal ensemble, without replicas. We discuss the incorporation of single-molecule measurements and dynamic observables such as relaxation parameters. The theoretical analysis of different Bayesian ensemble refinement approaches provides a basis for practical applications and a starting point for further investigations. [ABSTRACT FROM AUTHOR]

Binding of the spike protein of SARS-CoV-2 to the human angiotensinconverting enzyme 2 (ACE2) receptor triggers translocation of the virus into cells. Both the ACE2 receptor and the spike protein are heavily glycosylated, including at sites near their binding interface. We built fully glycosylated models of the ACE2 receptor bound to the receptor binding domain (RBD) of the SARS-CoV-2 spike protein. Using atomistic molecular dynamics (MD) simulations, we found that the glycosylation of the human ACE2 receptor contributes substantially to the binding of the virus. Interestingly, the glycans at two glycosylation sites, N90 and N322, have opposite effects on spike protein binding. The glycan at the N90 site partly covers the binding interface of the spike RBD. Therefore, this glycan can interfere with the binding of the spike protein and protect against docking of the virus to the cell. By contrast, the glycan at the N322 site interacts tightly with the RBD of the ACE2-bound spike protein and strengthens the complex. Remarkably, the N322 glycan binds to a conserved region of the spike protein identified previously as a cryptic epitope for a neutralizing antibody. By mapping the glycan binding sites, our MD simulations aid in the targeted development of neutralizing antibodies and SARS-CoV-2 fusion inhibitors. [ABSTRACT FROM AUTHOR]

We derive simple analytical expressions for the error and computational efficiency of simulated tempering (ST) simulations. The theory applies to the important case of systems whose dynamics at long times is dominated by the slow interconversion between two metastable states. An extension to the multistate case is described. We show that the relative gain in efficiency of ST simulations over regular molecular dynamics (MD) or Monte Carlo (MC) simulations is given by the ratio of their reactive fluxes, i.e., the number of transitions between the two states summed over all ST temperatures divided by the number of transitions at the single temperature of the MD or MC simulation. This relation for the efficiency is derived for the limit in which changes in the ST temperature are fast compared to the two-state transitions. In this limit, ST is most efficient. Our expression for the maximum efficiency gain of ST simulations is essentially identical to the corresponding expression derived by us for replica exchange MD and MC simulations [E. Rosta and G. Hummer, J. Chem. Phys. 131, 165102 (2009)] on a different route. We find quantitative agreement between predicted and observed efficiency gains in a test against ST and replica exchange MC simulations of a two-dimensional Ising model. Based on the efficiency formula, we provide recommendations for the optimal choice of ST simulation parameters, in particular, the range and number of temperatures, and the frequency of attempted temperature changes. [ABSTRACT FROM AUTHOR]

We derive simple analytical expressions for the error and computational efficiency of replica exchange molecular dynamics (REMD) simulations (and by analogy replica exchange Monte Carlo simulations). The theory applies to the important case of systems whose dynamics at long times is dominated by the slow interconversion between two metastable states. As a specific example, we consider the folding and unfolding of a protein. The efficiency is defined as the rate with which the error in an estimated equilibrium property, as measured by the variance of the estimator over repeated simulations, decreases with simulation time. For two-state systems, this rate is in general independent of the particular property. Our main result is that, with comparable computational resources used, the relative efficiency of REMD and molecular dynamics (MD) simulations is given by the ratio of the number of transitions between the two states averaged over all replicas at the different temperatures, and the number of transitions at the single temperature of the MD run. This formula applies if replica exchange is frequent, as compared to the transition times. High efficiency of REMD is thus achieved by including replica temperatures in which the frequency of transitions is higher than that at the temperature of interest. In tests of the expressions for the error in the estimator, computational efficiency, and the rate of equilibration we find quantitative agreement with the results both from kinetic models of REMD and from actual all-atom simulations of the folding of a peptide in water. [ABSTRACT FROM AUTHOR]

We describe a reverse integration approach for the exploration of low-dimensional effective potential landscapes. Coarse reverse integration initialized on a ring of coarse states enables efficient navigation on the landscape terrain: Escape from local effective potential wells, detection of saddle points, and identification of significant transition paths between wells. We consider several distinct ring evolution modes: Backward stepping in time, solution arc length, and effective potential. The performance of these approaches is illustrated for a deterministic problem where the energy landscape is known explicitly. Reverse ring integration is then applied to noisy problems where the ring integration routine serves as an outer wrapper around a forward-in-time inner simulator. Two versions of such inner simulators are considered: A Gillespie-type stochastic simulator and a molecular dynamics simulator. In these “equation-free” computational illustrations, estimation techniques are applied to the results of short bursts of inner simulation to obtain the unavailable (in closed-form) quantities (local drift and diffusion coefficient estimates) required for reverse ring integration; this naturally leads to approximations of the effective landscape. [ABSTRACT FROM AUTHOR]

The reactive flux between folded and unfolded states of a two-state protein, whose coarse-grained dynamics is described by a master equation, is expressed in terms of the commitment or splitting probabilities of the microstates in the bottleneck region. This allows one to determine how much each transition through a dividing surface contributes to the reactive flux. By repeating the analysis for a series of dividing surfaces or, alternatively, by partitioning the reactive flux into contributions of unidirectional pathways that connect reactants and products, insight can be gained into the mechanism of protein folding. Our results for the flux in a network with complex connectivity, obtained using the discrete counterpart of Kramers’ theory of activated rate processes, show that the number of reactive transitions is typically much smaller than the total number of transitions that cross a dividing surface at equilibrium. [ABSTRACT FROM AUTHOR]

We present a recently developed one-dimensional dipole lattice model that accurately captures the key properties of water in narrow nanopores. For this model, we derive three equivalent representations of the Hamiltonian that together yield a transparent physical picture of the energetics of the water chain and permit efficient computer simulations. In the charge representation, the Hamiltonian consists of nearest-neighbor interactions and Coulomb-like interactions of effective charges at the ends of dipole ordered segments. Approximations based on the charge picture shed light on the influence of the Coulomb-like interactions on the structure of nanopore water. We use Monte Carlo simulations to study the system behavior of the full Hamiltonian and its approximations as a function of chemical potential and system size and investigate the bimodal character of the density distribution occurring at small system sizes. [ABSTRACT FROM AUTHOR]

Molecular dynamics simulations on quantum energy surfaces are carried out to study the effects of perturbing electric fields on proton transport (PT) in protonated water chains. As an idealized model of a hydrophobic cavity in the interior of a protein the water molecules are confined into a carbon nanotube (CNT). The water chain connects a hydrated hydronium ion (H3O+) at one end of the CNT and an imidazole molecule at the other end. Without perturbing electric fields PT from the hydronium proton donor to the imidazole acceptor occurs on a picosecond time scale. External perturbations to PT are created by electric fields of varying intensities, normal to the CNT axis, generated by a neutral pair of charges on the nanotube wall. For fields above ∼0.5 V/Å, the hydronium ion is effectively trapped at the CNT center, and PT blocked. Fields of comparable strength are generated inside proteins by nearby polar/charged amino acids. At lower fields the system displays a rich dynamic behavior, where the excess charge shuttles back and forth along the water chain before reaching the acceptor group on the picosecond time scale. The effects of the perturbing field on the proton movement are analyzed in terms of structural and dynamic properties of the water chain. The implications of these observations on PT in biomolecular systems and its control by external perturbing fields are discussed. [ABSTRACT FROM AUTHOR]

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. © 2004 American Institute of Physics. [ABSTRACT FROM AUTHOR]

We present a “coarse molecular dynamics” approach and apply it to studying the kinetics and thermodynamics of a peptide fragment dissolved in water. Short bursts of appropriately initialized simulations are used to infer the deterministic and stochastic components of the peptide motion parametrized by an appropriate set of coarse variables. Techniques from traditional numerical analysis (Newton–Raphson, coarse projective integration) are thus enabled; these techniques help analyze important features of the free-energy landscape (coarse transition states, eigenvalues and eigenvectors, transition rates, etc.). Reverse integration of coarse variables backward in time can assist escape from free energy minima and trace low-dimensional free energy surfaces. To illustrate the coarse molecular dynamics approach, we combine multiple short (0.5 ps) replica simulations to map the free energy surface of the “alanine dipeptide” in water, and to determine the ∼1/(1000 ps) rate of interconversion between the two stable configurational basins corresponding to the α-helical and extended minima. © 2003 American Institute of Physics. [ABSTRACT FROM AUTHOR]

A class of simple expressions of increasing accuracy for the free-energy difference between two states is derived based on numerical thermodynamic integration. The implementation of these formulas requires simulations of the initial and final (and possibly a few intermediate) states. They involve higher free-energy derivatives at these states which are related to the moments of the probability distribution of the perturbation. Given a specified number of such derivatives, these integration formulas are optimal in the sense that they are exact to the highest possible order of free-energy perturbation theory. The utility of this approach is illustrated for the hydration free energy of water. This problem provides a quite stringent test because the free energy is a highly nonlinear function of the charge so that even fourth order perturbation theory gives a very poor estimate of the free-energy change. Our results should prove most useful for complex, computationally demanding problems where free-energy differences arise primarily from changes in the electrostatic interactions (e.g., electron transfer, charging of ions, protonation of amino acids in proteins). © 1996 American Institute of Physics. [ABSTRACT FROM AUTHOR]

The pair and triplet correlations of a restricted primitive model (RPM) description of aqueous NaCl solutions have been computed in the concentration regime 0.5–4 M via Monte Carlo simulations, using both the method of configuration statistics and the particle insertion (chemical potential) method. Comparison of the triplet data to Kirkwood’s superposition approximation and the Abe and convolution approximations shows that deviations occur only for distances close to contact. The entropy computed via the entropy expansion including terms up to third order and the insertion method is dominated essentially by packing effects. The pair correlations and thermodynamic properties derived from the simulations were compared to the hypernetted chain (HNC) approximation which was found to be very good corroborating the conclusions drawn from earlier RPM studies. [ABSTRACT FROM AUTHOR]

Examines the use of Ewald summation in calculating electrostatic interactions in computer simulations of polar and charged substances. Derivation of formula to be implemented in standard molecular dynamics; Difference between virial and thermodynamic pressures; Model used in investigating the difference between virial and thermodynamic pressures.

Performs Monte Carlo simulations of the tetramethylammonium ion hydrated by simple charge water molecules. Relation of the derivatives of the free energy with the fluctuations in the ion-water interaction energies; Dependence of the calculated hydration free energies on system size and the method used in calculating electrostatic interactions.

Among other improvements, the Martini 3 coarse-grained force field provides a more accurate description of the solvation of protein pockets and channels through the consistent use of various bead types and sizes. Here, we show that the representation of Na+ and Cl− ions as "tiny" (TQ5) beads limits the accessible time step to 25 fs. By contrast, with Martini 2, time steps of 30–40 fs were possible for lipid bilayer systems without proteins. This limitation is relevant for systems that require long equilibration times. We derive a quantitative kinetic model of time-integration instabilities in molecular dynamics (MD) as a function of the time step, ion concentration and mass, system size, and simulation time. We demonstrate that ion–water interactions are the main source of instability at physiological conditions, followed closely by ion–ion interactions. We show that increasing the ionic masses makes it possible to use time steps up to 40 fs with minimal impact on static equilibrium properties and dynamical quantities, such as lipid and solvent diffusion coefficients. Increasing the size of the bead representing the ions (and thus changing their hydration) also permits longer time steps. For a soluble protein, we find that increasing the mass of tiny beads also on the protein permits simulations with 30-fs time steps. The use of larger time steps in Martini 3 results in a more efficient exploration of configuration space. The kinetic model of MD simulation crashes can be used to determine the maximum allowed time step upfront for an efficient use of resources and whenever sampling efficiency is critical. [ABSTRACT FROM AUTHOR]

Channelrhodopsin-2 (ChR2) is a light-sensitive ion channel widely used in optogenetics. Photoactivation triggers a trans-to-cis isomerization of a covalently bound retinal. Ensuing conformational changes open a cation-selective channel. We explore the structural dynamics in the early photocycle leading to channel opening by classical (MM) and quantum mechanical (QM) molecular simulations. With QM/MM simulations, we generated a protein-adapted force field for the retinal chromophore, which we validated against absorption spectra. In a 4-µs MM simulation of a dark-adapted ChR2 dimer, water entered the vestibules of the closed channel. Retinal all-trans to 13-cis isomerization, simulated with metadynamics, triggered a major restructuring of the charge cluster forming the channel gate. On a microsecond time scale, water penetrated the gate to form a membrane-spanning preopen pore between helices H1, H2, H3, and H7. This influx of water into an ion-impermeable preopen pore is consistent with time-resolved infrared spectroscopy and electrophysiology experiments. In the retinal 13-cis state, D253 emerged as the proton acceptor of the Schiff base. Upon proton transfer from the Schiff base to D253, modeled by QM/MM simulations, we obtained an early-M/P2 390-like intermediate. Rapid rotation of the unprotonated Schiff base toward the cytosolic side effectively prevents its reprotonation from the extracellular side. From MM and QM simulations, we gained detailed insight into the mechanism of ChR2 photoactivation and early events in pore formation. By rearranging the network of charges and hydrogen bonds forming the gate, water emerges as a key player in light-driven ChR2 channel opening. [ABSTRACT FROM AUTHOR]

Likelihood-based analysis of single-particle electron microscopy images has contributed much to the recent improvements in resolution. By treating particle orientations and classes probabilistically, uncertainties in the reconstruction process are explicitly accounted for, and the risk of bias towards the initial model is diminished. As a result, the quality and reliability of the reconstructions have greatly improved at manageable computational cost. Likelihood-based analysis of electron microscopy images also offers a route to direct coordinate refinement for dynamic systems, as an alternative to 3D density reconstruction. Here, we review recent developments in the algorithms used for reconstructions of high-resolution maps, and in the integrative framework of combining likelihood methods with simulations to address conformational variability in cryo-electron microscopy. [ABSTRACT FROM AUTHOR]

A one-dimensional diffusion equation is derived for the time evolution of the orientational factor, κ2, in the Förster energy transfer rate. The κ2-dependent diffusion coefficient is obtained in three different ways: (1) by requiring the κ2 autocorrelation function, calculated using the κ2 diffusion equation, to be single-exponential with the exact characteristic time; (2) by projecting the multidimensional diffusion equation for the transition dipoles onto κ2 using the local equilibrium approximation; and (3) by requiring exact and approximate κ2 trajectories to be as close as possible using a Bayesian approach. Within the framework of this simple theory, the distance dependence of the fluorescence resonance energy transfer (FRET) efficiency can be calculated for all values of the ratio of the rotational correlation time of the transition dipoles to the lifetime of the donor excited state. The theoretical predictions are compared to the exact values obtained from Brownian dynamics simulations of the reorientation of the donor and acceptor transition dipoles. [ABSTRACT FROM AUTHOR]

Condensed-phase simulations commonly use periodic boundary conditions (PBCs) to represent the thermodynamic limit. For the vapor to liquid transfer of an ion, the gas/liquid boundary and its associated potential change are then missing. Furthermore, the electric potential and field at a given point are given by conditionally convergent infinite series, for which different summation schemes give different results. Nevertheless, standard simulation protocols can be used to compute experimental quantities unambiguously. In particular, using an auxiliary test particle and a multistep solvation path, we show that particle-based, Ewald, and common molecule-based summation schemes for the potential and field are all essentially equivalent. However, all methods require prior knowledge of the gas/liquid boundary potential to compute ionic solvation free energies using PBC protocols for both force-field and quantum-mechanical models. [ABSTRACT FROM AUTHOR]

Studies the free energies of ionic solvation calculated from computer simulations. System of size dependence exhibited by the free energies; Estimation of the Born ion size; Application of the finite-size correlation to systems.

We investigate the dependence of single-particle diffusion coefficients on the size and shape of the simulation box in molecular dynamics simulations of fluids and lipid membranes. We find that the diffusion coefficients of lipids and a carbon nanotube embedded in a lipid membrane diverge with the logarithm of the box width. For a neat Lennard-Jones fluid in flat rectangular boxes, diffusion becomes anisotropic, diverging logarithmically in all three directions with increasing box width. In elongated boxes, the diffusion coefficients normal to the long axis diverge linearly with the height-to-width ratio. For both lipid membranes and neat fluids, this behavior is predicted quantitatively by hydrodynamic theory. Mean-square displacements in the neat fluid exhibit intermediate regimes of anomalous diffusion, with t ln t and t 3/2 components in flat and elongated boxes, respectively. For membranes, the large finite-size effects, and the apparent inability to determine a well-defined lipid diffusion coefficient from simulation, rationalize difficulties in comparing simulation results to each other and to those from experiments. [ABSTRACT FROM AUTHOR]

Molecular modeling and molecular dynamics simulations play an important role in the functional characterization of complex I. With its large size and complicated function, linking quinone reduction to proton pumping across a membrane, complex I poses unique modeling challenges. Nonetheless, simulations have already helped in the identification of possible proton transfer pathways. Simulations have also shed light on the coupling between electron and proton transfer, thus pointing the way in the search for the mechanistic principles underlying the proton pump. In addition to reviewing what has already been achieved in complex I modeling, we aim here to identify pressing issues and to provide guidance for future research to harness the power of modeling in the functional characterization of complex I. This article is part of a Special Issue entitled Respiratory complex I, edited by Volker Zickermann and Ulrich Brandt. [ABSTRACT FROM AUTHOR]

All-atom molecular dynamics simulations now allow us to create movies of proteins folding and unfolding. However, it is difficult to assess the accuracy of the folding mechanisms observed because experiments cannot yet directly resolve events occurring along the transition paths between unfolded and folded states. Protein folding ϕ-values provide residue-resolved information about folding mechanisms by comparing the effects of mutations on folding rates and stability, but determining ϕ-values by separately simulating mutant proteins would be computationally demanding and prone to large statistical errors. Here we use transition path theory to develop a method for computing ϕ-values directly from the transition path ensemble, without the need for additional simulations. This path-based approach uses the full transition path information available from equilibrium folding and unfolding trajectories, or from transition path sampling, and does not require identification of folding transition states. Applying our approach to a set of simulations of 10 small proteins by Shaw and coworkers [Lindorff-Larsen K, Piana S, Dror RO, Shaw DE (2011) Science 334(6055):517-520; Piana S, Lindorff-Larsen K, Shaw DE (2011) Biophys J 100(9):L47-L49; and Piana S, Lindorff-Larsen K, Shaw DE (2013) Proc Natl Acad Sci USA 110(15):5915-5920], we find good agreement with experiments in most cases where data are available. We can further resolve the contributions to fractional ϕ-values coming from partial contact formation versus transition path heterogeneity. Although in some cases, there is substantial heterogeneity of folding mechanism, in others, such as Ubiquitin, the mechanism is strongly conserved. [ABSTRACT FROM AUTHOR]

In a typical single-molecule force spectroscopy experiment, the ends of the molecule of interest are connected by long polymer linkers to a pair of mesoscopic beads trapped in the focus of two laser beams. At constant force load, the total extension, i.e., the end-to-end distance of the molecule plus linkers, is measured as a function of time. In the simplest systems, the measured extension fluctuates about two values characteristic of folded and unfolded states, with occasional transitions between them. We have recently shown that molecular (un)folding rates can be recovered from such trajectories, with a small linker correction, as long as the characteristic time of the bead fluctuations is shorter than the residence time in the unfolded (folded) state. Here, we show that accurate measurements of the molecular transition path times require an even faster apparatus response. Transition paths, the trajectory segments in which the molecule (un)folds, are properly resolved only if the beads fluctuate more rapidly than the end-to-end distance of the molecule. Therefore, over a wide regime, the measured rates may be meaningful but not the transition path times. Analytic expressions for the measured mean transition path times are obtained for systems diffusing anisotropically on a two-dimensional free energy surface. The transition path times depend on the properties both of the molecule and of the pulling device. [ABSTRACT FROM AUTHOR]

We combine molecular simulations and mechanical modeling to explore the mechanism of energy conversion in the coupled rotary motors of FoF1-ATP synthase. A torsional viscoelastic model with frictional dissipation quantitatively reproduces the dynamics and energetics seen in atomistic molecular dynamics simulations of torque-driven γ-subunit rotation in the F1-ATPase rotary motor. The torsional elastic coefficients determined from the simulations agree with results from independent single-molecule experiments probing different segments of the γ-subunit, which resolves a long-lasting controversy. At steady rotational speeds of ~1 kHz corresponding to experimental turnover, the calculated frictional dissipation of less than kBT per rotation is consistent with the high thermodynamic efficiency of the fully reversible motor. Without load, the maximum rotational speed during transitions between dwells is reached at ~1 MHz. Energetic constraints dictate a unique pathway for the coupled rotations of the Fo and F1 rotary motors in ATP synthase, and explain the need for the finer stepping of the F1 motor in the mammalian system, as seen in recent experiments. Compensating for incommensurate eightfold and threefold rotational symmetries in Fo and F1, respectively, a significant fraction of the external mechanical work is transiently stored as elastic energy in the γ-subunit. The general framework developed here should be applicable to other molecular machines. [ABSTRACT FROM AUTHOR]

Translational diffusion coefficients are routinely estimated from molecular dynamics simulations. Linear fits to mean squared displacement (MSD) curves have become the de facto standard, from simple liquids to complex biomacromolecules. Nonlinearities in MSD curves at short times are handled with a wide variety of ad hoc practices, such as partial and piece-wise fitting of the data. Here, we present a rigorous framework to obtain reliable estimates of the self-diffusion coefficient and its statistical uncertainty. We also assess in a quantitative manner if the observed dynamics is, indeed, diffusive. By accounting for correlations between MSD values at different times, we reduce the statistical uncertainty of the estimator and, thereby, increase its efficiency. With a Kolmogorov–Smirnov test, we check for possible anomalous diffusion. We provide an easy-to-use Python data analysis script for the estimation of self-diffusion coefficients. As an illustration, we apply the formalism to molecular dynamics simulation data of pure TIP4P-D water and a single ubiquitin protein. In another paper [S. von Bülow, J. T. Bullerjahn, and G. Hummer, J. Chem. Phys. 153, 021101 (2020)], we demonstrate its ability to recognize deviations from regular diffusion caused by systematic errors in a common trajectory "unwrapping" scheme that is implemented in popular simulation and visualization software. [ABSTRACT FROM AUTHOR]

In molecular dynamics simulations under periodic boundary conditions, particle positions are typically wrapped into a reference box. For diffusion coefficient calculations using the Einstein relation, the particle positions need to be unwrapped. Here, we show that a widely used heuristic unwrapping scheme is not suitable for long simulations at constant pressure. Improper accounting for box-volume fluctuations creates, at long times, unphysical trajectories and, in turn, grossly exaggerated diffusion coefficients. We propose an alternative unwrapping scheme that resolves this issue. At each time step, we add the minimal displacement vector according to periodic boundary conditions for the instantaneous box geometry. Here and in another paper [J. T. Bullerjahn, S. von Bülow, and G. Hummer, J. Chem. Phys. 153, 024116 (2020)], we apply the new unwrapping scheme to extensive molecular dynamics and Brownian dynamics simulation data. We provide practitioners with a formula to assess if and by how much earlier results might have been affected by the widely used heuristic unwrapping scheme. [ABSTRACT FROM AUTHOR]

As shown by Jarzynski, free energy differences between equilibrium states can be expressed in terms of the statistics of work carried out on a system during non-equilibrium transformations. This exact result, as well as the related Crooks fluctuation theorem, provide the basis for the computation of free energy differences from fast switching molecular dynamics simulations, in which an external parameter is changed at a finite rate, driving the system away from equilibrium. In this article, we first briefly review the Jarzynski identity and the Crooks fluctuation theorem and then survey various algorithms building on these relations. We pay particular attention to the statistical efficiency of these methods and discuss practical issues arising in their implementation and the analysis of the results. [ABSTRACT FROM AUTHOR]

This special issue presents a series of papers that highlight a new and evolving view of water's molecular role in biological structure and dynamics. Increasingly water is appreciated for playing an active role in biological function rather than just serving as a spectator. The areas represented include molecular interfaces with water, biological self-assembly, conformational changes by macromolecules, and chemical reactivity [ABSTRACT FROM AUTHOR]

Abstract: We develop a method to extract structural information from electron microscopy (EM) images of dynamic and heterogeneous molecular assemblies. To overcome the challenge of disorder in the imaged structures, we analyze each image individually, avoiding information loss through clustering or averaging. The Bayesian inference of EM (BioEM) method uses a likelihood-based probabilistic measure to quantify the consistency between each EM image and given structural models. The likelihood function accounts for uncertainties in the molecular position and orientation, variations in the relative intensities and noise in the experimental images. The BioEM formalism is physically intuitive and mathematically simple. We show that for experimental GroEL images, BioEM correctly identifies structures according to the functional state. The top-ranked structure is the corresponding X-ray crystal structure, followed by an EM structure generated previously from a superset of the EM images used here. To analyze EM images of highly flexible molecules, we propose an ensemble refinement procedure, and validate it with synthetic EM maps of the ESCRT-I–II supercomplex. Both the size of the ensemble and its structural members are identified correctly. BioEM offers an alternative to 3D-reconstruction methods, extracting accurate population distributions for highly flexible structures and their assemblies. We discuss limitations of the method, and possible applications beyond ensemble refinement, including the cross-validation and unbiased post-assessment of model structures, and the structural characterization of systems where traditional approaches fail. Overall, our results suggest that the BioEM framework can be used to analyze EM images of both ordered and disordered molecular systems. [Copyright &y& Elsevier]

The recent availability of long equilibrium simulations of protein folding in atomistic detail for more than 10 proteins allows us to identify the key interactions driving folding. We find that the collective fraction of native amino acid contacts, Q, captures remarkably well the transition states for all the proteins with a folding free energy barrier. Going beyond this global picture, we devise two different measures to quantify the importance of individual interresidue contacts in the folding mechanism: (i) the log-ratio of lifetimes of contacts during folding transition paths and in the unfolded state and (ii) a Bayesian measure of how predictive the formation of each contact is for being on a transition path. Both of these measures indicate that native, or near-native, contacts are important for determining mechanism, as might be expected. More remarkably, however, we found that for almost all the proteins, with the designed protein α3D being a notable exception, nonnative contacts play no significant part in determining folding mechanisms. [ABSTRACT FROM AUTHOR]

F1-ATPase, the catalytic domain of ATP synthase, synthesizes most of the ATP in living organisms. Running in reverse powered by ATP hydrolysis, this hexameric ring-shaped molecular motor formed by three αβ-dimers creates torque on its central γ-subunit. This reverse operation enables detailed explorations of the mechanochemical coupling mechanisms in experiment and simulation. Here, we use molecular dynamics simulations to construct a first atomistic conformation of the intermediate state following the 40° substep of rotary motion, and to study the timing and molecular mechanism of inorganic phosphate (Pi) release coupled to the rotation. In response to torque-driven rotation of the γ-subunit in the hydrolysis direction, the nucleotide-free αβE interface forming the "empty" E site loosens and singly charged Pi readily escapes to the P loop. By contrast, the interface stays closed with doubly charged Pi. The γ-rotation tightens the ATP-bound αβTP interface, as required for hydrolysis. The calculated rate for the outward release of doubly charged Pi from the αβE interface 120° after ATP hydrolysis closely matches the ~1-ms functional timescale. Conversely, Pi release from the ADP-bound αβDP interface postulated in earlier models would occur through a kinetically infeasible inward-directed pathway. Our simulations help reconcile conflicting interpretations of single-molecule experiments and crystallographic studies by clarifying the timing of Pi exit, its pathway and kinetics, associated changes in Pi protonation, and changes of the F1-ATPase structure in the 40° substep. Important elements of the molecular mechanism of Pi release emerging from our simulations appear to be conserved in myosin despite the different functional motions. [ABSTRACT FROM AUTHOR]

The pair-distance distribution function (PDDF) contains all structural information probed in an elastic scattering experiment of macromolecular solutions. However, in small-angle x-ray scattering (SAXS) or small-angle neutron scattering (SANS) experiments only their Fourier transform is measured over a restricted range of scattering angles. We therefore developed a mathematically simple and computationally efficient method to calculate the PDDFs as well as accurate scattering intensities from molecular dynamics simulations. The calculated solution scattering intensities are in excellent agreement with SAXS and wide-angle x-ray scattering (WAXS) experiments for a series of proteins. The corresponding PDDFs are remarkably rich in features reporting on the detailed protein structure. Using an inverse Fourier transform method, most of these features can be recovered if scattering intensities are measured up to a momentum transfer of q ≈ 2-3 Å -1. Our results establish that high-precision solution scattering experiments utilizing x-ray free-electron lasers and third generation synchrotron sources can resolve subnanometer structural detail, well beyond size, shape, and fold. [ABSTRACT FROM AUTHOR]

We define an interface-interaction network (IIN) to capture the specificity and competition between protein-protein interactions (PPI). This new type of network represents interactions between individual interfaces used in functional protein binding and thereby contains the detail necessary to describe the competition and cooperation between any pair of binding partners. Here we establish a general framework for the construction of IINs that merges computational structure-based interface assignment with careful curation of available literature. To complement limited structural data, the inclusion of biochemical data is critical for achieving the accuracy and completeness necessary to analyze the specificity and competition between the protein interactions. Firstly, this procedure provides a means to clarify the information content of existing data on purported protein interactions and to remove indirect and spurious interactions. Secondly, the IIN we have constructed here for proteins involved in clathrin-mediated endocytosis (CME) exhibits distinctive topological properties. In contrast to PPI networks with their global and relatively dense connectivity, the fragmentation of the IIN into distinctive network modules suggests that different functional pressures act on the evolution of its topology. Large modules in the IIN are formed by interfaces sharing specificity for certain domain types, such as SH3 domains distributed across different proteins. The shared and distinct specificity of an interface is necessary for effective negative and positive design of highly selective binding targets. Lastly, the organization of detailed structural data in a network format allows one to identify pathways of specific binding interactions and thereby predict effects of mutations at specific surfaces on a protein and of specific binding inhibitors, as we explore in several examples. Overall, the endocytosis IIN is remarkably complex and rich in features masked in the coarser PPI, and collects relevant detail of protein association in a readily interpretable format. [ABSTRACT FROM AUTHOR]

We explore the theoretical foundation of different stringmethodsused to find dominant reaction pathways in high-dimensional configurationspaces. Pathways are assessed by the amount of reactive flux theycarry and by their orientation relative to the committor function.By examining the effects of transforming between different collectivecoordinates that span the same underlying space, we unmask artificialcoordinate dependences in strings optimized to follow the free energygradient. In contrast, strings optimized to follow the drift vectorproduce reaction pathways that are significantly less sensitive toreparameterizations of the collective coordinates. The differencesin these paths arise because the drift vector depends on both thefree energy gradient and the diffusion tensor of the coarse collectivevariables. Anisotropy and position dependence of diffusion tensorsarise commonly in spaces of coarse variables, whose generally slowdynamics are obtained by nonlinear projections of the strongly coupledatomic motions. We show here that transition paths constructed toaccount for dynamics by following the drift vector will (to a closeapproximation) carry the maximum reactive flux both in systems withisotropic position dependent diffusion and in systems with constantbut anisotropic diffusion. We derive a simple method for calculatingthe committor function along paths that follow the reactive flux.Lastly, we provide guidance for the practical implementation of thedynamic string method. [ABSTRACT FROM AUTHOR]

Abstract: Cytochrome c oxidase is an efficient energy transducer that reduces oxygen to water and converts the released chemical energy into an electrochemical membrane potential. As a true proton pump, cytochrome c oxidase translocates protons across the membrane against this potential. Based on a wealth of experiments and calculations, an increasingly detailed picture of the reaction intermediates in the redox cycle has emerged. However, the fundamental mechanism of proton pumping coupled to redox chemistry remains largely unresolved. Here we examine and extend a kinetic master-equation approach to gain insight into redox-coupled proton pumping in cytochrome c oxidase. Basic principles of the cytochrome c oxidase proton pump emerge from an analysis of the simplest kinetic models that retain essential elements of the experimentally determined structure, energetics, and kinetics, and that satisfy fundamental physical principles. The master-equation models allow us to address the question of how pumping can be achieved in a system in which all reaction steps are reversible. Whereas proton pumping does not require the direct modulation of microscopic reaction barriers, such kinetic gating greatly increases the pumping efficiency. Further efficiency gains can be achieved by partially decoupling the proton uptake pathway from the active-site region. Such a mechanism is consistent with the proposed Glu valve, in which the side chain of a key glutamic acid shuttles between the D channel and the active-site region. We also show that the models predict only small proton leaks even in the absence of turnover. The design principles identified here for cytochrome c oxidase provide a blueprint for novel biology-inspired fuel cells, and the master-equation formulation should prove useful also for other molecular machines. This article is part of a Special Issue entitled: Respiratory Oxidases. [Copyright &y& Elsevier]

Complex I (NADH-ubiquinone oxidoreductase) in the respiratory chain of mitochondria and several bacteria functions as a redoxdriven proton pump that contributes to the generation of the protonmotive force across the inner mitochondrial or bacterial membrane and thus to the aerobic synthesis of ATP. The stoichiometry of proton translocation is thought to be 4 H+ per NADH oxidized (2 e-). Here we show that a H+/2 e ratio of 3 appears more likely on the basis of the recently determined H+/ATP ratio of the mitochondrial F1F0-ATP synthase of animal mitochondria and of a set of carefully determined ATP/2 e ratios for different segments of the mitochondrial respiratory chain. This lower H+/2 e ratio of 3 is independently supported by thermodynamic analyses of experiments with both mitochondria and submitochondrial particles. A reduced H+/2 e stoichiometry of 3 has important mechanistic implications for this proton pump. In a rough mechanistic model, we suggest a concerted proton translocation mechanism in the three homologous and tightly packed antiporterlike subunits L, M, and N of the proton-translocating membrane domain of complex I. [ABSTRACT FROM AUTHOR]

The weighted histogram analysis method (WHAM) has become the standard technique for the analysis of umbrella sampling simulations. In this article, we address the challenges (1) of obtaining fast and accurate solutions of the coupled nonlinear WHAM equations, (2) of quantifying the statistical errors of the resulting free energies, (3) of diagnosing possible systematic errors, and (4) of optimally allocating of the computational resources. Traditionally, the WHAM equations are solved by a fixed-point direct iteration method, despite poor convergence and possible numerical inaccuracies in the solutions. Here, we instead solve the mathematically equivalent problem of maximizing a target likelihood function, by using superlinear numerical optimization algorithms with a significantly faster convergence rate. To estimate the statistical errors in one-dimensional free energy profiles obtained from WHAM, we note that for densely spaced umbrella windows with harmonic biasing potentials, the WHAM free energy profile can be approximated by a coarse-grained free energy obtained by integrating the mean restraining forces. The statistical errors of the coarse-grained free energies can be estimated straightforwardly and then used for the WHAM results. A generalization to multidimensional WHAM is described. We also propose two simple statistical criteria to test the consistency between the histograms of adjacent umbrella windows, which help identify inadequate sampling and hysteresis in the degrees of freedom orthogonal to the reaction coordinate. Together, the estimates of the statistical errors and the diagnostics of inconsistencies in the potentials of mean force provide a basis for the efficient allocation of computational resources in free energy simulations. © 2011 Wiley Periodicals, Inc. J Comput Chem, 2012 [ABSTRACT FROM AUTHOR]

Proton-coupled electron transfer (PCET) is an elementary chemical reaction crucial for biological oxidoreduction. We perform quantum chemical calculations to study the direct and water-mediated PCET between two stacked tyrosines, TyrO• + TyrOH → TyrOH + TyrO•, to mimic a key step in the catalytic reaction of class Ia ribonucleotide reductase (RNR). The energy surfaces of electronic ground and excited states are separated by a large gap of 20 kcal mol-1, indicative of an electronically adiabatic transfer mechanism. In response to chemical substitutions of the proton donor, the energy of the transition state for direct PCET shifts by exactly half of the change in energetic driving force, resulting in a linear free energy relation with a Brønsted slope of 1/2. In contrast, for water-mediated PCET, we observe integer Brønsted slopes of 1 and 0 for proton acceptor and donor modifications, respectively. Our calculations suggest that the π-stacking of the tyrosine dimer in RNR results in strong electronic coupling and adiabatic PCET. Water participation in the PCET can be identified perturbatively in a Brønsted analysis. [ABSTRACT FROM AUTHOR]

Multicellular organisms, from Caenorhabditis elegans to humans, have roughly the same number of protein encoding genes. We show that the need to prevent disease-causing nonspecific interactions between proteins provides a simple physical reason why organism complexity is not reflected in the number of distinct proteins. By collective evolution of the amino acid sequences of protein binding interfaces we estimate the degree of misbinding as a function of the number of distinct proteins. Protein interaction energies are calculated with an empirical, residue-specific. energy function tuned for protein binding. We show that the achievable energy gap favoring specific over nonspecific binding decreases with protein number in a power-law fashion. From the fraction of proteins involved in nonspecific complexes as a function of increasing protein number and decreasing energy gap, we predict the limits these binding requirements place on the number of different proteins that can function effectively in a given cellular compartment. Remarkably, the optimization of binding interfaces favors networks in which a few proteins have many partners, and most proteins have few partners, consistent with a scale-free network topology. We conclude that nonspecific binding adds to the evolutionary pressure to develop scale-free protein-protein interaction networks. [ABSTRACT FROM AUTHOR]

Nonequilibrium pulling experiments provide detailed information about the thermodynamic and kinetic properties of molecules. We show that unperturbed free energy profiles as a function of molecular extension can be obtained rigorously from such experiments without using work-weighted position histograms. An inverse Weierstrass transform is used to relate the system free energy obtained from the Jarzynski equality directly to the underlying molecular free energy surface. An accurate approximation for the free energy surface is obtained by using the method of steepest descent to evaluate the inverse transform. The formalism is applied to simulated data obtained from a kinetic model of RNA folding, in which the dynamics consists of jumping between linker-dominated folded and unfolded free energy surfaces. [ABSTRACT FROM AUTHOR]

Nerve signaling in humans and chemical sensing in bacteria both rely on the controlled opening and closing of the ion-conducting pore in pentameric ligand-gated ion channels. With the help of a multiscale simulation approach that combines mixed elastic network model calculations with molecular dynamics simulations, we study the opening and closing of the pore in Gloeobacter violaceus channel GLIC at atomic resolution. In oursimulations of the GLIC transmembrane domain, we first verify that the two endpoints of the transition are open and closed to sodium ion conduction, respectively. We then show that a two-stage tilting of the pore- lining helices induces cooperative drying and iris-like closing of the channel pore. From the free energy profile of the gating transi- tion and from unrestrained simulations, we conclude that the pore of the isolated GLIC transmembrane domain closes spontaneously. The mechanical work of opening the pore is performed primarily on the M2-M3 loop. Strong interactions of this short and conserved loop with the extracellular domain are therefore crucial to couple ligand binding to channel opening. [ABSTRACT FROM AUTHOR]

Diffusion on a low-dimensional free-energy surface is a remarkably successful model for the folding dynamics of small single-domain proteins. Complicating the interpretation of both simulations and experiments is the expectation that the effective diffusion coefficient D will in general depend on the position along the folding coordinate, and this dependence may vary for different coordinates. Here we explore the position dependence of D. its connection to protein internal friction, and the consequences for the interpretation of single-molecule experiments. We find a large decrease in D from unfolded to folded, for reaction coordinates that directly measure fluctuations in Cartesian configuration space, including those probed in single-molecule experiments. In contrast, D is almost independent of Q, the fraction of native amino acid contacts: Near the folded state, small fluctuations in position cause large fluctuations in Q, and vice versa for the unfolded state. In general, position-dependent free energies and diffusion coefficients for any two good reaction coordinates that separate reactant, product, and transition states, are related by a simple transformation, as we demonstrate. With this transformation, we obtain reaction coordinates with position-invariant D. The corresponding free-energy surfaces allow us to justify the assumptions used in estimating the speed limit for protein folding from experimental measurements of the reconfiguration time in the unfolded state, and also reveal intermediates hidden in the original free-energy projection. Lastly, we comment on the design of future single-molecule experiments that probe the position dependence of D directly. [ABSTRACT FROM AUTHOR]

The author discusses the adenosine triphosphate (ATP) synthases. Topics discussed by the author include the regeneration of molecular fuel for cellular processes, the efficieny of c-ring rotation in membrane proteins interfaces, and interaction between ATP synthase and cardiolipin (CL) establishes activity in mitochondria.