Improving the Force Field Description of Tyrosine-Choline Cation-π Interactions: QM Investigation of Phenol-N(Me)4+ Interactions

CHARMM at 45: Enhancements in Accessibility, Functionality, and Speed

Since its inception nearly a half century ago, CHARMM has been playing a central role in computational biochemistry and biophysics. Commensurate with the developments in experimental research and advances in computer hardware, the range of methods and applicability of CHARMM have also grown. This review summarizes major developments that occurred after 2009 when the last review of CHARMM was published. They include the following: new faster simulation engines, accessible user interfaces for convenient workflows, and a vast array of simulation and analysis methods that encompass quantum mechanical, atomistic, and coarse-grained levels, as well as extensive coverage of force fields. In addition to providing the current snapshot of the CHARMM development, this review may serve as a starting point for exploring relevant theories and computational methods for tackling contemporary and emerging problems in biomolecular systems. CHARMM is freely available for academic and nonprofit research at https://academiccharmm.org/program.

Model Mechanism for Lipid Uptake by the Human STARD2/PC-TP Phosphatidylcholine Transfer Protein

The human StAR-related lipid transfer domain protein 2 (STARD2), also known as phosphatidylcholine (PC) transfer protein, is a single-domain lipid transfer protein thought to transfer PC lipids between intracellular membranes. We performed extensive μs-long molecular dynamics simulations of STARD2 of its apo and holo forms in the presence or absence of complex lipid bilayers. The simulations in water reveal ligand-dependent conformational changes. In the 2 μs-long simulations of apo STARD2 in the presence of a lipid bilayer, we observed spontaneous reproducible PC lipid uptake into the protein hydrophobic cavity. We propose that the lipid extraction mechanism involves one to two metastable states stabilized by choline-tyrosine or choline-tryptophane cation-π interactions. Using free energy perturbation, we evaluate that PC-tyrosine cation-π interactions contribute 1.8 and 2.5 kcal/mol to the affinity of a PC-STARD2 metastable state, thus potentially providing a significant decrease of the energy barrier required for lipid desorption.

A Membrane-Assisted Mechanism for the Release of Ceramide from the CERT START Domain

Ceramide transfer protein CERT is the mediator of nonvesicular transfer of ceramide from the ER to Golgi. In CERT, START is the domain responsible for the binding and transport of ceramide. A wealth of structural data has revealed a helix-grip fold surrounding a large hydrophobic cavity holding the ceramide. Yet, little is known about the mechanisms by which START releases the ceramide through the polar region and into the packed environment of cellular membranes. As such events do not lend themselves easily to experimental investigations, we used multiple unbiased microsecond-long molecular simulations. We propose a membrane-assisted mechanism in which the membrane acts as an allosteric effector initiating the release of ceramide and where the passage of the ceramide acyl chains is facilitated by the intercalation of a single phosphatidylcholine lipid in the cavity, practically greasing the ceramide way out. We verify using free energy calculation and experimental lipidomics data that CERT forms stable complexes with phosphatidylcholine lipids, in addition to ceramide, thus providing validation for the proposed mechanism.

Membrane specificity of the human cholesterol transfer protein STARD4

STARD4 regulates cholesterol homeostasis by transferring cholesterol between the plasma membrane and endoplasmic reticulum. The STARD4 structure features a helix-grip fold surrounding a large hydrophobic cavity holding the sterol. Its access is controlled by a gate formed by the Ω1 and Ω4 loops and the C-terminal α-helix. Little is known about the mechanisms by which STARD4 binds to membranes and extracts/releases cholesterol. All available structures of STARD4 are without a bound sterol and display the same closed conformation of the gate. The cholesterol transfer activity of the mouse STARD4 is enhanced in the presence of anionic lipids, and in particular of phosphatidylinositol biphosphates (PIP2) for which two binding sites were proposed on the mouse STARD4 surface. Yet only one of these sites is conserved in human STARD4. We here report the results of a liposome microarray-based assay and microseconds-long molecular dynamics simulations of human STARD4 with complex lipid bilayers mimicking the composition of the donor and acceptor membranes. We show that the binding of apo form of human STARD4 is sensitive to the presence of PIP2 through two specific binding sites, one of which was not identified on mouse STARD4. We report two novel conformations of the gate in holo-STARD4: a yet-unobserved close conformation and an open conformation of Ω4 shedding light on the opening/closure mechanism needed for cholesterol uptake/release. Overall, the modulation of human STARD4 membrane-binding by lipid composition, and by the presence of the cargo supports the capacity of human STARD4 to achieve directed transfer between specific organelle membranes.

Standard Binding Free Energy and Membrane Desorption Mechanism for a Phospholipase C

Peripheral membrane proteins (PMPs) bind temporarily to cellular membranes and play important roles in signaling, lipid metabolism, and membrane trafficking. Obtaining accurate membrane-PMP affinities using experimental techniques is more challenging than for protein-ligand affinities in an aqueous solution. At the theoretical level, calculation of the standard protein-membrane binding free energy using molecular dynamics simulations remains a daunting challenge owing to the size of the biological objects at play, the slow lipid diffusion, and the large variation in configurational entropy that accompanies the binding process. To overcome these challenges, we used a computational framework relying on a series of potential-of-mean-force (PMF) calculations including a set of geometrical restraints on collective variables. This methodology allowed us to determine the standard binding free energy of a PMP to a phospholipid bilayer using an all-atom force field. Bacillus thuringiensis phosphatidylinositol-specific phospholipase C (BtPI-PLC) was chosen due to its importance as a virulence factor and owing to the host of experimental affinity data available. We computed a standard binding free energy of -8.2 ± 1.4 kcal/mol in reasonable agreement with the reported experimental values (-6.6 ± 0.2 kcal/mol). In light of the 2.3-μs separation PMF calculation, we investigated the mechanism whereby BtPI-PLC disengages from interactions with the lipid bilayer during separation. We describe how a short amphipathic helix engages in transitory interactions to ease the passage of its hydrophobes through the interfacial region upon desorption from the bilayer.

Specific interactions of peripheral membrane proteins with lipids: what can molecular simulations show us?

Peripheral membrane proteins (PMPs) can reversibly and specifically bind to biological membranes to carry out functions such as cell signalling, enzymatic activity, or membrane remodelling. Structures of these proteins and of their lipid-binding domains are typically solved in a soluble form, sometimes with a lipid or lipid headgroup at the binding site. To provide a detailed molecular view of PMP interactions with the membrane, computational methods such as molecular dynamics (MD) simulations can be applied. Here, we outline recent attempts to characterise these binding interactions, focusing on both intracellular proteins, such as phosphatidylinositol phosphate (PIP)-binding domains, and extracellular proteins such as glycolipid-binding bacterial exotoxins. We compare methods used to identify and analyse lipid-binding sites from simulation data and highlight recent work characterising the energetics of these interactions using free energy calculations. We describe how improvements in methodologies and computing power will help MD simulations to continue to contribute to this field in the future.

Demystifying DPP III Catalyzed Peptide Hydrolysis-Computational Study of the Complete Catalytic Cycle of Human DPP III Catalyzed Tynorphin Hydrolysis

Dipeptidyl peptides III (DPP III) is a dual-domain zinc exopeptidase that hydrolyzes peptides of varying sequence and size. Despite attempts to elucidate and narrow down the broad substrate-specificity of DPP III, there is no explanation as to why some of them, such as tynorphin (VVYPW), the truncated form of the endogenous heptapeptide spinorphin, are the slow-reacting substrates of DPP III compared to others, such as Leu-enkephalin. Using quantum molecular mechanics calculations followed by various molecular dynamics techniques, we describe for the first time the entire catalytic cycle of human DPP III, providing theoretical insight into the inhibitory mechanism of tynorphin. The chemical step of peptide bond hydrolysis and the substrate binding to the active site of the enzyme and release of the product were described for DPP III in complex with tynorphin and Leu-enkephalin and their products. We found that tynorphin is cleaved by the same reaction mechanism determined for Leu-enkephalin. More importantly, we showed that the product stabilization and regeneration of the enzyme, but not the nucleophilic attack of the catalytic water molecule and inversion at the nitrogen atom of the cleavable peptide bond, correspond to the rate-determining steps of the overall catalytic cycle of the enzyme.

Specificity of Loxosceles α clade phospholipase D enzymes for choline-containing lipids: Role of a conserved aromatic cage

Spider venom GDPD-like phospholipases D (SicTox) have been identified to be one of the major toxins in recluse spider venom. They are divided into two major clades: the α clade and the β clade. Most α clade toxins present high activity against lipids with choline head groups such as sphingomyelin, while activities in β clade toxins vary and include preference for substrates containing ethanolamine headgroups (Sicarius terrosus, St_βIB1). A structural comparison of available structures of phospholipases D (PLDs) reveals a conserved aromatic cage in the α clade. To test the potential influence of the aromatic cage on membrane-lipid specificity we performed molecular dynamics (MD) simulations of the binding of several PLDs onto lipid bilayers containing choline headgroups; two SicTox from the α clade, Loxosceles intermedia αIA1 (Li_αIA) and Loxosceles laeta αIII1 (Ll_αIII1), and one from the β clade, St_βIB1. The simulation results reveal that the aromatic cage captures a choline-headgroup and suggest that the cage plays a major role in lipid specificity. We also simulated an engineered St_βIB1, where we introduced the aromatic cage, and this led to binding with choline-containing lipids. Moreover, a multiple sequence alignment revealed the conservation of the aromatic cage among the α clade PLDs. Here, we confirmed that the i-face of α and β clade PLDs is involved in their binding to choline and ethanolamine-containing bilayers, respectively. Furthermore, our results suggest a major role in choline lipid recognition of the aromatic cage of the α clade PLDs. The MD simulation results are supported by in vitro liposome binding assay experiments.

Envenomation following bites from recluse spiders (Loxosceles) causes loxoscelism, a necrotic tissue breakdown in mammals, and leads to skin degeneration and systemic reactions in the worst case. Recluse spiders belong to the Sicariidae family which also includes six-eyed sand spiders in the genera Sicarius and Hexopthalma. While sicariid spiders are found natively on all continents except Australia, treatments of loxoscelism are typically antibody based and available in some regions of the Americas. Sphingomyelinase D/phospholipase D enzymes are one of the major toxins in venom of sicariid spiders, and have been divided in two clades called α and β. The activity of α and β clades toxins differs; most α clade toxins present high activity against lipids with choline headgroups (-N (CH₃)₃⁺) such as sphingomyelin, while activities in β clade toxins vary and include preference for substrates containing ethanolamine headgroups (-NH₃⁺). When comparing the structures of two α clade toxins and one β clade toxin, we noticed the presence in the α clade toxins only of a cage consisting of three aromatic amino acids. In this work we used numerical molecular simulations to probe the role of this cage in the preference of α clade toxins for choline head groups over ethanolamine head groups.

Physisorption of bio oil nitrogen compounds onto montmorillonite

We assess computationally the adsorption of a series of nitrogen containing heterocycles and fatty acid amides from bio-oil on a model clay surface, Na-montmorillonite. The adsorption energies and conformations predicted by atomistic detail molecular dynamics (MD) simulations are compared against density functional theory (DFT) based molecular electrostatic potentials (MEP) and Hirshfeld, AIM, Merz-Singh-Kollman, and ChelpG charges. MD predicts systematically adsorption via cation bridging with adsorption strength of the heterocycles following purine > pyridine > imidazole > pyrrole > indole > quinoline. The fatty acid amides adsorption strength follows the steric availability and bulkiness of the head group. A comparison against the DFT calculations shows that MEP predicts adsorption geometries and the MD simulations reproduce the conformations for single adsorption site species. However, the DFT derived charge distibutions show that MD force-fields with non-polarizable fixed partial charge representations parametrized for aqueous environments cannot be used in apolar solvent environments without careful accuracy considerations. The overall trends in adsorption energies are reproduced by the Charmm GenFF employed in the MD simulations but the adsorption energies are systematically overestimated in this apolar solvent environment. The work has significance both for revealing nitrogen compound adsorption trends in technologically relevant bio oil environments but also as a methodological assessment revealing the limits of state of the art biomolecular force-fields and simulation protocols in apolar bioenvironments.

PyAutoFEP: An Automated Free Energy Perturbation Workflow for GROMACS Integrating Enhanced Sampling Methods

Free energy perturbation (FEP) calculations are now routinely used in drug discovery to estimate the relative FEB (RFEB) of small molecules to a biomolecular target of interest. Using enhanced sampling can improve the correlation between predictions and experimental data, especially in systems with conformational changes. Due to the large number of perturbations required in drug discovery campaigns, the manual setup of FEP calculations is no longer viable. Here, we introduce PyAutoFEP, a flexible and open-source tool to aid the setup of RFEB FEP. PyAutoFEP is written in Python3, and automates the generation of perturbation maps, dual topologies, system building and molecular dynamics (MD), and analysis. PyAutoFEP supports multiple force fields, incorporates replica exchange with solute tempering (REST) and replica exchange with solute scaling (REST2) enhanced sampling methods, and allows flexible λ values along perturbation windows. To validate PyAutoFEP, it was applied to a set of 14 Farnesoid X receptor ligands, a system included in the drug design data resource grand challenge 2. An 88% mean correct sign prediction was achieved, and 75% of the predictions had an error below 1.5 kcal/mol. Results using Amber03/GAFF, CHARMM36m/CGenFF, and OPLS-AA/M/LigParGen had Pearson's r values of 0.71 ± 0.13, 0.30 ± 0.27, and 0.66 ± 0.20, respectively. The Amber03/GAFF and OPLS-AA/M/LigParGen results were on par with the top grand challenge 2 submissions. Applying REST2 improved the results using CHARMM36m/CGenFF (Pearson's r = 0.43 ± 0.21) but had little impact on the other force fields. CHARMM36-YF and CHARMM36-WYF modifications did not yield improved predictions compared to CHARMM36m. Finally, we estimated the probability of finding a molecule 1 pK_i better than a lead when using PyAutoFEP to screen 10 or 100 analogues. The probabilities, when compared to random sampling, increased up to sevenfold when 100 molecules were to be screened, suggesting that PyAutoFEP would likely be useful for lead optimization. PyAutoFEP is available on GitHub at https://github.com/lmmpf/PyAutoFEP.

Accuracy of Alternate Nonpolarizable Force Fields for the Determination of Protein-Ligand Binding Affinities Dominated by Cation-π Interactions

Modifying pair-specific Lennard-Jones parameters through the nonbonded FIX (NBFIX) feature of the CHARMM36 force field has proven cost-effective for improving the description of cation-π interactions in biological objects by means of pairwise additive potential energy functions. Here, two sets of newly optimized CHARMM36 force-field parameters including NBFIX corrections, coined CHARMM36m-NBF and CHARMM36-WYF, and the original force fields, namely CHARMM36m and Amber ff14SB, are used to determine the standard binding free energies of seven protein-ligand complexes containing cation-π interactions. Compared with precise experimental measurements, our results indicate that the uncorrected, original force fields significantly underestimate the binding free energies, with a mean error of 5.3 kcal/mol, while the mean errors of CHARMM36m-NBF and CHARMM36-WYF amount to 0.8 and 2.1 kcal/mol, respectively. The present study cogently demonstrates that the use of modified parameters jointly with NBFIX corrections dramatically increases the accuracy of the standard binding free energy of protein-ligand complexes dominated by cation-π interactions, most notably with CHARMM36m-NBF.

Automated Parameterization of Quantum Mechanically Derived Force Fields for Soft Materials and Complex Fluids: Development and Validation

The reliability of molecular dynamics (MD) simulations in predicting macroscopic properties of complex fluids and soft materials, such as liquid crystals, colloidal suspensions, or polymers, relies on the accuracy of the adopted force field (FF). We present an automated protocol to derive specific and accurate FFs, fully based on ab initio quantum mechanical (QM) data. The integration of the Joyce and Picky procedures, recently proposed by our group to provide an accurate description of simple liquids, is here extended to larger molecules, capable of exhibiting more complex fluid phases. While the standard Joyce protocol is employed to parameterize the intramolecular FF term, a new automated procedure is here proposed to handle the computational cost of the QM calculations required for the parameterization of the intermolecular FF term. The latter is thus obtained by integrating the old Picky procedure with a fragmentation reconstruction method (FRM) that allows for a reliable, yet computationally feasible sampling of the intermolecular energy surface at the QM level. The whole FF parameterization protocol is tested on a benchmark liquid crystal, and the performances of the resulting quantum mechanically derived (QMD) FF were compared with those delivered by a general-purpose, transferable one, and by the third, "hybrid" FF, where only the bonded terms were refined against QM data. Lengthy atomistic MD simulations are carried out with each FF on extended 5CB systems in both isotropic and nematic phases, eventually validating the proposed protocol by comparing the resulting macroscopic properties with other computational models and with experiments. The QMD-FF yields the best performances, reproducing both phases in the correct range of temperatures and well describing their structure, dynamics, and thermodynamic properties, thus providing a clear protocol that may be explored to predict such properties on other complex fluids or soft materials.

NADES-based surfactant-free microemulsions for solubilization and extraction of curcumin from Curcuma Longa

A choline chloride + lactic acid (1:1) natural deep eutectic solvent (NADES) is used as an adjuvant to ethanol/triacetin mixtures to solubilize and extract curcumin from Curcuma Longa. The obtained NADES/ethanol/triacetin mixtures are homogeneous, transparent and of low viscosity even in the absence of water. Dynamic light scattering revealed significant nanostructures, typical of surfactant-free microemulsions. A twofold increase of curcumin solubility and remarkable extraction power (yield of ~90%) can be achieved in the ternary system including the NADES, although curcumin is hydrophobic and the used NADES are very polar. Due to the elevated solubility of curcumin, more extraction cycles can be made than in the previously published aqueous systems with the same amount of solution. As a result, less solvent is required to achieve the same extraction yield.

Reentrant liquid condensate phase of proteins is stabilized by hydrophobic and non-ionic interactions

Liquid-liquid phase separation of proteins underpins the formation of membraneless compartments in living cells. Elucidating the molecular driving forces underlying protein phase transitions is therefore a key objective for understanding biological function and malfunction. Here we show that cellular proteins, which form condensates at low salt concentrations, including FUS, TDP-43, Brd4, Sox2, and Annexin A11, can reenter a phase-separated regime at high salt concentrations. By bringing together experiments and simulations, we demonstrate that this reentrant phase transition in the high-salt regime is driven by hydrophobic and non-ionic interactions, and is mechanistically distinct from the low-salt regime, where condensates are additionally stabilized by electrostatic forces. Our work thus sheds light on the cooperation of hydrophobic and non-ionic interactions as general driving forces in the condensation process, with important implications for aberrant function, druggability, and material properties of biomolecular condensates.

Reentrant liquid condensate phase of proteins is stabilized by hydrophobic and non-ionic interactions

Liquid-liquid phase separation of proteins underpins the formation of membraneless compartments in living cells. Elucidating the molecular driving forces underlying protein phase transitions is therefore a key objective for understanding biological function and malfunction. Here we show that cellular proteins, which form condensates at low salt concentrations, including FUS, TDP-43, Brd4, Sox2, and Annexin A11, can reenter a phase-separated regime at high salt concentrations. By bringing together experiments and simulations, we demonstrate that this reentrant phase transition in the high-salt regime is driven by hydrophobic and non-ionic interactions, and is mechanistically distinct from the low-salt regime, where condensates are additionally stabilized by electrostatic forces. Our work thus sheds light on the cooperation of hydrophobic and non-ionic interactions as general driving forces in the condensation process, with important implications for aberrant function, druggability, and material properties of biomolecular condensates.

Explicit Representation of Cation-π Interactions in Force Fields with 1/r4 Nonbonded Terms

The binding energies for cation-π complexation are underestimated by traditional fixed-charge force fields owing to their lack of explicit treatment of ion-induced dipole interactions. To address this deficiency, an explicit treatment of cation-π interactions has been introduced into the OPLS-AA force field. Following prior work with atomic cations, it is found that cation-π interactions can be handled efficiently by augmenting the usual 12-6 Lennard-Jones potentials with 1/r⁴ terms. Results are provided for prototypical complexes as well as protein-ligand systems of relevance for drug design. Alkali cation, ammonium, guanidinium, and tetramethylammonium were chosen for the representative cations, while benzene and six heteroaromatic molecules were used as the π systems. The required nonbonded parameters were fit to reproduce structure and interaction energies for gas-phase complexes from density functional theory (DFT) calculations at the ωB97X-D/6-311++G(d,p) level. The impact of the solvent was then examined by computing potentials of mean force (pmfs) in both aqueous and tetrahydrofuran (THF) solutions using the free-energy perturbation (FEP) theory. Further testing was carried out for two cases of strong and one case of weak cation-π interactions between druglike molecules and their protein hosts, namely, the JH2 domain of JAK2 kinase and macrophage migration inhibitory factor. FEP results reveal greater binding by 1.5-4.4 kcal/mol from the addition of the explicit cation-π contributions. Thus, in the absence of such treatment of cation-π interactions, errors for computed binding or inhibition constants of 10¹-10³ are expected.

Accurate Description of Cation-π Interactions in Proteins with a Nonpolarizable Force Field at No Additional Cost

Cation-π interactions play a significant role in a host of processes eminently relevant to biology. However, polarization effects arising from the interaction of cations with aromatic moieties have long been recognized to be inadequately described by pairwise additive force fields. In the present work, we address this longstanding shortcoming through the nonbonded FIX (NBFIX) feature of the CHARMM36 force field, modifying pair-specific Lennard-Jones (LJ) parameters, while circumventing the limitations of the Lorentz-Berthelot combination rules. The potentials of mean force (PMFs) characterizing prototypical cation-π interactions in aqueous solutions are first determined using a hybrid quantum mechanical/molecular mechanics (QM/MM) strategy in conjunction with an importance-sampling algorithm. The LJ parameters describing the cation-π pairs are then optimized to match the QM/MM PMFs. The standard binding free energies of nine cation-π complexes, i.e., toluene, para-cresol, and 3-methyl-indole interacting with either ammonium, guanidinium, or tetramethylammonium, determined with this new set of parameters agree well with the experimental measurements. Additional simulations were carried out on three different classes of biological objects featuring cation-π interactions, including five individual proteins, three protein-ligand complexes, and two protein-protein complexes. Our results indicate that the description of cation-π interactions is overall improved using NBFIX corrections, compared with the standard pairwise additive force field. Moreover, an accurate binding free energy calculation for a protein-ligand complex containing cation-π interactions (2BOK) shows that using the new parameters, the experimental binding affinity can be reproduced quantitatively. Put together, the present work suggests that the NBFIX parameters optimized here can be broadly utilized in the simulation of proteins in an aqueous solution to enhance the representation of cation-π interactions, at no additional computational cost.

A Spectroscopic Study of Xanthene Dyes on a Polystyrene Surface: an Investigation of Ion-π Interactions at Polymer Interfaces

Thin films of three differently charged xanthene dyes: rhodamine 6G (Rh6G, cationic), fluorescein 27 (F27, neutral), and disodium fluorescein (DSF, anionic) were cast onto a polystyrene (PS) coated glass substrate to investigate ion-π interactions. Absorbance spectroscopy was used to determine the aggregation state of the dyes on the PS surface. Deconvolution of the spectra of films revealed multiple peaks for all dyes assigned to isolated monomers, aggregates, and complexes between the dye and the polymer substrate. The shift of the low energy peak relative to the monomer peak was used as an indication of interaction strength of that species with the PS π system, which followed the trend of Rh6G > DSF > F27. Increase in the interaction energy is attributed to stabilization from ion-π interactions. Steady-state emission spectra and excited state lifetime measurements were performed on all films. The formation of a weakly emissive exciplex was found for Rh6G and DSF, consistent with ion-π interactions, but no evidence of an exciplex is found for the F27 films.

Insights into the Sensitivity of Arginine Concentration to Preserve the Folded Form of Insulin Monomer under Thermal Stress

Arginine, although popularly known as aggregation suppressor additive, has been found to quench proteins' structure and function by destabilizing their conformations. Driven by such controversial evidence, in this work we performed a series of atomistic molecular dynamics simulations of insulin monomer, a biologically active hormone protein, in arginine solution of varying concentrations (0.5, 1, and 2 M) at ambient and elevated temperature (400 K) to explore the arginine concentration driven structure-based stability of the protein. Our study reveals that the flexibility of the protein's structure is dependent on the arginine concentration, and among all the used solutions, 2 M arginine, a "neutral crowder" that mimics the cellular environment, can preserve the native folded form of the protein at ambient temperature in an excellent manner. Further, while the protein unfolds at 400 K in pure water, this solution worked satisfactorily to preserve the protein's folded conformation more firmly than the other solutions. The replica-exchange MD of insulin in 2 M arginine solution further supports the fact. In this aspect an important issue in molecular pharmacology is to identify and recognize the physical origin of the stability of a protein, i.e, in this case, how arginine directs the conformational flexibility of the protein and preserves its native folded form. We identified that the exclusion of arginine from the protein surface increases the local structuration of water around the protein, thereby preserving its "biological water" layer, and makes the protein more hydrated at 2 M concentration as compared to the other arginine solutions. Additionally, our microscopic investigation on the interactions of the protein-solvation layer revealed that the structural heterogeneity of the protein surface, arising from the differential physicochemical nature of the amino acid residues, controls the favorable formation of sluggish water-arginine mixed solvation layer at higher arginine concentration that helps the protein to maintain its structural rigidity. Importantly, apart from the protein-solvent hydrogen-bonding interactions, the anion-pi interactions, established between the carboxyl group of arginine and the aromatic amino acid residues of insulin, were recognized to facilitate the protein to maintain its native folded form at the experimental temperatures.

Capturing Choline-Aromatics Cation-π Interactions in the MARTINI Force Field

Cation−π interactions play an important role in biomolecular recognition, including interactions between membrane phosphatidylcholine lipids and aromatic amino acids of peripheral proteins. While molecular mechanics coarse grain (CG) force fields are particularly well suited to simulate membrane proteins in general, they are not parameterized to explicitly reproduce cation−π interactions. We here propose a modification of the polarizable MARTINI coarse grain (CG) model enabling it to model membrane binding events of peripheral proteins whose aromatic amino acid interactions with choline headgroups are crucial for their membrane binding. For this purpose, we first collected and curated a dataset of eight peripheral proteins from different families. We find that the MARTINI CG model expectedly underestimates aromatics–choline interactions and is unable to reproduce membrane binding of the peripheral proteins in our dataset. Adjustments of the relevant interactions in the polarizable MARTINI force field yield significant improvements in the observed binding events. The orientation of each membrane-bound protein is comparable to reference data from all-atom simulations and experimental binding data. We also use negative controls to ensure that choline–aromatics interactions are not overestimated. We finally check that membrane properties, transmembrane proteins, and membrane translocation potential of mean force (PMF) of aromatic amino acid side-chain analogues are not affected by the new parameter set. This new version “MARTINI 2.3P” is a significant improvement over its predecessors and is suitable for modeling membrane proteins including peripheral membrane binding of peptides and proteins.

Improved Modeling of Cation-π and Anion-Ring Interactions Using the Drude Polarizable Empirical Force Field for Proteins

Cation-π interactions are noncovalent interactions between a π-electron system and a positively charged ion that are regarded as a strong noncovalent interaction and are ubiquitous in biological systems. Similarly, though less studied, anion-ring interactions are present in proteins along with in-plane interactions of anions with aromatic rings. As these interactions are between a polarizing ion and a polarizable π system, the accuracy of the treatment of these interactions in molecular dynamics (MD) simulations using additive force fields (FFs) may be limited. In the present work, to allow for a better description of ion-π interactions in proteins in the Drude-2013 protein polarizable FF, we systematically optimized the parameters for these interactions targeting model compound quantum mechanical (QM) interaction energies with atom pair-specific Lennard-Jones parameters along with virtual particles as selected ring centroids introduced to target the QM interaction energies and geometries. Subsequently, MD simulations were performed on a series of protein structures where ion-π pairs occur to evaluate the optimized parameters in the context of the Drude-2013 FF. The resulting FF leads to a significant improvement in reproducing the ion-π pair distances observed in experimental protein structures, as well as a smaller root-mean-square differences and fluctuations of the overall protein structures from experimental structures. Accordingly, the optimized Drude-2013 protein polarizable FF is suggested for use in MD simulations of proteins where cation-π and anion-ring interactions are critical. © 2019 Wiley Periodicals, Inc.

Planarity and out-of-plane vibrational modes of tryptophan and tyrosine in biomolecular modeling

Tryptophan and tyrosine are amino acids that play significant roles in the folding processes of proteins at water-membrane interfaces because of their amphipathic heteroaromatic rings. Employing appropriate heteroaromatic molecular structures is essential for obtaining accurate dynamics and predictive capabilities in molecular simulations of these amino acids. In this study, molecular dynamics simulations that applied the most recent version of the CHARMM36 force field were conducted on aqueous solutions of tryptophan and of tyrosine. Geometric analysis and dynamics quantified how aromatic rings deviated from planar structures and exhibited out-of-plane fluctuations. Radial distribution functions showed possible biological significance because the extent of ring planarity slightly affected local water concentrations near aromatic rings. Instantaneous all-atom normal mode analysis (NMA) and Fourier transformation of time autocorrelation functions of out-of-plane displacements were applied to study out-of-plane vibrations of atoms in these rings. The NMA started with minimum energy configurations and then averaged over fluctuations in aqueous solution. The frequencies and frequency patterns that were obtained for tryptophan and tyrosine with CHARMM36 differed from literature reports of Raman spectra, infrared spectra, and frequencies calculated using quantum mechanics, with some out-of-plane modes found at higher frequencies. Effects of imposing improper torsion potentials and changing torsion angle force constants were investigated for all atoms in the rings of tryptophan and tyrosine. Results show that these coarse force field variations only affect planarity and out-of-plane vibrations of atoms within the rings, and not other vibrations. Although increasing improper torsion force constants reduced deviations from aromatic ring planarity significantly, it increased out-of-plane mode frequencies. Reducing torsion angle force constants (with and without improper torsions) shifted modes to lower frequencies. A combination of decreasing most torsion angle force constants for ring atoms in both amino acids and including improper torsion forces attained frequencies and frequency patterns for out-of-plane normal modes that were more similar to the literature spectra. These force field variations decreased the extents of out-of-plane vibrations within the heteroaromatic rings of tryptophan, especially around the nitrogen atom in the ring, but not within the heteroaromatic ring of tyrosine. Conclusions were unaffected by the peptide endgroup, water, or simulation ensemble.

Interfacial Aromatics Mediating Cation-π Interactions with Choline-Containing Lipids Can Contribute as Much to Peripheral Protein Affinity for Membranes as Aromatics Inserted below the Phosphates

Membrane-binding interfaces of peripheral proteins are restricted to a small part of their exposed surface, so the ability to engage in strong selective interactions with membrane lipids at various depths in the interface, both below and above the phosphates, is an advantage. Driven by their hydrophobicity, aromatic amino acids preferentially partition into membrane interfaces, often below the phosphates, yet enthalpically favorable interactions with the lipid headgroups, above the phosphate plane, are likely to further stabilize high interfacial positions. Using free-energy perturbation, we calculate the energetic cost of alanine substitution for 11 interfacial aromatic amino acids from 3 peripheral proteins. We show that the involvement in cation-π interactions with the headgroups (i) increases the ΔΔG_transfer as compared with insertion at the same depth without cation-π stabilization and (ii) can contribute at least as much as deeper insertion below the phosphates, highlighting the multiple roles of aromatics in peripheral membrane protein affinity.

Cation-π Interactions between Methylated Ammonium Groups and Tryptophan in the CHARMM36 Additive Force Field

Cation-π interactions between tryptophan and choline or trimethylated lysines are vital for many biological processes. The performance of the additive CHARMM36 force field against target quantum mechanical data is shown to reproduce QM equilibrium geometries but required modified Lennard-Jones potentials to accurately reproduce the QM interaction energies. The modified parameter set allows accurate modeling, including free energies, of cation-π indole-choline and indole-trimethylated lysines interactions relevant for protein-ligand, protein-membrane, and protein-protein interfaces.

Cooperative Nonbonded Forces Control Membrane Binding of the pH-Low Insertion Peptide pHLIP

Peptides with the ability to bind and insert into the cell membrane have immense potential in biomedical applications. pH (low) insertion peptide (pHLIP), a water-soluble polypeptide derived from helix C of bacteriorhodopsin, can insert into a membrane at acidic pH to form a stable transmembrane α-helix. The insertion process takes place in three stages: pHLIP is unstructured and soluble in water at neutral pH (state I), unstructured and bound to the surface of a membrane at neutral pH (state II), and inserted into the membrane as an α-helix at low pH (state III). Using molecular dynamics simulations, we have modeled state II of pHLIP and a fast-folding variant of pHLIP, in which each peptide is bound to a 1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine bilayer surface. Our results provide strong support for recently published spectroscopic studies, namely that pHLIP preferentially binds to the bilayer surface as a function of location of anionic amino acids and that backbone dehydration occurs upon binding. Unexpectedly, we also observed several instances of segments of pHLIP folding into a stable helical turn. Our results provide a molecular level of detail that is essential to providing new insights into pHLIP function and to facilitate design of variants with improved membrane-active capabilities.

Development and Validation of Quantum Mechanically Derived Force-Fields: Thermodynamic, Structural, and Vibrational Properties of Aromatic Heterocycles

A selection of several aromatic molecules, representative of the important class of heterocyclic compounds, has been considered for testing and validating an automated Force Field (FF) parametrization protocol, based only on Quantum Mechanical data. The parametrization is carried out separately for the intra- and intermolecular contributions, employing respectively the Joyce and Picky software packages, previously implemented and refined in our research group. The whole approach is here automated and integrated with a computationally effective yet accurate method, devised very recently ( J. Chem.

THEORY

Comput., 2018, 14, 543-556) to evaluate a large number of dimer interaction energies. The resulting quantum mechanically derived FFs are then used in extensive molecular dynamics simulations, in order to evaluate a number of thermodynamic, structural, and dynamic properties of the heterocycle's gas and liquid phases. The comparison with the available experimental data is good and furnishes a validation of the presented approach, which can be confidently exploited for the design of novel and more complex materials.

The Color of Cation-π Interactions: Subtleties of Amine-Tryptophan Interaction Energetics Allow for Radical-like Visible Absorbance and Fluorescence

Several peptides and a protein with an inter- or intramolecular cation-π interaction between tryptophan (Trp) and an amine cation are shown to absorb and fluoresce in the visible region of the spectrum. Titration of indole with sodium hydroxide or ammonium hydroxide yields an increasing visible fluorescence as well. Visible absorption and multipeaked fluorescence excitation spectra correlate with experimental absorption spectra and the vibrational modes of calculated absorption spectra for the neutral Trp radical. The radical character of the cation-indole interaction is predicted to stem from the electrostatic dislocation of indole highest occupied molecular orbital (HOMO) charge density toward the cation with a subsequent electronic transition from the HOMO-2 to the HOMO. Because this is a vertical transition, fluorescence is possible. Hydrogen bonding at the indole amine most likely stabilizes the radical-like state. These results provide new spectroscopic tools for the investigation of cation-π interactions in numerous biological systems, among them, proteins and their myriad ligands, and show that one, or at most, two, point mutations with natural amino acids are all that is required to impart visible fluorescence to proteins.

Exploring Covalent Allosteric Inhibition of Antigen 85C from Mycobacterium tuberculosis by Ebselen Derivatives

Previous studies identified ebselen as a potent in vitro and in vivo inhibitor of the Mycobacterium tuberculosis (Mtb) antigen 85 (Ag85) complex, comprising three homologous enzymes required for the biosynthesis of the mycobacterial cell wall. In this study, the Mtb Ag85C enzyme was cocrystallized with azido and adamantyl ebselen derivatives, resulting in two crystallographic structures of 2.01 and 1.30 Å resolution, respectively. Both structures displayed the anticipated covalent modification of the solvent accessible, noncatalytic Cys209 residue forming a selenenylsulfide bond. Continuous difference density for both thiol modifiers allowed for the assessment of interactions that influence ebselen binding and inhibitor orientation that were unobserved in previous Ag85C ebselen structures. The k_inact/K_I values for ebselen, adamantyl ebselen, and azido ebselen support the importance of observed constructive chemical interactions with Arg239 for increased in vitro efficacy toward Ag85C. To better understand the in vitro kinetic properties of these ebselen derivatives, the energetics of specific protein-inhibitor interactions and relative reaction free energies were calculated for ebselen and both derivatives using density functional theory. These studies further support the different in vitro properties of ebselen and two select ebselen derivatives from our previously published ebselen library with respect to kinetics and protein-inhibitor interactions. In both structures, the α9 helix was displaced farther from the enzyme active site than the previous Ag85C ebselen structure, resulting in the restructuring of a connecting loop and imparting a conformational change to residues believed to play a role in substrate binding specific to Ag85C. These notable structural changes directly affect protein stability, reducing the overall melting temperature by up to 14.5 °C, resulting in the unfolding of protein at physiological temperatures. Additionally, this structural rearrangement due to covalent allosteric modification creates a sizable solvent network that encompasses the active site and extends to the modified Cys209 residue. In all, this study outlines factors that influence enzyme inhibition by ebselen and its derivatives while further highlighting the effects of the covalent modification of Cys209 by said inhibitors on the structure and stability of Ag85C. Furthermore, the results suggest a strategy for developing new classes of Ag85 inhibitors with increased specificity and potency.

Structural and energetic study of cation–π–cation interactions in proteins

Statistical and energetic analysis of cation–π–cation motifs in protein structures suggests a potential stabilizing role in the protein fold.