8-Oxoguanine (8-oxoG), a mutagenic DNA lesion generated under oxidative stress, differs from its precursor guanine by only two substitutions (O8 and H7). Human 8-oxoguanine glycosylase 1 (OGG1) can locate and remove 8-oxoG through extrusion and excision. To date, it remains unclear how OGG1 efficiently distinguishes 8-oxoG from a large excess of undamaged DNA bases. We recently showed that formamidopyrimidine–DNA glycosylase (Fpg), a bacterial functional analog of OGG1, can selectively facilitate eversion of oxoG by stabilizing several intermediate states, and it is intriguing whether OGG1 also employs a similar mechanism in lesion recognition. Here, we use molecular dynamics simulations to explore the mechanism by which OGG1 discriminates between 8-oxoG and guanine along the base-eversion pathway. The MD results suggest an important role for kinking of the DNA by the glycosylase, which positions DNA phosphates in a way that assists lesion recognition during base eversion. The computational predictions were validated through experimental enzyme assays on phosphorothioate substrate analogs. Our simulations suggest that OGG1 distinguishes between 8-oxoG and G using their chemical dissimilarities not only at the active site but also at earlier stages during base eversion, and this mechanism is at least partially conserved in Fpg despite a lack of structural homology. The similarity also suggests that lesion recognition through multiple gating steps may be a common theme in DNA repair. Our results provide new insight into how enzymes can exploit kinetics and DNA conformational changes to probe the chemical modifications present in DNA lesions.

A general method is presented to characterize the helical properties of potentially irregular helices, such as those found in protein secondary and tertiary structures and nucleic acids. The method was validated using artificial helices with varying numbers of points, points per helical turn, pitch, and radius. The sensitivity of the method was validated by applying increasing amounts of random perturbation to the coordinates of these helices; 399 360 helices in total were evaluated. In addition, the helical parameters of protein secondary structure elements and nucleic acid helices were analyzed. Generally, at least seven points were required to recapitulate the parameters of a helix using our method. The method can also be used to calculate the helical parameters of nucleic acid-binding proteins, like TALE, enabling direct analysis of their helix complementarity to sequence-dependent DNA distortions.

The rational and predictable enhancement of protein stability is an important goal in protein design. Most efforts target the folded state, however stability is the free energy difference between the folded and unfolded states thus both are suitable targets. Strategies directed at the unfolded state usually seek to decrease chain entropy by introducing cross-links or by replacing glycines. Cross-linking has led to mixed results. Replacement of glycine with an l-amino acid, while reducing the entropy of the unfolded state, can introduce unfavorable steric interactions in the folded state, since glycine is often found in conformations that require a positive φ angle such as helical C-capping motifs or type I′ and II″ β-turns. l-Amino acids are strongly disfavored in these conformations, but d-amino acids are not. However, there are few reported examples and conflicting results have been obtained when glycines are replaced with d-Ala. We critically examine the effect of Gly-to-d-Ala substitutions on protein stability using experimental approaches together with molecular dynamics simulations and free energy calculations. The data, together with a survey of high resolution structures, show that the vast majority of proteins can be stabilized by substitution of C-capping glycines with d-Ala. Sites suitable for substitutions can be identified via sequence alignment with a high degree of success. Steric clashes in the native state due to the new side chain are rarely observed, but are likely responsible for the destabilizing or null effect observed for the small subset of Gly-to-d-Ala substitutions which are not stabilizing. Changes in backbone solvation play less of a role. Favorable candidates for d-Ala substitution can be identified using a rapid algorithm based on molecular mechanics.

The enoyl-ACP reductase (ENR) catalyzes the last reaction in the elongation cycle of the bacterial type II fatty acid biosynthesis (FAS-II) pathway. While the FabI ENR is a well-validated drug target in organisms such as Mycobacterium tuberculosis and Staphylococcus aureus, alternate ENR isoforms have been discovered in other pathogens, including the FabV enzyme that is the sole ENR in Yersinia pestis (ypFabV). Previously, we showed that the prototypical ENR inhibitor triclosan was a poor inhibitor of ypFabV and that inhibitors based on the 2-pyridone scaffold were more potent [Hirschbeck, M. (2012) Structure 20 (1), 89–100]. These studies were performed with the T276S FabV variant. In the work presented here, we describe a detailed examination of the mechanism and inhibition of wild-type ypFabV and the T276S variant. The T276S mutation significantly reduces the affinity of diphenyl ether inhibitors for ypFabV (20-fold → 100-fold). In addition, while T276S ypFabV generally displays an affinity for 2-pyridone inhibitors higher than that of the wild-type enzyme, the 4-pyridone scaffold yields compounds with similar affinity for both wild-type and T276S ypFabV. T276 is located at the N-terminus of the helical substrate-binding loop, and structural studies coupled with site-directed mutagenesis reveal that alterations in this residue modulate the size of the active site portal. Subsequently, we were able to probe the mechanism of time-dependent inhibition in this enzyme family by extending the inhibition studies to include P142W ypFabV, a mutation that results in a gain of slow-onset inhibition for the 4-pyridone PT156.

The Generalized Born (GB) implicit solvent model has undergone significant improvements in accuracy for modeling of proteins and small molecules. However, GB still remains a less widely explored option for nucleic acid simulations, in part because fast GB models are often unable to maintain stable nucleic acid structures or they introduce structural bias in proteins, leading to difficulty in application of GB models in simulations of protein–nucleic acid complexes. Recently, GB-neck2 was developed to improve the behavior of protein simulations. In an effort to create a more accurate model for nucleic acids, a similar procedure to the development of GB-neck2 is described here for nucleic acids. The resulting parameter set significantly reduces absolute and relative energy error relative to Poisson–Boltzmann for both nucleic acids and nucleic acid–protein complexes, when compared to its predecessor GB-neck model. This improvement in solvation energy calculation translates to increased structural stability for simulations of DNA and RNA duplexes, quadruplexes, and protein–nucleic acid complexes. The GB-neck2 model also enables successful folding of small DNA and RNA hairpins to near native structures as determined from comparison with experiment. The functional form and all required parameters are provided here and also implemented in the AMBER software.

Molecular mechanics is powerful for its speed in atomistic simulations, but an accurate force field is required. The Amber ff99SB force field improved protein secondary structure balance and dynamics from earlier force fields like ff99, but weaknesses in side chain rotamer and backbone secondary structure preferences have been identified. Here, we performed a complete refit of all amino acid side chain dihedral parameters, which had been carried over from ff94. The training set of conformations included multidimensional dihedral scans designed to improve transferability of the parameters. Improvement in all amino acids was obtained as compared to ff99SB. Parameters were also generated for alternate protonation states of ionizable side chains. Average errors in relative energies of pairs of conformations were under 1.0 kcal/mol as compared to QM, reduced 35% from ff99SB. We also took the opportunity to make empirical adjustments to the protein backbone dihedral parameters as compared to ff99SB. Multiple small adjustments of φ and ψ parameters were tested against NMR scalar coupling data and secondary structure content for short peptides. The best results were obtained from a physically motivated adjustment to the φ rotational profile that compensates for lack of ff99SB QM training data in the β-ppII transition region. Together, these backbone and side chain modifications (hereafter called ff14SB) not only better reproduced their benchmarks, but also improved secondary structure content in small peptides and reproduction of NMR χ1 scalar coupling measurements for proteins in solution. We also discuss the Amber ff12SB parameter set, a preliminary version of ff14SB that includes most of its improvements.

Slow-onset enzyme inhibitors are the subject of considerable interest as an approach to increasing the potency of pharmaceutical compounds by extending the residence time of the inhibitor on the target (the lifetime of the drug–receptor complex). However, rational modulation of residence time presents significant challenges because it requires additional mechanistic insight, such as the nature of the transition state for postbinding isomerization. Our previous work, based on X-ray crystallography, enzyme kinetics, and molecular dynamics simulation, suggested that the slow step in inhibition of the Mycobacterium tuberculosis enoyl-ACP reductase InhA involves a change in the conformation of the substrate binding loop from an open state in the initial enzyme–inhibitor complex to a closed state in the final enzyme–inhibitor complex. Here, we use multidimensional free energy landscapes for loop isomerization to obtain a computational model for the transition state. The results suggest that slow-onset inhibitors crowd key side chains on helices that slide past each other during isomerization, resulting in a steric clash. The landscapes become significantly flatter when residues involved in the steric clash are replaced with alanine. Importantly, this lower barrier can be increased by rational inhibitor redesign to restore the steric clash. Crystallographic studies and enzyme kinetics confirm the predicted effects on loop structure and flexibility, as well as inhibitor residence time. These loss and regain of function studies validate our mechanistic hypothesis for interactions controlling substrate binding loop isomerization, providing a platform for the future design of inhibitors with longer residence times and better in vivo potency. Similar opportunities for slow-onset inhibition via the same mechanism are identified in other pathogens.

The millisecond time scale needed for molecular dynamics simulations to approach the quantitative study of protein folding is not yet routine. One approach to extend the simulation time scale is to perform long simulations on specialized and expensive supercomputers such as Anton. Ideally, however, folding simulations would be more economical while retaining reasonable accuracy, and provide feedback on structure, stability and function rapidly enough if partnered directly with experiment. Approaches to this problem typically involve varied compromises between accuracy, precision, and cost; the goal here is to address whether simple implicit solvent models have become sufficiently accurate for their weaknesses to be offset by their ability to rapidly provide much more precise conformational data as compared to explicit solvent. We demonstrate that our recently developed physics-based model performs well on this challenge, enabling accurate all-atom simulated folding for 16 of 17 proteins with a variety of sizes, secondary structure, and topologies. The simulations were carried out using the Amber software on inexpensive GPUs, providing ∼1 μs/day per GPU, and >2.5 ms data presented here. We also show that native conformations are preferred over misfolded structures for 14 of the 17 proteins. For the other 3, misfolded structures are thermodynamically preferred, suggesting opportunities for further improvement.

Slow-onset enzyme inhibitors are of great interest for drug discovery programs since the slow dissociation of the inhibitor from the drug–target complex results in sustained target occupancy leading to improved pharmacodynamics. However, the structural basis for slow-onset inhibition is often not fully understood, hindering the development of structure-kinetic relationships and the rational optimization of drug-target residence time. Previously we demonstrated that slow-onset inhibition of the Mycobacterium tuberculosis enoyl-ACP reductase InhA correlated with motions of a substrate-binding loop (SBL) near the active site. In the present work, X-ray crystallography and molecular dynamics simulations have been used to map the structural and energetic changes of the SBL that occur upon enzyme inhibition. Helix-6 within the SBL adopts an open conformation when the inhibitor structure or binding kinetics is substrate-like. In contrast, slow-onset inhibition results in large-scale local refolding in which helix-6 adopts a closed conformation not normally populated during substrate turnover. The open and closed conformations of helix-6 are hypothesized to represent the EI and EI* states on the two-step induced-fit reaction coordinate for enzyme inhibition. These two states were used as the end points for nudged elastic band molecular dynamics simulations resulting in two-dimensional potential energy profiles that reveal the barrier between EI and EI*, thus rationalizing the binding kinetics observed with different inhibitors. Our findings indicate that the structural basis for slow-onset kinetics can be understood once the structures of both EI and EI* have been identified, thus providing a starting point for the rational control of enzyme–inhibitor binding kinetics.

The generalized Born (GB) model is one of the fastest implicit solvent models, and it has become widely adopted for Molecular Dynamics (MD) simulations. This speed comes with trade-offs, and many reports in the literature have pointed out weaknesses with GB models. Because the quality of a GB model is heavily affected by empirical parameters used in calculating solvation energy, in this work we have refit these parameters for GB-Neck, a recently developed GB model, in order to improve the accuracy of both the solvation energy and effective radii calculations. The data sets used for fitting are significantly larger than those used in the past. Comparing to other pairwise GB models like GB-OBC and the original GB-Neck, the new GB model (GB-Neck2) has better agreement with Poisson–Boltzmann (PB) in terms of reproducing solvation energies for a variety of systems ranging from peptides to proteins. Secondary structure preferences are also in much better agreement with those obtained from explicit solvent MD simulations. We also obtain near-quantitative reproduction of experimental structure and thermal stability profiles for several model peptides with varying secondary structure motifs. Extension to nonprotein systems will be explored in the future.

Base eversion is a fundamental process in the biochemistry of nucleic acids, allowing proteins engaged in DNA repair and epigenetic modifications to access target bases in DNA. Crystal structures reveal end points of these processes, but not the pathways involved in the dynamic process of base recognition. To elucidate the pathway taken by 8-oxoguanine during base excision repair by Fpg, we calculated free energy surfaces during eversion of the damaged base through the major and minor grooves. The minor groove pathway and free energy barrier (6–7 kcal/mol) are consistent with previously reported results (Qi, Y.; Spong, M. C.; Nam, K.; Banerjee, A.; Jiralerspong, S.; Karplus, M.; Verdine, G. L. Nature2009, 462, 762.) However, eversion of 8-oxoG through the major groove encounters a significantly lower barrier (3–4 kcal/mol) more consistent with experimentally determined rates of enzymatic sliding during lesion search (Blainey, P. C.; van Oijent, A. M.; Banerjee, A.; Verdine, G. L.; Xie, X. S. Proc. Natl. Acad. Sci. U.S.A.2006, 103, 5752.). Major groove eversion has been suggested for other glycosylases, suggesting that in addition to function, dynamics of base eversion may also be conserved.

Base flipping is a common strategy utilized by many enzymes to gain access to the functional groups of nucleic acid bases in duplex DNA which are otherwise protected by the DNA backbone and hydrogen bonding with their partner bases. Several X-ray crystallography studies have revealed flipped conformations of nucleotides bound to enzymes. However, little is known about the base-flipping process itself, even less about the role of the enzymes. Computational studies have used umbrella sampling to elicit the free energy profile of the base-flipping process using a pseudodihedral angle to represent the reaction coordinate. In this study, we have used an unrestrained trajectory in which a flipped base spontaneously reinserted into the helix in order to evaluate and improve the previously defined pseudodihedral angle. Our modified pseudodihedral angles use a new atom selection to improve the numerical stability of the restraints and also provide better correlation to the extent of flipping observed in simulations. Furthermore, on the basis of the comparison of potential of mean force (PMF) generated using different reaction coordinates, we observed that the shape of a flipping PMF profile is strongly dependent on the definition of the reaction coordinate, even for the same data set.

We report unrestrained, all-atom molecular dynamics simulations of HIV-1 protease that sample large conformational changes of the active site flaps. In particular, the unliganded protease undergoes multiple conversions between the “closed” and “semiopen” forms observed in crystal structures of inhibitor-bound and unliganded protease, respectively, including reversal of flap “handedness.” Simulations in the presence of a cyclic urea inhibitor yield stable closed flaps. Furthermore, we observe several events in which the flaps of the unliganded protease open to a much greater degree than observed in crystal structures and subsequently return to the semiopen state. Our data strongly support the hypothesis that the unliganded protease predominantly populates the semiopen conformation, with closed and fully open structures being a minor component of the overall ensemble. The results also provide a model for the flap opening and closing that is considered to be essential to enzyme function.

HIV-1 protease remains an important anti-AIDS drug target. Although it has been known that ligand binding induces large conformational changes in the protease, the dynamic aspects of binding have been largely ignored. Several computational models describing protease dynamics have been reported recently. These have reproduced experimental observations, and have also explained how ligands gain access to the binding site through dynamic behavior of the protease. Specifically, the transitions between three different conformations of the protein have been modeled in atomic detail. Two of these forms were determined by crystallography, and the third was implied by NMR experiments. Based on these computational models, it has been suggested that binding of inhibitors in allosteric sites might affect protease flexibility and disrupt its function.

Force-field validation is essential for the identification of weaknesses in current models and the development of more accurate models of biomolecules. NMR coupling and relaxation methods have been used to effectively diagnose the strengths and weaknesses of many existing force fields. Studies using the ff99SB force field have shown excellent agreement between experimental and calculated order parameters and residual dipolar calculations. However, recent studies have suggested that ff99SB demonstrates poor agreement with J-coupling constants for short polyalanines. We performed extensive replica-exchange molecular-dynamics simulations on Ala3 and Ala5 in TIP3P and TIP4P-Ew solvent models. Our results suggest that the performance of ff99SB is among the best of currently available models. In addition, scalar coupling constants derived from simulations in the TIP4P-Ew model show a slight improvement over those obtained using the TIP3P model. Despite the overall excellent agreement, the data suggest areas for possible improvement.