Kartograf: A Geometrically Accurate Atom Mapper for Hybrid-Topology Relative Free Energy Calculations.

Relative binding free energy (RBFE) calculations have emerged as a powerful tool that supports ligand optimization in drug discovery. Despite many successes, the use of RBFEs can often be limited by automation problems, in particular, the setup of such calculations. Atom mapping algorithms are an essential component in setting up automatic large-scale hybrid-topology RBFE calculation campaigns. Traditional algorithms typically employ a 2D subgraph isomorphism solver (SIS) in order to estimate the maximum common substructure. SIS-based approaches can be limited by time-intensive operations and issues with capturing geometry-linked chemical properties, potentially leading to suboptimal solutions. To overcome these limitations, we have developed Kartograf, a geometric-graph-based algorithm that uses primarily the 3D coordinates of atoms to find a mapping between two ligands. In free energy approaches, the ligand conformations are usually derived from docking or other previous modeling approaches, giving the coordinates a certain importance. By considering the spatial relationships between atoms related to the molecule coordinates, our algorithm bypasses the computationally complex subgraph matching of SIS-based approaches and reduces the problem to a much simpler bipartite graph matching problem. Moreover, Kartograf effectively circumvents typical mapping issues induced by molecule symmetry and stereoisomerism, making it a more robust approach for atom mapping from a geometric perspective. To validate our method, we calculated mappings with our novel approach using a diverse set of small molecules and used the mappings in relative hydration and binding free energy calculations. The comparison with two SIS-based algorithms showed that Kartograf offers a fast alternative approach. The code for Kartograf is freely available on GitHub (https://github.com/OpenFreeEnergy/kartograf). While developed for the OpenFE ecosystem, Kartograf can also be utilized as a standalone Python package.

Path Sampling Simulations Reveal How the Q61L Mutation Alters the Dynamics of KRas.

Flexibility is essential for many proteins to function, but can be difficult to characterize. Experiments lack resolution in space and time, while the time scales involved are prohibitively long for straightforward molecular dynamics simulations. In this work, we present a multiple state transition path sampling simulation study of a protein that has been notoriously difficult to characterize in its active state. The GTPase enzyme KRas is a signal transduction protein in pathways for cell differentiation, growth, and division. When active, KRas tightly binds guanosine triphosphate (GTP) in a rigid state. The protein-GTP complex can also visit more flexible states, in which it is not active. KRas mutations can affect the conversion between these rigid and flexible states, thus prolonging the activation of signal transduction pathways, which may result in tumor formation. In this work, we apply path sampling simulations to investigate the dynamic behavior of KRas-4B (wild type, WT) and the oncogenic mutant Q61L (Q61L). Our results show that KRas visits several conformational states, which are the same for WT and Q61L. The multiple state transition path sampling (MSTPS) method samples transitions between the different states in a single calculation. Tracking which transitions occur shows large differences between WT and Q61L. The MSTPS results further reveal that for Q61L, a route to a more flexible state is inaccessible, thus shifting the equilibrium to more rigid states. The methodology presented here enables a detailed characterization of protein flexibility on time scales not accessible with brute-force molecular dynamics simulations.

Atomistic insight into the kinetic pathways for Watson-Crick to Hoogsteen transitions in DNA.

DNA predominantly contains Watson-Crick (WC) base pairs, but a non-negligible fraction of base pairs are in the Hoogsteen (HG) hydrogen bonding motif at any time. In HG, the purine is rotated ∼180° relative to the WC motif. The transitions between WC and HG may play a role in recognition and replication, but are difficult to investigate experimentally because they occur quickly, but only rarely. To gain insight into the mechanisms for this process, we performed transition path sampling simulations on a model nucleotide sequence in which an AT pair changes from WC to HG. This transition can occur in two ways, both starting with loss of hydrogen bonds in the base pair, followed by rotation around the glycosidic bond. In one route the adenine base converts from WC to HG geometry while remaining entirely within the double helix. The other route involves the adenine leaving the confines of the double helix and interacting with water. Our results indicate that this outside route is more probable. We used transition interface sampling to compute rate constants and relative free energies for the transitions between WC and HG. Our results agree with experiments, and provide highly detailed insights into the mechanisms of this important process.

OpenPathSampling: A Python Framework for Path Sampling Simulations. 1. Basics.

Transition path sampling techniques allow molecular dynamics simulations of complex systems to focus on rare dynamical events, providing insight into mechanisms and the ability to calculate rates inaccessible by ordinary dynamics simulations. While path sampling algorithms are conceptually as simple as importance sampling Monte Carlo, the technical complexity of their implementation has kept these techniques out of reach of the broad community. Here, we introduce an easy-to-use Python framework called OpenPathSampling (OPS) that facilitates path sampling for (bio)molecular systems with minimal effort and yet is still extensible. Interfaces to OpenMM and an internal dynamics engine for simple models are provided in the initial release, but new molecular simulation packages can easily be added. Multiple ready-to-use transition path sampling methodologies are implemented, including standard transition path sampling (TPS) between reactant and product states and transition interface sampling (TIS) and its replica exchange variant (RETIS), as well as recent multistate and multiset extensions of transition interface sampling (MSTIS, MISTIS). In addition, tools are provided to facilitate the implementation of new path sampling schemes built on basic path sampling components. In this paper, we give an overview of the design of this framework and illustrate the simplicity of applying the available path sampling algorithms to a variety of benchmark problems.

OpenPathSampling: A Python Framework for Path Sampling Simulations. 2. Building and Customizing Path Ensembles and Sample Schemes.

The OpenPathSampling (OPS) package provides an easy-to-use framework to apply transition path sampling methodologies to complex molecular systems with a minimum of effort. Yet, the extensibility of OPS allows for the exploration of new path sampling algorithms by building on a variety of basic operations. In a companion paper [ Swenson et al. J. Chem. Theory Comput. 2018 , 10.1021/acs.jctc.8b00626 ] we introduced the basic concepts and the structure of the OPS package, and how it can be employed to perform standard transition path sampling and (replica exchange) transition interface sampling. In this paper, we elaborate on two theoretical developments that went into the design of OPS. The first development relates to the construction of path ensembles, the what is being sampled. We introduce a novel set-based notation for the path ensemble, which provides an alternative paradigm for constructing path ensembles and allows building arbitrarily complex path ensembles from fundamental ones. The second fundamental development is the structure for the customization of Monte Carlo procedures; how path ensembles are being sampled. We describe in detail the OPS objects that implement this approach to customization, the MoveScheme and the PathMover, and provide tools to create and manipulate these objects. We illustrate both the path ensemble building and sampling scheme customization with several examples. OPS thus facilitates both standard path sampling application in complex systems as well as the development of new path sampling methodology, beyond the default.

The opposing effects of isotropic and anisotropic attraction on association kinetics of proteins and colloids.

The association and dissociation of particles via specific anisotropic interactions is a fundamental process, both in biology (proteins) and in soft matter (colloidal patchy particles). The presence of alternative binding sites can lead to multiple productive states and also to non-productive "decoy" or intermediate states. Besides anisotropic interactions, particles can experience non-specific isotropic interactions. We employ single replica transition interface sampling to investigate how adding a non-productive binding site or a nonspecific isotropic interaction alters the dimerization kinetics of a generic patchy particle model. The addition of a decoy binding site reduces the association rate constant, independent of the site's position, while adding an isotropic interaction increases it due to an increased rebinding probability. Surprisingly, the association kinetics becomes non-monotonic for a tetramer complex formed by multivalent patchy particles. While seemingly identical to two-particle binding with a decoy state, the cooperativity of binding multiple particles leads to a kinetic optimum. Our results are relevant for the understanding and modeling of biochemical networks and self-assembly processes.

A replica exchange transition interface sampling method with multiple interface sets for investigating networks of rare events.

The multiple state transition interface sampling (TIS) framework in principle allows the simulation of a large network of complex rare event transitions, but in practice suffers from convergence problems. To improve convergence, we combine multiple state TIS [J. Rogal and P. G. Bolhuis, J. Chem. Phys. 129, 224107 (2008)] with replica exchange TIS [T. S. van Erp, Phys. Rev. Lett. 98, 268301 (2007)]. In addition, we introduce multiple interface sets, which allow more than one order parameter to be defined for each state. We illustrate the methodology on a model system of multiple independent dimers, each with two states. For reaction networks with up to 64 microstates, we determine the kinetics in the microcanonical ensemble, and discuss the convergence properties of the sampling scheme. For this model, we find that the kinetics depend on the instantaneous composition of the system. We explain this dependence in terms of the system's potential and kinetic energy.

A Cartesian quasi-classical model to nonequilibrium quantum transport: the Anderson impurity model.

We apply the recently proposed quasi-classical approach for a second quantized many-electron Hamiltonian in Cartesian coordinates [B. Li and W. H. Miller, J. Chem. Phys. 137, 154107 (2012)] to correlated nonequilibrium quantum transport. The approach provides accurate results for the resonant level model for a wide range of temperatures, bias, and gate voltages, correcting the flaws of our recently proposed mapping using action-angle variables. When electron-electron interactions are included, a Gaussian function scheme is required to map the two-electron integrals, leading to quantitative results for the Anderson impurity model. In particular, we show that the current mapping is capable of capturing quantitatively the Coulomb blockade effect and the temperature dependence of the current below and above the blockade.

Application of a semiclassical model for the second-quantized many-electron Hamiltonian to nonequilibrium quantum transport: the resonant level model.

A semiclassical approach is developed for nonequilibrium quantum transport in molecular junctions. Following the early work of Miller and White [J. Chem. Phys. 84, 5059 (1986)], the many-electron Hamiltonian in second quantization is mapped onto a classical model that preserves the fermionic character of electrons. The resulting classical electronic Hamiltonian allows for real-time molecular dynamics simulations of the many-body problem from an uncorrelated initial state to the steady state. Comparisons with exact results generated for the resonant level model reveal that a semiclassical treatment of transport provides a quantitative description of the dynamics at all relevant timescales for a wide range of bias and gate potentials, and for different temperatures. The approach opens a door to treating nontrivial quantum transport problems that remain far from the reach of fully quantum methodologies.

Remarkable features in the interactions of quadrupolar molecules.

Quadrupolar charge fields of molecules and of molecular fragments give rise to unique features in weakly interacting clusters and aggregations. Relative to dipole-dipole interactions, the interactions among quadrupolar molecules tend to allow for greater orientational distortions away from equilibrium. Potential surface regions have been found for several clusters that are attractive and yet very flat for certain directions. There is a notable slipperiness for the interactions in some of these cases. This implies significant vibrational excursions even in the ground state. Furthermore, the coupling of rotations among nearby molecules in pure clusters of quadrupolar molecules is different than for dipolar species, and it can lead to unexpectedly small internal rotation barriers. How these and other features develop and what they might imply for materials and biomolecular simulations are discussed here.