Validating Solution Ensembles from Molecular Dynamics Simulation by Wide-Angle X-ray Scattering Data

Wide-angle x-ray scattering (WAXS) experiments of biomolecules in solution have become increasingly popular because of technical advances in light sources and detectors. However, the structural interpretation of WAXS profiles is problematic, partly because accurate calculations of WAXS profiles from structural models have remained challenging. In this work, we present the calculation of WAXS profiles from explicit-solvent molecular dynamics (MD) simulations of five different proteins. Using only a single fitting parameter that accounts for experimental uncertainties because of the buffer subtraction and dark currents, we find excellent agreement to experimental profiles both at small and wide angles. Because explicit solvation eliminates free parameters associated with the solvation layer or the excluded solvent, which would require fitting to experimental data, we minimize the risk of overfitting. We further find that the influence from water models and protein force fields on calculated profiles are insignificant up to q ≈ 15 nm^?1. Using a series of simulations that allow increasing flexibility of the proteins, we show that incorporating thermal fluctuations into the calculations significantly improves agreement with experimental data, demonstrating the importance of protein dynamics in the interpretation of WAXS profiles. In addition, free MD simulations up to one microsecond suggest that the calculated profiles are highly sensitive with respect to minor conformational rearrangements of proteins, such as an increased flexibility of a loop or an increase of the radius of gyration by  <  1%. The present study suggests that quantitative comparison between MD simulations and experimental WAXS profiles emerges as an accurate tool to validate solution ensembles of biomolecules.     

The addition reactions of germylenoid H<Subscript>2</Subscript>GeAlCl<Subscript>3</Subscript> with ethylene: a theoretical investigation

Theoretical calculations using the M062X and QCISD methods were performed on the addition reactions of the aluminum germylenoid H₂GeAlCl₃ with ethylene. The most two stable structures of germylenoid H₂GeAlCl₃, i.e., the p-complex and three-membered ring structures, respectively, were employed as reactants. The calculated results indicate that, for the p-complex, H₂GeAlCl₃ there are two pathways, I and II, of which path I involves just one transition state, while path II involves two transition states between reactants and products. Comparing the reaction barrier heights of path I (44.6 kJ mol^?1) and II (37.6 kJ mol^?1), the two pathways are competitive, with similar barriers under the same conditions, while for the three-membered ring structure, another two pathways, III and IV, also exist. Path III has one transition state; however, in path IV, two transition states exist. By comparing their barrier heights, path III (barrier height 39.2 kJ mol^?1) could occur more easily than path IV (barrier height 92.8 kJ mol^?1). Considering solvent effects on these addition reactions, the PCM model and CH₂Cl₂ solvent were used in calculations, and the calculated results demonstrate that CH₂Cl₂ solvent is unfavorable for the reactions, except for path II. In CH₂Cl₂ solvent, paths II and III are more favorable than paths I and IV.     

Combined theoretical studies on solvation and hydrogen bond interactions in glycine tripeptide

The theoretical conformational analysis of glycine tripeptide (GT) has been carried out by molecular dynamics (MD) method in order to find minimum energy conformations. The MD studies on GT with water have been carried out for over 10 ns with a time step of 2 fs using fixed charge force field (AMBER ff03). By adding the solvation effect using water as a solvent, the GT conformers identified in this study exhibit α-helical conformation. Compared with the earlier reports, this MD study is able to identify the energetically favourable GT conformations. The obtained geometry of the five most stable GT conformations was optimised using the density functional theory method at B3LYP/6-311G** level of theory. Subsequently, the effects of solvation on the conformational characteristics of five most stable GT conformers with four water molecules (the number of water molecules in the first solvation shell of GT obtained from MD study) were investigated using the same method and the same level of theory. The effect of microsolvation on the fifth GT conformer has been studied with a cluster of 11 water molecules as the first hydration shell which generates folded structure. The interaction energies of all the complexes are calculated by correcting the basis set superposition error. The strong hydrogen bond mainly contributes to the interaction energies. The atoms in molecules theory and natural bond orbital analysis were used to study the origin of H-bonds. A good correlation between the structural parameters and the properties of charge density is found. NMR calculations show that the C = O carbons of the amine groups of the first and middle glycine fragments have maximum chemical shifts.     

Molecular dynamics simulations of void defects in the energetic material HMX

A molecular dynamics (MD) simulation was carried out to characterize the dynamic evolution of void defects in crystalline octahydro-1, 3, 5, 7-tetranitro-1, 3, 5, 7-tetrazocine (HMX). Different models were constructed with the same concentration of vacancies (10 %) to discuss the size effects of void. Energetic ground state properties were determined by annealing simulations. The void formation energy per molecule removed was found to be 55–63 kcal/mol^?1, and the average binding energy per molecule was between 32 and 34 kcal/mol^?1 according to the change in void size. Voids with larger size had lower formation energy. Local binding energies for molecules directly on the void surface decreased greatly compared to those in defect-free lattice, and then gradually increased until the distance away from the void surface was around 10 Å. Analysis of 1 ns MD simulations revealed that the larger the void size, the easier is void collapse. Mean square displacements (MSDs) showed that HMX molecules that had collapsed into void present liquid structure characteristics. Four unique low-energy conformers were found for HMX molecules in void: two whose conformational geometries corresponded closely to those found in HMX polymorphs and two, additional, lower energy conformers that were not seen in the crystalline phases. The ratio of different conformers changed with the simulated temperature, in that the ratio of α conformer increased with the increase in temperature.     

Theoretical study on the structures and properties of mixtures of urea and choline chloride

In this work, we investigated in detail the structural characteristics of mixtures of choline chloride and urea with different urea contents by performing molecular dynamic (MD) simulations, and offer possible explanations for the low melting point of the eutectic mixture of choline chloride and urea with a ratio of 1:2. The insertion of urea molecules was found to change the density distribution of cations and anions around the given cations significantly, disrupting the long-range ordered structure of choline chloride. Moreover, with increasing urea concentration, the hydrogen bond interactions between choline cations and Cl^? anions decreased, while those among urea molecules obviously increased. From the hydrogen bond lifetimes, it was found that a ratio of 1:2 between choline chloride and urea is necessary for a reasonable strength of hydrogen bond interaction to maintain the low melting point of the mixture of choline chloride with urea. In addition, it was also deduced from the interaction energies that a urea content of 67.7 % may make the interactions of cation–anion, cation–urea and anion–urea modest, and thus results in the lower melting point of the eutectic mixture of choline chloride and urea. The present results may offer assistance to some extent for understanding the physicochemical properties of the eutectic mixture of choline chloride and urea, and give valuable information for the further development and application of deep eutectic solvents.     

Calixarene-based artificial ionophores for chloride transport across natural liposomal bilayer: Synthesis,structure-function relationships,and computational study

An amphiphilic calix[6]arene, alone or complexed with an axle to form a pseudo-rotaxane, has been embedded into liposomes prepared from 1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine (POPC) and the permeability of the membrane-doped liposomes towards Cl^? ions has been evaluated by using lucigenin as the fluorescent probe. The pseudo-rotaxane promotes transmembrane transport of Cl^? ions more than calix[6]arene does. Surprisingly, the quenching of lucigenin was very fast for liposomes doped with the positively charged axle alone. Molecular dynamics (MD) simulations and quantum-chemical calculations were also carried out for providing a semi-quantitative support to the experimental results.     

Self-aggregation of trehalose in the mixed solvents of 1,3-dimethylimidazolium ionic liquid and water

Abstract

Conformational variations of solvated trehaloses in binary mixtures of 1,3-dialkylimidazolium ([dmim]Cl) ionic liquids and trehalose as well as ternary mixtures of trehalose, [dmim]Cl and water have been studied by molecular dynamics (MD) simulations with and without polarisable force fields. The interaction energy between anion Cl^? and water is stronger than that between water itself in the [dmim]Cl-water mixtures. Isolated water clusters were found in the binary [dmim]Cl-water mixtures with 60.0 and 75.0% mole fraction of water, but a continuous water network appears when the concentration of the mixture increases to 99.9%. In the case of binary mixtures of trehalose and [dmim]Cl, both non-polarisable and polarisable models demonstrated that the pyranose rings of trehalose displayed chair conformations. MD simulations with polarisation model could sample larger conformation space than that with non-polarisable model. A self-aggregation behaviour of trehalose was found in the ternary trehalose-[dmim]Cl-water mixtures, which can be rationalised by the stronger non-bonded interaction energy between trehalose molecules and anion Cl^? than that between trehalose molecules and water.     

Protonation and electronic structure of 2,6-dichlorophenolindophenolate during reduction. A theoretical study including explicit solvent

Protonation in the two-electron/two-proton reduction processes of 2,6-dichlorophenolindophenolate (DCIP) is investigated combining density functional theory (DFT) and molecular dynamics (MD) methods. DCIP (anion), DCIP^?– (radical anion), and DCIP^2? (dianion) are considered, including the electronic structure analysis from the prospective of quantum theory of atoms and molecules (QTAIM). It is shown that oxygen on the indophenolate moiety and nitrogen are the first and/or the second proton acceptor sites and their energetic order depends on the total charge of the system. MD simulations of differently charged species interacting with the solvent molecules have been performed for methanol, water, and oxonium cation (H₃O⁺). Methanol and water molecules are found to form only hydrogen bonds with the solute irrespective of its charge. The calculated pK_a values show that the imino group of DCIPH^? is a weaker acid than water. While in the case of DCIP (and DCIP^?–) plus oxonium cation, proton transfer from the solvent to the solute was evidenced for both aforementioned acceptor sites. In addition, MD simulations of bulks containing 15 and 43 molecules of water around the DCIP molecule have been performed, revealing the formation of 2–4 hydrogen bonds.

Open image in new window

Graphical Abstract 2,6-Dichlorophenolindophenolate interacts with solvent molecules (water, oxonium cation and methanol). Hydrogen transfer and electronic structure are studied by DFT and molecular dynamics methods

     

Theoretical study of the triangular bonding complex formed by carbon tetrabromide, a halide, and a solvent molecule in the gas phase

MP2(full)/aug-cc-pVDZ(-PP) computations predict that new triangular bonding complexes (where X^? is a halide and H–C refers to a protic solvent molecule) consist of one halogen bond and two hydrogen bonds in the gas phase. Carbon tetrabromide acts as the donor in the halogen bond, while it acts as an acceptor in the hydrogen bond. The halide (which commonly acts as an acceptor) can interact with both carbon tetrabromide and solvent molecule (CH₃CN, CH₂Cl₂, CHCl₃) to form a halogen bond and a hydrogen bond, respectively. The strength of the halogen bond obeys the order CBr₄???Cl^? > CBr₄???Br^? > CBr₄???I^?. For the hydrogen bonds formed between various halides and the same solvent molecule, the strength of the hydrogen bond obeys the order C-H???Cl^? > C-H???Br^? > C-H???I^?. For the hydrogen bonds formed between the same halide and various solvent molecules, the interaction strength is proportional to the acidity of the hydrogen in the solvent molecule. The diminutive effect is present between the hydrogen bonds and the halogen bond in chlorine and bromine triangular bonding complexes. Complexes containing iodide ion show weak cooperative effects.

Vibrational spectra of an RDX film over an aluminum substrate from molecular dynamics simulations and density functional theory

We report calculated vibrational spectra in the range of 0–3,500 cm^?1 of RDX (hexahydro-1,3,5-trinitro-1,3,5-triazine) molecules adsorbed on a model aluminum surface. A molecular film was modeled using two approaches: (1) density functional theory (DFT) was used to optimize a single RDX molecule interacting with its periodic images, and (2) a group of nine molecules extracted from the crystal structure was deposited on the surface and interacted with its periodic images via molecular dynamics (MD) simulations. In both cases, the molecule was initialized in the AAA conformer geometry having the three nitro groups in axial positions, and kept that conformation in the DFT examination, but some molecules were found to change to the AAE conformer (two nitro groups in axial and one in equatorial position) in the MD analysis. The vibrational spectra obtained from both methods are similar to each other, except in the regions where collective RDX intermolecular interactions (captured by MD simulations) are important, and compare fairly well with experimental findings.

Atomic insight into designed carbamate-based derivatives as acetylcholine esterase (AChE) inhibitors: a computational study by multiple molecular docking and molecular dynamics simulation

Over 100 variants have been designed and studied, using multiple docking methods such as Autodock Vina, ArgusLab, Molegro Virtual Docker, and Hex-Cuda, to study the effect of alteration in the structure of carbamate-based acetylcholyne esterase (AChE) inhibitors. Sixteen selected systems were then subjected to 14 ns molecular dynamics (MD) simulations. Results from all the docking methods are in agreement. Variants that involved biphenyl substituents possess the most negative binding energies in the ?37.64 to ?39.31 kJ mol^?1 range due to their π–π interactions with AChE aromatic residues. The root mean square deviation values showed that all of these components achieved equilibration after 6 ns. Gyration radius (R_g) and solvent accessibility surface area were calculated to further investigate the AChE conformational changes in the presence of these components. MD simulation results suggested that these components might interact with AChE, possibly with no major changes in AChE secondary and tertiary structures.     

Building a foundation for structure‐based cellulosome design for cellulosic ethanol: Insight into cohesin‐dockerin complexation from computer simulation

The organization and assembly of the cellulosome, an extracellular multienzyme complex produced by anaerobic bacteria, is mediated by the high‐affinity interaction of cohesin domains from scaffolding proteins with dockerins of cellulosomal enzymes. We have performed molecular dynamics simulations and free energy calculations on both the wild type (WT) and D39N mutant of the C. thermocellum Type I cohesin‐dockerin complex in aqueous solution. The D39N mutation has been experimentally demonstrated to disrupt cohesin‐dockerin binding. The present MD simulations indicate that the substitution triggers significant protein flexibility and causes a major change of the hydrogen‐bonding network in the recognition strips—the conserved loop regions previously proposed to be involved in binding—through electrostatic and salt‐bridge interactions between β‐strands 3 and 5 of the cohesin and α‐helix 3 of the dockerin. The mutation‐induced subtle disturbance in the local hydrogen‐bond network is accompanied by conformational rearrangements of the protein side chains and bound water molecules. Additional free energy perturbation calculations of the D39N mutation provide differences in the cohesin‐dockerin binding energy, thus offering a direct, quantitative comparison with experiments. The underlying molecular mechanism of cohesin‐dockerin complexation is further investigated through the free energy profile, that is, potential of mean force (PMF) calculations of WT cohesin‐dockerin complex. The PMF shows a high‐free energy barrier against the dissociation and reveals a stepwise pattern involving both the central β‐sheet interface and its adjacent solvent‐exposed loop/turn regions clustered at both ends of the β‐barrel structure.     

Chloride transport and compartmentation within main and lateral roots of two grapevine rootstocks differing in salt tolerance

Root Cl^? transport was investigated using ³⁶Cl^? flux analysis in two grapevine (Vitis sp.) rootstock hybrids differing in salt tolerance; 1103 Paulsen (salt-tolerant) and K 51–40 (salt sensitive). Initial ³⁶Cl^? influx to the root was greater in Paulsen than K 51–40. This flux, attributed to the Cl^? influx to the cytoplasm (Φ _oc) increased with increasing external concentrations of Cl^? for plants adapted to growth in 30 mM NaCl. The concentration kinetics in this high concentration range could be fit to a Michaeils–Menton equation. There was no significant difference between genotypes in Km (28.68 ± 15.76 and 24.27 ± 18.51 mM for Paulsen and K 51–40, respectively), but Paulsen had greater V _max (0.127 ± 0.042) compared to K 51–40 (0.059 ± 0.026 μm g^?1 FW min^?1). In Paulsen, the main root had greater contribution to ³⁶Cl^? uptake than lateral roots, there being no significant difference in lateral root influx between the genotypes. ³⁶Cl^? transport to the shoot of K 51–40 was greater than for Paulsen. It was estimated that efflux rate from the xylem parenchyma cells to the xylem vessels (Φ _cx) in K 51–40 was twice that of Paulsen. Compartmental analysis from ³⁶Cl^? efflux kinetics confirmed the larger Φ _oc and the higher ratio of main to lateral root Φ _oc for Paulsen. Efflux from the cytoplasm (Φ _co) was higher than 95 % of Φ _oc indicating a high degree of cycling across the plasma membrane in roots at these high external Cl^? concentrations. Paulsen appears to keep the cytoplasmic Cl^? concentration in roots lower than K 51–40 via greater efflux to the vacuole and to the outside medium. The difference in salt tolerance between the genotypes can be attributed to different Cl^? transport properties at the plasma membrane and tonoplast and particularly in Cl^? efflux to the xylem.     

Theoretical investigations on Azotobacter vinelandii ferredoxin I: effects of electron transfer on protein dynamics

Structural, energetic, and dynamical studies of Azotobacter vinelandii ferredoxin I are presented for native and mutant forms. The protein contains two iron-sulfur clusters, one of which ([3Fe-4S]) is believed to play a central role in the electron-coupled proton transfer. Different charge sets for the [3Fe-4S] cluster in its reduced and oxidized state are calculated with broken symmetry ab initio density functional theory methods and used in molecular dynamics (MD) simulations. The validity of the ab initio calculations is assessed by comparing partially optimized structures of the [3Fe-4S] clusters with x-ray structures. Possible proton transfer pathways between the protein and the iron-sulfur cluster are examined by both MD simulations and ab initio calculations. The MD simulations identify three main-chain hydrogen atoms--HN(13), HN(14), and HN(16)--that are within H-bonding distance of the [3Fe-4S] cluster throughout the MD simulations. They could thus play a role in the proton transfer from the protein to the iron-sulfur cluster. By contrast, the HD2(15) atom of the Asp-15 is seldom close enough to the [3Fe-4S] cluster to transfer a proton. Poisson-Boltzmann calculations indicate that there is a low, but nonzero probability, that Asp-15 is protonated at pH 7; this is a requirement for it to serve as a proton donor. Ab initio calculations with a fragment model for the protein find similar behavior for the transfer of a proton from the OH of the protonated side chain and the main-chain NH of Asp-15. The existence of a stable salt bridge between Asp-15 and Lys-84 in the D15E mutant, versus its absence in the wild-type, has been suggested as the cause of the difference in the rate of proton transfer. Extensive MD simulations were done to test this idea; the results do not support the proposal. The present findings, together with the available data, serve as the basis for an alternative proposal for the mechanism of the coupled electron-proton transfer reaction in ferredoxin I.     

Molecular dynamics and quantum chemical studies of solvent effects on cyclo glycylglycine and glycylalanine dipeptides

Hydrogen bond (H-bond) interactions between the two cyclo dipeptides, cyclo(glycyl-glycine) (CGG) and cyclo(glycyl-alanine) (CGA), and water have been studied using molecular dynamics (MD) and quantum chemical methods. The MD studies have been carried out on CGG and CGA in water using fixed charge force field (AMBER ff03) for over 10 ns with a MD time step of 2 fs. The results of this study show that the solvation pattern influences the conformations of the cyclo dipeptides. Following molecular simulations, post Hartree–Fock and density functional theory methods have been used to explore the molecular properties of the cyclo dipeptides in gaseous and aqueous phase environments. The self-consistent reaction field theory has been used to optimise the cyclopeptides in diethyl ether (? = 4.3) and water (? = 78.5), and the solvent effects have been analysed. A cluster of eight water molecules leads to the formation of first solvation shell of CGG and CGA and the strong H-bonding mainly contributes to the interaction energies. The H-bond interactions have been analysed by the calculation of electron density ρ(r) and its Laplacian ▽²ρ(r) at bond critical points using atoms in molecules theory. The natural bond orbital analysis was carried out to reveal the nature of H-bond interactions. In the solvated complexes, the keto carbons registered the maximum NMR chemical shifts.     

Binding of NDGA and morin with phospholipase A2: experimental and computational evidences

The effects of morin and nordihydroguaiaretic acid (NDGA), two plant secondary metabolites, on porcine pancreatic phospholipase A₂ (PLA₂) were investigated by isothermal titration calorimetry (ITC) and in silico docking analyses. The binding energies obtained for NDGA and morin from the ITC studies are ? 6.36 and ? 5.91 kcal mol^? 1, respectively. Similarly, the glide scores obtained for NDGA and morin towards PLA₂ were ? 7.32 and ? 7.23 kcal mol^? 1, respectively. Further the docked complexes were subjected to MD simulation in the presence of explicit water molecules to check the binding stability of the ligands in the active site of PLA₂. The bound ligands make hydrogen bonds with the active site residues of the enzyme and coordinate bonds with catalytically important Ca²⁺ ion. The binding of ligands at the active site of PLA₂ may also contribute to the reported anti-inflammatory properties of NDGA and morin.     

Dynamic properties of SOD1 mutants can predict survival time of patients carrying familial amyotrophic lateral sclerosis

One of the reasons for the death of motor neurons of the brain and spinal cord in patients with amyotrophic lateral sclerosis is known to be formation of subcellular protein aggregates that are caused by mutations in the SOD1 gene. Patient survival time was earlier shown to have limiting correlation with thermostability change of SOD1 mutant forms of patients’ carriers. We hypothesized that aggregation of mutant SOD1 may occur not only due to the protein destabilization, but through formation of novel interatomic bonds which stabilize “pathogenic” conformations of the mutant as well. To estimate these effects in the present paper, we performed statistical analysis of occupancy of intramolecular hydrogen bonds, hydrogen bonds between the protein and water molecules, and water bridges with use of molecular dynamics simulation for 38 mutant SOD1 forms. Multiple regression model based on these kinds of bonds demonstrated correlation with patient survival time significantly better (R = .9, p-value < 10^?11) than the thermostability of SOD1 mutants only. It was shown that the occupancy of intramolecular hydrogen bonds between amino acid residues is a key determinant (R = .89, p-value < 10^?10) in predicting patients’ survival time.     

The role of CS<Subscript>2</Subscript> in CS<Subscript>2</Subscript>/NMP mixed solvent in weakening the hydrogen bond of OH–N in coal: a DFT investigation

The interaction processes of trace amounts of N-methyl-2-pyrrolidinone (NMP), CS₂/NMP (1:1 by volume) and pure NMP solvent with the hydrogen bond of OH?N in coal were constructed and simulated by density functional theory methods. The distances and bond orders between the main related atoms, and the hydrogen bond energy of OH?N were calculated. The calculated results show that pure NMP solvent does not weaken the hydrogen bond of OH?N in coal. However, trace amounts of NMP and CS₂/NMP (1:1 by volume) have a strong capacity to weaken the hydrogen bond of OH?N in coal. The H2–N3 distances are elongated from 1.87 Å to 3.80 Å and 3.44 Å, the bond orders of H2–N3 all disappear, and the corresponding hydrogen bond energies of OH?N in coal decrease from 45.72 kJ mol^?1 to 7.06 and 11.24 kJ mol^?1, respectively. These results show that CS₂ added to pure NMP solvent plays an important role in releasing the original capacity of NMP to weaken the hydrogen bond of OH?N in coal, in agreement with experimental observations.