Free energy of a potassium ion in a model of the channel formed by an amphipathic leucine-serine peptide

We use molecular dynamics simulations to investigate the position-dependent free energy of a potassium ion in a model of an ion channel formed by the synthetic amphipathic leucine-serine peptide, LS3. The channel model is a parallel bundle of six LS3 helices around which are packed 146 methane-like spheres in order to mimic a membrane. At either end of and within the channel are 1051 water molecules, plus four ions (two potassium and two chloride). The free energy of a potassium ion in the channel was estimated using the weighted histogram analysis (WHAM) method. This is the first time to our knowledge that such a calculation has been carried out as a function of the position of an ion in three dimensions within a channel. The results indicate that for this channel, which is lined by hydrophilic serine sidechains, there is a relatively weak dependence of the free energy on the axial/off-axial position of the ion. There are some off-axis local minima, especially in the C-terminal half of the channel. Using the free energy results, a single channel current-voltage curve was estimated using a one-dimensional Nernst-Planck equation. Although reasonable agreement with experiment is achieved for K(+) ions flowing from the N-terminal to the C-terminal mouth, in the opposite direction the current is underestimated. This underestimation may be a consequence of under-sampling of the conformational dynamics of the channel. We suggest that our simulations may have captured, for example, a sub-conductance level (i.e. an incompletely open state) of the LS3 channel.     

The calculated free energy effects of 5-methyl cytosine on the B to Z transition in DNA

We have examined the free energy effects of 5-methylation of cytosine on the B in equilibrium Z conformational equilibrium in DNA. Free energy differences were calculated using the free energy perturbation approach, which uses an easily derived equation from classical statistical mechanics to relate the free energy difference between two states to the ensemble average of the potential energy difference between the states. Calculations were carried both in explicit solvent and (for comparison) in vacuo. The free energy values obtained for the explicit solvent systems are total free energies, with contributions from all parts of the system (solvent + solute), and so are relevant to the B in equilibrium Z transitions observed under real (physiological) conditions. We calculate that in solution, methylation makes the B in equilibrium Z transition more favorable by about -0.4 kcal/mole base pair (bp) in free energy. This value compares well with approximate experimentally derived values of about -0.3 kcal/mole-bp. We also discuss a method for determining the free energy difference between conformational states poorly maintained by a potential energy model. Finally, the effects of methylation on the melting temperature of DNA are examined.     

Calculation of conformational free energies by confinement simulations in explicit water with implicit desolvation

Abstract

The confinement method is a robust and conceptually simple free energy simulation method that allows the calculation of conformational free energy differences between highly dissimilar states. Application of the method to explicitly solvated systems requires a multi-stage simulation protocol for the calculation of desolvation free energies. Here we show that these desolvation free energies can be readily obtained from an implicit treatment, which is simpler and less costly. The accuracy and robustness of this protocol was shown by the calculation of conformational free energy differences of a series of explicitly solvated test systems. Given the accuracy and ease by which these free energy differences were obtained, the confinement method is promising for the treatment of conformational changes in large and complex systems.     

Theoretical study of large conformational transitions in DNA: the B<-->A conformational change in water and ethanol/water

We explore here the possibility of determining theoretically the free energy change associated with large conformational transitions in DNA, like the solvent-induced BA conformational change. We find that a combination of targeted molecular dynamics (tMD) and the weighted histogram analysis method (WHAM) can be used to trace this transition in both water and ethanol/water mixture. The pathway of the transition in the A→B direction mirrors the B→A pathway, and is dominated by two processes that occur somewhat independently: local changes in sugar puckering and global rearrangements (particularly twist and roll) in the structure. The B→A transition is found to be a quasi-harmonic process, which follows closely the first spontaneous deformation mode of B-DNA, showing that a physiologically-relevant deformation is in coded in the flexibility pattern of DNA.     

Efficient parameter estimation for RNA secondary structure prediction

MOTIVATION: Accurate prediction of RNA secondary structure from the base sequence is an unsolved computational challenge. The accuracy of predictions made by free energy minimization is limited by the quality of the energy parameters in the underlying free energy model. The most widely used model, the Turner99 model, has hundreds of parameters, and so a robust parameter estimation scheme should efficiently handle large data sets with thousands of structures. Moreover, the estimation scheme should also be trained using available experimental free energy data in addition to structural data. RESULTS: In this work, we present constraint generation (CG), the first computational approach to RNA free energy parameter estimation that can be efficiently trained on large sets of structural as well as thermodynamic data. Our CG approach employs a novel iterative scheme, whereby the energy values are first computed as the solution to a constrained optimization problem. Then the newly computed energy parameters are used to update the constraints on the optimization function, so as to better optimize the energy parameters in the next iteration. Using our method on biologically sound data, we obtain revised parameters for the Turner99 energy model. We show that by using our new parameters, we obtain significant improvements in prediction accuracy over current state of-the-art methods. AVAILABILITY: Our CG implementation is available at http://www.rnasoft.ca/CG/.     

A Dictionary Learning Approach with Overlap for the Low Dose Computed Tomography Reconstruction and Its Vectorial Application to Differential Phase Tomography

X-ray based Phase-Contrast Imaging (PCI) techniques have been demonstrated to enhance the visualization of soft tissues in comparison to conventional imaging methods. Nevertheless the delivered dose as reported in the literature of biomedical PCI applications often equals or exceeds the limits prescribed in clinical diagnostics. The optimization of new computed tomography strategies which include the development and implementation of advanced image reconstruction procedures is thus a key aspect. In this scenario, we implemented a dictionary learning method with a new form of convex functional. This functional contains in addition to the usual sparsity inducing and fidelity terms, a new term which forces similarity between overlapping patches in the superimposed regions. The functional depends on two free regularization parameters: a coefficient multiplying the sparsity-inducing norm of the patch basis functions coefficients, and a coefficient multiplying the norm of the differences between patches in the overlapping regions. The solution is found by applying the iterative proximal gradient descent method with FISTA acceleration. The gradient is computed by calculating projection of the solution and its error backprojection at each iterative step. We study the quality of the solution, as a function of the regularization parameters and noise, on synthetic data for which the solution is a-priori known. We apply the method on experimental data in the case of Differential Phase Tomography. For this case we use an original approach which consists in using vectorial patches, each patch having two components: one per each gradient component. The resulting algorithm, implemented in the European Synchrotron Radiation Facility tomography reconstruction code PyHST, has proven to be efficient and well-adapted to strongly reduce the required dose and the number of projections in medical tomography.     

Improved Ligand-Protein Binding Affinity Predictions Using Multiple Binding Modes

Accurate ligand-protein binding affinity prediction, for a set of similar binders, is a major challenge in the lead optimization stage in drug development. In general, docking and scoring functions perform unsatisfactorily in this application. Docking calculations, followed by molecular dynamics simulations and free energy calculations can be applied to improve the predictions. However, for targets with large, flexible binding sites, with no experimentally determined binding modes for a set of ligands, insufficient sampling can decrease the accuracy of the free energy calculations. Cytochrome P450s, a protein family of major importance for drug metabolism, is an example of a challenging target for binding affinity predictions. As a result, the choice of starting structure from the docking solutions becomes crucial. In this study, an iterative scheme is introduced that includes multiple independent molecular dynamics simulations to obtain weighted ensemble averages to be used in the linear interaction energy method. The proposed scheme makes the initial pose selection less crucial for further simulation, as it automatically calculates the relative weights of the various poses. It also properly takes into account the possibility that multiple binding modes contribute similarly to the overall affinity, or of similar compounds occupying very different poses. The method was applied to a set of 12 compounds binding to cytochrome P450 2C9 and it displayed a root mean-square error of 2.9 kJ/mol.     

On the interpretation of the observed linear free energy relationship in phosphate hydrolysis: a thorough computational study of phosphate diester hydrolysis in solution

The hydrolysis of phosphate esters is crucially important to biological systems, being involved in, among other things, signaling, energy transduction, biosynthesis, and the regulation of protein function. Despite this, there are many questions that remain unanswered in this important field, particularly with regard to the preferred mechanism of hydrolysis of phosphate esters, which can proceed through any of multiple pathways that are either associative or dissociative in nature. Previous comparisons of calculated and observed linear free energy relationships (LFERs) for phosphate monoester dianions with different leaving groups showed that the TS character gradually changes from associative to dissociative with the increasing acidity of the leaving group, while reproducing the experimental LFER. Here, we have generated ab initio potential energy surfaces for the hydrolysis of phosphate diesters in solution, with a variety of leaving groups. Once again, the reaction changes from a compact concerted pathway to one that is more expansive in character when the acidity of the leaving group increases. When such systems are examined in solution, it is essential to take into consideration the contribution of solute to the overall activation entropy, which remains a major computational challenge. The popular method of calculating the entropy using a quasi-harmonic approximation appears to markedly overestimate the configurational entropy for systems with multiple occupied energy wells. We introduce an improved restraint release approach for evaluating configurational entropies and apply this approach to our systems. We demonstrate that when this factor is taken into account, it is possible to reproduce the experimental LFER for this system with reasonable accuracy.     

A method to determine dielectric constants in nonhomogeneous systems: application to biological membranes

Continuum electrostatic models have had quantitative success in describing electrostatic-mediated phenomena on atomistic scales; however, there continues to be significant disagreement about how to assign dielectric constants in mixed, nonhomogeneous systems. We introduce a method for determining a position-dependent dielectric profile from molecular dynamics simulations. In this method, the free energy of introducing a test charge is computed two ways: from a free energy perturbation calculation and from a numerical solution to Poisson's Equation. The dielectric profile of the system is then determined by minimizing the discrepancy between these two calculations simultaneously for multiple positions of the test charge. We apply this method to determine the dielectric profile of a lipid bilayer surrounded by water. We find good agreement with dielectric models for lipid bilayers obtained by other approaches. The free energy of transferring an ion from bulk water to the lipid bilayer computed from the atomistic simulations indicates that large errors are introduced when the bilayer is represented as a single slab of low dielectric embedded in the higher-dielectric solvent. Significant improvement results from introducing an additional layer of intermediate dielectric (∼3) on each side of the low dielectric core extending from ∼12 Å to 18 Å. A small dip in transfer free energy just outside the lipid headgroups indicates the presence of a very high dielectric. These results have implications for the design of implicit membrane models and our understanding of protein-membrane interactions.     

An Adiabatic Molecular Dynamics Method for the Calculation of Free Energy Profiles

A new molecular dynamics method for calculating free energy profiles for rare events is presented. The new method is based on the creation of an adiabatic separation between the reaction coordinate subspace and the remaining degrees of freedom within a molecular dynamics run. This is achieved by associating with the reaction coordinate(s) a high temperature and large mass, thereby allowing the activated process to occur while permitting the remaining degrees of freedom to respond adiabatically. In this limit, by applying a formal multiple time scale Liouville operator factorization, it can be rigorously shown that the free energy profiles are obtained directly from the probability distribution of the reaction coordinate subspace and, therefore, require no postprocessing of the output data. The new method is applied to a variety of model problems and its performance tested against free energy calculations using the "bluemoon ensemble" approach. The comparison shows that free energy profiles can be calculated with greater ease and efficiency using the new method.     

Free energy changes associated with amino acid substitution in proteins

The estimation of free energy differences from computer simulation of macromolecular systems is important for rational strategies for drug design and for protein engineering. As an example of one mutation, we have studied the free energy change resulting from the conversion of a polar group (OH) to an apolar group (CH3) in aqueous solution. We have estimated the effect of various local environments on the magnitude of the free energy difference and find that significant environmental effects are found. We have also studied the reliability of the results in detail.     

Calculation of the concentrations of free cations and cation-ligand complexes in solutions containing multiple divalent cations and ligands.

The method described permits the computation of the concentrations of free ions and ion-ligand complexes in a solution containing arbitrary numbers of divalent cations and ligands. It is required that the pH be known, along with appropriate sets of ligand-hydrogen and ligand-divalent cation concentration binding constants. It is assumed that these sets of constants are chosen to be consistent with the ionic strength of the complete solution which contains the divalent cations and ligands. The technique is an iterative one which provides upper and lower bounds for the values of the unknowns. The method does not require initial guesses at the values of the unknowns, and it gives correct answers even when the concentrations involved are many orders of magnitude apart. The present formulation of the problem is restricted to the case where only one cation can bind to a given ligand at any one time. The method is applicable to large molecules with multiple "sub-ligands" provided these sub-ligands are independent in their function as ion-binding sites. These sub-ligands need not all have the same properties. It is also shown that a simple modification of the method permits the determination of the subset of total ion concentrations that are required in order to produce a specified subset of free ion concentrations. The modifications required to include monovalent cation binding are presented in outline form.     

Consensus folding of aligned sequences as a new measure for the detection of functional RNAs by comparative genomics

Facing the ever-growing list of newly discovered classes of functional RNAs, it can be expected that further types of functional RNAs are still hidden in recently completed genomes. The computational identification of such RNA genes is, therefore, of major importance. While most known functional RNAs have characteristic secondary structures, their free energies are generally not statistically significant enough to distinguish RNA genes from the genomic background. Additional information is required. Considering the wide availability of new genomic data of closely related species, comparative studies seem to be the most promising approach. Here, we show that prediction of consensus structures of aligned sequences can be a significant measure to detect functional RNAs. We report a new method to test multiple sequence alignments for the existence of an unusually structured and conserved fold. We show for alignments of six types of well-known functional RNA that an energy score consisting of free energy and a covariation term significantly improves sensitivity compared to single sequence predictions. We further test our method on a number of non-coding RNAs from Caenorhabditis elegans/Caenorhabditis briggsae and seven Saccharomyces species. Most RNAs can be detected with high significance. We provide a Perl implementation that can be used readily to score single alignments and discuss how the methods described here can be extended to allow for efficient genome-wide screens.     

Solvation effects in calculated electrostatic association free energies for the C3d‐CR2 complex and comparison with experimental data

The complement system is an integral part of the innate immune system that participates in the clearance of pathogens from the body. The association between complement protein fragment C3d and B or T cell‐receptor complement receptor (CR) 2 represents a crucial link between innate and adaptive immunities. The goal of this study is to predict association abilities of C3d and CR2 mutants by theoretically calculating electrostatic free energies of association and to assess the importance of solvation effects in the calculations. We demonstrate that calculated solvation free energy differences and Coulombic free energies of association are more sensitive than electrostatic free energies of association in solution and, thus, more accurate in predicting previously published experimental data for the association abilities (relative to the parent proteins) of specific C3d and CR2 mutants. We show that a proportional relationship exists between the predicted solvation free energy differences and the experimental data, while an inversely proportional relationship exists between the predicted Coulombic free energies of association and the experimental data. Our results yield new insights into the physicochemical properties underlying C3d‐CR2 association. We discuss the predictive validity of Coulombic, solvation, and solution electrostatic free energies of association and the generalization of our method for theoretical mutagenesis studies of other systems. This is a basic study, aimed toward improving our understanding of the theoretical basis of immune system regulation at the molecular level. Such insight can serve as the groundwork for the design of regulators with tailored properties, vaccines, and other biotechnology products. © 2010 Wiley Periodicals, Inc. Biopolymers 93:509–519, 2010. This article was originally published online as an accepted preprint. The “Published Online” date corresponds to the preprint version. You can request a copy of the preprint by emailing the Biopolymers editorial office at biopolymers@wiley.com     

Is the computed speciation of copper in a wide range of Chinese soils reliable?

Free Cu species in soils is a key issue to its bioavailability. However, predictive models for Cu speciation across a wide range of soils were still unavailable. In this study, Cu speciation in 34 contaminated soil samples were investigated via analytical technique and predictive models. The results showed that most of free Cu²⁺ was underestimated when using default log K_CuFA and 65% active fulvic acid as inputs in models of WHAM VI and NICA-Donnan. The best prediction was found when using either adjusted active fulvic acid from 10% to 125% for WHAM VI or from 15% to 65% for NICA-Donnan model with the RMSE < 0.32 and r² > 0.96. In contrast, NICA-Donnan demonstrated a slightly stronger binding for Cu than WHAM VI due to extra 26% of samples was underestimated. This work presents a comprehensive database of Cu speciation and an effective attempt of free Cu²⁺ prediction in a wide range of Chinese soils.     

Binding of phosphorus-containing inhibitors to thermolysin studied by the Poisson-Boltzmann method.

Zinc endopeptidase thermolysin can be inhibited by a series of phosphorus-containing peptide analogues, Cbz-Gly-psi (PO2)-X-Leu-Y-R (ZGp(X)L(y)R), where X = NH, O, or CH2; Y = NH or O; R = Leu, Ala, Gly, Phe, H, or CH3. The affinity correlation as well as an X-ray crystallography study suggest that these inhibitors bind to thermolysin in an identical mode. In this work, we calculate the electrostatic binding free energies for a series of 13 phosphorus-containing inhibitors with modifications at X, Y, and R moieties using finite difference solution to the Poisson-Boltzmann equation. A method has been developed to include the solvation entropy changes due to binding different ligands to a macromolecule. We demonstrate that the electrostatic energy and empirically derived solvation entropy can account for most of the binding energy differences in this series. By analyzing the binding contribution from individual residues, we show that the energy of a hydrogen bond is not confined to the donor and acceptor. In particular, the positive charges on Zn and Arg 203, which are not the acceptors, contribute significantly to the hydrogen bonds between two amides of ZGpLL and the thermolysin.     

基于堆积协变信息与最小自由能预测含伪结的RNA二级结构

RNA伪结预测是RNA研究的一个难点问题。文中提出一种基于堆积协变信息与最小自由能的RNA伪结预测方法。该方法使用已知结构的RNA比对序列(ClustalW比对和结构比对)测试此方法, 侧重考虑相邻碱基对之间相互作用形成的堆积协变信息, 并结合最小自由能方法对碱基配对综合评分, 通过逐步迭代求得含伪结的RNA二级结构。结果表明, 此方法能正确预测伪结, 其平均敏感性和特异性优于参考算法, 并且结构比对的预测性能比ClustalW比对的预测性能更加稳定。文中同时讨论了不同协变信息权重因子对预测性能的影响, 发现权重因子比值在l1: l2=5:1时, 预测性能达到最优。     

Free energy distributions in proteins.

Protein molecules in solution have a broad distribution of enthalpy states. A good approximation to the distribution function for enthalpy states can be calculated, using the maximum-entropy method, from the moments of the distribution that, in turn, are obtained from the experimental temperature dependence of the heat capacity. In the present paper, we show that the enthalpy probability distribution can then be formulated in terms of a free energy function that gives the free energy of the protein corresponding to a particular value of the enthalpy. By the location of the minima in this function, the free energy distribution graphically indicates the most probable values of the enthalpy for the protein. We find that the behavior of the free energy functions for proteins falls somewhere between two different cases: a two-state like function with two minima, the relative levels of the two states changing with temperature; and, a single-minimum function where the position of the minimum shifts to higher enthalpy values as the temperature is increased. We show that the temperature dependence of the free energy function can be expressed in terms of a central free energy distribution for a given, fixed temperature (which is most conveniently chosen as the temperature of the maximum in the heat capacity). The nature of this central free energy function for a given protein thus yields all of the thermodynamic behavior of that protein over the temperature range of the denaturation process.     

Solubility of sickle hemoglobin measured by a kinetic micromethod.

We have developed a photolytic method to determine the concentration of reactive hemes in a solution in the presence of a trace amount of CO. By measurement of the bimolecular rate of CO binding, and by calibration of the rate constant under equivalent conditions, the concentration of the reactive hemes can be determined. In a solution of sickle hemoglobin, the molecules in the gel contribute negligibly to the recombination rate, allowing the concentration of the molecules in the solution phase to be determined. To optimize signal to noise, modulated excitation methods were employed, although the method could also be used with pulse techniques and suitable signal averaging. Because the optical method employs a microspectrophotometer, only a few microliters of concentrated Hb solution is required to reproduce the entire temperature dependence of the solubility previously determined by centrifugation using milliliter quantities of solutions of the same concentration. This should be especially useful for studies of site-directed mutants, and we present results obtained on one such HbS in which Leu 88 beta has been replaced by Ala. The free energy difference in the polymerization of the Leu 88 beta double mutant is consistent with known differences in the amino acid hydrophobicities. The calibration required for these experiments also provides an excellent determination of the activation energy for binding the first CO to deoxy Hb.     

A Comparison of Seven Fast but Approximate Methods to Compute the Free Energy of Deprotonation for Amino Acids in Aqueous Solution

Abstract

In recent years, a variety of methods based on statistical mechanics have been successfully applied to calculate free energy differences of chemical reactions from molecular simulation. The accuracy and computational efficiency vary strongly between these methods. Seven approximate but fast methods to calculate free energy differences are compared in terms of accuracy and efficiency with the accurate but expensive thermodynamic integration method as reference, using 28 protonation and deprotonation reactions of aspartic acid in aqueous solution as test cases. At least two simulations are required to obtain an accurate free energy difference between two states of the system. Both, the averaged one-step perturbation method and the linear response method yield the most accurate results, while the latter method shows the fastest convergence.