 |
|
 |
 | We have recently changed our website. Looking for the previous site? |
Abstracts of Articles in 2015
To request any journal reprints, please contact us.- A Molecular Dynamics Investigation of Mycobacterium tuberculosis Prenyl Synthases: Conformational Flexibility, and Implications for Computer-Aided Drug Discovery
- Molecular dynamics simulations of lipid bilayers: structure, dynamics and interactions with antimicrobial peptides
- 18F-RB390: Innovative Ligand for Imaging the T877A Androgen Receptor Mutant in Prostate Cancer via Positron Emission Tomography (PET)
- In silico Screening for Plasmodium falciparum Enoyl-ACP-Reductase Inhibitors
- Computer-aided drug discovery approach finds calcium sensitizer of the cardiac troponin complex
- Effects of HCM cTnI Mutation R145G on Troponin Structure and Modulation by PKA Phosphorylation elucidated by Molecular Dynamics Simulations
- Poisson-Boltzmann vs. Size-modified Poisson-Boltzmann Electrostatics Applied to Lipid Bilayers
- Membranes Serve as Allosteric Activators of Phospholipase A2 Enabling it to Extract, Bind, and Hydrolyze Phospholipid Substrate
- Identification of Protein-Ligand Binding Sites by the Level-Set Variational Implicit Solvent Approach
- Improved cryoEM-Guided Iterative Molecular Dynamics–Rosetta Protein Structure Refinement Protocol for High Precision Protein Structure Prediction
- LS-VISM: A Software Package for Analysis of Biomolecular Solvation
- Allosteric Effects of Sodium Ion Binding on Activation of the M3 Muscarinic G-Protein Coupled Receptor
- Investigation of the conformational dynamics of the apo A2A adenosine receptor
- Accelerated Molecular Dynamics Simulations of Ligand Binding to a Muscarinic G-protein Coupled Receptor
- How to Deal with Multiple Binding Poses in Alchemical Relative Protein-Ligand Binding Free Energy Calculations
- Enhanced Ligand Sampling for Relative Protein-Ligand Binding Free Energy Calculations
- Conformational Dynamics and Binding Free Energies of Inhibitors of BACE-1: From the Perspective of Protonation Equilibria.
- Accelerated Molecular Dynamics Simulations of Protein Folding.
- Broken-Symmetry DFT Computations for the Reaction Pathway of IspH, an Iron−Sulfur Enzyme in Pathogenic Bacteria
- Gaussian Accelerated Molecular Dynamics: Unconstrained Enhanced Sampling and Free Energy Calculation
- Enzyme localization, crowding and buffers collectively modulate diffusion influenced signal transduction: insights from continuum diffusion modeling
- Electrostatic steering enhances the rate of cAMP binding to phosphodiesterase: Brownian dynamics modeling
- A Self-Consistent Phase-Field Approach to Implicit Solvation of Charged Molecules with Poisson-Boltzmann Electrostatics
- Troponin I Mutations R146G and R21C Alter Cardiac Troponin Function, Contractile Properties and Modulation by PKA-mediated Phosphorylation
- Anti-infectives targeting enzymes, and the proton motive force
A Molecular Dynamics Investigation of Mycobacterium tuberculosis Prenyl Synthases: Conformational Flexibility, and Implications for Computer-Aided Drug DiscoveryChem. Biol. Drug Disc. 85, 756-769 (2015)
With the rise in antibiotic resistance, there is interest in discovering new drugs active against new targets. Here, we investigate the dynamic structures of three isoprenoid synthases from Mycobacterium tuberculosis using molecular dynamics (MD) methods with a view to discovering new drug leads. Two of the enzymes, cis-farnesyl diphosphate synthase (cis-FPPS) and cis-decaprenyl diphosphate synthase (cis-DPPS), are involved in bacterial cell wall biosynthesis, while the third, tuberculosinyl adenosine synthase (Rv3378c), is involved in virulence factor formation. The MD results for these three enzymes were then compared with previous results on undecaprenyl diphosphate synthase (UPPS) by means of active site volume fluctuation and principal component analyses. In addition, an analysis of the binding of prenyl diphosphates to cis-FPPS, cis-DPPS, and UPPS utilizing the new MD results is reported. We also screened libraries of inhibitors against cis-DPPS, finding ~1 ¿m inhibitors, and used the receiver operating characteristic-area under the curve (ROC-AUC) method to test the predictive power of X-ray and MD-derived cis-DPPS receptors. We found that one compound with potent M. tuberculosis cell growth inhibition activity was an IC50 ~0.5- to 20-¿m inhibitor (depending on substrate) of cis-DPPS, a ~660-nm inhibitor of Rv3378c as well as a 4.8-¿m inhibitor of cis-FPPS, opening up the possibility of multitarget inhibition involving both cell wall biosynthesis and virulence factor formation.
With the rise in antibiotic resistance, there is interest in discovering new drugs active against new targets. Here, we investigate the dynamic structures of three isoprenoid synthases from Mycobacterium tuberculosis using molecular dynamics (MD) methods with a view to discovering new drug leads. Two of the enzymes, cis-farnesyl diphosphate synthase (cis-FPPS) and cis-decaprenyl diphosphate synthase (cis-DPPS), are involved in bacterial cell wall biosynthesis, while the third, tuberculosinyl adenosine synthase (Rv3378c), is involved in virulence factor formation. The MD results for these three enzymes were then compared with previous results on undecaprenyl diphosphate synthase (UPPS) by means of active site volume fluctuation and principal component analyses. In addition, an analysis of the binding of prenyl diphosphates to cis-FPPS, cis-DPPS, and UPPS utilizing the new MD results is reported. We also screened libraries of inhibitors against cis-DPPS, finding ~1 ¿m inhibitors, and used the receiver operating characteristic-area under the curve (ROC-AUC) method to test the predictive power of X-ray and MD-derived cis-DPPS receptors. We found that one compound with potent M. tuberculosis cell growth inhibition activity was an IC50 ~0.5- to 20-¿m inhibitor (depending on substrate) of cis-DPPS, a ~660-nm inhibitor of Rv3378c as well as a 4.8-¿m inhibitor of cis-FPPS, opening up the possibility of multitarget inhibition involving both cell wall biosynthesis and virulence factor formation.
Molecular dynamics simulations of lipid bilayers: structure, dynamics and interactions with antimicrobial peptidesIn "Molecular Modeling at the Atomic Scale: Methods and Applications in Quantitative Biology" (R. Zhou, Ed.), CRC Press, Taylor & Francis, pp. 191-213 (2015)
18F-RB390: Innovative Ligand for Imaging the T877A Androgen Receptor Mutant in Prostate Cancer via Positron Emission Tomography (PET)The Prostate 75, 348–359 (2015)
BACKGROUND: Detecting prostate cancer before spreading or predicting a favorable therapy are challenging issues for impacting patient's survival. Presently, 2-[(18) F]-fluoro-2-deoxy-D-glucose ((18) F-FDG) and/or (18) F-fluorocholine ((18) F-FCH) are the generally used PET-tracers in oncology yet do not emphasize the T877A androgen receptor (AR) mutation being exclusively present in cancerous tissue and escaping androgen deprivation treatment.
BACKGROUND: Detecting prostate cancer before spreading or predicting a favorable therapy are challenging issues for impacting patient's survival. Presently, 2-[(18) F]-fluoro-2-deoxy-D-glucose ((18) F-FDG) and/or (18) F-fluorocholine ((18) F-FCH) are the generally used PET-tracers in oncology yet do not emphasize the T877A androgen receptor (AR) mutation being exclusively present in cancerous tissue and escaping androgen deprivation treatment.
METHODS: We designed and synthesized fluorinated 5α-dihydrotestosterone (DHT) derivatives to target T877A-AR. We performed binding assays to select suitable candidates using COS-7 cells transfected with wild-type or T877A AR (WT-AR, T877A-AR) expressing plasmids and investigated cellular uptake of candidate (18) F-RB390. Stability, biodistribution analyses and PET-Imaging were assessed by injecting (18) F-RB390 (10MBq), with and without co-injection of an excess of unlabeled DHT in C4-2 and PC-3 tumor bearing male SCID mice (n = 12).
RESULTS: RB390 presented a higher relative binding affinity (RBA) (28.1%, IC50 = 32 nM) for T877A-AR than for WT-AR (1.7%, IC50 = 357 nM) related to DHT (RBA = 100%). A small fraction of (18) F-RB390 was metabolized when incubated with murine liver homogenate or human blood for 3 hr. The metabolite of RB390, 3-hydroxysteroid RB448, presented similar binding characteristics as RB390. (18) F-RB390 but not (18) F-FDG or (18) F-FCH accumulated 2.5× more in COS-7 cells transfected with pSG5AR-T877A than with control plasmid. Accumulation was reduced with an excess of DHT. PET/CT imaging and biodistribution studies revealed a significantly higher uptake of (18) F-RB390 in T877A mutation positive xenografts compared to PC-3 control tumors. This effect was blunted with DHT.
CONCLUSION: Given the differential binding capacity and the favorable radioactivity pattern, (18) F-RB390 represents the portrayal of the first imaging ligand with predictive potential for mutant T877A-AR in prostate cancer for guiding therapy. Prostate 75:348-359, 2015. © 2014 Wiley Periodicals, Inc.
In silico Screening for Plasmodium falciparum Enoyl-ACP-Reductase InhibitorsJ Comp Aid Mol Des 29, 79–87 (2015)
The need for novel therapeutics against Plasmodium falciparum is urgent due to recent emergence of multi-drug resistant malaria parasites. Since fatty acids are essential for both the liver and blood stages of the malarial parasite, targeting fatty acid biosynthesis is a promising strategy for combatting P. falciparum. We present a combined computational and experimental study to identify novel inhibitors of enoyl-acyl carrier protein reductase (PfENR) in the fatty acid biosynthesis pathway. A smallmolecule database from ChemBridge was docked into three distinct PfENR crystal structures that provide multiple receptor conformations. Two different docking algorithms were used to generate a consensus score in order to rank possible small molecule hits. Our studies led to the identification of five low-micromolar pyrimidine dione inhibitors of PfENR.
The need for novel therapeutics against Plasmodium falciparum is urgent due to recent emergence of multi-drug resistant malaria parasites. Since fatty acids are essential for both the liver and blood stages of the malarial parasite, targeting fatty acid biosynthesis is a promising strategy for combatting P. falciparum. We present a combined computational and experimental study to identify novel inhibitors of enoyl-acyl carrier protein reductase (PfENR) in the fatty acid biosynthesis pathway. A smallmolecule database from ChemBridge was docked into three distinct PfENR crystal structures that provide multiple receptor conformations. Two different docking algorithms were used to generate a consensus score in order to rank possible small molecule hits. Our studies led to the identification of five low-micromolar pyrimidine dione inhibitors of PfENR.
Computer-aided drug discovery approach finds calcium sensitizer of the cardiac troponin complexChem Biol Drug Des 85, 99–106 (2015)
In the fight against heart failure, therapeutics that have the ability to increase the contractile power of the heart are urgently needed. One possible route of action to improve heart contractile power is increasing the calcium sensitivity of the thin filament. From a pharmaceutical standpoint, calcium sensitizers have the distinct advantage of not altering cardiomyocyte calcium levels and thus have lower potential for side-effects. Small chemical molecules have been shown to bind to the interface between cTnC and the cTnI switch peptide and exhibit calcium-sensitizing properties, possibly by stabilizing cTnC in an open conformation. Building on existing structural data of a known calcium sensitizer bound to cardiac troponin, we combined computational structure-based virtual screening drug discovery methods and solution NMR titration assays to identify a novel calcium sensitizer 4-(4-(2,5-dimethylphenyl)-1-piperazinyl)-3-pyridinamine (NSC147866) which binds to cTnC and the cTnC-cTnI147-163 complex. Its presence increases the affinity of switch peptide to cTnC by approximately a factor of two. This action is comparable to that of known levosimendan analogues.
In the fight against heart failure, therapeutics that have the ability to increase the contractile power of the heart are urgently needed. One possible route of action to improve heart contractile power is increasing the calcium sensitivity of the thin filament. From a pharmaceutical standpoint, calcium sensitizers have the distinct advantage of not altering cardiomyocyte calcium levels and thus have lower potential for side-effects. Small chemical molecules have been shown to bind to the interface between cTnC and the cTnI switch peptide and exhibit calcium-sensitizing properties, possibly by stabilizing cTnC in an open conformation. Building on existing structural data of a known calcium sensitizer bound to cardiac troponin, we combined computational structure-based virtual screening drug discovery methods and solution NMR titration assays to identify a novel calcium sensitizer 4-(4-(2,5-dimethylphenyl)-1-piperazinyl)-3-pyridinamine (NSC147866) which binds to cTnC and the cTnC-cTnI147-163 complex. Its presence increases the affinity of switch peptide to cTnC by approximately a factor of two. This action is comparable to that of known levosimendan analogues.
Effects of HCM cTnI Mutation R145G on Troponin Structure and Modulation by PKA Phosphorylation elucidated by Molecular Dynamics SimulationsBiophysical Journal, Volume 108, Issue 2, 20 January 2015, Pages 395–407
Cardiac troponin (cTn) is a key molecule in the regulation of human cardiac muscle contraction. The N-terminal cardiac-specific peptide of the inhibitory subunit of troponin, cTnI (cTnI1-39), is a target for phosphorylation by protein kinase A (PKA) during β-adrenergic stimulation. We recently presented evidence indicating that this peptide interacts with the inhibitory peptide (cTnl137-147) when S23 and S24 are phosphorylated. The inhibitory peptide is also the target of the point mutation cTnI-R145G, which is associated with hypertrophic cardiomyopathy (HCM), a disease associated with sudden death in apparently healthy young adults. It has been shown that both phosphorylation and this mutation alter the cTnC-cTnI (C-I) interaction, which plays a crucial role in modulating contractile activation. However, little is known about the molecular-level events underlying this modulation. Here, we computationally investigated the effects of the cTnI-R145G mutation on the dynamics of cTn, cTnC Ca(2+) handling, and the C-I interaction. Comparisons were made with the cTnI-R145G/S23D/S24D phosphomimic mutation, which has been used both experimentally and computationally to study the cTnI N-terminal specific effects of PKA phosphorylation. Additional comparisons between the phosphomimic mutations and the real phosphorylations were made. For this purpose, we ran triplicate 150 ns molecular dynamics simulations of cTnI-R145G Ca(2+)-bound cTnC1-161-cTnI1-172-cTnT236-285, cTnI-R145G/S23D/S24D Ca(2+)-bound cTnC1-161-cTnI1-172-cTnT236-285, and cTnI-R145G/PS23/PS24 Ca(2+)-bound cTnC1-161-cTnI1-172-cTnT236-285, respectively. We found that the cTnI-R145G mutation did not impact the overall dynamics of cTn, but stabilized crucial Ca(2+)-coordinating interactions. However, the phosphomimic mutations increased overall cTn fluctuations and destabilized Ca(2+) coordination. Interestingly, cTnI-R145G blunted the intrasubunit interactions between the cTnI N-terminal extension and the cTnI inhibitory peptide, which have been suggested to play a crucial role in modulating troponin function during β-adrenergic stimulation. These findings offer a molecular-level explanation for how the HCM mutation cTnI-R145G reduces the modulation of cTn by phosphorylation of S23/S24 during β-adrenergic stimulation.
Cardiac troponin (cTn) is a key molecule in the regulation of human cardiac muscle contraction. The N-terminal cardiac-specific peptide of the inhibitory subunit of troponin, cTnI (cTnI1-39), is a target for phosphorylation by protein kinase A (PKA) during β-adrenergic stimulation. We recently presented evidence indicating that this peptide interacts with the inhibitory peptide (cTnl137-147) when S23 and S24 are phosphorylated. The inhibitory peptide is also the target of the point mutation cTnI-R145G, which is associated with hypertrophic cardiomyopathy (HCM), a disease associated with sudden death in apparently healthy young adults. It has been shown that both phosphorylation and this mutation alter the cTnC-cTnI (C-I) interaction, which plays a crucial role in modulating contractile activation. However, little is known about the molecular-level events underlying this modulation. Here, we computationally investigated the effects of the cTnI-R145G mutation on the dynamics of cTn, cTnC Ca(2+) handling, and the C-I interaction. Comparisons were made with the cTnI-R145G/S23D/S24D phosphomimic mutation, which has been used both experimentally and computationally to study the cTnI N-terminal specific effects of PKA phosphorylation. Additional comparisons between the phosphomimic mutations and the real phosphorylations were made. For this purpose, we ran triplicate 150 ns molecular dynamics simulations of cTnI-R145G Ca(2+)-bound cTnC1-161-cTnI1-172-cTnT236-285, cTnI-R145G/S23D/S24D Ca(2+)-bound cTnC1-161-cTnI1-172-cTnT236-285, and cTnI-R145G/PS23/PS24 Ca(2+)-bound cTnC1-161-cTnI1-172-cTnT236-285, respectively. We found that the cTnI-R145G mutation did not impact the overall dynamics of cTn, but stabilized crucial Ca(2+)-coordinating interactions. However, the phosphomimic mutations increased overall cTn fluctuations and destabilized Ca(2+) coordination. Interestingly, cTnI-R145G blunted the intrasubunit interactions between the cTnI N-terminal extension and the cTnI inhibitory peptide, which have been suggested to play a crucial role in modulating troponin function during β-adrenergic stimulation. These findings offer a molecular-level explanation for how the HCM mutation cTnI-R145G reduces the modulation of cTn by phosphorylation of S23/S24 during β-adrenergic stimulation.
Poisson-Boltzmann vs. Size-modified Poisson-Boltzmann Electrostatics Applied to Lipid BilayersJ. Phys. Chem. B, 2014, 118 (51), pp 14827–14832
Mean-field methods, such as the Poisson–Boltzmann equation (PBE), are often used to calculate the electrostatic properties of molecular systems. In the past two decades, an enhancement of the PBE, the size-modified Poisson–Boltzmann equation (SMPBE), has been reported. Here, the PBE and the SMPBE are reevaluated for realistic molecular systems, namely, lipid bilayers, under eight different sets of input parameters. The SMPBE appears to reproduce the molecular dynamics simulation results better than the PBE only under specific parameter sets, but in general, it performs no better than the Stern layer correction of the PBE. These results emphasize the need for careful discussions of the accuracy of mean-field calculations on realistic systems with respect to the choice of parameters and call for reconsideration of the cost-efficiency and the significance of the current SMPBE formulation.
Mean-field methods, such as the Poisson–Boltzmann equation (PBE), are often used to calculate the electrostatic properties of molecular systems. In the past two decades, an enhancement of the PBE, the size-modified Poisson–Boltzmann equation (SMPBE), has been reported. Here, the PBE and the SMPBE are reevaluated for realistic molecular systems, namely, lipid bilayers, under eight different sets of input parameters. The SMPBE appears to reproduce the molecular dynamics simulation results better than the PBE only under specific parameter sets, but in general, it performs no better than the Stern layer correction of the PBE. These results emphasize the need for careful discussions of the accuracy of mean-field calculations on realistic systems with respect to the choice of parameters and call for reconsideration of the cost-efficiency and the significance of the current SMPBE formulation.
Membranes Serve as Allosteric Activators of Phospholipase A2 Enabling it to Extract, Bind, and Hydrolyze Phospholipid SubstrateProc. Natl. Acad. Sci. USA in press (2015)
Defining the molecular details and consequences of the association of water-soluble proteins with membranes is fundamental to understanding protein–lipid interactions and membrane functioning. Phospholipase A2 (PLA2) enzymes, which catalyze the hydrolysis of phospholipid substrates that compose the membrane bilayers, provide the ideal system for studying protein–lipid interactions. Our study focuses on understanding the catalytic cycle of two different human PLA2s: the cytosolic Group IVA cPLA2 and calcium-independent Group VIA iPLA2. Computer-aided techniques guided by deuterium exchange mass spectrometry data, were used to create structural complexes of each enzyme with a single phospholipid substrate molecule, whereas the substrate extraction process was studied using steered molecular dynamics simulations. Molecular dynamic simulations of the enzyme–substrate–membrane systems revealed important information about the mechanisms by which these enzymes associate with the membrane and then extract and bind their phospholipid substrate. Our data support the hypothesis that the membrane acts as an allosteric ligand that binds at the allosteric site of the enzyme’s interfacial surface, shifting its conformation from a closed (inactive) state in water to an open (active) state at the membrane interface.
Defining the molecular details and consequences of the association of water-soluble proteins with membranes is fundamental to understanding protein–lipid interactions and membrane functioning. Phospholipase A2 (PLA2) enzymes, which catalyze the hydrolysis of phospholipid substrates that compose the membrane bilayers, provide the ideal system for studying protein–lipid interactions. Our study focuses on understanding the catalytic cycle of two different human PLA2s: the cytosolic Group IVA cPLA2 and calcium-independent Group VIA iPLA2. Computer-aided techniques guided by deuterium exchange mass spectrometry data, were used to create structural complexes of each enzyme with a single phospholipid substrate molecule, whereas the substrate extraction process was studied using steered molecular dynamics simulations. Molecular dynamic simulations of the enzyme–substrate–membrane systems revealed important information about the mechanisms by which these enzymes associate with the membrane and then extract and bind their phospholipid substrate. Our data support the hypothesis that the membrane acts as an allosteric ligand that binds at the allosteric site of the enzyme’s interfacial surface, shifting its conformation from a closed (inactive) state in water to an open (active) state at the membrane interface.
Identification of Protein-Ligand Binding Sites by the Level-Set Variational Implicit Solvent ApproachJ. Chem. Theory Comp. pp 753–765), 11 (2), (2015)
Protein–ligand binding is a key biological process at the molecular level. The identification and characterization of small-molecule binding sites on therapeutically relevant proteins have tremendous implications for target evaluation and rational drug design. In this work, we used the recently developed level-set variational implicit-solvent model (VISM) with the Coulomb field approximation (CFA) to locate and characterize potential protein–small-molecule binding sites. We applied our method to a data set of 515 protein–ligand complexes and found that 96.9% of the cocrystallized ligands bind to the VISM-CFA-identified pockets and that 71.8% of the identified pockets are occupied by cocrystallized ligands. For 228 tight-binding protein–ligand complexes (i.e, complexes with experimental pKd values larger than 6), 99.1% of the cocrystallized ligands are in the VISM-CFA-identified pockets. In addition, it was found that the ligand binding orientations are consistent with the hydrophilic and hydrophobic descriptions provided by VISM. Quantitative characterization of binding pockets with topological and physicochemical parameters was used to assess the “ligandability” of the pockets. The results illustrate the key interactions between ligands and receptors and can be very informative for rational drug design.
Protein–ligand binding is a key biological process at the molecular level. The identification and characterization of small-molecule binding sites on therapeutically relevant proteins have tremendous implications for target evaluation and rational drug design. In this work, we used the recently developed level-set variational implicit-solvent model (VISM) with the Coulomb field approximation (CFA) to locate and characterize potential protein–small-molecule binding sites. We applied our method to a data set of 515 protein–ligand complexes and found that 96.9% of the cocrystallized ligands bind to the VISM-CFA-identified pockets and that 71.8% of the identified pockets are occupied by cocrystallized ligands. For 228 tight-binding protein–ligand complexes (i.e, complexes with experimental pKd values larger than 6), 99.1% of the cocrystallized ligands are in the VISM-CFA-identified pockets. In addition, it was found that the ligand binding orientations are consistent with the hydrophilic and hydrophobic descriptions provided by VISM. Quantitative characterization of binding pockets with topological and physicochemical parameters was used to assess the “ligandability” of the pockets. The results illustrate the key interactions between ligands and receptors and can be very informative for rational drug design.
Improved cryoEM-Guided Iterative Molecular Dynamics–Rosetta Protein Structure Refinement Protocol for High Precision Protein Structure PredictionJ. Chem. Theory Comput., 11 (3), pp 1337–1346, (2015)
Many excellent methods exist that incorporate cryo-electron microscopy (cryoEM) data to constrain computational protein structure prediction and refinement. Previously, it was shown that iteration of two such orthogonal sampling and scoring methods – Rosetta and molecular dynamics (MD) simulations – facilitated exploration of conformational space in principle. Here, we go beyond a proof-of-concept study and address significant remaining limitations of the iterative MD – Rosetta protein structure refinement protocol. Specifically, all parts of the iterative refinement protocol are now guided by medium-resolution cryoEM density maps and previous knowledge about the native structure of the protein is no longer necessary. Models are identified solely based on score or simulation time. All four benchmark proteins showed substantial improvement through three rounds of the iterative refinement protocol. The best-scoring final models of two proteins had sub-Angstrom RMSD to the native structure over residues in secondary structure elements. Molecular dynamics was most efficient in refining secondary structure elements, and was thus highly complementary to the Rosetta refinement which is most powerful in refining side chains and loop regions.
Many excellent methods exist that incorporate cryo-electron microscopy (cryoEM) data to constrain computational protein structure prediction and refinement. Previously, it was shown that iteration of two such orthogonal sampling and scoring methods – Rosetta and molecular dynamics (MD) simulations – facilitated exploration of conformational space in principle. Here, we go beyond a proof-of-concept study and address significant remaining limitations of the iterative MD – Rosetta protein structure refinement protocol. Specifically, all parts of the iterative refinement protocol are now guided by medium-resolution cryoEM density maps and previous knowledge about the native structure of the protein is no longer necessary. Models are identified solely based on score or simulation time. All four benchmark proteins showed substantial improvement through three rounds of the iterative refinement protocol. The best-scoring final models of two proteins had sub-Angstrom RMSD to the native structure over residues in secondary structure elements. Molecular dynamics was most efficient in refining secondary structure elements, and was thus highly complementary to the Rosetta refinement which is most powerful in refining side chains and loop regions.
LS-VISM: A Software Package for Analysis of Biomolecular SolvationJ. Comp. Chem.,Volume 36, (14) pp 1047-1059, (2015)
We introduce a software package for the analysis of biomolecular solvation. The package collects computer codes that implement numerical methods for a variational implicit-solvent model (VISM). The input of the package includes the atomic data of biomolecules under consideration and the macroscopic parameters such as solute-solvent surface tension, bulk solvent density and ionic concentrations, and the dielectric coefficients. The output includes estimated solvation free energies and optimal macroscopic solute-solvent interfaces that are obtained by minimizing the VISM solvation free-energy functional among all possible solute-solvent interfaces enclosing the solute atoms. We review the VISM with various descriptions of electrostatics. We also review our numerical methods that consist mainly of the level-set method for relaxing the VISM free-energy functional and a compact coupling interface method for the dielectric Poisson-Boltzmann equation. Such numerical methods and algorithms constitute the central modules of the software package. We detail the structure of the package, format of input and output files, workflow of the codes, and the postprocessing of output data. Our demo application to a host-guest system illustrates how to use the package to perform solvation analysis for biomolecules, including ligand-receptor binding systems. The package is simple and flexible with respect to minimum adjustable parameters and a wide range of applications. Future extensions of the package use can include the efficient identification of ligand binding pockets on protein surfaces.
We introduce a software package for the analysis of biomolecular solvation. The package collects computer codes that implement numerical methods for a variational implicit-solvent model (VISM). The input of the package includes the atomic data of biomolecules under consideration and the macroscopic parameters such as solute-solvent surface tension, bulk solvent density and ionic concentrations, and the dielectric coefficients. The output includes estimated solvation free energies and optimal macroscopic solute-solvent interfaces that are obtained by minimizing the VISM solvation free-energy functional among all possible solute-solvent interfaces enclosing the solute atoms. We review the VISM with various descriptions of electrostatics. We also review our numerical methods that consist mainly of the level-set method for relaxing the VISM free-energy functional and a compact coupling interface method for the dielectric Poisson-Boltzmann equation. Such numerical methods and algorithms constitute the central modules of the software package. We detail the structure of the package, format of input and output files, workflow of the codes, and the postprocessing of output data. Our demo application to a host-guest system illustrates how to use the package to perform solvation analysis for biomolecules, including ligand-receptor binding systems. The package is simple and flexible with respect to minimum adjustable parameters and a wide range of applications. Future extensions of the package use can include the efficient identification of ligand binding pockets on protein surfaces.
Allosteric Effects of Sodium Ion Binding on Activation of the M3 Muscarinic G-Protein Coupled ReceptorBiophys. J.,Volume 8, (7) pp 1796-1806, (2015)
G-protein coupled receptors (GPCRs) are important membrane proteins that mediate cellular signaling and represent primary targets for about one third of currently marketed drugs. Recent X-ray crystallographic studies identified distinct conformations of GPCRs in the active and inactive states. An allosteric sodium ion was found bound to a highly conserved D2.50 residue in inactive GPCRs, while the D2.50 allosteric pocket became collapsed in active GPCR structures. However, the dynamic mechanisms of these observations remain elusive. In this study, we aim to understand the mechanistic effects of sodium ion binding on dynamic activation of the M3 muscarinic GPCR through long-timescale accelerated molecular dynamics (aMD) simulations. Results showed that with the D2.50 residue deprotonated, the M3 receptor is bound by an allosteric sodium ion and confined mostly in the inactive state with remarkably reduced flexibility. In contrast, the D2.50-protonated receptor does not exhibit sodium ion binding to the D2.50 allosteric site and samples significantly larger conformational space. The receptor activation is captured and characterized by large-scale structural rearrangements of the transmembrane helices via dynamic hydrogen bond and salt bridge interactions. The residue motions are highly correlated during receptor activation. Further network analysis revealed that the allosteric signaling between residue D2.50 and key residues in the intracellular, extracellular and orthosteric pockets is significantly weakened upon sodium ion binding.
G-protein coupled receptors (GPCRs) are important membrane proteins that mediate cellular signaling and represent primary targets for about one third of currently marketed drugs. Recent X-ray crystallographic studies identified distinct conformations of GPCRs in the active and inactive states. An allosteric sodium ion was found bound to a highly conserved D2.50 residue in inactive GPCRs, while the D2.50 allosteric pocket became collapsed in active GPCR structures. However, the dynamic mechanisms of these observations remain elusive. In this study, we aim to understand the mechanistic effects of sodium ion binding on dynamic activation of the M3 muscarinic GPCR through long-timescale accelerated molecular dynamics (aMD) simulations. Results showed that with the D2.50 residue deprotonated, the M3 receptor is bound by an allosteric sodium ion and confined mostly in the inactive state with remarkably reduced flexibility. In contrast, the D2.50-protonated receptor does not exhibit sodium ion binding to the D2.50 allosteric site and samples significantly larger conformational space. The receptor activation is captured and characterized by large-scale structural rearrangements of the transmembrane helices via dynamic hydrogen bond and salt bridge interactions. The residue motions are highly correlated during receptor activation. Further network analysis revealed that the allosteric signaling between residue D2.50 and key residues in the intracellular, extracellular and orthosteric pockets is significantly weakened upon sodium ion binding.
Investigation of the conformational dynamics of the apo A2A adenosine receptorProtein Sci. 24 (6) 1004-1012 (2015)
The activation/deactivation processes for G-protein coupled receptors (GPCRs) have been computationally studied for several different classes, including rhodopsin, the β2 adrenergic receptor, and the M2 muscarinic receptor. Despite determined co-crystal structures of the adenosine A2A receptor (A2A AR) in complex with antagonists, agonists and an antibody, the deactivation process of this GPCR is not completely understood. In this study, we investigate the convergence of two apo simulations, one starting with an agonist-bound conformation (PDB: 3QAK) and the other starting with an antagonist-bound conformation (PDB: 3EML). Despite the two simulations not completely converging, we were able to identify distinct intermediate steps of the deactivation process characterized by the movement of Y2887.53 in the NPxxY motif. We find that Y2887.53 contributes to the process by forming hydrogen bonds to residues in transmembrane helices 2 and 7 and losing these interactions upon full deactivation. Y1975.58 also plays a role in the process by forming a hydrogen bond only once the side chain moves from the lipid interface to the middle of the helical bundle.
The activation/deactivation processes for G-protein coupled receptors (GPCRs) have been computationally studied for several different classes, including rhodopsin, the β2 adrenergic receptor, and the M2 muscarinic receptor. Despite determined co-crystal structures of the adenosine A2A receptor (A2A AR) in complex with antagonists, agonists and an antibody, the deactivation process of this GPCR is not completely understood. In this study, we investigate the convergence of two apo simulations, one starting with an agonist-bound conformation (PDB: 3QAK) and the other starting with an antagonist-bound conformation (PDB: 3EML). Despite the two simulations not completely converging, we were able to identify distinct intermediate steps of the deactivation process characterized by the movement of Y2887.53 in the NPxxY motif. We find that Y2887.53 contributes to the process by forming hydrogen bonds to residues in transmembrane helices 2 and 7 and losing these interactions upon full deactivation. Y1975.58 also plays a role in the process by forming a hydrogen bond only once the side chain moves from the lipid interface to the middle of the helical bundle.
Accelerated Molecular Dynamics Simulations of Ligand Binding to a Muscarinic G-protein Coupled ReceptorQ Rev Biophysics, 48, (4) 479-487 (2015)
Elucidating the detailed process of ligand binding to a receptor is pharmaceutically important for identifying druggable binding sites. With the ability to provide atomistic detail, computational methods are well poised to study these processes. Here, accelerated molecular dynamics (aMD) is proposed to simulate processes of ligand binding to a G-protein-coupled receptor (GPCR), in this case the M3 muscarinic receptor, which is a target for treating many human diseases, including cancer, diabetes and obesity. Long-timescale aMD simulations were performed to observe the binding of three chemically diverse ligand molecules: antagonist tiotropium (TTP), partial agonist arecoline (ARc) and full agonist acetylcholine (ACh). In comparison with earlier microsecond-timescale conventional MD simulations, aMD greatly accelerated the binding of ACh to the receptor orthosteric ligand-binding site and the binding of TTP to an extracellular vestibule. Further aMD simulations also captured binding of ARc to the receptor orthosteric site. Additionally, all three ligands were observed to bind in the extracellular vestibule during their binding pathways, suggesting that it is a metastable binding site. This study demonstrates the applicability of aMD to protein–ligand binding, especially the drug recognition of GPCRs.
Elucidating the detailed process of ligand binding to a receptor is pharmaceutically important for identifying druggable binding sites. With the ability to provide atomistic detail, computational methods are well poised to study these processes. Here, accelerated molecular dynamics (aMD) is proposed to simulate processes of ligand binding to a G-protein-coupled receptor (GPCR), in this case the M3 muscarinic receptor, which is a target for treating many human diseases, including cancer, diabetes and obesity. Long-timescale aMD simulations were performed to observe the binding of three chemically diverse ligand molecules: antagonist tiotropium (TTP), partial agonist arecoline (ARc) and full agonist acetylcholine (ACh). In comparison with earlier microsecond-timescale conventional MD simulations, aMD greatly accelerated the binding of ACh to the receptor orthosteric ligand-binding site and the binding of TTP to an extracellular vestibule. Further aMD simulations also captured binding of ARc to the receptor orthosteric site. Additionally, all three ligands were observed to bind in the extracellular vestibule during their binding pathways, suggesting that it is a metastable binding site. This study demonstrates the applicability of aMD to protein–ligand binding, especially the drug recognition of GPCRs.
How to Deal with Multiple Binding Poses in Alchemical Relative Protein-Ligand Binding Free Energy CalculationsJ. Chem. Theory Comput., 11, (6) pp 2670-2679 (2015)
Recent advances in improved force fields and sampling methods have made it possible for the accurate calculation of protein–ligand binding free energies. Alchemical free energy perturbation (FEP) using an explicit solvent model is one of the most rigorous methods to calculate relative binding free energies. However, for cases where there are high energy barriers separating the relevant conformations that are important for ligand binding, the calculated free energy may depend on the initial conformation used in the simulation due to the lack of complete sampling of all the important regions in phase space. This is particularly true for ligands with multiple possible binding modes separated by high energy barriers, making it difficult to sample all relevant binding modes even with modern enhanced sampling methods. In this paper, we apply a previously developed method that provides a corrected binding free energy for ligands with multiple binding modes by combining the free energy results from multiple alchemical FEP calculations starting from all enumerated poses, and the results are compared with Glide docking and MM-GBSA calculations. From these calculations, the dominant ligand binding mode can also be predicted. We apply this method to a series of ligands that bind to c-Jun N-terminal kinase-1 (JNK1) and obtain improved free energy results. The dominant ligand binding modes predicted by this method agree with the available crystallography, while both Glide docking and MM-GBSA calculations incorrectly predict the binding modes for some ligands. The method also helps separate the force field error from the ligand sampling error, such that deviations in the predicted binding free energy from the experimental values likely indicate possible inaccuracies in the force field. An error in the force field for a subset of the ligands studied was identified using this method, and improved free energy results were obtained by correcting the partial charges assigned to the ligands. This improved the root-mean-square error (RMSE) for the predicted binding free energy from 1.9 kcal/mol with the original partial charges to 1.3 kcal/mol with the corrected partial charges.
Recent advances in improved force fields and sampling methods have made it possible for the accurate calculation of protein–ligand binding free energies. Alchemical free energy perturbation (FEP) using an explicit solvent model is one of the most rigorous methods to calculate relative binding free energies. However, for cases where there are high energy barriers separating the relevant conformations that are important for ligand binding, the calculated free energy may depend on the initial conformation used in the simulation due to the lack of complete sampling of all the important regions in phase space. This is particularly true for ligands with multiple possible binding modes separated by high energy barriers, making it difficult to sample all relevant binding modes even with modern enhanced sampling methods. In this paper, we apply a previously developed method that provides a corrected binding free energy for ligands with multiple binding modes by combining the free energy results from multiple alchemical FEP calculations starting from all enumerated poses, and the results are compared with Glide docking and MM-GBSA calculations. From these calculations, the dominant ligand binding mode can also be predicted. We apply this method to a series of ligands that bind to c-Jun N-terminal kinase-1 (JNK1) and obtain improved free energy results. The dominant ligand binding modes predicted by this method agree with the available crystallography, while both Glide docking and MM-GBSA calculations incorrectly predict the binding modes for some ligands. The method also helps separate the force field error from the ligand sampling error, such that deviations in the predicted binding free energy from the experimental values likely indicate possible inaccuracies in the force field. An error in the force field for a subset of the ligands studied was identified using this method, and improved free energy results were obtained by correcting the partial charges assigned to the ligands. This improved the root-mean-square error (RMSE) for the predicted binding free energy from 1.9 kcal/mol with the original partial charges to 1.3 kcal/mol with the corrected partial charges.
Enhanced Ligand Sampling for Relative Protein-Ligand Binding Free Energy CalculationsJ. Phys. Chem. B, 119(20), pp 6190-6197 (2015)
Free energy calculations are used to study how strongly potential drug molecules interact with their target receptors. The accuracy of these calculations depends on the accuracy of the molecular dynamics (MD) force field as well as proper sampling of the major conformations of each molecule. However, proper sampling of ligand conformations can be difficult when there are large barriers separating the major ligand conformations. An example of this is for ligands with an asymmetrically substituted phenyl ring, where the presence of protein loops hinders the proper sampling of the different ring conformations. These ring conformations become more difficult to sample when the size of the functional groups attached to the ring increases. The Adaptive Integration Method (AIM) has been developed, which adaptively changes the alchemical coupling parameter λ during the MD simulation so that conformations sampled at one λ can aid sampling at the other λ values. The Accelerated Adaptive Integration Method (AcclAIM) builds on AIM by lowering potential barriers for specific degrees of freedom at intermediate λ values. However, these methods may not work when there are very large barriers separating the major ligand conformations. In this work, we describe a modification to AIM that improves sampling of the different ring conformations, even when there is a very large barrier between them. This method combines AIM with conformational Monte Carlo sampling, giving improved convergence of ring populations and the resulting free energy. This method, called AIM/MC, is applied to study the relative binding free energy for a pair of ligands that bind to thrombin and a different pair of ligands that bind to aspartyl protease β-APP cleaving enzyme 1 (BACE1). These protein-ligand binding free energy calculations illustrate the improvements in conformational sampling and the convergence of the free energy compared to both AIM and AcclAIM.
Free energy calculations are used to study how strongly potential drug molecules interact with their target receptors. The accuracy of these calculations depends on the accuracy of the molecular dynamics (MD) force field as well as proper sampling of the major conformations of each molecule. However, proper sampling of ligand conformations can be difficult when there are large barriers separating the major ligand conformations. An example of this is for ligands with an asymmetrically substituted phenyl ring, where the presence of protein loops hinders the proper sampling of the different ring conformations. These ring conformations become more difficult to sample when the size of the functional groups attached to the ring increases. The Adaptive Integration Method (AIM) has been developed, which adaptively changes the alchemical coupling parameter λ during the MD simulation so that conformations sampled at one λ can aid sampling at the other λ values. The Accelerated Adaptive Integration Method (AcclAIM) builds on AIM by lowering potential barriers for specific degrees of freedom at intermediate λ values. However, these methods may not work when there are very large barriers separating the major ligand conformations. In this work, we describe a modification to AIM that improves sampling of the different ring conformations, even when there is a very large barrier between them. This method combines AIM with conformational Monte Carlo sampling, giving improved convergence of ring populations and the resulting free energy. This method, called AIM/MC, is applied to study the relative binding free energy for a pair of ligands that bind to thrombin and a different pair of ligands that bind to aspartyl protease β-APP cleaving enzyme 1 (BACE1). These protein-ligand binding free energy calculations illustrate the improvements in conformational sampling and the convergence of the free energy compared to both AIM and AcclAIM.
Conformational Dynamics and Binding Free Energies of Inhibitors of BACE-1: From the Perspective of Protonation Equilibria.PLoS Comp. Biol. in press (2015).
BACE-1 is the b-secretase responsible for the initial amyloidogenesis in Alzheimer’s disease, catalyzing hydrolytic cleavage of substrate in a pH-sensitive manner. The catalytic mechanism of BACE-1 requires water-mediated proton transfer from aspartyl dyad to the substrate, as well as structural flexibility in the flap region. Thus, the coupling of protonation and conformational equilibria is essential to a full in silico characterization of BACE-1. In this work, we perform constant pH replica exchange molecular dynamics simulations on both apo BACE-1 and five BACE-1-inhibitor complexes to examine the effect of pH on dynamics and inhibitor binding properties of BACE-1. In our simulations, we find that solution pH controls the conformational flexibility of apo BACE-1, whereas bound inhibitors largely limit the motions of the holo enzyme at all levels of pH. The microscopic pKa values of titratable residues in BACE-1 including its aspartyl dyad are computed and compared between apo and inhibitor-bound states. Changes in protonation between the apo and holo forms suggest a thermodynamic linkage between binding of inhibitors and protons localized at the dyad. Utilizing our recently developed computational protocol applying the binding polynomial formalism to the constant pH molecular dynamics (CpHMD) framework, we are able to obtain the pH-dependent binding free energy profiles for various BACE-1-inhibitor complexes. Our results highlight the importance of correctly addressing the binding-induced protonation changes in protein-ligand systems where binding accompanies a net proton transfer. This work comprises the first application of our CpHMD-based free energy computational method to protein-ligand complexes and illustrates the value of CpHMD as an all-purpose tool for obtaining pH-dependent dynamics and binding free energies of biological systems.
BACE-1 is the b-secretase responsible for the initial amyloidogenesis in Alzheimer’s disease, catalyzing hydrolytic cleavage of substrate in a pH-sensitive manner. The catalytic mechanism of BACE-1 requires water-mediated proton transfer from aspartyl dyad to the substrate, as well as structural flexibility in the flap region. Thus, the coupling of protonation and conformational equilibria is essential to a full in silico characterization of BACE-1. In this work, we perform constant pH replica exchange molecular dynamics simulations on both apo BACE-1 and five BACE-1-inhibitor complexes to examine the effect of pH on dynamics and inhibitor binding properties of BACE-1. In our simulations, we find that solution pH controls the conformational flexibility of apo BACE-1, whereas bound inhibitors largely limit the motions of the holo enzyme at all levels of pH. The microscopic pKa values of titratable residues in BACE-1 including its aspartyl dyad are computed and compared between apo and inhibitor-bound states. Changes in protonation between the apo and holo forms suggest a thermodynamic linkage between binding of inhibitors and protons localized at the dyad. Utilizing our recently developed computational protocol applying the binding polynomial formalism to the constant pH molecular dynamics (CpHMD) framework, we are able to obtain the pH-dependent binding free energy profiles for various BACE-1-inhibitor complexes. Our results highlight the importance of correctly addressing the binding-induced protonation changes in protein-ligand systems where binding accompanies a net proton transfer. This work comprises the first application of our CpHMD-based free energy computational method to protein-ligand complexes and illustrates the value of CpHMD as an all-purpose tool for obtaining pH-dependent dynamics and binding free energies of biological systems.
Accelerated Molecular Dynamics Simulations of Protein Folding.J. Comp. Chem., Volume 36, Issue 20, pp 1536-1549 (2015)
Folding of four fast-folding proteins, including chignolin, Trp-cage, villin headpiece and WW domain, was simulated via accelerated molecular dynamics (aMD). In comparison with hundred-of-microsecond timescale conventional molecular dynamics (cMD) simulations performed on the Anton supercomputer, aMD captured complete folding of the four proteins in significantly shorter simulation time. The folded protein conformations were found within 0.2–2.1 Å of the native NMR or X-ray crystal structures. Free energy profiles calculated through improved reweighting of the aMD simulations using cumulant expansion to the second-order are in good agreement with those obtained from cMD simulations. This allows us to identify distinct conformational states (e.g., unfolded and intermediate) other than the native structure and the protein folding energy barriers. Detailed analysis of protein secondary structures and local key residue interactions provided important insights into the protein folding pathways. Furthermore, the selections of force fields and aMD simulation parameters are discussed in detail. Our work shows usefulness and accuracy of aMD in studying protein folding, providing basic references in using aMD in future protein-folding studies. © 2015 Wiley Periodicals, Inc.
Folding of four fast-folding proteins, including chignolin, Trp-cage, villin headpiece and WW domain, was simulated via accelerated molecular dynamics (aMD). In comparison with hundred-of-microsecond timescale conventional molecular dynamics (cMD) simulations performed on the Anton supercomputer, aMD captured complete folding of the four proteins in significantly shorter simulation time. The folded protein conformations were found within 0.2–2.1 Å of the native NMR or X-ray crystal structures. Free energy profiles calculated through improved reweighting of the aMD simulations using cumulant expansion to the second-order are in good agreement with those obtained from cMD simulations. This allows us to identify distinct conformational states (e.g., unfolded and intermediate) other than the native structure and the protein folding energy barriers. Detailed analysis of protein secondary structures and local key residue interactions provided important insights into the protein folding pathways. Furthermore, the selections of force fields and aMD simulation parameters are discussed in detail. Our work shows usefulness and accuracy of aMD in studying protein folding, providing basic references in using aMD in future protein-folding studies. © 2015 Wiley Periodicals, Inc.
Broken-Symmetry DFT Computations for the Reaction Pathway of IspH, an Iron−Sulfur Enzyme in Pathogenic BacteriaInorg. Chem. 54, 6439−6461 (2015)
The recently discovered methylerythritol phosphate (MEP) pathway provides new targets for the development of antibacterial and antimalarial drugs. In the final step of the MEP pathway, the [4Fe-4S] IspH protein catalyzes the 2e(-)/2H(+) reductive dehydroxylation of (E)-4-hydroxy-3-methyl-but-2-enyl diphosphate (HMBPP) to afford the isoprenoid precursors isopentenyl pyrophosphate (IPP) and dimethylallyl pyrophosphate (DMAPP). Recent experiments have attempted to elucidate the IspH catalytic mechanism to drive inhibitor development. Two competing mechanisms have recently emerged, differentiated by their proposed HMBPP binding modes upon 1e(-) reduction of the [4Fe-4S] cluster: (1) a Birch reduction mechanism, in which HMBPP remains bound to the [4Fe-4S] cluster through its terminal C4-OH group (ROH-bound) until the -OH is cleaved as water; and (2) an organometallic mechanism, in which the C4-OH group rotates away from the [4Fe-4S] cluster, allowing the HMBPP olefin group to form a metallacycle complex with the apical iron (η(2)-bound). We perform broken-symmetry density functional theory computations to assess the energies and reduction potentials associated with the ROH- and η(2)-bound states implicated by these competing mechanisms. Reduction potentials obtained for ROH-bound states are more negative (-1.4 to -1.0 V) than what is typically expected of [4Fe-4S] ferredoxin proteins. Instead, we find that η(2)-bound states are lower in energy than ROH-bound states when the [4Fe-4S] cluster is 1e(-) reduced. Furthermore, η(2)-bound states can already be generated in the oxidized state, yielding reduction potentials of ca. -700 mV when electron addition occurs after rotation of the HMBPP C4-OH group. We demonstrate that such η(2)-bound states are kinetically accessible both when the IspH [4Fe-4S] cluster is oxidized and 1e(-) reduced. The energetically preferred pathway gives 1e(-) reduction of the cluster after substrate conformational change, generating the 1e(-) reduced intermediate proposed in the organometallic mechanism.
The recently discovered methylerythritol phosphate (MEP) pathway provides new targets for the development of antibacterial and antimalarial drugs. In the final step of the MEP pathway, the [4Fe-4S] IspH protein catalyzes the 2e(-)/2H(+) reductive dehydroxylation of (E)-4-hydroxy-3-methyl-but-2-enyl diphosphate (HMBPP) to afford the isoprenoid precursors isopentenyl pyrophosphate (IPP) and dimethylallyl pyrophosphate (DMAPP). Recent experiments have attempted to elucidate the IspH catalytic mechanism to drive inhibitor development. Two competing mechanisms have recently emerged, differentiated by their proposed HMBPP binding modes upon 1e(-) reduction of the [4Fe-4S] cluster: (1) a Birch reduction mechanism, in which HMBPP remains bound to the [4Fe-4S] cluster through its terminal C4-OH group (ROH-bound) until the -OH is cleaved as water; and (2) an organometallic mechanism, in which the C4-OH group rotates away from the [4Fe-4S] cluster, allowing the HMBPP olefin group to form a metallacycle complex with the apical iron (η(2)-bound). We perform broken-symmetry density functional theory computations to assess the energies and reduction potentials associated with the ROH- and η(2)-bound states implicated by these competing mechanisms. Reduction potentials obtained for ROH-bound states are more negative (-1.4 to -1.0 V) than what is typically expected of [4Fe-4S] ferredoxin proteins. Instead, we find that η(2)-bound states are lower in energy than ROH-bound states when the [4Fe-4S] cluster is 1e(-) reduced. Furthermore, η(2)-bound states can already be generated in the oxidized state, yielding reduction potentials of ca. -700 mV when electron addition occurs after rotation of the HMBPP C4-OH group. We demonstrate that such η(2)-bound states are kinetically accessible both when the IspH [4Fe-4S] cluster is oxidized and 1e(-) reduced. The energetically preferred pathway gives 1e(-) reduction of the cluster after substrate conformational change, generating the 1e(-) reduced intermediate proposed in the organometallic mechanism.
Gaussian Accelerated Molecular Dynamics: Unconstrained Enhanced Sampling and Free Energy CalculationJ. Chem. Theory Comput, 11 (8), pp 3584-3595 (2015)
A Gaussian accelerated molecular dynamics (GaMD) approach for simultaneous enhanced sampling and free energy calculation of biomolecules is presented. By constructing a boost potential that follows Gaussian distribution, accurate reweighting of the GaMD simulations is achieved using cumulant expansion to the second order. Here, GaMD is demonstrated on three biomolecular model systems: alanine dipeptide, chignolin folding, and ligand binding to the T4-lysozyme. Without the need to set predefined reaction coordinates, GaMD enables unconstrained enhanced sampling of these biomolecules. Furthermore, the free energy profiles obtained from reweighting of the GaMD simulations allow us to identify distinct low-energy states of the biomolecules and characterize the protein-folding and ligand-binding pathways quantitatively.
A Gaussian accelerated molecular dynamics (GaMD) approach for simultaneous enhanced sampling and free energy calculation of biomolecules is presented. By constructing a boost potential that follows Gaussian distribution, accurate reweighting of the GaMD simulations is achieved using cumulant expansion to the second order. Here, GaMD is demonstrated on three biomolecular model systems: alanine dipeptide, chignolin folding, and ligand binding to the T4-lysozyme. Without the need to set predefined reaction coordinates, GaMD enables unconstrained enhanced sampling of these biomolecules. Furthermore, the free energy profiles obtained from reweighting of the GaMD simulations allow us to identify distinct low-energy states of the biomolecules and characterize the protein-folding and ligand-binding pathways quantitatively.
Enzyme localization, crowding and buffers collectively modulate diffusion influenced signal transduction: insights from continuum diffusion modelingJ. Chem. Phys. 143, 094103 (2015).
Biochemical reaction networks consisting of coupled enzymes connect substrate signaling events with biological function. Substrates involved in these reactions can be strongly influenced by diffusion “barriers” arising from impenetrable cellular structures and macromolecules, as well as interactions with biomolecules, especially within crowded environments. For diffusion-influenced reactions, the spatial organization of diffusion barriers arising from intracellular structures, non-specific crowders, and specific-binders (buffers) strongly controls the temporal and spatial reaction kinetics. In this study, we use two prototypical biochemical reactions, a Goodwin oscillator, and a reaction with a periodic source/sink term to examine how a diffusion barrier that partitions substrates controls reaction behavior. Namely, we examine how conditions representative of a densely packed cytosol, including reduced accessible volume fraction, non-specific interactions, and buffers, impede diffusion over nanometer length-scales. We find that diffusion barriers can modulate the frequencies and amplitudes of coupled diffusion-influenced reaction networks, as well as give rise to “compartments” of decoupled reactant populations. These effects appear to be intensified in the presence of buffers localized to the diffusion barrier. These findings have strong implications for the role of the cellular environment in tuning the dynamics of signaling pathways.
Biochemical reaction networks consisting of coupled enzymes connect substrate signaling events with biological function. Substrates involved in these reactions can be strongly influenced by diffusion “barriers” arising from impenetrable cellular structures and macromolecules, as well as interactions with biomolecules, especially within crowded environments. For diffusion-influenced reactions, the spatial organization of diffusion barriers arising from intracellular structures, non-specific crowders, and specific-binders (buffers) strongly controls the temporal and spatial reaction kinetics. In this study, we use two prototypical biochemical reactions, a Goodwin oscillator, and a reaction with a periodic source/sink term to examine how a diffusion barrier that partitions substrates controls reaction behavior. Namely, we examine how conditions representative of a densely packed cytosol, including reduced accessible volume fraction, non-specific interactions, and buffers, impede diffusion over nanometer length-scales. We find that diffusion barriers can modulate the frequencies and amplitudes of coupled diffusion-influenced reaction networks, as well as give rise to “compartments” of decoupled reactant populations. These effects appear to be intensified in the presence of buffers localized to the diffusion barrier. These findings have strong implications for the role of the cellular environment in tuning the dynamics of signaling pathways.
Electrostatic steering enhances the rate of cAMP binding to phosphodiesterase: Brownian dynamics modelingProt. Sci. in press (2015).
Signaling in cells often involves co-localization of the signaling molecules. Most experimental evidence has shown that intracellular compartmentalization restricts the range of action of the second messenger, 3'-5'-cyclic adenosine monophosphate (cAMP), which is degraded by phosphodiesterases (PDEs). The objective of this study is to understand the details of molecular encounter that may play a role in efficient operation of the cAMP signaling apparatus. The results from electrostatic potential calculations and Brownian dynamics simulations suggest that positive potential of the active site from PDE enhances capture of diffusing cAMP molecules. This electrostatic steering between cAMP and the active site of a PDE plays a major role in the enzyme-substrate encounter, an effect that may be of significance in sequestering cAMP released from a nearby binding site or in attracting more freely diffusing cAMP molecules.
Signaling in cells often involves co-localization of the signaling molecules. Most experimental evidence has shown that intracellular compartmentalization restricts the range of action of the second messenger, 3'-5'-cyclic adenosine monophosphate (cAMP), which is degraded by phosphodiesterases (PDEs). The objective of this study is to understand the details of molecular encounter that may play a role in efficient operation of the cAMP signaling apparatus. The results from electrostatic potential calculations and Brownian dynamics simulations suggest that positive potential of the active site from PDE enhances capture of diffusing cAMP molecules. This electrostatic steering between cAMP and the active site of a PDE plays a major role in the enzyme-substrate encounter, an effect that may be of significance in sequestering cAMP released from a nearby binding site or in attracting more freely diffusing cAMP molecules.
A Self-Consistent Phase-Field Approach to Implicit Solvation of Charged Molecules with Poisson-Boltzmann ElectrostaticsJ. Chem Phys. 143, 243110 (2015)
Dielectric boundary based implicit-solvent models provide efficient descriptions of coarse-grained effects, particularly the electrostatic effect, of aqueous solvent. Recent years have seen the initial success of a new such model, variational implicit-solvent model (VISM) [Dzubiella, Swanson, and McCammon Phys. Rev. Lett. 96, 087802 (2006) and J. Chem. Phys. 124, 084905 (2006)], in capturing multiple dry and wet hydration states, describing the subtle electrostatic effect in hydrophobic interactions, and providing qualitatively good estimates of solvation free energies. Here, we develop a phase-field VISM to the solvation of charged molecules in aqueous solvent to include more flexibility. In this approach, a stable equilibrium molecular system is described by a phase field that takes one constant value in the solute region and a different constant value in the solvent region, and smoothly changes its value on a thin transition layer representing a smeared solute-solvent interface or dielectric boundary. Such a phase field minimizes an effective solvation free-energy functional that consists of the solute-solvent interfacial energy, solute-solvent van der Waals interaction energy, and electrostatic free energy described by the Poisson–Boltzmann theory. We apply our model and methods to the solvation of single ions, two parallel plates, and protein complexes BphC and p53/MDM2 to demonstrate the capability and efficiency of our approach at different levels. With a diffuse dielectric boundary, our new approach can describe the dielectric asymmetry in the solute-solvent interfacial region. Our theory is developed based on rigorous mathematical studies and is also connected to the Lum–Chandler–Weeks theory (1999). We discuss these connections and possible extensions of our theory and methods.
Dielectric boundary based implicit-solvent models provide efficient descriptions of coarse-grained effects, particularly the electrostatic effect, of aqueous solvent. Recent years have seen the initial success of a new such model, variational implicit-solvent model (VISM) [Dzubiella, Swanson, and McCammon Phys. Rev. Lett. 96, 087802 (2006) and J. Chem. Phys. 124, 084905 (2006)], in capturing multiple dry and wet hydration states, describing the subtle electrostatic effect in hydrophobic interactions, and providing qualitatively good estimates of solvation free energies. Here, we develop a phase-field VISM to the solvation of charged molecules in aqueous solvent to include more flexibility. In this approach, a stable equilibrium molecular system is described by a phase field that takes one constant value in the solute region and a different constant value in the solvent region, and smoothly changes its value on a thin transition layer representing a smeared solute-solvent interface or dielectric boundary. Such a phase field minimizes an effective solvation free-energy functional that consists of the solute-solvent interfacial energy, solute-solvent van der Waals interaction energy, and electrostatic free energy described by the Poisson–Boltzmann theory. We apply our model and methods to the solvation of single ions, two parallel plates, and protein complexes BphC and p53/MDM2 to demonstrate the capability and efficiency of our approach at different levels. With a diffuse dielectric boundary, our new approach can describe the dielectric asymmetry in the solute-solvent interfacial region. Our theory is developed based on rigorous mathematical studies and is also connected to the Lum–Chandler–Weeks theory (1999). We discuss these connections and possible extensions of our theory and methods.
Troponin I Mutations R146G and R21C Alter Cardiac Troponin Function, Contractile Properties and Modulation by PKA-mediated PhosphorylationJ. Biol. Chem. 290, 27749-27766 (2015)
Two hypertrophic cardiomyopathy-associated cardiac troponin I (cTnI) mutations, R146G and R21C, are located in different regions of cTnI, the inhibitory peptide and the cardiac-specific N terminus. We recently reported that these regions may interact when Ser-23/Ser-24 are phosphorylated, weakening the interaction of cTnI with cardiac TnC. Little is known about how these mutations influence the affinity of cardiac TnC for cTnI (KC-I) or contractile kinetics during β-adrenergic stimulation. Here, we tested how cTnI(R146G) or cTnI(R21C) influences contractile activation and relaxation and their response to protein kinase A (PKA). Both mutations significantly increased Ca(2+) binding affinity to cTn (KCa) and KC-I. PKA phosphorylation resulted in a similar reduction of KCa for all complexes, but KC-I was reduced only with cTnI(WT). cTnI(WT), cTnI(R146G), and cTnI(R21C) were complexed into cardiac troponin and exchanged into rat ventricular myofibrils, and contraction/relaxation kinetics were measured ± PKA phosphorylation. Maximal tension (Tmax) was maintained for cTnI(R146G)- and cTnI(R21C)-exchanged myofibrils, and Ca(2+) sensitivity of tension (pCa50) was increased. PKA phosphorylation decreased pCa50 for cTnI(WT)-exchanged myofibrils but not for either mutation. PKA phosphorylation accelerated the early slow phase relaxation for cTnI(WT) myofibrils, especially at Ca(2+) levels that the heart operates in vivo. Importantly, this effect was blunted for cTnI(R146G)- and cTnI(R21C)-exchanged myofibrils. Molecular dynamics simulations suggest both mutations inhibit formation of intra-subunit contacts between the N terminus and the inhibitory peptide of cTnI that is normally seen with WT-cTn upon PKA phosphorylation. Together, our results suggest that cTnI(R146G) and cTnI(R21C) blunt PKA modulation of activation and relaxation kinetics by prohibiting cardiac-specific N-terminal interaction with the cTnI inhibitory peptide.
Two hypertrophic cardiomyopathy-associated cardiac troponin I (cTnI) mutations, R146G and R21C, are located in different regions of cTnI, the inhibitory peptide and the cardiac-specific N terminus. We recently reported that these regions may interact when Ser-23/Ser-24 are phosphorylated, weakening the interaction of cTnI with cardiac TnC. Little is known about how these mutations influence the affinity of cardiac TnC for cTnI (KC-I) or contractile kinetics during β-adrenergic stimulation. Here, we tested how cTnI(R146G) or cTnI(R21C) influences contractile activation and relaxation and their response to protein kinase A (PKA). Both mutations significantly increased Ca(2+) binding affinity to cTn (KCa) and KC-I. PKA phosphorylation resulted in a similar reduction of KCa for all complexes, but KC-I was reduced only with cTnI(WT). cTnI(WT), cTnI(R146G), and cTnI(R21C) were complexed into cardiac troponin and exchanged into rat ventricular myofibrils, and contraction/relaxation kinetics were measured ± PKA phosphorylation. Maximal tension (Tmax) was maintained for cTnI(R146G)- and cTnI(R21C)-exchanged myofibrils, and Ca(2+) sensitivity of tension (pCa50) was increased. PKA phosphorylation decreased pCa50 for cTnI(WT)-exchanged myofibrils but not for either mutation. PKA phosphorylation accelerated the early slow phase relaxation for cTnI(WT) myofibrils, especially at Ca(2+) levels that the heart operates in vivo. Importantly, this effect was blunted for cTnI(R146G)- and cTnI(R21C)-exchanged myofibrils. Molecular dynamics simulations suggest both mutations inhibit formation of intra-subunit contacts between the N terminus and the inhibitory peptide of cTnI that is normally seen with WT-cTn upon PKA phosphorylation. Together, our results suggest that cTnI(R146G) and cTnI(R21C) blunt PKA modulation of activation and relaxation kinetics by prohibiting cardiac-specific N-terminal interaction with the cTnI inhibitory peptide.
Anti-infectives targeting enzymes, and the proton motive forceProc. Natl. Acad. Sci. USA in press (2015).
There is a growing need for new antibiotics. Compounds that target the proton motive force (PMF), uncouplers, represent one possible class of compounds that might be developed because they are already used to treat parasitic infections, and there is interest in their use for the treatment of other diseases, such as diabetes. Here, we tested a series of compounds, most with known antiinfective activity, for uncoupler activity. Many cationic amphiphiles tested positive, and some targeted isoprenoid biosynthesis or affected lipid bilayer structure. As an example, we found that clomiphene, a recently discovered undecaprenyl diphosphate synthase inhibitor active against Staphylococcus aureus, is an uncoupler. Using in silico screening, we then found that the anti-glioblastoma multiforme drug lead vacquinol is an inhibitor of Mycobacterium tuberculosis tuberculosinyl adenosine synthase, as well as being an uncoupler. Because vacquinol is also an inhibitor of M. tuberculosis cell growth, we used similarity searches based on the vacquinol structure, finding analogs with potent (∼0.5-2 μg/mL) activity against M. tuberculosis and S. aureus. Our results give a logical explanation of the observation that most new tuberculosis drug leads discovered by phenotypic screens and genome sequencing are highly lipophilic (logP ∼5.7) bases with membrane targets because such species are expected to partition into hydrophobic membranes, inhibiting membrane proteins, in addition to collapsing the PMF. This multiple targeting is expected to be of importance in overcoming the development of drug resistance because targeting membrane physical properties is expected to be less susceptible to the development of resistance.
There is a growing need for new antibiotics. Compounds that target the proton motive force (PMF), uncouplers, represent one possible class of compounds that might be developed because they are already used to treat parasitic infections, and there is interest in their use for the treatment of other diseases, such as diabetes. Here, we tested a series of compounds, most with known antiinfective activity, for uncoupler activity. Many cationic amphiphiles tested positive, and some targeted isoprenoid biosynthesis or affected lipid bilayer structure. As an example, we found that clomiphene, a recently discovered undecaprenyl diphosphate synthase inhibitor active against Staphylococcus aureus, is an uncoupler. Using in silico screening, we then found that the anti-glioblastoma multiforme drug lead vacquinol is an inhibitor of Mycobacterium tuberculosis tuberculosinyl adenosine synthase, as well as being an uncoupler. Because vacquinol is also an inhibitor of M. tuberculosis cell growth, we used similarity searches based on the vacquinol structure, finding analogs with potent (∼0.5-2 μg/mL) activity against M. tuberculosis and S. aureus. Our results give a logical explanation of the observation that most new tuberculosis drug leads discovered by phenotypic screens and genome sequencing are highly lipophilic (logP ∼5.7) bases with membrane targets because such species are expected to partition into hydrophobic membranes, inhibiting membrane proteins, in addition to collapsing the PMF. This multiple targeting is expected to be of importance in overcoming the development of drug resistance because targeting membrane physical properties is expected to be less susceptible to the development of resistance.
