Bayesian Mechanistic Inference, Statistical Mechanics, and a New Era for Monte Carlo.

On the one hand, much of computational chemistry is concerned with "bottom-up" calculations which elucidate observable behavior starting from exact or approximated physical laws, a paradigm exemplified by typical quantum mechanical calculations and molecular dynamics simulations. On the other hand, "top down" computations aiming to formulate mathematical models consistent with observed data, e.g., parametrizing force fields, binding or kinetic models, have been of interest for decades but recently have grown in sophistication with the use of Bayesian inference (BI). Standard BI provides an estimation of parameter values, uncertainties, and correlations among parameters. Used for "model selection," BI can also distinguish between model structures such as the presence or absence of individual states and transitions. Fortunately for physical scientists, BI can be formulated within a statistical mechanics framework, and indeed, BI has led to a resurgence of interest in Monte Carlo (MC) algorithms, many of which have been directly adapted from or inspired by physical strategies. Certain MC algorithms─notably procedures using an "infinite temperature" reference state─can be successful in a 5-20 parameter BI context which would be unworkable in molecular spaces of 103 coordinates and more. This Review provides a pedagogical introduction to BI and reviews key aspects of BI through a physical lens, setting the computations in terms of energy landscapes and free energy calculations and describing promising sampling algorithms. Statistical mechanics and basic probability theory also provide a reference for understanding intrinsic limitations of Bayesian inference with regard to model selection and the choice of priors.

From Average Transient Transporter Currents to Microscopic Mechanism─A Bayesian Analysis.

Electrophysiology studies of secondary active transporters have revealed quantitative mechanistic insights over many decades of research. However, the emergence of new experimental and analytical approaches calls for investigation of the capabilities and limitations of the newer methods. We examine the ability of solid-supported membrane electrophysiology (SSME) to characterize discrete-state kinetic models with >10 rate constants. We use a Bayesian framework applied to synthetic data for three tasks: to quantify and check (i) the precision of parameter estimates under different assumptions, (ii) the ability of computation to guide the selection of experimental conditions, and (iii) the ability of our approach to distinguish among mechanisms based on SSME data. When the general mechanism, i.e., event order, is known in advance, we show that a subset of kinetic parameters can be "practically identified" within ∼1 order of magnitude, based on SSME current traces that visually appear to exhibit simple exponential behavior. This remains true even when accounting for systematic measurement bias and realistic uncertainties in experimental inputs (concentrations) are incorporated into the analysis. When experimental conditions are optimized or different experiments are combined, the number of practically identifiable parameters can be increased substantially. Some parameters remain intrinsically difficult to estimate through SSME data alone, suggesting that additional experiments are required to fully characterize parameters. We also demonstrate the ability to perform model selection and determine the order of events when that is not known in advance, comparing Bayesian and maximum-likelihood approaches. Finally, our studies elucidate good practices for the increasingly popular but subtly challenging Bayesian calculations for structural and systems biology.

Single-cell morphodynamical trajectories enable prediction of gene expression accompanying cell state change.

Extracellular signals induce changes to molecular programs that modulate multiple cellular phenotypes, including proliferation, motility, and differentiation status. The connection between dynamically adapting phenotypic states and the molecular programs that define them is not well understood. Here we develop data-driven models of single-cell phenotypic responses to extracellular stimuli by linking gene transcription levels to "morphodynamics" - changes in cell morphology and motility observable in time-lapse image data. We adopt a dynamics-first view of cell state by grouping single-cell trajectories into states with shared morphodynamic responses. The single-cell trajectories enable development of a first-of-its-kind computational approach to map live-cell dynamics to snapshot gene transcript levels, which we term MMIST, Molecular and Morphodynamics-Integrated Single-cell Trajectories. The key conceptual advance of MMIST is that cell behavior can be quantified based on dynamically defined states and that extracellular signals alter the overall distribution of cell states by altering rates of switching between states. We find a cell state landscape that is bound by epithelial and mesenchymal endpoints, with distinct sequences of epithelial to mesenchymal transition (EMT) and mesenchymal to epithelial transition (MET) intermediates. The analysis yields predictions for gene expression changes consistent with curated EMT gene sets and provides a prediction of thousands of RNA transcripts through extracellular signal-induced EMT and MET with near-continuous time resolution. The MMIST framework leverages true single-cell dynamical behavior to generate molecular-level omics inferences and is broadly applicable to other biological domains, time-lapse imaging approaches and molecular snapshot data.

Morphodynamical cell state description via live-cell imaging trajectory embedding.

Time-lapse imaging is a powerful approach to gain insight into the dynamic responses of cells, but the quantitative analysis of morphological changes over time remains challenging. Here, we exploit the concept of "trajectory embedding" to analyze cellular behavior using morphological feature trajectory histories-that is, multiple time points simultaneously, rather than the more common practice of examining morphological feature time courses in single timepoint (snapshot) morphological features. We apply this approach to analyze live-cell images of MCF10A mammary epithelial cells after treatment with a panel of microenvironmental perturbagens that strongly modulate cell motility, morphology, and cell cycle behavior. Our morphodynamical trajectory embedding analysis constructs a shared cell state landscape revealing ligand-specific regulation of cell state transitions and enables quantitative and descriptive models of single-cell trajectories. Additionally, we show that incorporation of trajectories into single-cell morphological analysis enables (i) systematic characterization of cell state trajectories, (ii) better separation of phenotypes, and (iii) more descriptive models of ligand-induced differences as compared to snapshot-based analysis. This morphodynamical trajectory embedding is broadly applicable to the quantitative analysis of cell responses via live-cell imaging across many biological and biomedical applications.

A Suite of Tutorials for the WESTPA 2.0 Rare-Events Sampling Software [Article v2.0].

The weighted ensemble (WE) strategy has been demonstrated to be highly efficient in generating pathways and rate constants for rare events such as protein folding and protein binding using atomistic molecular dynamics simulations. Here we present two sets of tutorials instructing users in the best practices for preparing, carrying out, and analyzing WE simulations for various applications using the WESTPA software. The first set of more basic tutorials describes a range of simulation types, from a molecular association process in explicit solvent to more complex processes such as host-guest association, peptide conformational sampling, and protein folding. The second set ecompasses six advanced tutorials instructing users in the best practices of using key new features and plugins/extensions of the WESTPA 2.0 software package, which consists of major upgrades for larger systems and/or slower processes. The advanced tutorials demonstrate the use of the following key features: (i) a generalized resampler module for the creation of "binless" schemes, (ii) a minimal adaptive binning scheme for more efficient surmounting of free energy barriers, (iii) streamlined handling of large simulation datasets using an HDF5 framework, (iv) two different schemes for more efficient rate-constant estimation, (v) a Python API for simplified analysis of WE simulations, and (vi) plugins/extensions for Markovian Weighted Ensemble Milestoning and WE rule-based modeling for systems biology models. Applications of the advanced tutorials include atomistic and non-spatial models, and consist of complex processes such as protein folding and the membrane permeability of a drug-like molecule. Users are expected to already have significant experience with running conventional molecular dynamics or systems biology simulations.

Quantifying cooperative multisite binding in the hub protein LC8 through Bayesian inference.

Multistep protein-protein interactions underlie most biological processes, but their characterization through methods such as isothermal titration calorimetry (ITC) is largely confined to simple models that provide little information on the intermediate, individual steps. In this study, we primarily examine the essential hub protein LC8, a small dimer that binds disordered regions of 100+ client proteins in two symmetrical grooves at the dimer interface. Mechanistic details of LC8 binding have remained elusive, hampered in part by ITC data analyses employing simple models that treat bivalent binding as a single event with a single binding affinity. We build on existing Bayesian ITC approaches to quantify thermodynamic parameters for multi-site binding interactions impacted by significant uncertainty in protein concentration. Using a two-site binding model, we identify positive cooperativity with high confidence for LC8 binding to multiple client peptides. In contrast, application of an identical model to the two-site binding between the coiled-coil NudE dimer and the intermediate chain of dynein reveals little evidence of cooperativity. We propose that cooperativity in the LC8 system drives the formation of saturated induced-dimer structures, the functional units of most LC8 complexes. In addition to these system-specific findings, our work advances general ITC analysis in two ways. First, we describe a previously unrecognized mathematical ambiguity in concentrations in standard binding models and clarify how it impacts the precision with which binding parameters are determinable in cases of high uncertainty in analyte concentrations. Second, building on observations in the LC8 system, we develop a system-agnostic heat map of practical parameter identifiability calculated from synthetic data which demonstrates that the ability to determine microscopic binding parameters is strongly dependent on both the parameters themselves and experimental conditions. The work serves as a foundation for determination of multi-step binding interactions, and we outline best practices for Bayesian analysis of ITC experiments.

#COVIDisAirborne: AI-enabled multiscale computational microscopy of delta SARS-CoV-2 in a respiratory aerosol.

We seek to completely revise current models of airborne transmission of respiratory viruses by providing never-before-seen atomic-level views of the SARS-CoV-2 virus within a respiratory aerosol. Our work dramatically extends the capabilities of multiscale computational microscopy to address the significant gaps that exist in current experimental methods, which are limited in their ability to interrogate aerosols at the atomic/molecular level and thus obscure our understanding of airborne transmission. We demonstrate how our integrated data-driven platform provides a new way of exploring the composition, structure, and dynamics of aerosols and aerosolized viruses, while driving simulation method development along several important axes. We present a series of initial scientific discoveries for the SARS-CoV-2 Delta variant, noting that the full scientific impact of this work has yet to be realized.

Ligand-mediated Structural Dynamics of a Mammalian Pancreatic KATP Channel.

Regulation of pancreatic KATP channels involves orchestrated interactions of their subunits, Kir6.2 and SUR1, and ligands. Previously we reported KATP channel cryo-EM structures in the presence and absence of pharmacological inhibitors and ATP, focusing on the mechanisms by which inhibitors act as pharmacological chaperones of KATP channels (Martin et al., 2019). Here we analyzed the same cryo-EM datasets with a focus on channel conformational dynamics to elucidate structural correlates pertinent to ligand interactions and channel gating. We found pharmacological inhibitors and ATP enrich a channel conformation in which the Kir6.2 cytoplasmic domain is closely associated with the transmembrane domain, while depleting one where the Kir6.2 cytoplasmic domain is extended away into the cytoplasm. This conformational change remodels a network of intra- and inter-subunit interactions as well as the ATP and PIP2 binding pockets. The structures resolved key contacts between the distal N-terminus of Kir6.2 and SUR1's ABC module involving residues implicated in channel function and showed a SUR1 residue, K134, participates in PIP2 binding. Molecular dynamics simulations revealed two Kir6.2 residues, K39 and R54, that mediate both ATP and PIP2 binding, suggesting a mechanism for competitive gating by ATP and PIP2.

Continuum dynamics and statistical correction of compositional heterogeneity in multivalent IDP oligomers resolved by single-particle EM.

Multivalent intrinsically disordered protein (IDP) complexes are prevalent in biology and act in regulation of diverse processes, including transcription, signaling events, and the assembly and disassembly of complex macromolecular architectures. These systems pose significant challenges to structural investigation, due to continuum dynamics imparted by the IDP and compositional heterogeneity resulting from characteristic low-affinity interactions. Here, we developed a modular pipeline for automated single-particle electron microscopy (EM) distribution analysis of common but relatively understudied semi-ordered systems: 'beads-on-a-string' assemblies, composed of IDPs bound at multivalent sites to the ubiquitous ∼20 kDa cross-linking hub protein LC8. This approach quantifies conformational geometries and compositional heterogeneity on a single-particle basis, and statistically corrects spurious observations arising from random proximity of bound and unbound LC8. The statistical correction is generically applicable to oligomer characterization and not specific to our pipeline. Following validation, the approach was applied to the nuclear pore IDP Nup159 and the transcription factor ASCIZ. This analysis unveiled significant compositional and conformational diversity in both systems that could not be obtained from ensemble single particle EM class-averaging strategies, and new insights for exploring how these architectural properties might contribute to their physiological roles in supramolecular assembly and transcriptional regulation. We expect that this approach may be adopted to many other intrinsically disordered systems that have evaded traditional methods of structural characterization.

WESTPA 2.0: High-Performance Upgrades for Weighted Ensemble Simulations and Analysis of Longer-Timescale Applications.

The weighted ensemble (WE) family of methods is one of several statistical mechanics-based path sampling strategies that can provide estimates of key observables (rate constants and pathways) using a fraction of the time required by direct simulation methods such as molecular dynamics or discrete-state stochastic algorithms. WE methods oversee numerous parallel trajectories using intermittent overhead operations at fixed time intervals, enabling facile interoperability with any dynamics engine. Here, we report on the major upgrades to the WESTPA software package, an open-source, high-performance framework that implements both basic and recently developed WE methods. These upgrades offer substantial improvements over traditional WE methods. The key features of the new WESTPA 2.0 software enhance the efficiency and ease of use: an adaptive binning scheme for more efficient surmounting of large free energy barriers, streamlined handling of large simulation data sets, exponentially improved analysis of kinetics, and developer-friendly tools for creating new WE methods, including a Python API and resampler module for implementing both binned and "binless" WE strategies.

Vascular KATP channel structural dynamics reveal regulatory mechanism by Mg-nucleotides.

Vascular tone is dependent on smooth muscle KATP channels comprising pore-forming Kir6.1 and regulatory SUR2B subunits, in which mutations cause Cantú syndrome. Unique among KATP isoforms, they lack spontaneous activity and require Mg-nucleotides for activation. Structural mechanisms underlying these properties are unknown. Here, we determined cryogenic electron microscopy structures of vascular KATP channels bound to inhibitory ATP and glibenclamide, which differ informatively from similarly determined pancreatic KATP channel isoform (Kir6.2/SUR1). Unlike SUR1, SUR2B subunits adopt distinct rotational "propeller" and "quatrefoil" geometries surrounding their Kir6.1 core. The glutamate/aspartate-rich linker connecting the two halves of the SUR-ABC core is observed in a quatrefoil-like conformation. Molecular dynamics simulations reveal MgADP-dependent dynamic tripartite interactions between this linker, SUR2B, and Kir6.1. The structures captured implicate a progression of intermediate states between MgADP-free inactivated, and MgADP-bound activated conformations wherein the glutamate/aspartate-rich linker participates as mobile autoinhibitory domain, suggesting a conformational pathway toward KATP channel activation.

What Markov State Models Can and Cannot Do: Correlation versus Path-Based Observables in Protein-Folding Models.

Markov state models (MSMs) have been widely applied to study the kinetics and pathways of protein conformational dynamics based on statistical analysis of molecular dynamics (MD) simulations. These MSMs coarse-grain both configuration space and time in ways that limit what kinds of observables they can reproduce with high fidelity over different spatial and temporal resolutions. Despite their popularity, there is still limited understanding of which biophysical observables can be computed from these MSMs in a robust and unbiased manner, and which suffer from the space-time coarse-graining intrinsic in the MSM model. Most theoretical arguments and practical validity tests for MSMs rely on long-time equilibrium kinetics, such as the slowest relaxation time scales and experimentally observable time-correlation functions. Here, we perform an extensive assessment of the ability of well-validated protein folding MSMs to accurately reproduce path-based observable such as mean first-passage times (MFPTs) and transition path mechanisms compared to a direct trajectory analysis. We also assess a recently proposed class of history-augmented MSMs (haMSMs) that exploit additional information not accounted for in standard MSMs. We conclude with some practical guidance on the use of MSMs to study various problems in conformational dynamics of biomolecules. In brief, MSMs can accurately reproduce correlation functions slower than the lag time, but path-based observables can only be reliably reproduced if the lifetimes of states exceed the lag time, which is a much stricter requirement. Even in the presence of short-lived states, we find that haMSMs reproduce path-based observables more reliably.

Connexin 46 and connexin 50 gap junction channel properties are shaped by structural and dynamic features of their N-terminal domains.

KEY POINTS: Gap junctions formed by different connexins are expressed throughout the body and harbour unique channel properties that have not been fully defined mechanistically. Recent structural studies by cryo-electron microscopy have produced high-resolution models of the related but functionally distinct lens connexins (Cx50 and Cx46) captured in a stable open state, opening the door for structure-function comparison. Here, we conducted comparative molecular dynamics simulation and electrophysiology studies to dissect the isoform-specific differences in Cx46 and Cx50 intercellular channel function. We show that key determinants Cx46 and Cx50 gap junction channel open stability and unitary conductance are shaped by structural and dynamic features of their N-terminal domains, in particular the residue at the 9th position and differences in hydrophobic anchoring sites. The results of this study establish the open state Cx46/50 structural models as archetypes for structure-function studies targeted at elucidating the mechanism of gap junction channels and the molecular basis of disease-causing variants. ABSTRACT: Connexins form intercellular communication channels, known as gap junctions (GJs), that facilitate diverse physiological roles, from long-range electrical and chemical coupling to coordinating development and nutrient exchange. GJs formed by different connexin isoforms harbour unique channel properties that have not been fully defined mechanistically. Recent structural studies on Cx46 and Cx50 defined a novel and stable open state and implicated the amino-terminal (NT) domain as a major contributor for isoform-specific functional differences between these closely related lens connexins. To better understand these differences, we constructed models corresponding to wildtype Cx50 and Cx46 GJs, NT domain swapped chimeras, and point variants at the 9th residue for comparative molecular dynamics (MD) simulation and electrophysiology studies. All constructs formed functional GJ channels, except the chimeric Cx46-50NT variant, which correlated with an introduced steric clash and increased dynamical behaviour (instability) of the NT domain observed by MD simulation. Single channel conductance correlated well with free-energy landscapes predicted by MD, but resulted in a surprisingly greater degree of effect. Additionally, we observed significant effects on transjunctional voltage-dependent gating (Vj gating) and/or open state dwell times induced by the designed NT domain variants. Together, these studies indicate intra- and inter-subunit interactions involving both hydrophobic and charged residues within the NT domains of Cx46 and Cx50 play important roles in defining GJ open state stability and single channel conductance, and establish the open state Cx46/50 structural models as archetypes for structure-function studies targeted at elucidating GJ channel mechanisms and the molecular basis of cataract-linked connexin variants.

A gentle introduction to the non-equilibrium physics of trajectories: Theory, algorithms, and biomolecular applications.

Despite the importance of non-equilibrium statistical mechanics in modern physics and related fields, the topic is often omitted from undergraduate and core-graduate curricula. Key aspects of non-equilibrium physics, however, can be understood with a minimum of formalism based on a rigorous trajectory picture. The fundamental object is the ensemble of trajectories, a set of independent time-evolving systems, which easily can be visualized or simulated (e.g., for protein folding) and which can be analyzed rigorously in analogy to an ensemble of static system configurations. The trajectory picture provides a straightforward basis for understanding first-passage times, "mechanisms" in complex systems, and fundamental constraints on the apparent reversibility of complex processes. Trajectories make concrete the physics underlying the diffusion and Fokker-Planck partial differential equations. Last but not least, trajectory ensembles underpin some of the most important algorithms that have provided significant advances in biomolecular studies of protein conformational and binding processes.

Should Virus Capsids Assemble Perfectly? Theory and Observation of Defects.

Although published structural models of viral capsids generally exhibit a high degree of regularity or symmetry, structural defects might be expected because of the fluctuating environment in which capsids assemble and the requirement of some capsids for disassembly before genome delivery. Defective structures are observed in computer simulations, and are evident in single-particle cryoelectron microscopy studies. Here, we quantify the conditions under which defects might be expected, using a statistical mechanics model allowing for ideal, defective, and vacant sites. The model displays a threshold in affinity parameters below which there is an appreciable population of defective capsids. Even when defective sites are not allowed, there is generally some population of vacancies. Analysis of single particles in cryoelectron microscopy micrographs yields a confirmatory ≳15% of defective particles. Our findings suggest structural heterogeneity in virus capsids may be under-appreciated, and also points to a nontraditional strategy for assembly inhibition.

Accelerated Estimation of Long-Timescale Kinetics from Weighted Ensemble Simulation via Non-Markovian "Microbin" Analysis.

The weighted ensemble (WE) simulation strategy provides unbiased sampling of nonequilibrium processes, such as molecular folding or binding, but the extraction of rate constants relies on characterizing steady-state behavior. Unfortunately, WE simulations of sufficiently complex systems will not relax to steady state on observed simulation times. Here, we show that a postsimulation clustering of molecular configurations into "microbins" using methods developed in the Markov State Model (MSM) community can yield unbiased kinetics from WE data before steady-state convergence of the WE simulation itself. Because WE trajectories are directional and not equilibrium distributed, the history-augmented MSM (haMSM) formulation can be used, which yields the mean first-passage time (MFPT) without bias for arbitrarily small lag times. Accurate kinetics can be obtained while bypassing the often prohibitive convergence requirements of the nonequilibrium weighted ensemble. We validate the method in a simple diffusive process on a two-dimensional (2D) random energy landscape and then analyze atomistic protein folding simulations using WE molecular dynamics. We report significant progress toward the unbiased estimation of protein folding times and pathways, though key challenges remain.

Connexin-46/50 in a dynamic lipid environment resolved by CryoEM at 1.9 Å.

Gap junctions establish direct pathways for cells to transfer metabolic and electrical messages. The local lipid environment is known to affect the structure, stability and intercellular channel activity of gap junctions; however, the molecular basis for these effects remains unknown. Here, we incorporate native connexin-46/50 (Cx46/50) intercellular channels into a dual lipid nanodisc system, mimicking a native cell-to-cell junction. Structural characterization by CryoEM reveals a lipid-induced stabilization to the channel, resulting in a 3D reconstruction at 1.9 Å resolution. Together with all-atom molecular dynamics simulations, it is shown that Cx46/50 in turn imparts long-range stabilization to the dynamic local lipid environment that is specific to the extracellular lipid leaflet. In addition, ~400 water molecules are resolved in the CryoEM map, localized throughout the intercellular permeation pathway and contributing to the channel architecture. These results illustrate how the aqueous-lipid environment is integrated with the architectural stability, structure and function of gap junction communication channels.

A systems-biology approach to molecular machines: Exploration of alternative transporter mechanisms.

Motivated by growing evidence for pathway heterogeneity and alternative functions of molecular machines, we demonstrate a computational approach for investigating two questions: (1) Are there multiple mechanisms (state-space pathways) by which a machine can perform a given function, such as cotransport across a membrane? (2) How can additional functionality, such as proofreading/error-correction, be built into machine function using standard biochemical processes? Answers to these questions will aid both the understanding of molecular-scale cell biology and the design of synthetic machines. Focusing on transport in this initial study, we sample a variety of mechanisms by employing Metropolis Markov chain Monte Carlo. Trial moves adjust transition rates among an automatically generated set of conformational and binding states while maintaining fidelity to thermodynamic principles and a user-supplied fitness/functionality goal. Each accepted move generates a new model. The simulations yield both single and mixed reaction pathways for cotransport in a simple environment with a single substrate along with a driving ion. In a "competitive" environment including an additional decoy substrate, several qualitatively distinct reaction pathways are found which are capable of extremely high discrimination coupled to a leak of the driving ion, akin to proofreading. The array of functional models would be difficult to find by intuition alone in the complex state-spaces of interest.

A kinetic mechanism for enhanced selectivity of membrane transport.

Membrane transport is generally thought to occur via an alternating access mechanism in which the transporter adopts at least two states, accessible from two different sides of the membrane to exchange substrates from the extracellular environment and the cytoplasm or from the cytoplasm and the intracellular matrix of the organelles (only in eukaryotes). In recent years, a number of high resolution structures have supported this general framework for a wide class of transport molecules, although additional states along the transport pathway are emerging as critically important. Given that substrate binding is often weak in order to enhance overall transport rates, there exists the distinct possibility that transporters may transport the incorrect substrate. This is certainly the case for many pharmaceutical compounds that are absorbed in the gut or cross the blood brain barrier through endogenous transporters. Docking studies on the bacterial sugar transporter vSGLT reveal that many highly toxic compounds are compatible with binding to the orthosteric site, further motivating the selective pressure for additional modes of selectivity. Motivated by recent work in which we observed failed substrate delivery in a molecular dynamics simulation where the energized ion still goes down its concentration gradient, we hypothesize that some transporters evolved to harness this 'slip' mechanism to increase substrate selectivity and reduce the uptake of toxic molecules. Here, we test this idea by constructing and exploring a kinetic transport model that includes a slip pathway. While slip reduces the overall productive flux, when coupled with a second toxic molecule that is more prone to slippage, the overall substrate selectivity dramatically increases, suppressing the accumulation of the incorrect compound. We show that the mathematical framework for increased substrate selectivity in our model is analogous to the classic proofreading mechanism originally proposed for tRNA synthase; however, because the transport cycle is reversible we identified conditions in which the selectivity is essentially infinite and incorrect substrates are exported from the cell in a 'detoxification' mode. The cellular consequences of proofreading and membrane slippage are discussed as well as the impact on future drug development.