Deep Free Energy Research Articles

In silico screening of LRRK2 WDR domain inhibitors using deep docking and free energy simulations.

The Critical Assessment of Computational Hit-Finding Experiments (CACHE) Challenge series is focused on identifying small molecule inhibitors of protein targets using computational methods. Each challenge contains two phases, hit-finding and follow-up optimization, each of which is followed by experimental validation of the computational predictions. For the CACHE Challenge #1, the Leucine-Rich Repeat Kinase 2 (LRRK2) WD40 Repeat (WDR) domain was selected as the target for in silico hit-finding and optimization. Mutations in LRRK2 are the most common genetic cause of the familial form of Parkinson's disease. The LRRK2 WDR domain is an understudied drug target with no known molecular inhibitors. Herein we detail the first phase of our winning submission to the CACHE Challenge #1. We developed a framework for the high-throughput structure-based virtual screening of a chemically diverse small molecule space. Hit identification was performed using the large-scale Deep Docking (DD) protocol followed by absolute binding free energy (ABFE) simulations. ABFEs were computed using an automated molecular dynamics (MD)-based thermodynamic integration (TI) approach. 4.1 billion ligands from Enamine REAL were screened with DD followed by ABFEs computed by MD TI for 793 ligands. 76 ligands were prioritized for experimental validation, with 59 compounds successfully synthesized and 5 compounds identified as hits, yielding a 8.5% hit rate. Our results demonstrate the efficacy of the combined DD and ABFE approaches for hit identification for a target with no previously known hits. This approach is widely applicable for the efficient screening of ultra-large chemical libraries as well as rigorous protein-ligand binding affinity estimation leveraging modern computational resources.

Deep Variational Free Energy Approach to Dense Hydrogen.

We developed a deep generative model-based variational free energy approach to the equations of state of dense hydrogen. We employ a normalizing flow network to model the proton Boltzmann distribution and a fermionic neural network to model the electron wave function at given proton positions. By jointly optimizing the two neural networks we reached a comparable variational free energy to the previous coupled electron-ion MonteCarlo calculation. The predicted equation of state of dense hydrogen under planetary conditions is denser than the findings of abinitio molecular dynamics calculation and empirical chemical model. Moreover, direct access to the entropy and free energy of dense hydrogen opens new opportunities in planetary modeling and high-pressure physics research.

Substituted Oligosaccharides as Protein Mimics: Deep Learning Free Energy Landscapes.

Protein-protein complexes power the majority of cellular processes. Interfering with the formation of such complexes using well-designed mimics is a difficult, yet actively pursued, research endeavor. Due to the limited availability of results on the conformational preferences of oligosaccharides compared to polypeptides, the former have been much less explored than the latter as protein mimics, despite interesting ADMET characteristics. In this work, the conformational landscapes of a series of 956 substituted glucopyranose oligomers of lengths 3 to 12 designed as protein interface mimics are revealed using microsecond-time-scale, enhanced-sampling molecular dynamics simulations. Deep convolutional networks are trained on these large conformational ensembles, to predict the stability of longer oligosaccharide structures from those of their constituent trimer motifs. Deep generative adversarial networks are then designed to suggest plausible conformations for oligosaccharide mimics of arbitrary length and substituent sequences that can subsequently be used as input to docking simulations. Analyzing the performance of the neural networks also yields insights into the intricate collective effects that dominate oligosaccharide conformational dynamics.

Comparing α-Synuclein Fibrils Formed in the Absence and Presence of a Model Lipid Membrane: A Small and Wide-Angle X-Ray Scattering Study

Amyloid fibrils are associated with a number of different neurodegenerative diseases. Detailed knowledge of the fibril structure will be of importance in the search of therapy and may guide experiments to understand amyloid formation. In this paper we investigate the morphology of α-synuclein amyloid fibrils, associated with Parkinson’s disease, formed under different conditions. In particular, we study, by means of small and wide-angle X-ray scattering, whether the presence of model lipid membranes affect the overall structure of the fibrils formed, motivated by the fact that amyloid fibrils in vivo are formed in a highly lipid-rich environment. Comparing fibrils formed in the presence of lipid with fibrils formed in their absence, show that the presence of lipids has no detectable effect on the fibril cross-section radius and that the characteristic β-strand repeat distance of 4.7 Å of the extended intermolecular β-sheets remains unaffected. We also show that the observed fibril radius is consistent with a fibril structure composed of two protofilaments. This indicates overall that the particular fibril structure, with their stacks of two-dimensionally folded α-synuclein molecules, represent a deep free energy minimum, not largely affected by the co-aggregation with lipids.

A combined approach to RNA secondary structure prediction based on deep learning and minimum free energy model

A combined approach to RNA secondary structure prediction based on deep learning and minimum free energy model

Conductance selectivity of Na+ across the K+ channel via Na+ trapped in a tortuous trajectory

Ion selectivity of the potassium channel is crucial for regulating electrical activity in living cells; however, the mechanism underlying the potassium channel selectivity that favors large K+ over small Na+ remains unclear. Generally, Na+ is not completely excluded from permeation through potassium channels. Herein, the distinct nature of Na+ conduction through the prototypical KcsA potassium channel was examined. Single-channel current recordings revealed that, at a high Na+ concentration (200 mM), the channel was blocked by Na+, and this blocking was relieved at high membrane potentials, suggesting the passage of Na+ across the channel. At a 2,000 mM Na+ concentration, single-channel Na+ conductance was measured as one-eightieth of the K+ conductance, indicating that the selectivity filter allows substantial conduit of Na+ Molecular dynamics simulations revealed unprecedented atomic trajectories of Na+ permeation. In the selectivity filter having a series of carbonyl oxygen rings, a smaller Na+ was distributed off-center in eight carbonyl oxygen-coordinated sites as well as on-center in four carbonyl oxygen-coordinated sites. This amphipathic nature of Na+ coordination yielded a continuous but tortuous path along the filter. Trapping of Na+ in many deep free energy wells in the filter caused slow elution. Conversely, K+ is conducted via a straight path, and as the number of occupied K+ ions increased to three, the concerted conduction was accelerated dramatically, generating the conductance selectivity ratio of up to 80. The selectivity filter allows accommodation of different ion species, but the ion coordination and interactions between ions render contrast conduction rates, constituting the potassium channel conductance selectivity.

Sampling efficiency of the counting method for permeability calculations estimated with the inhomogeneous solubility-diffusion model.

Permeability is a key property in various fields such as membrane technology for chemical separation and transport of substances through cell membranes. At the molecular scale, the counting method uses the number of membrane crossings in a conventional unbiased molecular dynamics simulation to predict the permeability. This contribution investigates under which conditions the counting method has insufficient statistics. An equation is derived for a compartmental model based on the inhomogeneous solubility-diffusion (Smoluchowski) model, giving insight into how the flux correlates with the solubility of permeants. This equation shows that a membrane crossing is a rare event not only when the membrane forms a large free energy barrier but also when the membrane forms a deep free energy well that traps permeants. Such a permeant trap has a high permeability; yet, the counting method suffers from poor statistics. To illustrate this, coarse-grained MD was run for 16 systems of dipalmitoylphosphatidylcholine bilayer membranes with different permeant types. The composition rule for permeability is shown to also hold for fluxes, and it is highlighted that the considered thickness of the membrane causes uncertainty in the permeability calculation of highly permeable membranes. In conclusion, a high permeability in itself is not an effective indicator of the sampling efficiency of the counting method, and caution should be taken for permeants whose solubility varies greatly over the simulation box. A practical consequence relevant in, e.g., drug design is that a drug with high membrane permeability might get trapped by membranes thus reducing its efficacy.

Movement of Arginine through OprD: The Energetics of Permeation and the Role of Lipopolysaccharide in Directing Arginine to the Protein.

The outer membrane channel OprD from Pseudomonas aeruginosa transports basic amino acids and clinically relevant carbapenem antibiotics. Understanding the molecular basis of substrate permeation across this channel will therefore lead to better therapeutic designs to treat infections. Using umbrella sampling simulations, we calculated the potential of mean force (PMF) for the arginine permeation pathway through OprD. The PMF reveals a deep free energy well of ∼6 kT around the putative substrate binding site followed by a shallower well of ∼4 kT close to the most constricted region of the pore. Despite becoming partially dehydrated during translocation, some water molecules are retained to shield the guanidinium side chain of arginine from the ladder of basic residues in the protein. Sugars of the lipopolysaccharide headgroups form contacts with arginine and could potentially play an important role in transferring substrate from the external medium to OprD. The PMF through bulk membrane shows a large energetic barrier of ∼45 kT within the hydrophobic core of the membrane, suggesting that spontaneous translocation without OprD is highly unlikely. This significant energetic penalty is likely caused by the extensive distortion of the lower leaflet of the outer membrane as phospholipid headgroups sink inwards to interact with charged groups of arginine. Our results provide quantitative insights into solute permeation across the bacterial outer membrane.

Deep Free Energy Principle for Adaptive Behavior

Deep Free Energy Principle for Adaptive Behavior

Combined Experimental and Molecular Simulation Investigation of the Individual Effects of Corexit Surfactants on the Aerosolization of Oil Spill Matter.

We report laboratory aerosolization experiments and classical molecular dynamics (MD) simulations, with the objective of investigating the individual effects of the two Corexit surfactants Span 80 (nonionic) and dioctyl sodium sulfosuccinate (DOSS, ionic), on the aerosolization of oil spill matter to the atmosphere. Our simulation results show that Span 80, DOSS, and the oil alkanes n-pentadecane (C15) and n-triacontane (C30) exhibit deep free energy minima at the air/seawater interface. C15 and C30 exhibit deeper free energy minima at the interface when Span 80 is present, as compared to the situation when DOSS or no surfactants are at the interface. These results suggest that Span 80 makes these oil hydrocarbons more likely to be adsorbed at the surface of seawater droplets and carried out to the atmosphere, relative to DOSS or to the situation where no surfactants are present. These simulation trends are in qualitative agreement with our experimental observations in a bubble-column setup, where larger amounts of oil hydrocarbons are ejected when Span 80 is mixed with oil and injected into the column, as compared to when DOSS is used. Our simulations also indicate that Span 80 has a larger thermodynamic incentive than DOSS to move from the seawater phase and into the air/seawater interface. This observation is also in qualitative agreement with our experimental measurements, which indicate that Span 80 is ejected in larger quantities than DOSS. Our simulations also suggest that DOSS predominantly adopts a perpendicular orientation with respect to the air/seawater interface at a dispersant to oil ratio (DOR) of 1:20, but has a slight preference to lie parallel to the interfaces at a DOR = 1:5; in both cases, DOSS molecules have their tails wide open and stretched. In contrast, Span 80 has a slight preference to align parallel to the interfaces with a coiled conformation at both DOR values.

Interlayer water regulates the bio-nano interface of a β-sheet protein stacking on graphene.

Using molecular dynamics simulations, we investigated an integrated bio-nano interface consisting of a β-sheet protein stacked onto graphene. We found that the stacking assembly of the model protein on graphene could be controlled by water molecules. The interlayer water filled within interstices of the bio-nano interface could suppress the molecular vibration of surface groups on protein, and could impair the CH···π interaction driving the attraction of the protein and graphene. The intermolecular coupling of interlayer water would be relaxed by the relative motion of protein upon graphene due to the interaction between water and protein surface. This effect reduced the hindrance of the interlayer water against the assembly of protein on graphene, resulting an appropriate adsorption status of protein on graphene with a deep free energy trap. Thereby, the confinement and the relative sliding between protein and graphene, the coupling of protein and water, and the interaction between graphene and water all have involved in the modulation of behaviors of water molecules within the bio-nano interface, governing the hindrance of interlayer water against the protein assembly on hydrophobic graphene. These results provide a deep insight into the fundamental mechanism of protein adsorption onto graphene surface in water.

Bubble bursting as an aerosol generation mechanism during an oil spill in the deep-sea environment: molecular dynamics simulations of oil alkanes and dispersants in atmospheric air/salt water interfaces

Potential of mean force (PMF) calculations and molecular dynamics (MD) simulations were performed to investigate the properties of oil n-alkanes [i.e., n-pentadecane (C15), n-icosane (C20) and n-triacontane (C30)], as well as several surfactant species [i.e., the standard anionic surfactant sodium dodecyl sulfate (SDS), and three model dispersants similar to the Tween and Span species present in Corexit 9500A] at air/salt water interfaces. This study was motivated by the 2010 Deepwater Horizon (DWH) oil spill, and our simulation results show that, from the thermodynamic point of view, the n-alkanes and the model dispersants have a strong preference to remain at the air/salt water interface, as indicated by the presence of deep free energy minima at these interfaces. The free energy minimum of these n-alkanes becomes deeper as their chain length increases, and as the concentration of surfactant species at the interface increases. The n-alkanes tend to adopt a flat orientation and form aggregates at the bare air/salt water interface. When this interface is coated with surfactants, the n-alkanes tend to adopt more tilted orientations with respect to the vector normal to the interface. These simulation results are consistent with the experimental findings reported in the accompanying paper [Ehrenhauser et al., Environ. Sci.: Processes Impacts 2013, in press, (DOI: 10.1039/c3em00390f)]. The fact that these long-chain n-alkanes show a strong thermodynamic preference to remain at the air/salt water interfaces, especially if these interfaces are coated with surfactants, makes these species very likely to adsorb at the surface of bubbles or droplets and be ejected to the atmosphere by sea surface processes such as whitecaps (breaking waves) and bubble bursting. Finally, the experimental finding that more oil hydrocarbons are ejected when Corexit 9500A is present in the system is consistent with the deeper free energy minima observed for the n-alkanes at the air/salt water interface at increasing concentrations of surfactant species.

Molecular Modeling of the Green Leaf Volatile Methyl Salicylate on Atmospheric Air/Water Interfaces

Methyl salicylate (MeSA) is a green leaf volatile (GLV) compound that is emitted in significant amounts by plants, especially when they are under stress conditions. GLVs can then undergo chemical reactions with atmospheric oxidants, yielding compounds that contribute to the formation of secondary organic aerosols (SOAs). We investigated the adsorption of MeSA on atmospheric air/water interfaces at 298 K using thermodynamic integration (TI), potential of mean force (PMF) calculations, and classical molecular dynamics (MD) simulations. Our molecular models can reproduce experimental results of the 1-octanol/water partition coefficient of MeSA. A deep free energy minimum was found for MeSA at the air/water interface, which is mainly driven by energetic interactions between MeSA and water. At the interface, the oxygenated groups in MeSA tend to point toward the water side of the interface, with the aromatic group of MeSA lying farther away from water. Increases in the concentrations of MeSA lead to reductions in the height of the peaks in the MeSA-MeSA g(r) functions, a slowing down of the dynamics of both MeSA and water at the interface, and a reduction in the interfacial surface tension. Our results indicate that MeSA has a strong thermodynamic preference to remain at the air/water interface, and thus chemical reactions with atmospheric oxidants are more likely to take place at this interface, rather than in the water phase of atmospheric water droplets or in the gas phase.

Ice Growth from Supercooled Aqueous Solutions of Benzene, Naphthalene, and Phenanthrene

Classical molecular dynamics (MD) were performed to investigate the growth of ice from supercooled aqueous solutions of benzene, naphthalene, or phenanthrene. The main objective of this study is to explore the fate of those aromatic molecules after freezing of the supercooled aqueous solutions, i.e., if these molecules become trapped inside the ice lattice or if they are displaced to the QLL or to the interface with air. Ice growth from supercooled aqueous solutions of benzene, naphthalene, or phenanthrene result in the formation of quasi-liquid layers (QLLs) at the air/ice interface that are thicker than those observed when pure supercooled water freezes. Naphthalene and phenanthrene molecules in the supercooled aqueous solutions are displaced to the air/ice interface during the freezing process at both 270 and 260 K; no incorporation of these aromatics into the ice lattice is observed throughout the freezing process. Similar trends were observed during freezing of supercooled aqueous solutions of benzene at 270 K. In contrast, a fraction of the benzene molecules become trapped inside the ice lattice during the freezing process at 260 K, with the rest of the benzene molecules being displaced to the air/ice interface. These results suggest that the size of the aromatic molecule in the supercooled aqueous solution is an important parameter in determining whether these molecules become trapped inside the ice crystals. Finally, we also report potential of mean force (PMF) calculations aimed at studying the adsorption of gas-phase benzene and phenanthrene on atmospheric air/ice interfaces. Our PMF calculations indicate the presence of deep free energy minima for both benzene and phenanthrene at the air/ice interface, with these molecules adopting a flat orientation at the air/ice interface.

Simulations of anion transport through OprP reveal the molecular basis for high affinity and selectivity for phosphate

The outer membrane protein OprP from Pseudomonas aeruginosa forms a phosphate selective pore. To understand the mechanism of phosphate permeation and selectivity, we used three simulation techniques [equilibrium molecular dynamics simulations, steered molecular dynamics, and calculation of a potential of mean force (PMF)]. The PMF for phosphate reveals a deep free energy well midway along the OprP channel. Two adjacent phosphate-binding sites (W1 and W2), each with a well depth of approximately 8 kT, are identified close to the L3 loop in the most constricted region of the pore. A dissociation constant for phosphate of 6 microM is computed from the PMF, within the range of reported experimental values. The transfer of phosphate between sites W1 and W2 is correlated with changes in conformation of the sidechain of K121, which serves as a "charged brush" to facilitate phosphate passage between the two subsites. OprP also binds chloride, but less strongly than phosphate, as calculated from a Cl(-) PMF. The difference in affinity and hence selectivity is due to the "Lys-cluster" motif, the positive charges of which interact strongly with a partially dehydrated phosphate ion but are shielded from a Cl(-) by the hydration shell of the smaller ion. Our simulations suggest that OprP does not conform to the conventional picture of a channel with relatively flat energy landscape for permeant ions, but rather resembles a membrane-inserted binding protein with a high specificity that allows access to a centrally located binding site from both the extracellular and the periplasmic spaces.

Academic Life: The Whole Package

When I received the unexpected but welcome invitation to write a Reflections article, I wondered which aspects of my life in science are sufficiently distinctive to make the article worth reading. One aspect may be my survival as a department chair for nearly 25 years while struggling to contribute as well in teaching and research, in other words, embracing the “whole package” of academic life: teaching, research, and service. Another part may be my experience as coauthor of a biochemistry textbook. To deal with both, I will begin with my introduction to biochemistry.

Adsorption of polycyclic aromatic hydrocarbons at the air–water interface: Molecular dynamics simulations and experimental atmospheric observations

Adsorption of benzene, naphthalene, anthracene, and phenanthrene at the aqueous surface is investigated by means of molecular dynamics simulations. Potentials of mean force, i.e., free energy profiles obtained when moving the studied molecules across an aqueous slab were evaluated. In all cases, deep surface free energy minima, corresponding to orders of magnitude of surface enhancement of the aromatic molecule, were located. This enhancement, which increases with the size of the solute, points to the importance of the aqueous surface for the chemistry of polycyclic aromatic hydrocarbons (PAHs). Supporting evidence in the atmospheric environment related to the heterogeneous chemistry of PAHs on water droplets and planar surfaces is summarized. There is good agreement between the hydration free energies computed from MD calculations and the experimentally determined values. Data pertaining to the importance of the air-water interface in the adsorption and transport of PAHs on micron sized water droplets are described. The relevant data on adsorption and reaction (ozonation and photochemical) at the air-water interface of planar surfaces and droplets are also summarized.

Use of Strain in a Stereospecific Catalytic Mechanism: Crystal Structures of Escherichia coli Thymidylate Synthase Bound to FdUMP and Methylenetetrahydrofolate,

Two crystal structures for E. coli thymidylate synthase (TS) bound to the mechanism-based inhibitor 5-fluoro-dUMP (FdUMP) and methylenetetrahydrofolate (CH2THF) have been determined to 2.6 and 2.2 A nominal resolutions, with crystallographic R factors of 0.180 and 0.178, respectively. The inhibitor and cofactor are well ordered in both structures and display covalent links to each other and to Cys 146 in the TS active site. The structures are in general agreement with a previous report for this complex (D. A. Matthews et al. (1990) J. Mol. Biol. 214, 937-948), but differ in two key respects: (i) the methylene bridge linking FdUMP and CH2THF is rotated about 60 degrees to a different position and (ii) the electron density for C6 of FdUMP, which is covalently linked to Cys 146, is more diffuse than for the other atoms in the pyrimidine ring. The ligand arrangement observed in the previous structure led the authors to propose that a large conformational change in ligand geometry must occur in order to facilitate catalysis and yield the correct chirality in the methyl of product dTMP. The new structures suggest a different mechanism for product formation that does not require ligands to greatly alter their conformations during catalysis and which makes use of instability in the nucleotide-Cys 146 thiol adduct to avoid a deep free energy well and assist in proton abstraction from dUMP. All intermediates in the proposed mechanism were modeled and energy minimized in the TS active site, and all can be accommodated in the present structures. The role of ligand-induced conformational change in the TS mechanism and the possibility of Tyr 94 acting as a base during catalysis are also discussed.