Computational investigation of the HIV-1 Rev multimerization using molecular dynamics simulations and binding free energy calculations

The HIV Rev protein mediates the nuclear export of viral mRNA, and is thereby essential for the production of late viral proteins in the replication cycle. Rev forms a large organized multimeric protein-protein complex for proper functioning. Recently, the three-dimensional structures of a Rev dimer and tetramer have been resolved and provide the basis for a thorough structural analysis of the binding interaction. Here, molecular dynamics (MD) and binding free energy calculations were performed to elucidate the forces thriving dimerization and higher order multimerization of the Rev protein. It is found that despite the structural differences between each crystal structure, both display a similar behavior according to our calculations. Our analysis based on a molecular mechanics-generalized Born surface area (MM/GBSA) and a configurational entropy approach demonstrates that the higher order multimerization site is much weaker than the dimerization site. In addition, a quantitative hot spot analysis combined with a mutational analysis reveals the most contributing amino acid residues for protein interactions in agreement with experimental results. Additional residues were found in each interface, which are important for the protein interaction. The investigation of the thermodynamics of the Rev multimerization interactions performed here could be a further step in the development of novel antiretrovirals using structure based drug design. Moreover, the variability of the angle between each Rev monomer as measured during the MD simulations suggests a role of the Rev protein in allowing flexibility of the arginine rich domain (ARM) to accommodate RNA binding.     

Molecular dynamics and free energy studies on the Drosophila melanogaster and Leptinotarsa decemlineata ecdysone receptor complexed with agonists: Mechanism for binding and selectivity

The ecdysone receptor is a nuclear hormone receptor that plays a pivotal role in the insect metamorphosis and development. To address the molecular mechanisms of binding and selectivity, the interactions of two typical agonists Ponasterone A and 20-Hydroxyecdysone with Drosophila melanogaster (DME) and Leptinotarsa decemlineata ecdysone (LDE) receptors were investigated by homology modeling, molecular docking, molecular dynamic simulation, and thermodynamic analysis. We discover that 1) the L5-loop, L11-loop, and H12 helix for DME, L7-loop, and L11-loop for LDE are more flexible, which affect the global dynamics of the ligand-binding pocket, thus facilitating the ligand recognition of ecdysone receptor; 2) several key residues (Thr55/Thr37, Phe109/Phe91, Arg95/Arg77, Arg99/Arg81, Phe108/Leu90, and Ala110/Val92) are responsible for the binding of the proteins; 3) the binding-free energy is mainly contributed by the van der Waals forces as well as the electrostatic interactions of ligand and receptor; 4) the computed binding-free energy difference between DME-C1 and LDE-C1 is –4.65 kcal/mol, explains that C1 can form many more interactions with the DME; 5) residues Phe108/Leu90 and Ala110/Val92 have relatively position and orientation difference in the two receptors, accounting most likely for the ligand selectivity of ecdysone receptor from different orders of insects. This study underscores the expectation that different insect pests should be able to discriminate among compounds from different as yet undiscovered compounds, and the results firstly show a structural and functional relay between the agonists and receptors (DME and LDE), which can provide an avenue for the development of target-specific insecticides.

Communicated by Ramaswamy H. Sarma     

Study of the insulin dimerization: binding free energy calculations and per-residue free energy decomposition

A calculation of the binding free energy for the dimerization of insulin has been performed using the molecular mechanics-generalized Born surface area approach. The calculated absolute binding free energy is -11.9 kcal/mol, in approximate agreement with the experimental value of -7.2 kcal/mol. The results show that the dimerization is mainly due to nonpolar interactions. The role of the hydrogen bonds between the 2 monomers appears to give the direction of the interactions. A per-atom decomposition of the binding free energy has been performed to identify the residues contributing most to the self association free energy. Residues B24-B26 are found to make the largest favorable contributions to the dimerization. Other residues situated at the interface between the 2 monomers were found to make favorable but smaller contributions to the dimerization: Tyr B16, Val B12, and Pro B28, and to an even lesser extent, Gly B23. The energy decomposition on a per-residue basis is in agreement with experimental alanine scanning data. The results obtained from a single trajectory (i.e., the dimer trajectory is also used for the monomer analysis) and 2 trajectories (i.e., separate trajectories are used for the monomer and dimer) are similar.     

Absolute binding free energy calculations: On the accuracy of computational scoring of protein–ligand interactions

Calculating the absolute binding free energies is a challenging task. Reliable estimates of binding free energies should provide a guide for rational drug design. It should also provide us with deeper understanding of the correlation between protein structure and its function. Further applications may include identifying novel molecular scaffolds and optimizing lead compounds in computer‐aided drug design. Available options to evaluate the absolute binding free energies range from the rigorous but expensive free energy perturbation to the microscopic linear response approximation (LRA/β version) and related approaches including the linear interaction energy (LIE) to the more approximated and considerably faster scaled protein dipoles Langevin dipoles (PDLD/S‐LRA version) as well as the less rigorous molecular mechanics Poisson–Boltzmann/surface area (MM/PBSA) and generalized born/surface area (MM/GBSA) to the less accurate scoring functions. There is a need for an assessment of the performance of different approaches in terms of computer time and reliability. We present a comparative study of the LRA/β, the LIE, the PDLD/S‐LRA/β, and the more widely used MM/PBSA and assess their abilities to estimate the absolute binding energies. The LRA and LIE methods perform reasonably well but require specialized parameterization for the nonelectrostatic term. The PDLD/S‐LRA/β performs effectively without the need of reparameterization. Our assessment of the MM/PBSA is less optimistic. This approach appears to provide erroneous estimates of the absolute binding energies because of its incorrect entropies and the problematic treatment of electrostatic energies. Overall, the PDLD/S‐LRA/β appears to offer an appealing option for the final stages of massive screening approaches. Proteins 2010. © 2010 Wiley‐Liss, Inc.     

Position-resolved free energy of solvation for amino acids in lipid membranes from molecular dynamics simulations

Studies of insertion and interactions of amino acids in lipid membranes are pivotal to our understanding of membrane protein structure and function. Calculating the insertion cost as a function of transmembrane helix sequence is thus an important step towards improved membrane protein prediction and eventually drug design. Here, we present position-dependent free energies of solvation for all amino acid analogs along the membrane normal. The profiles cover the entire region from bulk water to hydrophobic core, and were produced from all-atom molecular dynamics simulations. Experimental differences corresponding to mutations and costs for entire segments match experimental data well, and in addition the profiles provide the spatial resolution currently not available from experiments. Polar side-chains largely maintain their hydration and assume quite ordered conformations, which indicates the solvation cost is mainly entropic. The cost of solvating charged side-chains is not only significantly lower than for implicit solvation models, but also close to experiments, meaning these could well maintain their protonation states inside the membrane. The single notable exception to the experimental agreement is proline, which is quite expensive to introduce in vivo despite its hydrophobicity--a difference possibly explained by kinks making it harder to insert helices in the translocon.     

Exploring the mechanism how AF9 recognizes and binds H3K9ac by molecular dynamics simulations and free energy calculations

Histone acetylation is a very important regulatory mechanism in gene expression in the chromatin context. A new protein family‐YEATS domains have been found as a novel histone acetylation reader, which could specific recognize the histone lysine acetylation. AF9 is an important one in the YEATS family. Focused on the AF9‐H3K9ac (K9 acetylation) complex (ALY) (PDB code: 4TMP) and a serials of mutants, MUT (the acetyllsine of H3K9ac was mutated to lysine), F59A, G77A, and D103A, we applied molecular dynamics simulation and molecular mechanics Poisson?Boltzmann (MM‐PBSA) free energy calculations to examine the role of AF9 protein in recognition interaction. The simulation results and analysis indicate that some residues of the protein have significant influence on recognition and binding to H3K9ac peptides and hydrophobic surface show the hydrophobic interactions play an important role in the binding. Our work can give important information to understand how the protein AF9 recognizes the peptides H3K9ac. © 2016 Wiley Periodicals, Inc. Biopolymers 105: 779–786, 2016.     

Ferric uptake regulator protein: binding free energy calculations and per-residue free energy decomposition

Iron homeostasis is, in many bacterial species, mediated by the ferric uptake regulator (Fur). A regulatory site able to bind iron to activate Fur for DNA binding has been described, and a structural zinc site essential for the dimerization has also been proposed. They have been localized and named site 1 and site 2, respectively, from the crystal structure of a zinc-substituted Pseudomonas aeruginosa Fur (PA-Fur). Notwithstanding the studies on Fur proteins from various species, both the precise site of iron binding and the effect on DNA binding affinity are still controversial. These issues were investigated here by molecular dynamics simulations and free energy calculations. Simulations were performed for eight molecular systems represented by the three forms of Fur, that is, apo Fur, metal-substituted Fur, and Fur complexed with DNA. Because of the lack of a Fur-DNA complex crystal structure, the recently published model based on mass spectrometry experiments on Escherichia coli Fur (EC-Fur), and the crystal structure of PA-Fur, was used, after adjustment to adopt a symmetric conformation. The simulation results suggest that the formerly proposed site 2 is, in fact, the regulatory iron-sensing site. The calculations also predict that Fe(2+) at site 2 is hexacoordinated having an octahedral environment with only nitrogen and oxygen atoms, which is in accordance with previous spectroscopic characterizations. Energy decomposition pinpoints H87 as an additional amino acid that defines the regulatory metal site. Finally, free energy decomposition analysis reveals a number of amino acids potentially important in dimerization and in DNA binding.     

Docking and free energy simulations to predict conformational domains involved in hCG-LH receptor interactions using recombinant antibodies

Single chain fragment variables (ScFvs) have been extensively employed in studying the protein-protein interactions. ScFvs derived from phage display libraries have an additional advantage of being generated against a native antigen, circumventing loss of information on conformational epitopes. In the present study, an attempt has been made to elucidate human chorionic gonadotropin (hCG)-luteinizing hormone (LH) receptor interactions by using a neutral and two inhibitory ScFvs against hCG. The objective was to dock a computationally derived model of these ScFvs onto the crystal structure of hCG and understand the differential roles of the mapped epitopes in hCG-LH receptor interactions. An anti-hCG ScFv, whose epitope was mapped previously using biochemical tools, served as the positive control for assessing the quality of docking analysis. To evaluate the role of specific side chains at the hCG-ScFv interface, binding free energy as well as residue interaction energies of complexes in solution were calculated using molecular mechanics Poisson-Boltzmann/surface area method after performing the molecular dynamic simulations on the selected hCG-ScFv models and validated using biochemical and SPR analysis. The robustness of these calculations was demonstrated by comparing the theoretically determined binding energies with the experimentally obtained kinetic parameters for hCG-ScFv complexes. Superimposition of hCG-ScFv model onto a model of hCG complexed with the 51-266 residues of LH receptor revealed importance of the residues previously thought to be unimportant for hormone binding and response. This analysis provides an alternate tool for understanding the structure-function analysis of ligand-receptor interactions.     

Helix-capping interaction in λ cro protein: A free energy simulation analysis

The stability mutant Tyr-26 → Asp was studied in the Cro protein from bacteriophage λ using free energy molecular dynamics simulations. The mutant was calculated to be more stable than the wild type by 3.0 ± 1.7 kcal/mol/monomer, in reasonable agreement with experiment (1.4 kcal/mol/monomer). Moreover, the aspartic acid in the mutant was found to form a capping interation with the amino terminus of the third α-helix of Cro. The simulations were analyzed to understand better the source of the stability of this helix-capping interaction and to examine the results in light of previous explanations of stabilizing helix caps-namely, a model of local unsatisfied hydrogen bonds at the helix termini and the helix macro dipole model. Analysis of the simulations shows that the stabilizing effect of this charged helical cap is due both to favorable hydrogen bonds with backbone NH groups at the helix terminus and to favorable electrostatic interactions (but not hydrogen bonds) with their carbonyls (effectively the next row of local dipoles in the helix). However, electrostatic interactions are weak or negligible with backbone dipolar groups in the helix further away from the terminus. Moreover, the importance of other local electrostatic interactions with polar side chains near the helix terminus, which are neglected in most treatments of this effect, are shown to be important. Thus, the results support a model that is intermediate between the two previous explanations: both unsatisfied hydrogen bonds at the helix terminus and other, local preoriented dipolar groups stabilize the helix cap. These findings suggest that similar interactions with preoriented dipolar groups may be important for cooperativity in other charge–dipole interactions and may be employed to advantage for molecular design. © 1994 Wiley-Liss, Inc.     

The effect of structural parameters and positive charge distance on the interaction free energy of antimicrobial peptides with membrane surface

Many attempts have been made to find hints explaining the relationship between physicochemical and structural properties of antimicrobial peptides (AMPs) which are relevant to their antimicrobial activities. We here found that there is a difference in the percentages of hydrophobic, hydrophilic, and charged residues between AMPs killing both bacteria and fungi (Group A) and AMPs that only kill bacteria (Group B). The percentage of charged residues in Group A AMPs is highly elevated, while in Group B the percentage of hydrophobic residues is increased. This result suggests a sequence-based mechanism of selectivity for AMPs. Moreover, we examined how the distance between basic residues affects the interaction free energy of AMPs with the membrane surface, since most of the known AMPs act by membrane perturbation. We measured the average distance between basic residues throughout the 3D structure of AMPs by defining D_pr parameter and calculated the interaction free energy for 10 AMPs that interacted with the DPPC membrane using molecular dynamics simulation. We found that the changes of the interaction free energy correlates with the change of D_pr by a linear regression coefficient of r^2?=?.47 and a cubic regression coefficient of r^2?=?.70.     

Energetic evaluation of binding modes in the C3d and Factor H (CCP 19-20) complex

As a part of innate immunity, the complement system relies on activation of the alternative pathway (AP). While feed-forward amplification generates an immune response towards foreign surfaces, the process requires regulation to prevent an immune response on the surface of host cells. Factor H (FH) is a complement protein secreted by native cells to negatively regulate the AP. In terms of structure, FH is composed of 20 complement-control protein (CCP) modules that are structurally homologous but vary in composition and function. Mutations in these CCPs have been linked to states of autoimmunity. In particular, several mutations in CCP 19-20 are correlated to atypical hemolytic uremic syndrome (aHUS). From crystallographic structures there are three putative binding sites of CCP 19-20 on C3d. Since there has been some controversy over the primary mode of binding from experimental studies, we approach characterization of binding using computational methods. Specifically, we compare each binding mode in terms of electrostatic character, structural stability, dissociative and associative properties, and predicted free energy of binding. After a detailed investigation, we found two of the three binding sites to be similarly stable while varying in the number of contacts to C3d and in the energetic barrier to complex dissociation. These sites are likely physiologically relevant and may facilitate multivalent binding of FH CCP 19-20 to C3b and either C3d or host glycosaminoglycans. We propose thermodynamically stable binding with modules 19 and 20, the latter driven by electrostatics, acting synergistically to increase the apparent affinity of FH for host surfaces.     

Affinity and specificity of protein U1A-RNA complex formation based on an additive component free energy model

An MM-GBSA computational protocol was used to investigate wild-type U1A-RNA and F56 U1A mutant experimental binding free energies. The trend in mutant binding free energies compared to wild-type is well-reproduced. Following application of a linear-response-like equation to scale the various energy components, the binding free energies agree quantitatively with observed experimental values. Conformational adaptation contributes to the binding free energy for both the protein and the RNA in these systems. Small differences in DeltaGs are the result of different and sometimes quite large relative contributions from various energetic components. Residual free energy decomposition indicates differences not only at the site of mutation, but throughout the entire protein. MM-GBSA and ab initio calculations performed on model systems suggest that stacking interactions may nearly, but not completely, account for observed differences in mutant binding affinities. This study indicates that there may be different underlying causes of ostensibly similar experimentally observed binding affinities of different mutants, and thus recommends caution in the interpretation of binding affinities and specificities purely by inspection.     

Protein-protein recognition: an experimental and computational study of the R89K mutation in Raf and its effect on Ras binding

Binding of the protein Raf to the active form of Ras promotes activation of the MAP kinase signaling pathway, triggering cell growth and differentiation. Raf/Arg89 in the center of the binding interface plays an important role determining Ras-Raf binding affinity. We have investigated experimentally and computationally the Raf-R89K mutation, which abolishes signaling in vivo. The binding to [gamma-35S]GTP-Ras of a fusion protein between the Raf-binding domain (RBD) of Raf and GST was reduced at least 175-fold by the mutation, corresponding to a standard binding free energy decrease of at least 3.0 kcal/mol. To compute this free energy and obtain insights into the microscopic interactions favoring binding, we performed alchemical simulations of the RBD, both complexed to Ras and free in solution, in which residue 89 is gradually mutated from Arg into Lys. The simulations give a standard binding free energy decrease of 2.9+/-1.9 kcal/mol, in agreement with experiment. The use of numerous runs with three different force fields allows insights into the sources of uncertainty in the free energy and its components. The binding decreases partly because of a 7 kcal/mol higher cost to desolvate Lys upon binding, compared to Arg, due to better solvent interactions with the more concentrated Lys charge in the unbound state. This effect is expected to be general, contributing to the lower propensity of Lys to participate in protein-protein interfaces. Large contributions to the free energy change also arise from electrostatic interactions with groups up to 8 A away, namely residues 37-41 in the conserved effector domain of Ras (including 4 kcal/mol from Ser39 which loses a bifurcated hydrogen bond to Arg89), the conserved Lys84 and Lys87 of Raf, and 2-3 specific water molecules. This analysis will provide insights into the large experimental database of Ras-Raf mutations.     

Engineering Strontium Binding Affinity in an EF-hand Motif: A Quantum Chemical and Molecular Dynamics Study

Abstract

Proteins with the ability to specifically bind strontium would potentially be of great use in the field of nuclear waste management. Unfortunately, no such peptides or proteins are known—indeed, it is uncertain whether they exist under natural conditions due to low environmental concentrations of strontium. To investigate the possibility of devising such molecules, one of us (CV), in a previous experimental study [J. Biol. Inorg. Chem. 8, 33440 (2003)], proposed starting from an EF-hand motif of the protein calmodulin and mutating some residues to change the motif's specificity for calcium into one for strontium. In this paper, which represents a theoretical complement to the experimental work, we analyzed small-molecule crystallographic structures and performed quantum chemical calculations to identify possible mutations. We then constructed seven mutant sequences of the EF-hand motif and analyzed their dynamical and binding behaviors using molecular dynamics simulations and free-energy calculations (using the MM/PBSA method). As a result of these analyzes we were able to isolate some characteristics that could lead to mutant peptides with enhanced strontium affinity.     

Zinc finger protein binding to DNA: an energy perspective using molecular dynamics simulation and free energy calculations on mutants of both zinc finger domains and their specific DNA bases

Energy calculations based on MM-GBSA were employed to study various zinc finger protein (ZF) motifs binding to DNA. Mutants of both the DNA bound to their specific amino acids were studied. Calculated energies gave evidence for a relationship between binding energy and affinity of ZF motifs to their sites on DNA. ΔG values were ?15.82(12), ?3.66(12), and ?12.14(11.6) kcal/mol for finger one, finger two, and finger three, respectively. The mutations in the DNA bases reduced the value of the negative energies of binding (maximum value for ΔΔG = 42Kcal/mol for F1 when GCG mutated to GGG, and ΔΔG = 22 kcal/mol for F2, the loss in total energy of binding originated in the loss in electrostatic energies upon mutation (r = .98). The mutations in key amino acids in the ZF motif in positions-1, 2, 3, and 6 showed reduced binding energies to DNA with correlation coefficients between total free energy and electrostatic was .99 and with Van der Waal was .93. Results agree with experimentally found selectivity which showed that Arginine in position-1 is specific to G, while Aspartic acid (D) in position 2 plays a complicated role in binding. There is a correlation between the MD calculated free energies of binding and those obtained experimentally for prepared ZF motifs bound to triplet bases in other reports (), our results may help in the design of ZF motifs based on the established recognition codes based on energies and contributing energies to the total energy.     

Molecular dynamics simulations and free energy analyses on the dimer formation of an amyloidogenic heptapeptide from human beta2-microglobulin: implication for the protofibril structure

Amyloid formation is associated with many neurodegenerative diseases. Recent findings suggest that early oligomeric aggregates could be major sources of toxicity. We present a computational investigation of the first step of amyloid initiation-dimer formation of a seven residue peptide (NHVTLSQ) from human beta2-microglobulin at pH 2.0, which renders +2.0 units charges to each peptide. A total of over 1.2 micros of simulations with explicit solvent and 1.0 micros of simulations with implicit solvent were conducted. Main-chain conformational restraint was applied to facilitate the formation of ordered dimers. An antiparallel beta-sheet with six main-chain hydrogen bonds was dominant in the implicit solvent simulations. In contrast, no stable dimers were observed in the two negative controls, the mouse heptapeptide (KHDSMAE, +3.0 units charges) and the scrambled human heptapeptide (QVLHTSN). Explicit solvent simulations presented a more complex scenario. The wild-type human heptapeptide formed predominantly antiparallel beta-sheets ( approximately 38%) although parallel ones ( approximately 12%) were also observed. Hydrophobic contacts preceded hydrogen bond saturation in the majority of the association events in the explicit solvent simulations, highlighting the important role of hydrophobic interaction in amyloid initiation. The fact that the mouse dimer dissociated immediately after the removal of conformational restraint suggests that the higher conformational entropy barrier, along with the stronger charge repulsion and weaker hydrophobic interaction, contributed to its inability to form amyloid fibril. The closeness of positive charge pairs in the dimers of the scrambled human heptapeptide may prohibit further beta-sheet extension and fibril growth. Combining the results from simulations and free energy analyses, we propose that the building block for this amyloid fibril is an antiparallel dimer with a two-residue register shift and six main-chain hydrogen bonds. A double-layer protofibril structure is also proposed in which two antiparallel beta-sheets face each other and are held together by hydrophobic staples and hydrogen bonds of the polar side-chains.     

Conservative and nonconservative mutations in proteins: anomalous mutations in a transport receptor analyzed by free energy and quantum chemical calculations.

Experimental studies on a bacterial sulfate receptor have indicated anomalous relative binding affinities for the mutations Ser130-->Cys,Ser130-->Gly, and Ser130-->Ala. The loss of affinity for sulfate in the former mutation was previously attributed to a greater steric effect on the part of the Cys side chain relative to the Ser side chain, whereas the relatively small loss of binding affinity for the latter two mutations was attributed to the loss of a single hydrogen bond. In this report we present quantum chemical and statistical thermodynamic studies of these mutations. Qualitative results from these studies indicate that for the Ser130-->Cys mutation the large decrease in binding affinity is in part caused by steric effects, but also significantly by the differential work required to polarize the Cys thiol group relative to the Ser hydroxyl group. The Gly mutant cobinds a water molecule in the same location as the Ser side chain resulting in a relatively small decrease in binding affinity. Results for the Ala mutant are in disagreement with experimental results but are likely to be limited by insufficient sampling of configuration space due to physical constraints applied during the simulation.     

Role of the electrostatic energy between the photo-products and ionic groups on the photoinduced electron transfer rates from aromatic amino acids to the excited flavin in five single-point substitution isoforms of the charged amino acid residue-13 in the FMN-binding protein

Effect of the charge (negative, positive or neutral) of amino acid residue-13 on the photoinduced electron transfer (ET) from Trp32, Tyr35 and Trp106 to the excited isoalloxazine was evaluated in the flavin mononucleotide-binding protein from Desulfovibrio vulgaris isolate Miyazaki F (DvFBP). The protein structures of the wild type and the four isoforms where glutamic acid-13 is replaced with lysine (E13K), arginine (E13R), threonine (E13T) and glutamine (E13Q) in aqueous solution were obtained by molecular dynamics simulation. The distances between the amino acid residue-13 and isoalloxazine (Iso), and between the amino acid residue-13 and the ET donors were longer than 1 nm. The ET rates were evaluated with the Kakitani and Mataga model (KM theory) from their ultrafast fluorescence dynamics by means of a non-linear least squares method. Electrostatic (ES) energies between the photo-products and other ionic groups in the proteins markedly varied among ET donors and among the DvFBP isoforms, while the other physical quantities related to the ET rates, the solvent reorganisation and ES energies between the Iso anion and the donor cations did not vary much between the proteins and donors. A plot of the logarithmic ET rates versus either the total free energy gaps or the net ES energies between the photo-products and the other ionic groups both displayed a parabolic function and so the net ES energies are an important influential factor upon the ET rate, in addition to the donor–acceptor distance.