GPU-specific algorithms for improved solute sampling in grand canonical Monte Carlo simulations.

The Grand Canonical Monte Carlo (GCMC) ensemble defined by the excess chemical potential, µex , volume, and temperature, in the context of molecular simulations allows for variations in the number of particles in the system. In practice, GCMC simulations have been widely applied for the sampling of rare gasses and water, but limited in the context of larger molecules. To overcome this limitation, the oscillating µex GCMC method was introduced and shown to be of utility for sampling small solutes, such as formamide, propane, and benzene, as well as for ionic species such as monocations, acetate, and methylammonium. However, the acceptance of GCMC insertions is low, and the method is computationally demanding. In the present study, we improved the sampling efficiency of the GCMC method using known cavity-bias and configurational-bias algorithms in the context of GPU architecture. Specifically, for GCMC simulations of aqueous solution systems, the configurational-bias algorithm was extended by applying system partitioning in conjunction with a random interval extraction algorithm, thereby improving the efficiency in a highly parallel computing environment. The method is parallelized on the GPU using CUDA and OpenCL, allowing for the code to run on both Nvidia and AMD GPUs, respectively. Notably, the method is particularly well suited for GPU computing as the large number of threads allows for simultaneous sampling of a large number of configurations during insertion attempts without additional computational overhead. In addition, the partitioning scheme allows for simultaneous insertion attempts for large systems, offering considerable efficiency. Calculations on the BK Channel, a transporter, including a lipid bilayer with over 760,000 atoms, show a speed up of ~53-fold through the use of system partitioning. The improved algorithm is then combined with an enhanced µex oscillation protocol and shown to be of utility in the context of the site-identification by ligand competitive saturation (SILCS) co-solvent sampling approach as illustrated through application to the protein CDK2.

Extension of the CHARMM Classical Drude Polarizable Force Field to N- and O-Linked Glycopeptides and Glycoproteins.

Molecular dynamic simulations are an effective tool to study complex molecular systems and are contingent upon the availability of an accurate and reliable molecular mechanics force field. The Drude polarizable force field, which allows for the explicit treatment of electronic polarization in a computationally efficient fashion, has been shown to reproduce experimental properties that were difficult or impossible to reproduce with the CHARMM additive force field, including peptide folding cooperativity, RNA hairpin structures, and DNA base flipping. Glycoproteins are essential components of glycoconjugate vaccines, antibodies, and many pharmaceutically important molecules, and an accurate polarizable force field that includes compatibility between the protein and carbohydrate aspect of the force field is essential to study these types of systems. In this work, we present an extension of the Drude polarizable force field to glycoproteins, including both N- and O-linked species. Parameter optimization focused on the dihedral terms using a reweighting protocol targeting NMR solution J-coupling data for model glycopeptides. Validation of the model include eight model glycopeptides and four glycoproteins with multiple N- and O-linked glycosylations. The new glycoprotein carbohydrate force field can be used in conjunction with the remainder of Drude polarizable force field through a variety of MD simulation programs including GROMACS, OPENMM, NAMD, and CHARMM and may be accessed through the Drude Prepper module in the CHARMM-GUI.

SILCS-RNA: Toward a Structure-Based Drug Design Approach for Targeting RNAs with Small Molecules.

RNA molecules can act as potential drug targets in different diseases, as their dysregulated expression or misfolding can alter various cellular processes. Noncoding RNAs account for ∼70% of the human genome, and these molecules can have complex tertiary structures that present a great opportunity for targeting by small molecules. In the present study, the site identification by ligand competitive saturation (SILCS) computational approach is extended to target RNA, termed SILCS-RNA. Extensions to the method include an enhanced oscillating excess chemical potential protocol for the grand canonical Monte Carlo calculations and individual simulations of the neutral and charged solutes from which the SILCS functional group affinity maps (FragMaps) are calculated for subsequent binding site identification and docking calculations. The method is developed and evaluated against seven RNA targets and their reported small molecule ligands. SILCS-RNA provides a detailed characterization of the functional group affinity pattern in the small molecule binding sites, recapitulating the types of functional groups present in the ligands. The developed method is also shown to be useful for identification of new potential binding sites and identifying ligand moieties that contribute to binding, granular information that can facilitate ligand design. However, limitations in the method are evident including the ability to map the regions of binding sites occupied by ligand phosphate moieties and to fully account for the wide range of conformational heterogeneity in RNA associated with binding of different small molecules, emphasizing inherent challenges associated with applying computer-aided drug design methods to RNA. While limitations are present, the current study indicates how the SILCS-RNA approach may enhance drug discovery efforts targeting RNAs with small molecules.

CHARMM-GUI Drude prepper for molecular dynamics simulation using the classical Drude polarizable force field.

Explicit treatment of electronic polarizability in empirical force fields (FFs) represents an extension over a traditional additive or pairwise FF and provides a more realistic model of the variations in electronic structure in condensed phase, macromolecular simulations. To facilitate utilization of the polarizable FF based on the classical Drude oscillator model, Drude Prepper has been developed in CHARMM-GUI. Drude Prepper ingests additive CHARMM protein structures file (PSF) and pre-equilibrated coordinates in CHARMM, PDB, or NAMD format, from which the molecular components of the system are identified. These include all residues and patches connecting those residues along with water, ions, and other solute molecules. This information is then used to construct the Drude FF-based PSF using molecular generation capabilities in CHARMM, followed by minimization and equilibration. In addition, inputs are generated for molecular dynamics (MD) simulations using CHARMM, GROMACS, NAMD, and OpenMM. Validation of the Drude Prepper protocol and inputs is performed through conversion and MD simulations of various heterogeneous systems that include proteins, nucleic acids, lipids, polysaccharides, and atomic ions using the aforementioned simulation packages. Stable simulations are obtained in all studied systems, including 5 µs simulation of ubiquitin, verifying the integrity of the generated Drude PSFs. In addition, the ability of the Drude FF to model variations in electronic structure is shown through dipole moment analysis in selected systems. The capabilities and availability of Drude Prepper in CHARMM-GUI is anticipated to greatly facilitate the application of the Drude FF to a range of condensed phase, macromolecular systems.

Contributions and competition of Mg2+ and K+ in folding and stabilization of the Twister ribozyme.

Native folded and compact intermediate states of RNA typically involve tertiary structures in the presence of divalent ions such as Mg2+ in a background of monovalent ions. In a recent study, we have shown how the presence of Mg2+ impacts the transition from partially unfolded to folded states through a "push-pull" mechanism where the ion both favors and disfavors the sampling of specific phosphate-phosphate interactions. To further understand the ion atmosphere of RNA in folded and partially folded states results from atomistic umbrella sampling and oscillating chemical potential grand canonical Monte Carlo/molecular dynamics (GCMC/MD) simulations are used to obtain atomic-level details of the distributions of Mg2+ and K+ ions around Twister RNA. Results show the presence of 100 mM Mg2+ to lead to increased charge neutralization over that predicted by counterion condensation theory. Upon going from partially unfolded to folded states, overall charge neutralization increases at all studied ion concentrations that, while associated with an increase in the number of direct ion-phosphate interactions, is fully accounted for by the monovalent K+ ions. Furthermore, K+ preferentially interacts with purine N7 atoms of helical regions in partially unfolded states, thereby potentially stabilizing the helical regions. Thus, both secondary helical structures and formation of tertiary structures leads to increased counterion condensation, thereby stabilizing those structural features of Twister. Notably, it is shown that K+ can act as a surrogate for Mg2+ by participating in specific interactions with nonsequential phosphate pairs that occur in the folded state, explaining the ability of Twister to self-cleave at submillimolar Mg2+ concentrations.

Balanced polarizable Drude force field parameters for molecular anions: phosphates, sulfates, sulfamates, and oxides.

Polarizable force fields are emerging as a more accurate alternative to additive force fields in terms of modeling and simulations of a variety of chemicals including biomolecules. Explicit treatment of induced polarization in charged species such as phosphates and sulfates offers the potential for achieving an improved atomistic understanding of the physical forces driving their interactions with their environments. To help achieve this, in this study we present balanced Drude polarizable force field parameters for molecular ions including phosphates, sulfates, sulfamates, and oxides. Better balance was primarily achieved in the relative values of minimum interaction energies and distances of the anionic model compounds with water at the Drude and quantum mechanical (QM) model chemistries. Parametrization involved reoptimizing available parameters as well as extending the force field to new molecules with the goal of achieving self-consistency with respect to the Lennard-Jones and electrostatic parameters targeting QM and experimental hydration free energies. The resulting force field parameters achieve consistent treatment across the studied anions, facilitating more balanced simulations of biomolecules and small organic molecules in the context of the classical Drude polarizable force field. Graphical abstract.

Mg2+ Impacts the Twister Ribozyme through Push-Pull Stabilization of Nonsequential Phosphate Pairs.

RNA molecules perform a variety of biological functions for which the correct three-dimensional structure is essential, including as ribozymes where they catalyze chemical reactions. Metal ions, especially Mg2+, neutralize these negatively charged nucleic acids and specifically stabilize RNA tertiary structures as well as impact the folding landscape of RNAs as they assume their tertiary structures. Specific binding sites of Mg2+ in folded conformations of RNA have been studied extensively; however, the full range of interactions of the ion with compact intermediates and unfolded states of RNA is challenging to investigate, and the atomic details of the mechanism by which the ion facilitates tertiary structure formation is not fully known. Here, umbrella sampling combined with oscillating chemical potential Grand Canonical Monte Carlo/molecular dynamics simulations are used to capture the energetics and atomic-level details of Mg2+-RNA interactions that occur along an unfolding pathway of the Twister ribozyme. The free energy profiles reveal stabilization of partially unfolded states by Mg2+, as observed in unfolding experiments, with this stabilization being due to increased sampling of simultaneous interactions of Mg2+ with two or more nonsequential phosphate groups. Notably, these results indicate a push-pull mechanism in which the Mg2+-RNA interactions actually lead to destabilization of specific nonsequential phosphate-phosphate interactions (i.e., pushed apart), whereas other interactions are stabilized (i.e., pulled together), a balance that stabilizes unfolded states and facilitates the folding of Twister, including the formation of hydrogen bonds associated with the tertiary structure. This study establishes a better understanding of how Mg2+-ion interactions contribute to RNA structural properties and stability.

Predicting Partition Coefficients of Neutral and Charged Solutes in the Mixed SLES-Fatty Acid Micellar System.

Sodium laureth sulfate (SLES) and fatty acids are common ingredients in many cosmetic products. Understanding how neutral and charged fatty acid compounds partition between micellar and water phases is crucial to achieve the optimal design of the product formulation. In this paper, we first study the formation of mixed SLES and fatty acid micelles using molecular dynamics (MD) simulations. Micelle/water partition coefficients of neutral and charged fatty acids are then calculated using COSMOmic as well as a MD approach based on the potential of mean force (PMF) calculations performed using umbrella sampling (US). The combined US/PMF approach was performed with both the additive, non-polarizable CHARMM general force field (CGenFF) and the classical Drude polarizable force field. The partition coefficients for the neutral solutes are shown to be accurately calculated with the COSMOmic and additive CGenFF US/PMF approaches, while only the US/PMF approach with the Drude polarizable force field accurately calculated the experimental partition coefficient of the charged solute. These results indicate the utility of the Drude polarizable force field as a tool for the rational development of mixed micelles.

Cellulose-specific Type B carbohydrate binding modules: understanding oligomeric and non-crystalline substrate recognition mechanisms.

BACKGROUND: Effective enzymatic degradation of crystalline polysaccharides requires a synergistic cocktail of hydrolytic enzymes tailored to the wide-ranging degree of substrate crystallinity. To accomplish this type of targeted carbohydrate recognition, nature produces multi-modular enzymes, having at least one catalytic domain appended to one or more carbohydrate binding modules (CBMs). The Type B CBM categorization encompasses several families (i.e., protein folds) of CBMs that are generally thought to selectively bind oligomeric polysaccharides; however, a subset of cellulose-specific CBM families (17 and 28) appear to bind non-crystalline cellulose more tightly than oligomers and in a manner that discriminates between surface topology. RESULTS: To provide insight into this unexplained phenomenon, we investigated the molecular-level origins of oligomeric and non-crystalline carbohydrate recognition in cellulose-specific Type B CBMs using molecular dynamics (MD) simulation and free energy calculations. Examining two CBMs from three different families (4, 17, and 28), we describe how protein-ligand dynamics contribute to observed variations in binding affinity of oligomers within the same CBM family. Comparisons across the three CBM families identified factors leading to modified functionality prohibiting competitive binding, despite similarity in sequence and specificity. Using free energy perturbation with Hamiltonian replica exchange MD, we also examined the hypothesis that the open topology of the binding grooves in families 17 and 28 necessitates tight binding of an oligomer, while the more confined family 4 binding groove does not require the same degree of tight binding. Finally, we elucidated the mechanisms of non-crystalline carbohydrate recognition by modeling CBMs complexed with a partially decrystallized cellulose substrate. Molecular simulation provided structural and dynamic data for direct comparison to oligomeric modes of carbohydrate recognition, and umbrella sampling MD was used to determine ligand binding free energy. Comparing both protein-carbohydrate interactions and ligand binding free energies, which were in good agreement with experimental values, we confirmed the hypothesis that family 17 and 28 CBMs bind non-crystalline cellulose and oligomers with different affinities (i.e., high and low). CONCLUSIONS: Our study provides an unprecedented level of insight into the complex solid and soluble carbohydrate substrate recognition mechanisms of Type B CBMs, the findings of which hold considerable promise for enhancing lignocellulosic biomass conversion technology and development of plant cell wall probes.

Structural and molecular dynamics studies of a C1-oxidizing lytic polysaccharide monooxygenase from Heterobasidion irregulare reveal amino acids important for substrate recognition.

Lytic polysaccharide monooxygenases (LPMOs) are a group of recently discovered enzymes that play important roles in the decomposition of recalcitrant polysaccharides. Here, we report the biochemical, structural, and computational characterization of an LPMO from the white-rot fungus Heterobasidion irregulare (HiLPMO9B). This enzyme oxidizes cellulose at the C1 carbon of glycosidic linkages. The crystal structure of HiLPMO9B was determined at 2.1 Å resolution using X-ray crystallography. Unlike the majority of the currently available C1-specific LPMO structures, the HiLPMO9B structure contains an extended L2 loop, connecting ß-strands ß2 and ß3 of the ß-sandwich structure. Molecular dynamics (MD) simulations suggest roles for both aromatic and acidic residues in the substrate binding of HiLPMO9B, with the main contribution from the residues located on the extended region of the L2 loop (Tyr20) and the LC loop (Asp205, Tyr207, and Glu210). Asp205 and Glu210 were found to be involved in the hydrogen bonding with the hydroxyl group of the C6 carbon of glucose moieties directly or via a water molecule. Two different binding orientations were observed over the course of the MD simulations. In each orientation, the active-site copper of this LPMO preferentially skewed toward the pyranose C1 of the glycosidic linkage over the targeted glycosidic bond. This study provides additional insight into cellulose binding by C1-specific LPMOs, giving a molecular-level picture of active site substrate interactions. DATABASE: The atomic coordinates and structure factors for HiLPMO9B have been deposited in the Protein Data Bank with accession code 5NNS.

Inhibition of Mammalian Glycoprotein YKL-40: IDENTIFICATION OF THE PHYSIOLOGICAL LIGAND.

YKL-40 is a mammalian glycoprotein associated with progression, severity, and prognosis of chronic inflammatory diseases and a multitude of cancers. Despite this well documented association, identification of the lectin's physiological ligand and, accordingly, biological function has proven experimentally difficult. YKL-40 has been shown to bind chito-oligosaccharides; however, the production of chitin by the human body has not yet been documented. Possible alternative ligands include proteoglycans, polysaccharides, and fibers like collagen, all of which makeup the extracellular matrix. It is likely that YKL-40 is interacting with these alternative polysaccharides or proteins within the body, extending its function to cell biological roles such as mediating cellular receptors and cell adhesion and migration. Here, we consider the feasibility of polysaccharides, including cello-oligosaccharides, hyaluronan, heparan sulfate, heparin, and chondroitin sulfate, and collagen-like peptides as physiological ligands for YKL-40. We use molecular dynamics simulations to resolve the molecular level recognition mechanisms and calculate the free energy of binding the hypothesized ligands to YKL-40, addressing thermodynamic preference relative to chito-oligosaccharides. Our results suggest that chitohexaose and hyaluronan preferentially bind to YKL-40 over collagen, and hyaluronan is likely the preferred physiological ligand, because the negatively charged hyaluronan shows enhanced affinity for YKL-40 over neutral chitohexaose. Collagen binds in two locations at the YKL-40 surface, potentially related to a role in fibrillar formation. Finally, heparin non-specifically binds at the YKL-40 surface, as predicted from structural studies. Overall, YKL-40 likely binds many natural ligands in vivo, but its concurrence with physical maladies may be related to associated increases in hyaluronan.

Structural and Functional Characterization of a Lytic Polysaccharide Monooxygenase with Broad Substrate Specificity.

The recently discovered lytic polysaccharide monooxygenases (LPMOs) carry out oxidative cleavage of polysaccharides and are of major importance for efficient processing of biomass. NcLPMO9C from Neurospora crassa acts both on cellulose and on non-cellulose ß-glucans, including cellodextrins and xyloglucan. The crystal structure of the catalytic domain of NcLPMO9C revealed an extended, highly polar substrate-binding surface well suited to interact with a variety of sugar substrates. The ability of NcLPMO9C to act on soluble substrates was exploited to study enzyme-substrate interactions. EPR studies demonstrated that the Cu(2+) center environment is altered upon substrate binding, whereas isothermal titration calorimetry studies revealed binding affinities in the low micromolar range for polymeric substrates that are due in part to the presence of a carbohydrate-binding module (CBM1). Importantly, the novel structure of NcLPMO9C enabled a comparative study, revealing that the oxidative regioselectivity of LPMO9s (C1, C4, or both) correlates with distinct structural features of the copper coordination sphere. In strictly C1-oxidizing LPMO9s, access to the solvent-facing axial coordination position is restricted by a conserved tyrosine residue, whereas access to this same position seems unrestricted in C4-oxidizing LPMO9s. LPMO9s known to produce a mixture of C1- and C4-oxidized products show an intermediate situation.

Cello-oligomer-binding dynamics and directionality in family 4 carbohydrate-binding modules.

Carbohydrate-binding modules (CBMs) play significant roles in modulating the function of cellulases, and understanding the protein-carbohydrate recognition mechanisms by which CBMs selectively bind substrate is critical to development of enhanced biomass conversion technology. CBMs exhibit a limited range of specificity and appear to bind polysaccharides in a directional fashion dictated by the position of the ring oxygen relative to the protein fold. The two family 4 CBMs of Cellulomonas fimi Cel9B (CfCBM4) are reported to preferentially bind cellulosic substrates. However, experimental evidence suggests that these CBMs may not exhibit a thermodynamic preference for a particular orientation. We use molecular dynamics (MD) and free energy calculations to investigate protein-carbohydrate recognition mechanisms in CfCBM4-1 and CfCBM4-2 and to elucidate preferential ligand-binding orientation. We evaluate four cellopentaose orientations including that of the crystal structure and three others suggested by nuclear magnetic resonance (NMR). These four orientations differ based on position of the ligand reducing end (RE) and pyranose ring orientations relative to the protein core. MD simulations indicate that the plausible orientations reduce to two conformations. Calculated ligand-binding free energy discerns each of the orientations is equally favorable. The calculated free energies are in excellent agreement with isothermal titration calorimetry measurements from the literature. MD simulations further reveal the approximate structural symmetry of the oligosaccharides relative to the amino acids along the binding cleft plays a role in the promiscuity of ligand binding. A survey of ligand-bound structures suggests this phenomenon may be characteristic of the broader class of proteins belonging to the ß-sandwich fold.