3D Reference Interaction Site Model Research Articles

Modeling the aqueous ionic environment of biomolecules with the 3D-reference interaction site model of molecular solvation

Modeling the aqueous ionic environment of biomolecules with the 3D-reference interaction site model of molecular solvation

Effect of Main and Side Chains on the Folding Mechanism of the Trp-Cage Miniprotein.

Proteins that do not fold into their functional native state have been linked to diseases. In this study, the influence of the main and side chains of individual amino acids on the folding of the tryptophan cage (Trp-cage), a designed 20-residue miniprotein, was analyzed. For this purpose, we calculated the solvation free energy (SFE) contributions of individual atoms by using the 3D-reference interaction site model with the atomic decomposition method. The mechanism by which the Trp-cage is stabilized during the folding process was examined by calculating the total energy, which is the sum of the conformational energy and SFE. The folding process of the Trp-cage resulted in a stable native state, with a total energy that was 62.4 kcal/mol lower than that of the unfolded state. The solvation entropy, which is considered to be responsible for the hydrophobic effect, contributed 31.3 kcal/mol to structural stabilization. In other words, the contribution of the solvation entropy accounted for approximately half of the total contribution to Trp-cage folding. The hydrophobic core centered on Trp6 contributed 15.6 kcal/mol to the total energy, whereas the solvation entropy contribution was 6.3 kcal/mol. The salt bridge formed by the hydrophilic side chains of Asp9 and Arg16 contributed 10.9 and 5.0 kcal/mol, respectively. This indicates that not only the hydrophobic core but also the salt bridge of the hydrophilic side chains gain solvation entropy and contribute to stabilizing the native structure of the Trp-cage.

Evaluating performance of the approximate 3D-RISM-KH molecular solvation theory for solvation free energies in alkanes and alkane-water partition coefficients

The performance of the 3D-reference interaction site model (3D-RISM) with the Kovalenko-Hirata (KH) closure is benchmarked for liquid alkanes. The liquids are modelled using the united atom formalism. The intermolecular distribution functions obtained from the computations are in good agreement with the experimental reports. The solvation free energy (SFE) computed in the alkanes for a large set of solutes also agrees well with the experimental data. The molecular partition coefficients for alkane-water partition of the compounds were calculated for cyclohexane-water, alkane-water, and hexadecane-water partitions of a total of eight hundred and forty-three organic compounds with reliable accuracy. A calibration scale is proposed for alkane-water partition coefficient calculations based on the 3D-RISM-KH theory.

EPISOL: A software package with expanded functions to perform 3D-RISM calculations for the solvation of chemical and biological molecules.

Integral equation theory (IET) provides an effective solvation model for chemical and biological systems that balances computational efficiency and accuracy. We present a new software package, the expanded package for IET-based solvation (EPISOL), that performs 3D-reference interaction site model (3D-RISM) calculations to obtain the solvation structure and free energies of solute molecules in different solvents. In EPISOL, we have implemented 22 different closures, multiple free energy functionals, and new variations of 3D-RISM theory, including the recent hydrophobicity-induced density inhomogeneity (HI) theory for hydrophobic solutes and ion-dipole correction (IDC) theory for negatively charged solutes. To speed up the convergence and enhance the stability of the self-consistent iterations, we have introduced several numerical schemes in EPISOL, including a newly developed dynamic mixing approach. We show that these schemes have significantly reduced the failure rate of 3D-RISM calculations compared to AMBER-RISM software. EPISOL consists of both a user-friendly graphic interface and a kernel library that allows users to call its routines and adapt them to other programs. EPISOL is compatible with the force-field and coordinate files from both AMBER and GROMACS simulation packages. Moreover, EPISOL is equipped with an internal memory control to efficiently manage the use of physical memory, making it suitable for performing calculations on large biomolecules. We demonstrate that EPISOL can efficiently and accurately calculate solvation density distributions around various solute molecules (including a protein chaperone consisting of 120,715 atoms) and obtain solvent free energy for a wide range of organic compounds. We expect that EPISOL can be widely applied as a solvation model for chemical and biological systems. EPISOL is available at https://github.com/EPISOLrelease/EPISOL.

Identifying Systematic Force Field Errors Using a 3D-RISM Element Counting Correction.

Hydration free energies of small molecules are commonly used as benchmarks for solvation models. However, errors in predicting hydration free energies are partially due to the force fields used and not just the solvation model. To address this, we have used the 3D reference interaction site model (3D-RISM) of molecular solvation and existing benchmark explicit solvent calculations with a simple element count correction (ECC) to identify problems with the non-bond parameters in the general AMBER force field (GAFF). 3D-RISM was used to calculate hydration free energies of all 642 molecules in the FreeSolv database, and a partial molar volume correction (PMVC), ECC, and their combination (PMVECC) were applied to the results. The PMVECC produced a mean unsigned error of 1.01±0.04kcal/mol and root mean squared error of 1.44±0.07kcal/mol, better than the benchmark explicit solvent calculations from FreeSolv, and required less than 15 s of computing time per molecule on a single CPU core. Importantly, parameters for PMVECC showed systematic errors for molecules containing Cl, Br, I, and P. Applying ECC to the explicit solvent hydration free energies found the same systematic errors. The results strongly suggest that some small adjustments to the Lennard-Jones parameters for GAFF will lead to improved hydration free energy calculations for all solvent models.

Liquid structures of chloroethenes from molecular simulations and electronic structure calculations

The structures of six liquid chloroethenes viz., chloroethene, 1,1-Dichloroethene, cis- and trans-1,2-Dichloroethene, trichloro- and tetrachloroethene, are explored using the all-atom molecular dynamics simulations, molecular simulations with the 3D-reference interaction site model and the Kovalenko-Hirata closure relation (3D-RISM-KH), and dispersion corrected density functional theory. The liquid structures have different local orientations stabilized by possible π-stacking, π⋅⋅⋅Cl, and H⋅⋅⋅Cl intermolecular interactions. The RISM computations and MD simulations are in good agreement with each other, and provided important structural details. The electronic structure calculations support the presence of dimeric structures in the liquid form as well as the findings from the molecular simulations.

The Ion-Dipole Correction of the 3DRISM Solvation Model to Accurately Compute Water Distributions around Negatively Charged Biomolecules.

The 3D reference interaction site model (3DRISM) provides an efficient grid-based solvation model to compute the structural and thermodynamic properties of biomolecules in aqueous solutions. However, it remains challenging for existing 3DRISM methods to correctly predict water distributions around negatively charged solute molecules. In this paper, we first show that this challenge is mainly due to the orientation of water molecules in the first solvation shell of the negatively charged solute molecules. To properly consider this orientational preference, position-dependent two-body intramolecular correlations of solvent need to be included in the 3DRISM theory, but direct evaluations of these position-dependent two-body intramolecular correlations remain numerically intractable. To address this challenge, we introduce the Ion-Dipole Correction (IDC) to the 3DRISM theory, in which we incorporate the orientation preference of water molecules via an additional solute-solvent interaction term (i.e., the ion-dipole interaction) while keeping the formulism of the 3DRISM equation unchanged. We prove that this newly introduced IDC term is equivalent to an effective direct correlation function which can effectively consider the orientation effect that arises from position dependent two-body correlations. We first quantitatively validate our 3DRISM-IDC theory combined with the PSE3 closure on Cl-, [ClO]- (a two-site anion), and [NO2]- (a three-site anion). For all three anions, we show that our 3DRISM-IDC theory significantly outperforms the 3DRISM theory in accurately predicting the solvation structures in comparison to MD simulations, including RDFs and 3D water distributions. Furthermore, we have also demonstrated that the 3DRISM-IDC can improve the accuracy of hydration free-energy calculation for Cl-. We further demonstrate that our 3DRISM-IDC theory yields significant improvements over the 3DRISM theory when applied to compute the solvation structures for various negatively charged solute molecules, including adenosine triphosphate (ATP), a short peptide containing 19 residues, a DNA hairpin containing 24 nucleotides, and a riboswitch RNA molecule with 77 nucleotides. We expect that our 3DRISM-IDC-PSE3 solvation model holds great promise to be widely applied to study solvation properties for nucleic acids and other biomolecules containing negatively charged functional groups.

Understanding the Liquid States of Cyclic Hydrocarbons Containing N, O, and S Atoms via the 3D-RISM-KH Molecular Solvation Theory.

The 3D-reference interaction site model (3D-RISM) molecular solvation theory in combination with the Kovalenko–Hirata (KH) closure is extended to seven heterocyclic liquids to understand their liquid states and to test the performance of the theory in solvation free energy (SFE) calculations of solutes in select solvents. The computed solvent site distribution profiles were compared with the all-atom molecular dynamics (MD) simulations, showing comparable performances. The computational results were compared against the structural parameters for liquids, whenever available, as well as against the experimental SFEs. The liquids are found to have local ordered structures held together via weak interactions in both the RISM and MD simulations. The 3D-RISM-KH computed SFEs are in good agreement with the benchmark values for the tetrahydrothiophene-S,S-dioxide, and showed comparatively larger deviations in the case of the SFEs in the tetrahydrofuran continuum.

Accelerating the 3D reference interaction site model theory of molecular solvation with treecode summation and cut-offs.

The 3D reference interaction site model (3D-RISM) of molecular solvation is a powerful tool for computing the equilibrium thermodynamics and density distributions of solvents, such as water and co-ions, around solute molecules. However, 3D-RISM solutions can be expensive to calculate, especially for proteins and other large molecules where calculating the potential energy between solute and solvent requires more than half the computation time. To address this problem, we have developed and implemented treecode summation for long-range interactions and analytically corrected cut-offs for short-range interactions to accelerate the potential energy and long-range asymptotics calculations in non-periodic 3D-RISM in the AmberTools molecular modeling suite. For the largest single protein considered in this work, tubulin, the total computation time was reduced by a factor of 4. In addition, parallel calculations with these new methods scale almost linearly and the iterative solver remains the largest impediment to parallel scaling. To demonstrate the utility of our approach for large systems, we used 3D-RISM to calculate the solvation thermodynamics and density distribution of 7-ring microtubule, consisting of 910 tubulin dimers, over 1.2 million atoms.

Implementation and Optimization of the Embedded Cluster Reference Interaction Site Model with Atomic Charges.

In this work, we implemented the embedded cluster reference interaction site model (EC-RISM) originally developed by Kloss, Heil, and Kast (J. Phys. Chem. B2008, 112, 4337–4343). This method combines quantum mechanical calculations with the 3D reference interaction site model (3D-RISM). Numerous options, such as buffer, grid space, basis set, charge model, water model, closure relation, and so forth, were investigated to find the best settings. Additionally, the small point charges, which are derived from the solvent distribution from the 3D-RISM solution to represent the solvent in the QM calculation, were neglected to reduce the overhead without the loss of accuracy. On the MNSOL[a], MNSOL, and FreeSolv databases, our implemented and optimized method provides solvation free energies in water with 5.70, 6.32, and 6.44 kJ/mol root-mean-square deviations, respectively, but with different settings, 5.22, 6.08, and 6.63 kJ/mol can also be achieved. Only solvent models containing fitting parameters, like COSMO-RS and EC-RISM with universal correction and directly used electrostatic potential, perform better than our EC-RISM implementation with atomic charges.

Integral equation models for solvent in macromolecular crystals.

The solvent can occupy up to ∼70% of macromolecular crystals, and hence, having models that predict solvent distributions in periodic systems could improve the interpretation of crystallographic data. Yet, there are few implicit solvent models applicable to periodic solutes, and crystallographic structures are commonly solved assuming a flat solvent model. Here, we present a newly developed periodic version of the 3D-reference interaction site model (RISM) integral equation method that is able to solve efficiently and describe accurately water and ion distributions in periodic systems; the code can compute accurate gradients that can be used in minimizations or molecular dynamics simulations. The new method includes an extension of the Ornstein–Zernike equation needed to yield charge neutrality for charged solutes, which requires an additional contribution to the excess chemical potential that has not been previously identified; this is an important consideration for nucleic acids or any other charged system where most or all the counter- and co-ions are part of the “disordered” solvent. We present several calculations of proteins, RNAs, and small molecule crystals to show that x-ray scattering intensities and the solvent structure predicted by the periodic 3D-RISM solvent model are in closer agreement with the experiment than are intensities computed using the default flat solvent model in the refmac5 or phenix refinement programs, with the greatest improvement in the 2 to 4 Å range. Prospects for incorporating integral equation models into crystallographic refinement are discussed.

Asymmetric Interplay Between K+ and Blocker and Atomistic Parameters From Physiological Experiments Quantify K+ Channel Blocker Release.

Modulating the activity of ion channels by blockers yields information on both the mode of drug action and on the biophysics of ion transport. Here we investigate the interplay between ions in the selectivity filter (SF) of K+ channels and the release kinetics of the blocker tetrapropylammonium in the model channel KcvNTS. A quantitative expression calculates blocker release rate constants directly from voltage-dependent ion occupation probabilities in the SF. The latter are obtained by a kinetic model of single-channel currents recorded in the absence of the blocker. The resulting model contains only two adjustable parameters of ion-blocker interaction and holds for both symmetric and asymmetric ionic conditions. This data-derived model is corroborated by 3D reference interaction site model (3D RISM) calculations on several model systems, which show that the K+ occupation probability is unaffected by the blocker, a direct consequence of the strength of the ion-carbonyl attraction in the SF, independent of the specific protein background. Hence, KcvNTS channel blocker release kinetics can be reduced to a small number of system-specific parameters. The pore-independent asymmetric interplay between K+ and blocker ions potentially allows for generalizing these results to similar potassium channels.

A molecular solvation theory simulation of liquid alkyl esters of acetic acid with the 3D Reference Interaction Site Model

The liquid states of three aliphatic esters, viz. methyl-, ethyl-, and n-butyl-acetate are modeled using the three-dimensional reference interaction site model with the Kovalenko-Hirata closure in combination with the united atom AMBER force field. The distributions of different solvent sites were compared with those obtained from all-atom molecular dynamics simulations. The possibility of existence of dimeric and tetrameric clusters in liquids were explored using the dispersion corrected density functional theory. The molecular solvation theory showed good accuracy in solvation free energy calculations and in ester to water partition coefficient calculations for a large set of organic compounds.

Structural basis of mutants of PET-degrading enzyme from Saccharomonospora viridis AHK190 with high activity and thermal stability.

The cutinase-like enzyme from the thermophile Saccharomonospora viridis AHK190, Cut190, is a good candidate to depolymerize polyethylene terephthalate (PET) efficiently. We previously developed a mutant of Cut190 (S226P/R228S), which we designated as Cut190* that has both increased activity and stability and solved its crystal structure. Recently, we showed that mutation of D250C/E296C on one of the Ca2+ -binding sites resulted in a higher thermal stability while retaining its polyesterase activity. In this study, we solved the crystal structures of Cut190* mutants, Q138A/D250C-E296C/Q123H/N202H, designated as Cut190*SS, and its inactive S176A mutant, Cut190*SS_S176A, at high resolution. The overall structures were similar to those of Cut190* and Cut190*S176A reported previously. As expected, Cys250 and Cys296 were closely located to form a disulfide bond, which would assuredly contribute to increase the stability. Isothermal titration calorimetry experiments and 3D Reference Interaction Site Model calculations showed that the metal-binding properties of the Cut190*SS series were different from those of the Cut190* series. However, our results show that binding of Ca2+ to the weak binding site, site 1, would be retained, enabling Cut190*SS to keep its ability to use Ca2+ to accelerate the conformational change from the closed (inactive) to the open (active) form. While increasing the thermal stability, Cut190*SS could still express its enzymatic function. Even after incubation at 70°C, which corresponds to the glass transition temperature of PET, the enzyme retained its activity well, implying a high applicability for industrial PET depolymerization using Cut190*SS.

An online repository of solvation thermodynamic and structural maps of SARS-CoV-2 targets.

SARS-CoV-2 recently jumped species and rapidly spread via human-to-human transmission to cause a global outbreak of COVID-19. The lack of effective vaccine combined with the severity of the disease necessitates attempts to develop small molecule drugs to combat the virus. COVID19_GIST_HSA is a freely available online repository to provide solvation thermodynamic maps of COVID-19-related protein small molecule drug targets. Grid inhomogeneous solvation theory maps were generated using AmberTools cpptraj-GIST, 3D reference interaction site model maps were created with AmberTools rism3d.snglpnt and hydration site analysis maps were created using SSTMap code. The resultant data can be applied to drug design efforts: scoring solvent displacement for docking, rational lead modification, prioritization of ligand- and protein- based pharmacophore elements, and creation of water-based pharmacophores. Herein, we demonstrate the use of the solvation thermodynamic mapping data. It is hoped that this freely provided data will aid in small molecule drug discovery efforts to defeat SARS-CoV-2.

Application of the Approximate 3D-Reference Interaction Site Model (RISM) Molecular Solvation Theory to Acetonitrile as Solvent.

The integral equation formalism based reference interaction site model (RISM) molecular solvation theory is applied to pure liquid acetonitrile and water-acetonitrile binary mixtures of different compositions. Solvate formation of d- and f-block ions by ACN is also calculated to check applicability of the RISM theory. The generalized Amber force field (GAFF) and the universal force field (UFF) parameters were found to be suitable for applications with the RISM theory for acetonitrile solvent. The presence of local microsolvated clusters for water-acetonitrile mixtures is validated in the RISM calculations.

A molecular reconstruction approach to site-based 3D-RISM and comparison to GIST hydration thermodynamic maps in an enzyme active site.

Computed, high-resolution, spatial distributions of solvation energy and entropy can provide detailed information about the role of water in molecular recognition. While grid inhomogeneous solvation theory (GIST) provides rigorous, detailed thermodynamic information from explicit solvent molecular dynamics simulations, recent developments in the 3D reference interaction site model (3D-RISM) theory allow many of the same quantities to be calculated in a fraction of the time. However, 3D-RISM produces atomic-site, rather than molecular, density distributions, which are difficult to extract physical meaning from. To overcome this difficulty, we introduce a method to reconstruct molecular density distributions from atomic-site density distributions. Furthermore, we assess the quality of the resulting solvation thermodynamics density distributions by analyzing the binding site of coagulation Factor Xa with both GIST and 3D-RISM. We find good qualitative agreement between the methods for oxygen and hydrogen densities as well as direct solute-solvent energetic interactions. However, 3D-RISM predicts lower energetic and entropic penalties for moving water from the bulk to the binding site.

Recent Developments in Integral Equation Theory for Solvation to Treat Density Inhomogeneity at Solute–Solvent Interface

AbstractThe integration equation theory (IET) provides highly efficient tools for the calculation of structural and thermodynamic properties of molecular liquids. In recent years, the 3D reference interaction site model (3DRISM), the most developed IET for solvation, has been widely applied to study protein solvation, aggregation, and drug‐receptor binding. However, hydrophobic solutes with sufficient size (>nm) can induce water density depletion at the solute–solvent interface. This density depletion is not considered in the original 3DRISM theory. The authors here review the recent developments of 3DRISM at hydrophobic surfaces and related theories to address this challenge. At hydrophobic surfaces, an additional hydrophobicity‐induced density inhomogeneity equation is introduced to 3DRISM theory to consider this density depletion. Accordingly, several new closures equations including D2 closure and D2MSA closures are developed to enable stable numerical solutions of 3DRISM equations. These newly developed theories hold great promise for an accurate and rapid calculation of the solvation effect for complex molecular systems such as proteins. At the end of the report, the authors also provide a perspective on other challenges of the IETs as an efficient solvation model.

Metal binding to cutinase-like enzyme from Saccharomonospora viridis AHK190 and its effects on enzyme activity and stability.

A cutinase from Saccharomonospora viridis AHK190, Cut190, can hydrolyze polyethylene terephthalate and has a unique feature that the activity and stability are regulated by Ca2+ binding. Our recent structural and functional analyses showed three Ca2+ binding sites and their respective roles. Here, we analysed the binding thermodynamics of Mn2+, Zn2+ and Mg2+ to Cut190 and their effects on the catalytic activity and thermal stability. The binding affinities of Mn2+ and Zn2+ were higher than that of Mg2+ and are all entropy driven with a binding stoichiometry of three, one and one for Zn2+, Mn2+ and Mg2+, respectively. The catalytic activity was measured in the presence of the respective metals, where the activity of 0.25 mM Mn2+ was comparable to that of 2.5 mM Ca2+. Our 3D Reference Interaction Site Model calculations suggested that all the ions exhibited a high occupancy rate for Site 2. Thus, Mn2+ and Mg2+ would most likely bind to Site 2 (contributes to stability) with high affinity, while to Sites 1 and 3 (contributes to activity) with low affinity. We elucidate the metal-dependent structural and functional properties of Cut190 and show the subtle balance on structure stability and flexibility is controlled by specific metal ions.

Role of Mg2+ Ions in DNA Hydrolysis by EcoRV, Studied by the 3D-Reference Interaction Site Model and Molecular Dynamics.

The role of Mg2+ ions during precursor formation in DNA hydrolysis by the homodimeric restriction enzyme EcoRV was elucidated based on the 3D-reference interaction site model (RISM) theory and the molecular dynamics (MD) simulation. From an analysis of the spatial distribution of Mg2+ in an active site using 3D-RISM, we identified a new position for Mg2+ in the X-ray EcoRV-DNA complex structure ( 1rvb ). We refer to the position as site IV†. Site IV† is almost at the same position as that of a Ca2+ ion in the superimposed X-ray crystallographic active-site structure of the PvuII-DNA complex ( 1f0o ). 3D-RISM was also used to locate the position of water molecules, including the water nucleophile at the active site. MD simulations were carried out with the initial structure having two Mg2+ ions at site IV† and at site I*, experimentally identified by Horton et al., to find a stable complex structure in which the DNA fragment was rearranged to orient the scissile bond direction toward the water nucleophile. The equilibrium active-site structure of the EcoRV-DNA complex obtained from the MD simulation was similar to the superimposed X-ray crystallographic structure of the BamHI-DNA complex ( 2bam ). In the active-site structure, two metal ions have almost the same position (≤1.0 Å) as that of 2bam , and the scissile phosphate is twisted to orient the scissile bond toward the water nucleophile, as is the case in 2bam . We propose the equilibrium active-site structure obtained in this study as a precursor for the hydrolysis reaction of EcoRV.