The Amber ff99 Force Field Predicts Relative Free Energy Changes for RNA Helix Formation

Exploring the neurotoxicity of chiral dinotefuran towards nicotinic acetylcholine receptors: Enantioselective insights into species selectivity

Dinotefuran is a chiral neonicotinoid that is widely distributed in environmental matrices, but its health risks to different organisms are poorly understood. This study investigated the neurotoxic responses of honeybee/cotton aphid nicotinic acetylcholine receptors (nAChRs) to chiral dinotefuran at the enantiomeric scale and demonstrated the microscopic mechanism of species selectivity in nAChR-mediated enantioselective neurotoxicity. The findings indicated that (S)-dinotefuran had a higher affinity for honeybee nAChR than (R)-dinotefuran whereas both enantiomers exhibited similar bioactivity toward cotton aphid nAChR. The results of dynamic neurotoxic processes indicated the association of conformational changes induced by chiral dinotefuran with its macroscopic neurotoxicity, and (R)-dinotefuran, which exhibit low toxicity to honeybee, was found to induce significant conformational changes in the enantioselective neurotoxic reaction, as supported by the average root-mean-square fluctuation (0.35 nm). Energy decomposition results indicated that electrostatic contribution (ΔGele) is the critical energy term that leads to substantial enantioselectivity, and both Trp-51 (－2.57 kcal mol-1) and Arg-75 (－4.86 kcal mol-1), which form a hydrogen-bond network, are crucial residues in mediating the species selectivity for enantioselective neurotoxic responses. Clearly, this study provides experimental evidence for a comprehensive assessment of the health hazards of chiral dinotefuran.

Interactions between Water and a Hydrophobic Polymer

To investigate the mechanisms of interactions between a hydrophobic polymer and water, molecular dynamics calculations and Raman spectroscopic measurements of cis-1,4-polyisoprene (PI)-water systems were performed. The results show that PI in water undergoes a coil-globule transition at around 248 K. The transition is attributed to changes in the density and diffusivity of water. The volume expansion of the supercooled liquid water induces the coil structure of PI. The phase separation of PI from water with an increase in the self-diffusion coefficient of water molecules results in the globule structure of PI. The self-diffusion coefficient of free water with PI is larger than that of pure water because PI has an effect to decrease the hydrogen-bonding strength of water. The result suggests that the effects of the coexisting water are important factors governing the physical and chemical properties of hydrophobic polymers.

QSAR modelling, molecular docking, molecular dynamic and ADMET prediction of pyrrolopyrimidine derivatives as novel Bruton's tyrosine kinase (BTK) inhibitors

In recent years, there has been a focus on developing and discovering novel Bruton's tyrosine kinase (BTK) inhibitors, as they offer an effective treatment strategy for B-cell malignancies. BTK plays a crucial role in B cell receptor (BCR)-mediated activation and proliferation by regulating downstream factors such as the NF-κB and MAP kinase pathways. To address this challenge and propose potential therapeutic options for B-cell lymphomas, researchers conducted 2D-QSAR and ADMET studies on pyrrolopyrimidine derivatives that act as inhibitors of the BCR site in cytochrome b. These studies aim to improve and identify new compounds that could serve as more potent potential BTK inhibitors, which would lead to the identification of new drug candidates in this field. In our study, we used 2D-QSAR (multiple linear regression, multiple nonlinear regression, and artificial neural networks), molecular docking, molecular dynamics, and ADMET properties to investigate the potential of 35 pyrrolopyrimidine derivatives as BTK inhibitors. A molecular docking study and molecular dynamics simulations of molecule 13 over 10 ns revealed that it establishes multiple hydrogen bonds with several residues and exhibits frequent stability throughout the simulation period. Based on the results obtained by molecular modeling, we proposed six new compounds (Pred1, Pred2, Pred3, Pred4, Pred5, and Pred6) with highly significant predicted activity by MLR models. A study based on the in silico evaluation of the predicted ADMET properties of the new candidate molecules is strongly recommended to classify these molecules as promising candidates for new anticancer agents specifically designed to target Bruton's tyrosine kinase (BTK) inhibition.

Combination of Resonance and Non-Resonance Chiral Raman Scattering in a Cobalt(III) Complex

Resonance Raman optical activity (RROA) spectra with high sensitivity reveal details on molecular structure, chirality, and excited electronic properties. Despite the difficulty of the measurements, the recorded data for the Co(III) complex with S,S-N,N-ethylenediaminedisuccinic acid are of exceptional quality and, coupled with the theory, spectacularly document the molecular behavior in resonance. This includes a huge enhancement of the chiral scattering, contribution of the antisymmetric polarizabilities to the signal, and the Herzberg-Teller effect significantly shaping the spectra. The chiral component is by about one order of magnitude bigger than for an analogous aluminum complex. The band assignment and intensity profile were confirmed by simulations based on density functional and vibronic theories. The resonance was attributed to the S0 →S3 transition, with the strongest signal enhancement of Raman and ROA spectral bands below about 800 cm-1 . For higher wavenumbers, other excited electronic states contribute to the scattering in a less resonant way. RROA spectroscopy thus appears as a unique tool to study the structure and electronic states of absorbing molecules in analytical chemistry, biology, and material science.

Identification of potential mutational hotspots in serratiopeptidase to address its poor pH tolerance issue

Serratiopeptidase is the multifunctionality metalloendopeptidase extensively employed in biopharmaceutical and industrial biotechnology. Despite its poor pH tolerance, serratiopeptidase must withstand the highly acidic environment of the gastrointestinal tract to be used as a potent anti-inflammatory and analgesic medication. In earlier studies, post-translational deamination related mutations showed alteration in the net charge of protein's surface. Therefore, the current study aimed to enhance the acid resistance of serratiopeptidase via implementing computational interventions to screen out the most stable mutational hotspot. The methodology used in this study is as follows: (a) Higher accessibility to surface (b) 4 Å away from active site region to avoid interference with its proteolytic activity, and (c) By converting non-conserved amide residues to acidic residues. A docking study has been conducted to establish the substrate specificity and binding affinity to native and mutant proteins. The docking outcomes were then validated using molecular dynamic simulations to clarify each mutant's molecular stability and conformation while preserving their activity. The results showed that N412D is the best-screened mutant with negative electrostatic potential that can alter the overall charge on the protein's surface with increased H+ ions. Alteration in overall charge leads to protein surface more acidic that causes a common ion effect in stomach pH and act as a buffer which could stabilize the serratiopeptidase amid extreme pH.Communicated by Ramaswamy H. Sarma.

Structure and Diffusive Properties of Water in Polymer Hydrogels

To investigate the diffusive properties of water in hydrogels, ab initio molecular orbital and molecular dynamics calculations of poly-N,N-dimethylacrylamide-water systems were performed. The results show that the mean diffusion coefficient of water drastically decreases in the middle dehydration stage. In this stage, the mobilities of water are restrained because part of the water forms hydrogen bonds to bind polymer chains due to a glass transition. In addition to the three well-known types of water (i.e., bound, intermediate, and free water) around hydrophilic polymers in hydrogels, our results suggest that the intermediate water can be further classified into two types: first and second intermediate water. The bound and first intermediate water acts as a cross-linker between polymer chains, even if the polymer does not form intra- or intermolecular bonds itself. The diffusive properties of water might have important implications for the interpretation of properties of polymer hydrogels.

A Test and Refinement of Folding Free Energy Nearest Neighbor Parameters for RNA Including N6-Methyladenosine

RNA folding free energy change parameters are widely used to predict RNA secondary structure and to design RNA sequences. These parameters include terms for the folding free energies of helices and loops. Although the full set of parameters has only been traditionally available for the four common bases and backbone, it is well known that covalent modifications of nucleotides are widespread in natural RNAs. Covalent modifications are also widely used in engineered sequences. We recently derived a full set of nearest neighbor terms for RNA that includes N6-methyladenosine (m6A). In this work, we test the model using 98 optical melting experiments, matching duplexes with or without N6-methylation of A. Most experiments place RRACH, the consensus site of N6-methylation, in a variety of contexts, including helices, bulge loops, internal loops, dangling ends, and terminal mismatches. For matched sets of experiments that include either A or m6A in the same context, we find that the parameters for m6A are as accurate as those for A. Across all experiments, the root mean squared deviation between estimated and experimental free energy changes is 0.67 kcal/mol. We used the new experimental data to refine the set of nearest neighbor parameter terms for m6A. These parameters enable prediction of RNA secondary structures including m6A, which can be used to model how N6-methylation of A affects RNA structure.

Secondary structure prediction for RNA sequences including N6-methyladenosine

There is increasing interest in the roles of covalently modified nucleotides in RNA. There has been, however, an inability to account for modifications in secondary structure prediction because of a lack of software and thermodynamic parameters. We report the solution for these issues for N⁶-methyladenosine (m⁶A), allowing secondary structure prediction for an alphabet of A, C, G, U, and m⁶A. The RNAstructure software now works with user-defined nucleotide alphabets of any size. We also report a set of nearest neighbor parameters for helices and loops containing m⁶A, using experiments. Interestingly, N⁶-methylation decreases folding stability for adenosines in the middle of a helix, has little effect on folding stability for adenosines at the ends of helices, and increases folding stability for unpaired adenosines stacked on a helix. We demonstrate predictions for an N⁶-methylation-activated protein recognition site from MALAT1 and human transcriptome-wide effects of N⁶-methylation on the probability of adenosine being buried in a helix.

RNA folding free energy nearest neighbor parameters were determined for sequences with the nucleotide m⁶A. The RNAstructure software package can accommodate modified nucleotides, enabling secondary structure prediction of sequences with m⁶A.

The development of an Amber-compatible organosilane force field for drug-like small molecules

As members of the group IVA elements, silicon and carbon have long been thought of as isosteres of each other in drug design. However, the lack of silicon parameters in current main stream force fields hinders the computational study of this important element in drug discovery. Thus, in this study, we attempted to supplement the parameters of organosilanes in the General Amber Force Field (GAFF2). The parameters have been designed following the principles of GAFF2 to make it compatible with the Amber force field family. The accuracy of the parameters was discussed by comparing the pair interaction energy, the liquid properties, and the structures and alchemical binding free energy differences for a set of protein-ligand complexes.

Acetylcholinesterase inhibitory activities and bioguided fractionation of the Ocotea percoriacea extracts: HPLC-DAD-MS/MS characterization and molecular modeling of their alkaloids in the active fraction

In vitro acetylcholinesterase activities of the hexane, dichloromethane, ethyl acetate, n-butanol and aqueous extracts of leaves of Ocotea percoriacea Kosterm. (Lauraceae) were evaluated. The bioguided fractionation of the most active extract (dichloromethane) using silica gel open-column chromatography led to an active alkaloidal fraction composed of isocorydine N-oxide, isocorydine N-oxide derivative, palmatine, roemerine and roemerine N-Oxide. The identification of the chemical structure of these compounds was carried out with high-performance liquid chromatography coupled to electrospray ionization multiple-stage mass spectrometry (HPLC-ESI-MS/MS). Aiming to understand their inhibitory activities, these alkaloids were docked into a 3D model of Electrophorus electricus Acetylcholinesterase (EelAChE) built in the Modeller 9.18 employing homology modeling approach. The results suggest that the alkaloids had the same binding mode and, possibly, the inhibition mechanism of classic drugs (ex. tacrine and donepezil). The structural difference of these compounds opens a new opportunity for the optimization of leading compounds.

Chemically Accurate Relative Folding Stability of RNA Hairpins from Molecular Simulations

To benchmark RNA force fields, we compared the folding stabilities of three 12-nucleotide hairpin stem loops estimated by simulation to stabilities determined by experiment. We used umbrella sampling and a reaction coordinate of end-to-end (5' to 3' hydroxyl oxygen) distance to estimate the free energy change of the transition from the native conformation to a fully extended conformation with no hydrogen bonds between non-neighboring bases. Each simulation was performed four times using the AMBER FF99+bsc0+χ_OL3 force field, and each window, spaced at 1 Å intervals, was sampled for 1 μs, for a total of 552 μs of simulation. We compared differences in the simulated free energy changes to analogous differences in free energies from optical melting experiments using thermodynamic cycles where the free energy change between stretched and random coil sequences is assumed to be sequence-independent. The differences between experimental and simulated ΔΔ G° are, on average, 0.98 ± 0.66 kcal/mol, which is chemically accurate and suggests that analogous simulations could be used predictively. We also report a novel method to identify where replica free energies diverge along a reaction coordinate, thus indicating where additional sampling would most improve convergence. We conclude by discussing methods to more economically perform these simulations.

Thermodynamics of the fourU RNA thermal switch derived from molecular dynamics simulations and spectroscopic techniques

Bacterial 5' untranslated regions of mRNA, termed thermal switches or thermometers, change their structure in response to temperature change. This structural change provides for the regulation of gene expression. One of such thermal switches, called fourU, is present in the Salmonella species. Mutations of fourU were found to abrogate its regulatory properties. We investigated the thermodynamics of the fourU fragment responsible for its structural changes. All-atom molecular dynamics simulations at various temperatures and spectroscopic experiments in solution were performed for the wild-type fourU and its mutants. We found that the U11C and A8C mutations stabilize the fourU structure in comparison to the wild-type fourU, and the double-point G14A/C25U mutant has the most destabilizing effect on the fourU hairpin 2 responsible for temperature sensing. The G14A/C25U mutant is also the easiest to strand-invade by a complementary oligonucleotide as indicated by fluorescence spectroscopy experiments.

Flavonones from Penthorum chinense Ameliorate Hepatic Steatosis by Activating the SIRT1/AMPK Pathway in HepG2 Cells

Pinocembrin-7-O-β-d-glucoside (PCBG), pinocembrin (PCB), and 5-methoxy-pinocembrin-7-O-β-d-glucoside (MPG) are three flavonones isolated from Penthorum chinense Pursh (P. chinense). The effects of the three flavonones on hepatic steatosis and their molecular mechanisms in HepG2 cells were investigated in this study for the first time. A model of hepatic steatosis in HepG2 cells was induced by free fatty acid (FFA), and co-treated with the three flavonones as mentioned. Intracellular lipid droplets were detected by Oil Red O staining. PCB, PCBG, and MPG suppressed oxidative stress by decreasing malondialdehyde (MDA) levels and increasing superoxide dismutase (SOD) and glutathione peroxidase (GSH-Px) activities. The levels of aspartate aminotransferase (AST) and alanine aminotransferase (ALT) were ameliorated. Moreover, these flavonones enhanced the phosphorylation of AMP-activated protein kinase (AMPK) and the expression of silent mating type information regulation 2 homolog 1 (SIRT1) and peroxisome proliferator-activated receptor α (PPARα), and reduced the expression of sterol regulatory element binding protein-1c (SREBP1c) and the downstream targets fatty acid synthase (FAS), acetyl-CoA carboxylase (ACC), and stearoyl-CoA desaturase 1 (SCD1). Molecular docking was used to predict the interaction and combination patterns between the three flavonones and the enzymes above. The results revealed that the SIRT1/AMPK pathway is involved in the functions of the three flavonones, and the most effective flavonone against hepatic steatosis might be PCBG, followed by MPG and PCB. Therefore, the three flavonones from P. chinense were found to exert preventive effects against hepatic steatosis by regulating the SIRT1/AMPK pathway.

Molecular dynamics correctly models the unusual major conformation of the GAGU RNA internal loop and with NMR reveals an unusual minor conformation

The RNA "GAGU" duplex, (5'GACGAGUGUCA)2, contains the internal loop (5'-GAGU-3')2 , which has two conformations in solution as determined by NMR spectroscopy. The major conformation has a loop structure consisting of trans-Watson-Crick/Hoogsteen GG pairs, A residues stacked on each other, U residues bulged outside the helix, and all sugars with a C2'-endo conformation. This differs markedly from the internal loops, (5'-GAGC-3')2, (5'-AAGU-3')2, and (5'-UAGG-3')2, which all have cis-Watson-Crick/Watson-Crick AG "imino" pairs flanked by cis-Watson-Crick/Watson-Crick canonical pairs resulting in maximal hydrogen bonding. Here, molecular dynamics was used to test whether the Amber force field (ff99 + bsc0 + OL3) approximates molecular interactions well enough to keep stable the unexpected conformation of the GAGU major duplex structure and the NMR structures of the duplexes containing (5'-GAGC-3')2, (5'-AAGU-3')2, and (5'-UAGG-3')2 internal loops. One-microsecond simulations were repeated four times for each of the duplexes starting in their NMR conformations. With the exception of (5'-UAGG-3')2, equivalent simulations were also run starting with alternative conformations. Results indicate that the Amber force field keeps the NMR conformations of the duplexes stable for at least 1 µsec. They also demonstrate an unexpected minor conformation for the (5'-GAGU-3')2 loop that is consistent with newly measured NMR spectra of duplexes with natural and modified nucleotides. Thus, unrestrained simulations led to the determination of the previously unknown minor conformation. The stability of the native (5'-GAGU-3')2 internal loop as compared to other loops can be explained by changes in hydrogen bonding and stacking as the flanking bases are changed.

Effect of PDGF-B aptamer on PDGFRβ/PDGF-B interaction: Molecular dynamics study

PDGFRβ/PDGF-B interaction plays a role in angiogenesis, and is mandatory in wound healing and cancer treatment. It has been reported that the PDGF-B aptamer was able to bind to PDGF-B, thus regulating the angiogenesis. However, the binding interaction between the aptamer and the growth factor, including the binding sites, has not been well investigated. This study applied a molecular dynamics (MD) simulation to investigate the aptamer-growth factor interaction in the presence or absence of a receptor (PDGFRβ). Characterization of the structure of an aptamer-growth factor complex revealed binding sites from each section in the complex. Upon the complex formation, PDGF-B and its aptamer exhibited less flexibility in their molecular movement, as indicated by the minimum values of RMSD, RMSF, loop-to-loop distance, and the summation of PCA eigenvalues. Our study of residue pairwise interaction demonstrated that the binding interaction was mainly contributed by electrostatic interaction between the positively-charged amino acid and the negatively-charged phosphate backbone. The role of the PDGF-B aptamer in PDGFRβ/PDGF-B interaction was also investigated. We demonstrated that the stability of the Apt-PDGF-B complex could prevent the presence of a competitor, of PDGFRβ, interrupting the binding process. Because the aptamer was capable of binding with PDGF-B, and blocking the growth factor from the PDGFRβ, it could down regulate the consequent signaling pathway. We provide evidence that the PDGF-BB aptamer is a promising molecule for regulation of angiogenesis. The MD study provides a molecular understanding to modification of the aptamer binding interaction, which could be used in a number of medical applications.

RNA Structural Dynamics As Captured by Molecular Simulations: A Comprehensive Overview

With both catalytic and genetic functions, ribonucleic acid (RNA) is perhaps the most pluripotent chemical species in molecular biology, and its functions are intimately linked to its structure and dynamics. Computer simulations, and in particular atomistic molecular dynamics (MD), allow structural dynamics of biomolecular systems to be investigated with unprecedented temporal and spatial resolution. We here provide a comprehensive overview of the fast-developing field of MD simulations of RNA molecules. We begin with an in-depth, evaluatory coverage of the most fundamental methodological challenges that set the basis for the future development of the field, in particular, the current developments and inherent physical limitations of the atomistic force fields and the recent advances in a broad spectrum of enhanced sampling methods. We also survey the closely related field of coarse-grained modeling of RNA systems. After dealing with the methodological aspects, we provide an exhaustive overview of the available RNA simulation literature, ranging from studies of the smallest RNA oligonucleotides to investigations of the entire ribosome. Our review encompasses tetranucleotides, tetraloops, a number of small RNA motifs, A-helix RNA, kissing-loop complexes, the TAR RNA element, the decoding center and other important regions of the ribosome, as well as assorted others systems. Extended sections are devoted to RNA-ion interactions, ribozymes, riboswitches, and protein/RNA complexes. Our overview is written for as broad of an audience as possible, aiming to provide a much-needed interdisciplinary bridge between computation and experiment, together with a perspective on the future of the field.

Erythromycin leads to differential protein expression through differences in electrostatic and dispersion interactions with nascent proteins

The antibiotic activity of erythromycin, which reversibly binds to a site within the bacterial ribosome exit tunnel, against many gram positive microorganisms indicates that it effectively inhibits the production of proteins. Similar to other macrolides, the activity of erythromycin is far from universal, as some peptides can bypass the macrolide-obstructed exit tunnel and become partially or fully synthesized. It is unclear why, at the molecular level, some proteins can be synthesized while others cannot. Here, we use steered molecular dynamics simulations to examine how erythromycin inhibits synthesis of the peptide ErmCL but not the peptide H-NS. By pulling these peptides through the exit tunnel of the E.coli ribosome with and without erythromycin present, we find that erythromycin directly interacts with both nascent peptides, but the force required for ErmCL to bypass erythromycin is greater than that of H-NS. The largest forces arise three to six residues from their N-terminus as they start to bypass Erythromycin. Decomposing the interaction energies between erythromycin and the peptides at this point, we find that there are stronger electrostatic and dispersion interactions with the more C-terminal residues of ErmCL than with H-NS. These results suggest that erythromycin slows or stalls synthesis of ErmCL compared to H-NS due to stronger interactions with particular residue positions along the nascent protein.

Physics-based all-atom modeling of RNA energetics and structure

The database of RNA sequences is exploding, but knowledge of energetics, structures, and dynamics lags behind. All-atom computational methods, such as molecular dynamics, hold promise for closing this gap. New algorithms and faster computers have accelerated progress in improving the reliability and accuracy of predictions. Currently, the methods can facilitate refinement of experimentally determined nuclear magnetic resonance and x-ray structures, but are 'unreliable' for predictions based only on sequence. Much remains to be discovered, however, about the many molecular interactions driving RNA folding and the best way to approximate them quantitatively. The large number of parameters required means that a wide variety of experimental results will be required to benchmark force fields and different approaches. As computational methods become more reliable and accessible, they will be used by an increasing number of biologists, much as x-ray crystallography has expanded. Thus, many fundamental physical principles underlying the computational methods are described. This review presents a summary of the current state of molecular dynamics as applied to RNA. It is designed to be helpful to students, postdoctoral fellows, and faculty who are considering or starting computational studies of RNA. WIREs RNA 2017, 8:e1422. doi: 10.1002/wrna.1422.

Conformational stabilization of FOX-DNA complex architecture to sensitize prostate cancer chemotherapy

The forkhead box (FOX) transcription factor is a family of tumor suppressors that negatively regulates the tumorigenesis activity of prostate cancer; stabilization of FOX-DNA complex architecture has been recognized as a new and promising strategy for sensitizing cancer chemotherapy. Here, we described a systematic method that combined in silico analysis and in vitro assay to investigate the intermolecular interaction between FOX DNA-binding domain (DBD) and its cognate DNA partner. The structural and energetic information harvested from the molecular investigation were used to guide high-throughput virtual screening against a structurally diverse, nonredundant library of natural product compounds, aiming at discovery of novel small-molecule medicines that can conformationally stabilize and promote FOX-DNA recognition and interaction. The screening identified a number of theoretically promising hits, which were then examined by using fluorescence anisotropy assay to determine their binding potency for FOX DBD domain. The antitumor activity of identified high-affinity compounds was also tested at cellular level. Structural dynamics analysis found that the small-molecule stabilizers can shift the conformational equilibrium of FOX DBD to DNA-bound state, thus promoting the protein domain to bind tightly with its DNA partner.

Energetics and structural characterization of the “DFG-flip” conformational transition of B-RAF kinase: a SITS molecular dynamics study

A combination of a homology modeling technique and an enhanced sampling molecular dynamics simulation implemented using the SITS method is employed to compute a detailed map of the free-energy landscape and explore the conformational transition pathway of B-RAF kinase.

Blind tests of RNA nearest-neighbor energy prediction

The predictive modeling and design of biologically active RNA molecules requires understanding the energetic balance among their basic components. Rapid developments in computer simulation promise increasingly accurate recovery of RNA's nearest-neighbor (NN) free-energy parameters, but these methods have not been tested in predictive trials or on nonstandard nucleotides. Here, we present, to our knowledge, the first such tests through a RECCES-Rosetta (reweighting of energy-function collection with conformational ensemble sampling in Rosetta) framework that rigorously models conformational entropy, predicts previously unmeasured NN parameters, and estimates these values' systematic uncertainties. RECCES-Rosetta recovers the 10 NN parameters for Watson-Crick stacked base pairs and 32 single-nucleotide dangling-end parameters with unprecedented accuracies: rmsd of 0.28 kcal/mol and 0.41 kcal/mol, respectively. For set-aside test sets, RECCES-Rosetta gives rmsd values of 0.32 kcal/mol on eight stacked pairs involving G-U wobble pairs and 0.99 kcal/mol on seven stacked pairs involving nonstandard isocytidine-isoguanosine pairs. To more rigorously assess RECCES-Rosetta, we carried out four blind predictions for stacked pairs involving 2,6-diaminopurine-U pairs, which achieved 0.64 kcal/mol rmsd accuracy when tested by subsequent experiments. Overall, these results establish that computational methods can now blindly predict energetics of basic RNA motifs, including chemically modified variants, with consistently better than 1 kcal/mol accuracy. Systematic tests indicate that resolving the remaining discrepancies will require energy function improvements beyond simply reweighting component terms, and we propose further blind trials to test such efforts.

Are Waters around RNA More than Just a Solvent? - An Insight from Molecular Dynamics Simulations

Hydrating water molecules are believed to be an inherent part of the RNA structure and have a considerable impact on RNA conformation. However, the magnitude and mechanism of the interplay between water molecules and the RNA structure are still poorly understood. In principle, such hydration effects can be studied by molecular dynamics (MD) simulations. In our recent MD studies, we observed that the choice of water model has a visible impact on the predicted structure and structural dynamics of RNA and, in particular, has a larger effect than type, parametrization, and concentration of the ions. Furthermore, the water model effect is sequence dependent and modulates the sequence dependence of A-RNA helical parameters. Clearly, the sensitivity of A-RNA structural dynamics to the water model parametrization is a rather spurious effect that complicates MD studies of RNA molecules. These results nevertheless suggest that the sequence dependence of the A-RNA structure, usually attributed to base stacking, might be driven by the structural dynamics of specific hydration. Here, we present a systematic MD study that aimed to (i) clarify the atomistic mechanism of the water model sensitivity and (ii) discover whether and to what extent specific hydration modulates the A-RNA structural variability. We carried out an extended set of MD simulations of canonical A-RNA duplexes with TIP3P, TIP4P/2005, TIP5P, and SPC/E explicit water models and found that different water models provided a different extent of water bridging between 2'-OH groups across the minor groove, which in turn influences their distance and consequently also inclination, roll, and slide parameters. Minor groove hydration is also responsible for the sequence dependence of these helical parameters. Our simulations suggest that TIP5P is not optimal for RNA simulations.

Evaluation of the Gibbs Free Energy Changes and Melting Temperatures of DNA/DNA Duplexes Using Hybridization Enthalpy Calculated by Molecular Dynamics Simulation

A molecular dynamics simulation approach was applied for the prediction of the thermal stability of oligonucleotide duplexes. It was shown that the enthalpy of the DNA/DNA complex formation could be calculated using this approach. We have studied the influence of various simulation parameters on the secondary structure and the hybridization enthalpy value of Dickerson-Drew dodecamer. The optimal simulation parameters for the most reliable prediction of the enthalpy values were determined. The thermodynamic parameters (enthalpy and entropy changes) of a duplex formation were obtained experimentally for 305 oligonucleotides of various lengths and GC-content. The resulting database was studied with molecular dynamics (MD) simulation using the optimized simulation parameters. Gibbs free energy changes and the melting temperatures were evaluated using the experimental correlation between enthalpy and entropy changes of the duplex formation and the enthalpy values calculated by the MD simulation. The average errors in the predictions of enthalpy, the Gibbs free energy change, and the melting temperature of oligonucleotide complexes were 11%, 10%, and 4.4 °C, respectively. We have shown that the molecular dynamics simulation gives a possibility to calculate the thermal stability of native DNA/DNA complexes a priori with an unexpectedly high accuracy.

Predicting RNA Duplex Dimerization Free-Energy Changes upon Mutations Using Molecular Dynamics Simulations

The dimerization free energies of RNA-RNA duplexes are fundamental values that represent the structural stability of RNA complexes. We report a comparative analysis of RNA-RNA duplex dimerization free-energy changes upon mutations, estimated from a molecular dynamics simulation and experiments. A linear regression for nine pairs of double-stranded RNA sequences, six base pairs each, yielded a mean absolute deviation of 0.55 kcal/mol and an R(2) value of 0.97, indicating quantitative agreement between simulations and experimental data. The observed accuracy indicates that the molecular dynamics simulation with the current molecular force field is capable of estimating the thermodynamic properties of RNA molecules.

Stacking in RNA: NMR of Four Tetramers Benchmark Molecular Dynamics

Molecular dynamics (MD) simulations for RNA tetramers r(AAAA), r(CAAU), r(GACC), and r(UUUU) are benchmarked against 1H-1H NOESY distances and 3J scalar couplings to test effects of RNA torsion parametrizations. Four different starting structures were used for r(AAAA), r(CAAU), and r(GACC), while five starting structures were used for r(UUUU). On the basis of X-ray structures, criteria are reported for quantifying stacking. The force fields, AMBER ff99, parmbsc0, parm99χ_Yil, ff10, and parmTor, all predict experimentally unobserved stacks and intercalations, e.g., base 1 stacked between bases 3 and 4, and incorrect χ, ϵ, and sugar pucker populations. The intercalated structures are particularly stable, often lasting several microseconds. Parmbsc0, parm99χ_Yil, and ff10 give similar agreement with NMR, but the best agreement is only 46%. Experimentally unobserved intercalations typically are associated with reduced solvent accessible surface area along with amino and hydroxyl hydrogen bonds to phosphate nonbridging oxygens. Results from an extensive set of MD simulations suggest that recent force field parametrizations improve predictions, but further improvements are necessary to provide reasonable agreement with NMR. In particular, intramolecular stacking and hydrogen bonding interactions may not be well balanced with the TIP3P water model. NMR data and the scoring method presented here provide rigorous benchmarks for future changes in force fields and MD methods.

Stacking Free Energies of All DNA and RNA Nucleoside Pairs and Dinucleoside-Monophosphates Computed Using Recently Revised AMBER Parameters and Compared with Experiment

We report the results of a series of 1-μs-long explicit-solvent molecular dynamics (MD) simulations performed to compare the free energies of stacking (ΔGstack) of all possible combinations of DNA and RNA nucleoside (NS) pairs and dinucleoside-monophosphates (DNMPs). For both NS pairs and DNMPs, we show that the computed stacking free energies are in reasonable qualitative agreement with experimental measurements and appear to provide the closest correspondence with experimental data yet found among computational studies; in all cases, however, the computed stacking free energies are too favorable relative to experimental data. Comparisons of NS-pair systems indicate that stacking interactions are very similar in RNA and DNA systems except when a thymine or uracil base is involved: the presence of a thymine base favors stacking by ∼0.3 kcal/mol relative to a uracil base. One exception is found in the self-stacking of cytidines, which are found to be significantly more favorable for the DNA form; an analysis of the rotational orientations sampled during stacking events suggests that this is likely to be due to more favorable sugar-sugar interactions in stacked complexes of deoxycytidines. Comparisons of the DNMP systems indicate that stacking interactions are more favorable in RNA than in DNA except, again, when thymine or uracil bases are involved. Finally, additional simulations performed using a previous generation of the AMBER force field-in which the description of glycosidic bond rotations was less than optimal-produce computed stacking free energies that are in poorer agreement with experimental data. Overall, the simulations provide a comprehensive view of stacking thermodynamics in NS pairs and in DNMPs as predicted by a state-of-the-art MD force field.

Modified Amber Force Field Correctly Models the Conformational Preference for Tandem GA pairs in RNA

Molecular mechanics with all-atom models was used to understand the conformational preference of tandem guanine-adenine (GA) noncanonical pairs in RNA. These tandem GA pairs play important roles in determining stability, flexibility, and structural dynamics of RNA tertiary structures. Previous solution structures showed that these tandem GA pairs adopt either imino (cis Watson-Crick/Watson-Crick A-G) or sheared (trans Hoogsteen/sugar edge A-G) conformations depending on the sequence and orientation of the adjacent closing base pairs. The solution structures (GCGGACGC)2 [Biochemistry, 1996, 35, 9677-9689] and (GCGGAUGC)2 [Biochemistry, 2007, 46, 1511-1522] demonstrate imino and sheared conformations for the two central GA pairs, respectively. These systems were studied using molecular dynamics and free energy change calculations for conformational changes, using umbrella sampling. For the structures to maintain their native conformations during molecular dynamics simulations, a modification to the standard Amber ff10 force field was required, which allowed the amino group of guanine to leave the plane of the base [J. Chem. Theory Comput., 2009, 5, 2088-2100] and form out-of-plane hydrogen bonds with a cross-strand cytosine or uracil. The requirement for this modification suggests the importance of out-of-plane hydrogen bonds in stabilizing the native structures. Free energy change calculations for each sequence demonstrated the correct conformational preference when the force field modification was used, but the extent of the preference is underestimated.