Structural and kinetic evidence of aging after organophosphate inhibition of human Cathepsin A.

Human Cathepsin A (CatA) is a lysosomal serine carboxypeptidase of the renin-angiotensin system (RAS) and is structurally similar to acetylcholinesterase (AChE). CatA can remove the C-terminal amino acids of endothelin I, angiotensin I, Substance P, oxytocin, and bradykinin, and can deamidate neurokinin A. Proteomic studies identified CatA and its homologue, SCPEP1, as potential targets of organophosphates (OP). CatA could be stably inhibited by low µM to high nM concentrations of racemic sarin (GB), soman (GD), cyclosarin (GF), VX, and VR within minutes to hours at pH 7. Cyclosarin was the most potent with a kinetically measured dissociation constant (KI) of 2 µM followed by VR (KI = 2.8 µM). Bimolecular rate constants for inhibition by cyclosarin and VR were 1.3 × 103 M-1sec-1 and 1.2 × 103 M-1sec-1, respectively, and were approximately 3-orders of magnitude lower than those of human AChE indicating slower reactivity. Notably, both AChE and CatA bound diisopropylfluorophosphate (DFP) comparably and had KIDFP = 13 µM and 11 µM, respectively. At low pH, greater than 85% of the enzyme spontaneously reactivated after OP inhibition, conditions under which OP-adducts of cholinesterases irreversibly age. At pH 6.5 CatA remained stably inhibited by GB and GF and <10% of the enzyme spontaneously reactivated after 200 h. A crystal structure of DFP-inhibited CatA was determined and contained an aged adduct. Similar to AChE, CatA appears to have a "backdoor" for product release. CatA has not been shown previously to age. These results may have implications for: OP-associated inflammation; cardiovascular effects; and the dysregulation of RAS enzymes by OP.

Paired Carboxylic Acids in Enzymes and Their Role in Selective Substrate Binding, Catalysis, and Unusually Shifted pKa Values.

Cathepsin A (CatA, EC 3.4.16.5, UniProtKB P10619 ) is a human lysosomal carboxypeptidase. Counterintuitively, crystal structures of CatA and its homologues show a cluster of Glu and Asp residues binding the C-terminal carboxylic acid of the product or inhibitor. Each of these enzymes functions in an acidic environment and contains a highly conserved pair of Glu residues with side chain carboxyl group oxygens that are approximately 2.3-2.6 Šapart. In small molecules, carboxyl groups separated by ∼3 Šcan overcome the repulsive interaction by protonation of one of the two groups. The pKa of one group increases (pKa ∼ 11) and can be as much as ∼6 pH units higher than the paired group. Consequently, at low and neutral pH, one carboxylate can carry a net negative charge while the other can remain protonated and neutral. In CatA, E69 and E149 form a Glu pair that is important to catalysis as evidenced by the 56-fold decrease in kcat/Km in the E69Q/E149Q variant. Here, we have measured the pH dependencies of log(kcat), log(Km), and log(kcat/Km) for wild type CatA and its variants and have compared the measured pKa with calculated values. We propose a substrate-assisted mechanism in which the high pKa of E149 (>8.5) favors the binding of the carboxylate form of the substrate and promotes the abstraction of the proton from H429 of the catalytic triad effectively decreasing its pKa in a low-pH environment. We also identify a similar motif consisting of a pair of histidines in S-formylglutathione hydrolase.

Quantitative Analysis of Repertoire-Scale Immunoglobulin Properties in Vaccine-Induced B-Cell Responses.

Recent advances in the next-generation sequencing of B-cell receptors (BCRs) enable the characterization of humoral responses at a repertoire-wide scale and provide the capability for identifying unique features of immune repertoires in response to disease, vaccination, or infection. Immunosequencing now readily generates 103-105 sequences per sample; however, statistical analysis of these repertoires is challenging because of the high genetic diversity of BCRs and the elaborate clonal relationships among them. To date, most immunosequencing analyses have focused on reporting qualitative trends in immunoglobulin (Ig) properties, such as usage or somatic hypermutation (SHM) percentage of the Ig heavy chain variable (IGHV) gene segment family, and on reducing complex Ig property distributions to simple summary statistics. However, because Ig properties are typically not normally distributed, any approach that fails to assess the distribution as a whole may be inadequate in (1) properly assessing the statistical significance of repertoire differences, (2) identifying how two repertoires differ, and (3) determining appropriate confidence intervals for assessing the size of the differences and their potential biological relevance. To address these issues, we have developed a technique that uses Wilcox' robust statistics toolbox to identify statistically significant vaccine-specific differences between Ig repertoire properties. The advantage of this technique is that it can determine not only whether but also where the distributions differ, even when the Ig repertoire properties are non-normally distributed. We used this technique to characterize murine germinal center (GC) B-cell repertoires in response to a complex Ebola virus-like particle (eVLP) vaccine candidate with known protective efficacy. The eVLP-mediated GC B-cell responses were highly diverse, consisting of thousands of clonotypes. Despite this staggering diversity, we identified statistically significant differences between non-immunized, vaccine only, and vaccine-plus-adjuvant groups in terms of Ig properties, including IGHV-family usage, SHM percentage, and characteristics of the BCR complementarity-determining region. Most notably, our analyses identified a robust eVLP-specific feature-enhanced IGHV8-family usage in B-cell repertoires. These findings demonstrate the utility of our technique in identifying statistically significant BCR repertoire differences following vaccination. More generally, our approach is potentially applicable to a wide range of studies in infection, vaccination, auto-immunity, and cancer.

Dengue virus antibody database: Systematically linking serotype-specificity with epitope mapping in dengue virus.

BACKGROUND: A majority infections caused by dengue virus (DENV) are asymptomatic, but a higher incidence of severe illness, such as dengue hemorrhagic fever, is associated with secondary infections, suggesting that pre-existing immunity plays a central role in dengue pathogenesis. Primary infections are typically associated with a largely serotype-specific antibody response, while secondary infections show a shift to a broadly cross-reactive antibody response. METHODS/PRINCIPAL FINDINGS: We hypothesized that the basis for the shift in serotype-specificity between primary and secondary infections can be found in a change in the antibody fine-specificity. To investigate the link between epitope- and serotype-specificity, we assembled the Dengue Virus Antibody Database, an online repository containing over 400 DENV-specific mAbs, each annotated with information on 1) its origin, including the immunogen, host immune history, and selection methods, 2) binding/neutralization data against all four DENV serotypes, and 3) epitope mapping at the domain or residue level to the DENV E protein. We combined epitope mapping and activity information to determine a residue-level index of epitope propensity and cross-reactivity and generated detailed composite epitope maps of primary and secondary antibody responses. We found differing patterns of epitope-specificity between primary and secondary infections, where secondary responses target a distinct subset of epitopes found in the primary response. We found that secondary infections were marked with an enhanced response to cross-reactive epitopes, such as the fusion-loop and E-dimer region, as well as increased cross-reactivity in what are typically more serotype-specific epitope regions, such as the domain I-II interface and domain III. CONCLUSIONS/SIGNIFICANCE: Our results support the theory that pre-existing cross-reactive memory B cells form the basis for the secondary antibody response, resulting in a broadening of the response in terms of cross-reactivity, and a focusing of the response to a subset of epitopes, including some, such as the fusion-loop region, that are implicated in poor neutralization and antibody-dependent enhancement of infection.

Modeling the Role of Epitope Arrangement on Antibody Binding Stoichiometry in Flaviviruses.

Cryo-electron-microscopy (cryo-EM) structures of flaviviruses reveal significant variation in epitope occupancy across different monoclonal antibodies that have largely been attributed to epitope-level differences in conformation or accessibility that affect antibody binding. The consequences of these variations for macroscopic properties such as antibody binding and neutralization are the results of the law of mass action-a stochastic process of innumerable binding and unbinding events between antibodies and the multiple binding sites on the flavivirus in equilibrium-that cannot be directly imputed from structure alone. We carried out coarse-grained spatial stochastic binding simulations for nine flavivirus antibodies with epitopes defined by cryo-EM or x-ray crystallography to assess the role of epitope spatial arrangement on antibody-binding stoichiometry, occupancy, and neutralization. In our simulations, all epitopes were equally competent for binding, representing the upper limit of binding stoichiometry that results from epitope spatial arrangement alone. Surprisingly, our simulations closely reproduced the relative occupancy and binding stoichiometry observed in cryo-EM, without having to account for differences in epitope accessibility or conformation, suggesting that epitope spatial arrangement alone may be sufficient to explain differences in binding occupancy and stoichiometry between antibodies. Furthermore, we found that there was significant heterogeneity in binding configurations even at saturating antibody concentrations, and that bivalent antibody binding may be more common than previously thought. Finally, we propose a structure-based explanation for the stoichiometric threshold model of neutralization.

BRILIA: Integrated Tool for High-Throughput Annotation and Lineage Tree Assembly of B-Cell Repertoires.

The somatic diversity of antigen-recognizing B-cell receptors (BCRs) arises from Variable (V), Diversity (D), and Joining (J) (VDJ) recombination and somatic hypermutation (SHM) during B-cell development and affinity maturation. The VDJ junction of the BCR heavy chain forms the highly variable complementarity determining region 3 (CDR3), which plays a critical role in antigen specificity and binding affinity. Tracking the selection and mutation of the CDR3 can be useful in characterizing humoral responses to infection and vaccination. Although tens to hundreds of thousands of unique BCR genes within an expressed B-cell repertoire can now be resolved with high-throughput sequencing, tracking SHMs is still challenging because existing annotation methods are often limited by poor annotation coverage, inconsistent SHM identification across the VDJ junction, or lack of B-cell lineage data. Here, we present B-cell repertoire inductive lineage and immunosequence annotator (BRILIA), an algorithm that leverages repertoire-wide sequencing data to globally improve the VDJ annotation coverage, lineage tree assembly, and SHM identification. On benchmark tests against simulated human and mouse BCR repertoires, BRILIA correctly annotated germline and clonally expanded sequences with 94 and 70% accuracy, respectively, and it has a 90% SHM-positive prediction rate in the CDR3 of heavily mutated sequences; these are substantial improvements over existing methods. We used BRILIA to process BCR sequences obtained from splenic germinal center B cells extracted from C57BL/6 mice. BRILIA returned robust B-cell lineage trees and yielded SHM patterns that are consistent across the VDJ junction and agree with known biological mechanisms of SHM. By contrast, existing BCR annotation tools, which do not account for repertoire-wide clonal relationships, systematically underestimated both the size of clonally related B-cell clusters and yielded inconsistent SHM frequencies. We demonstrate BRILIA's utility in B-cell repertoire studies related to VDJ gene usage, mechanisms for adenosine mutations, and SHM hot spot motifs. Furthermore, we show that the complete gene usage annotation and SHM identification across the entire CDR3 are essential for studying the B-cell affinity maturation process through immunosequencing methods.

A reaction path study of the catalysis and inhibition of the Bacillus anthracis CapD γ-glutamyl transpeptidase.

The CapD enzyme of Bacillus anthracis is a γ-glutamyl transpeptidase from the N-terminal nucleophile hydrolase superfamily that covalently anchors the poly-γ-D-glutamic acid (pDGA) capsule to the peptidoglycan. The capsule hinders phagocytosis of B. anthracis by host cells and is essential for virulence. The role CapD plays in capsule anchoring and remodeling makes the enzyme a promising target for anthrax medical countermeasures. Although the structure of CapD is known, and a covalent inhibitor, capsidin, has been identified, the mechanisms of CapD catalysis and inhibition are poorly understood. Here, we used a computational approach to map out the reaction steps involved in CapD catalysis and inhibition. We found that the rate-limiting step of either CapD catalysis or inhibition was a concerted asynchronous formation of the tetrahedral intermediate with a barrier of 22-23 kcal/mol. However, the mechanisms of these reactions differed for the two amides. The formation of the tetrahedral intermediate with pDGA was substrate-assisted with two proton transfers. In contrast, capsidin formed the tetrahedral intermediate in a conventional way with one proton transfer. Interestingly, capsidin coupled a conformational change in the catalytic residue of the tetrahedral intermediate to stretching of the scissile amide bond. Furthermore, capsidin took advantage of iminol-amide tautomerism of its diacetamide moiety to convert the tetrahedral intermediate to the acetylated CapD. As evidence of the promiscuous nature of CapD, the enzyme cleaved the amide bond of capsidin by attacking it on the opposite side compared to pDGA.

Exploring chemical reaction mechanisms through harmonic Fourier beads path optimization.

Here, we apply the harmonic Fourier beads (HFB) path optimization method to study chemical reactions involving covalent bond breaking and forming on quantum mechanical (QM) and hybrid QM∕molecular mechanical (QM∕MM) potential energy surfaces. To improve efficiency of the path optimization on such computationally demanding potentials, we combined HFB with conjugate gradient (CG) optimization. The combined CG-HFB method was used to study two biologically relevant reactions, namely, L- to D-alanine amino acid inversion and alcohol acylation by amides. The optimized paths revealed several unexpected reaction steps in the gas phase. For example, on the B3LYP∕6-31G(d,p) potential, we found that alanine inversion proceeded via previously unknown intermediates, 2-iminopropane-1,1-diol and 3-amino-3-methyloxiran-2-ol. The CG-HFB method accurately located transition states, aiding in the interpretation of complex reaction mechanisms. Thus, on the B3LYP∕6-31G(d,p) potential, the gas phase activation barriers for the inversion and acylation reactions were 50.5 and 39.9 kcal∕mol, respectively. These barriers determine the spontaneous loss of amino acid chirality and cleavage of peptide bonds in proteins. We conclude that the combined CG-HFB method further advances QM and QM∕MM studies of reaction mechanisms.

Quantitative predictions of binding free energy changes in drug-resistant influenza neuraminidase.

Quantitatively predicting changes in drug sensitivity associated with residue mutations is a major challenge in structural biology. By expanding the limits of free energy calculations, we successfully identified mutations in influenza neuraminidase (NA) that confer drug resistance to two antiviral drugs, zanamivir and oseltamivir. We augmented molecular dynamics (MD) with Hamiltonian Replica Exchange and calculated binding free energy changes for H274Y, N294S, and Y252H mutants. Based on experimental data, our calculations achieved high accuracy and precision compared with results from established computational methods. Analysis of 15 micros of aggregated MD trajectories provided insights into the molecular mechanisms underlying drug resistance that are at odds with current interpretations of the crystallographic data. Contrary to the notion that resistance is caused by mutant-induced changes in hydrophobicity of the binding pocket, our simulations showed that drug resistance mutations in NA led to subtle rearrangements in the protein structure and its dynamics that together alter the active-site electrostatic environment and modulate inhibitor binding. Importantly, different mutations confer resistance through different conformational changes, suggesting that a generalized mechanism for NA drug resistance is unlikely.

Probing the donor and acceptor substrate specificity of the γ-glutamyl transpeptidase.

γ-Glutamyl transpeptidase (GGT) is a two-substrate enzyme that plays a central role in glutathione metabolism and is a potential target for drug design. GGT catalyzes the cleavage of γ-glutamyl donor substrates and the transfer of the γ-glutamyl moiety to an amine of an acceptor substrate or water. Although structures of bacterial GGT have revealed details of the protein-ligand interactions at the donor site, the acceptor substrate site is relatively undefined. The recent identification of a species-specific acceptor site inhibitor, OU749, suggests that these inhibitors may be less toxic than glutamine analogues. Here we investigated the donor and acceptor substrate preferences of Bacillus anthracis GGT (CapD) and applied computational approaches in combination with kinetics to probe the structural basis of the enzyme's substrate and inhibitor binding specificities and compare them with human GGT. Site-directed mutagenesis studies showed that the R432A and R520S variants exhibited 6- and 95-fold decreases in hydrolase activity, respectively, and that their activity was not stimulated by the addition of the l-Cys acceptor substrate, suggesting an additional role in acceptor binding and/or catalysis of transpeptidation. Rat GGT (and presumably HuGGT) has strict stereospecificity for L-amino acid acceptor substrates, while CapD can utilize both L- and D-acceptor substrates comparably. Modeling and kinetic analysis suggest that R520 and R432 allow two alternate acceptor substrate binding modes for L- and D-acceptors. R432 is conserved in Francisella tularensis, Yersinia pestis, Burkholderia mallei, Helicobacter pylori and Escherichia coli, but not in human GGT. Docking and MD simulations point toward key residues that contribute to inhibitor and acceptor substrate binding, providing a guide to designing novel and specific GGT inhibitors.

Fungal bis-Naphthopyrones as Inhibitors of Botulinum Neurotoxin Serotype A.

An in silico screen of the NIH Molecular Library Small Molecule Repository (MLSMR) of ∼350000 compounds and confirmatory bioassays led to identification of chaetochromin A (1) as an inhibitor of botulinum neurotoxin serotype A (BoNT A). Subsequent acquisition and testing of analogues of 1 uncovered two compounds, talaroderxines A (2) and B (3), with improved activity. These are the first fungal metabolites reported to exhibit BoNT/A inhibitory activity.

Improved Binding Free Energy Predictions from Single-Reference Thermodynamic Integration Augmented with Hamiltonian Replica Exchange.

Reliable predictions of relative binding free energies are essential in drug discovery, where chemists modify promising compounds with the aim of increasing binding affinity. Conventional Thermodynamic Integration (TI) approaches can estimate corresponding changes in binding free energies, but suffer from inadequate sampling due to ruggedness of the molecular energy surfaces. Here, we present an improved TI strategy for computing relative binding free energies of congeneric ligands. This strategy employs a specific, unphysical single-reference (SR) state and Hamiltonian replica exchange (HREX) to locally enhance sampling. We then apply this strategy to compute relative binding free energies of twelve ligands in the L99A mutant of T4 Lysozyme. Besides the ligands, our approach enhances hindered rotations of the important V111, as well as V87 and L118 sidechains. Concurrently, we devise practical strategies to monitor and improve HREX-SRTI efficiency. Overall, the HREX-SRTI results agree well (R(2) = 0.76, RMSE = 0.3 kcal/mol) with available experimental data. When optimized for efficiency, the HREX-SRTI precision matches that of experimental measurements.

Separation of Betti Reaction Product Enantiomers: Absolute Configuration and Inhibition of Botulinum Neurotoxin A.

The racemic product of the Betti reaction of 5-chloro-8-hydroxyquinoline, benzaldehyde and 2-aminopyridine was separated by chiral HPLC to determine which enantiomer inhibited botulinum neurotoxin serotype A. When the enantiomers unexpectedly proved to have comparable activity, the absolute structures of (+)-(R)-1 and (-)-(S)-1 were determined by comparison of calculated and observed circular dichroism spectra. Molecular modeling studies were undertaken in an effort to understand the observed bioactivity and revealed different ensembles of binding modes, with roughly equal binding energies, for the two enantiomers.

Computing Relative Free Energies of Solvation using Single Reference Thermodynamic Integration Augmented with Hamiltonian Replica Exchange.

This paper introduces an efficient single-topology variant of Thermodynamic Integration (TI) for computing relative transformation free energies in a series of molecules with respect to a single reference state. The presented TI variant that we refer to as Single-Reference TI (SR-TI) combines well-established molecular simulation methodologies into a practical computational tool. Augmented with Hamiltonian Replica Exchange (HREX), the SR-TI variant can deliver enhanced sampling in select degrees of freedom. The utility of the SR-TI variant is demonstrated in calculations of relative solvation free energies for a series of benzene derivatives with increasing complexity. Noteworthy, the SR-TI variant with the HREX option provides converged results in a challenging case of an amide molecule with a high (13-15 kcal/mol) barrier for internal cis/trans interconversion using simulation times of only 1 to 4 ns.

A transition path ensemble study reveals a linchpin role for Mg(2+) during rate-limiting ADP release from protein kinase A.

Protein kinases are key regulators of diverse signaling networks critical for growth and development. Protein kinase A (PKA) is an important kinase prototype that phosphorylates protein targets at Ser and Thr residues by converting ATP to ADP. Mg(2+) ions play a crucial role in regulating phosphoryl transfer and can limit overall enzyme turnover by affecting ADP release. However, the mechanism by which Mg(2+) participates in ADP release is poorly understood. Here we use a novel transition path ensemble technique, the harmonic Fourier beads method, to explore the atomic and energetic details of the Mg(2+)-dependent ADP binding and release. Our studies demonstrate that adenine-driven ADP binding to PKA creates three ion-binding sites at the ADP/PKA interface that are absent otherwise. Two of these sites bind the previously characterized Mg(2+) ions, whereas the third site binds a monovalent cation with high affinity. This third site can bind the P-3 residue of substrate proteins and may serve as a reporter of the active site occupation. Binding of Mg(2+) ions restricts mobility of the Gly-rich loop that closes over the active site. We find that simultaneous release of ADP with Mg(2+) ions from the active site is unfeasible. Thus, we conclude that Mg(2+) ions act as a linchpin and that at least one ion must be removed prior to pyrophosphate-driven ADP release. The results of the present study enhance understanding of Mg(2+)-dependent association of nucleotides with protein kinases.

Free energy for the permeation of Na(+) and Cl(-) ions and their ion-pair through a zwitterionic dimyristoyl phosphatidylcholine lipid bilayer by umbrella integration with harmonic fourier beads.

Understanding the mechanism of ion permeation across lipid bilayers is key to controlling osmotic pressure and developing new ways of delivering charged, drug-like molecules inside cells. Recent reports suggest ion-pairing as the mechanism to lower the free energy barrier for the ion permeation in disagreement with predictions from the simple electrostatic models. In this paper we quantify the effect of ion-pairing or charge quenching on the permeation of Na(+) and Cl(-) ions across DMPC lipid bilayer by computing the corresponding potentials of mean force (PMFs) using fully atomistic molecular dynamics simulations. We find that the free energy barrier to permeation reduces in the order Na(+)-Cl(-) ion-pair (27.6 kcal/mol) > Cl(-) (23.6 kcal/mol) > Na(+) (21.9 kcal/mol). Furthermore, with the help of these PMFs we derive the change in the binding free energy between the Na(+) and Cl(-) with respect to that in water as a function of the bilayer permeation depth. Despite the fact that the bilayer boosts the Na(+)-Cl(-) ion binding free energy by as high as 17.9 kcal/mol near its center, ion-pairing between such hydrophilic ions as Na(+) and Cl(-) does not assist their permeation. However, based on a simple thermodynamic cycle, we suggest that ion-pairing between ions of opposite charge and solvent philicity could enhance ion permeation. Comparison of the computed permeation barriers for Na(+) and Cl(-) ions with available experimental data supports this notion. This work establishes general computational methodology to address ion-pairing in fluid anisotropic media and details the ion permeation mechanism on atomic level.

Computing accurate potentials of mean force in electrolyte solutions with the generalized gradient-augmented harmonic Fourier beads method.

We establish the accuracy of the novel generalized gradient-augmented harmonic Fourier beads (ggaHFB) method in computing free-energy profiles or potentials of mean force (PMFs) through comparison with two independent conventional techniques. In particular, we employ umbrella sampling with one dimensional weighted histogram analysis method (WHAM) and free molecular dynamics simulation of radial distribution functions to compute the PMF for the Na(+)-Cl(-) ion-pair separation to 16 A in 1.0M NaCl solution in water. The corresponding ggaHFB free-energy profile in six dimensional Cartesian space is in excellent agreement with the conventional benchmarks. We then explore changes in the PMF in response to lowering the NaCl concentration to physiological 0.3 and 0.1M, and dilute 0.0M concentrations. Finally, to expand the scope of the ggaHFB method, we formally develop the free-energy gradient approximation in arbitrary nonlinear coordinates. This formal development underscores the importance of the logarithmic Jacobian correction to reconstruct true PMFs from umbrella sampling simulations with either WHAM or ggaHFB techniques when nonlinear coordinate restraints are used with Cartesian propagators. The ability to employ nonlinear coordinates and high accuracy of the computed free-energy profiles further advocate the use of the ggaHFB method in studies of rare events in complex systems.

Intrinsic Free Energy of the Conformational Transition of the KcsA Signature Peptide from Conducting to Nonconducting State.

We explore a conformational transition of the TATTVGYG signature peptide of the KcsA ion selectivity filter and its GYG to AYA mutant from the conducting α-strand state into the nonconducting pII-like state using a novel technique for multidimensional optimization of transition path ensembles and free energy calculations. We find that the wild type peptide, unlike the mutant, intrinsically favors the conducting state due to G77 backbone propensities and additional hydrophobic interaction between the V76 and Y78 side chains in water. The molecular mechanical free energy profiles in explicit water are in very good agreement with the corresponding adiabatic energies from the Generalized Born Molecular Volume (GBMV) implicit solvent model. However comparisons of the energies to higher level B3LYP/6-31G(d) Density Functional Theory calculations with Polarizable Continuum Model (PCM) suggest that the nonconducting state might be more favorable than predicted by molecular mechanics simulations. By extrapolating the single peptide results to the tetrameric channel, we propose a novel hypothesis for the ion selectivity mechanism.

Generalized gradient-augmented harmonic Fourier beads method with multiple atomic and/or center-of-mass positional restraints.

We describe a generalization of the gradient-augmented harmonic Fourier beads method for finding minimum free-energy transition path ensembles and similarly minimum potential energy paths to allow positional restraints on the centers of mass of selected atoms. The generalized gradient-augmented harmonic Fourier beads (ggaHFB) method further extends the scope of the HFB methodology to studying molecule transport across various mobile phases such as lipid membranes. Furthermore, the new implementation improves the applicability of the HFB method to studies of ligand binding, protein folding, and enzyme catalysis as well as modeling equilibrium pulling experiments. Like its predecessor, the ggaHFB method provides accurate energy profiles along the specified paths and in certain simple cases avoids the need for path optimization. The utility of the ggaHFB method is demonstrated with an application to the water permeation through a single-wall (5,5) carbon nanotube with a diameter of 6.78 A and length of 16.0 A. We provide a simple rationale as to why water enters the hydrophobic nanotube and why it does so in pulses and in wire assembly.

Conformational change of the methionine 20 loop of Escherichia coli dihydrofolate reductase modulates pKa of the bound dihydrofolate.

We evaluate the pK(a) of dihydrofolate (H(2)F) at the N(5) position in three ternary complexes with Escherichia coli dihydrofolate reductase (ecDHFR), namely ecDHFR(NADP(+):H(2)F) in the closed form (1), and the Michaelis complexes ecDHFR(NADPH:H(2)F) in the closed (2) and occluded (3) forms, by performing free energy perturbation with molecular dynamics simulations (FEP/MD). Our simulations suggest that in the Michaelis complex the pK(a) is modulated by the Met20 loop fluctuations, providing the largest pK(a) shift in substates with a "tightly closed" loop conformation; in the "partially closed/open" substates, the pK(a) is similar to that in the occluded complex. Conducive to the protonation, tightly closing the Met20 loop enhances the interactions of the cofactor and the substrate with the Met20 side chain and aligns the nicotinamide ring of the cofactor coplanar with the pterin ring of the substrate. Overall, the present study favors the hypothesis that N(5) is protonated directly from solution and provides further insights into the mechanism of the substrate protonation.