Direct Calculation of Protein Fitness Landscapes through Computational Protein Design.

Naturally selected amino-acid sequences or experimentally derived ones are often the basis for understanding how protein three-dimensional conformation and function are determined by primary structure. Such sequences for a protein family comprise only a small fraction of all possible variants, however, representing the fitness landscape with limited scope. Explicitly sampling and characterizing alternative, unexplored protein sequences would directly identify fundamental reasons for sequence robustness (or variability), and we demonstrate that computational methods offer an efficient mechanism toward this end, on a large scale. The dead-end elimination and A(∗) search algorithms were used here to find all low-energy single mutant variants, and corresponding structures of a G-protein heterotrimer, to measure changes in structural stability and binding interactions to define a protein fitness landscape. We established consistency between these algorithms with known biophysical and evolutionary trends for amino-acid substitutions, and could thus recapitulate known protein side-chain interactions and predict novel ones.

Rational modification of protein stability by targeting surface sites leads to complicated results.

The rational modification of protein stability is an important goal of protein design. Protein surface electrostatic interactions are not evolutionarily optimized for stability and are an attractive target for the rational redesign of proteins. We show that surface charge mutants can exert stabilizing effects in distinct and unanticipated ways, including ones that are not predicted by existing methods, even when only solvent-exposed sites are targeted. Individual mutation of three solvent-exposed lysines in the villin headpiece subdomain significantly stabilizes the protein, but the mechanism of stabilization is very different in each case. One mutation destabilizes native-state electrostatic interactions but has a larger destabilizing effect on the denatured state, a second removes the desolvation penalty paid by the charged residue, whereas the third introduces unanticipated native-state interactions but does not alter electrostatics. Our results show that even seemingly intuitive mutations can exert their effects through unforeseen and complex interactions.

Biophysical and functional analyses suggest that adenovirus E4-ORF3 protein requires higher-order multimerization to function against promyelocytic leukemia protein nuclear bodies.

The early region 4 open reading frame 3 protein (E4-ORF3; UniProt ID P04489) is the most highly conserved of all adenovirus-encoded gene products at the amino acid level. A conserved attribute of the E4-ORF3 proteins of different human adenoviruses is the ability to disrupt PML nuclear bodies from their normally punctate appearance into heterogeneous filamentous structures. This E4-ORF3 activity correlates with the inhibition of PML-mediated antiviral activity. The mechanism of E4-ORF3-mediated reorganization of PML nuclear bodies is unknown. Biophysical analysis of the purified WT E4-ORF3 protein revealed an ordered secondary/tertiary structure and the ability to form heterogeneous higher-order multimers in solution. Importantly, a nonfunctional E4-ORF3 mutant protein, L103A, forms a stable dimer with WT secondary structure content. Because the L103A mutant is incapable of PML reorganization, this result suggests that higher-order multimerization of E4-ORF3 may be required for the activity of the protein. In support of this hypothesis, we demonstrate that the E4-ORF3 L103A mutant protein acts as a dominant-negative effector when coexpressed with the WT E4-ORF3 in mammalian cells. It prevents WT E4-ORF3-mediated PML track formation presumably by binding to the WT protein and inhibiting the formation of higher-order multimers. In vitro protein binding studies support this conclusion as demonstrated by copurification of coexpressed WT and L103A proteins in Escherichia coli and coimmunoprecipitation of WT·L103A E4-ORF3 complexes in mammalian cells. These results provide new insight into the properties of the Ad E4-ORF3 protein and suggest that higher-order protein multimerization is essential for E4-ORF3 activity.

Ionic strength effects on amyloid formation by amylin are a complicated interplay among Debye screening, ion selectivity, and Hofmeister effects.

Amyloid formation plays a role in a wide range of human diseases. The rate and extent of amyloid formation depend on solution conditions, including pH and ionic strength. Amyloid fibrils often adopt structures with parallel, in-register ß-sheets, which generate quasi-infinite arrays of aligned side chains. These arrangements can lead to significant electrostatic interactions between adjacent polypeptide chains. The effect of ionic strength and ion composition on the kinetics of amyloid formation by islet amyloid polypeptide (IAPP) is examined. IAPP is a basic 37-residue polypeptide responsible for islet amyloid formation in type 2 diabetes. Poisson-Boltzmann calculations revealed significant electrostatic repulsion in a model of the IAPP fibrillar state. The kinetics of IAPP amyloid formation are strongly dependent on ionic strength, varying by a factor of >10 over the range of 20-600 mM NaCl at pH 8.0, but the effect is not entirely due to Debye screening. At low ionic strengths, the rate depends strongly on the identity of the anion, varying by a factor of nearly 4, and scales with the electroselectivity series, implicating anion binding. At high ionic strengths, the rate varies by only 8% and scales with the Hofmeister series. At intermediate ionic strengths, no clear trend is detected, likely because of the convolution of different effects. The effects of salts on the growth phase and lag phase of IAPP amyloid formation are strongly correlated. At pH 5.5, where the net charge on IAPP is higher, the effect of different anions scales with the electroselectivity series at all salt concentrations.

Carbohydrate recognition by the antiviral lectin cyanovirin-N.

Cyanovirin-N (CVN) is a cyanobacterial lectin with potent antiviral activity and has been the focus of extensive preclinical investigation as a potential prophylactic for the prevention of the sexual transmission of the human immunodeficiency virus (HIV). Here we present a detailed analysis of carbohydrate recognition by this important protein, using a combination of computational methods, including extensive molecular dynamics simulations and molecular mechanics/Poisson-Boltzmann surface area (MM/PBSA) energetic analysis. The simulation results strongly suggest that the observed tendency of wild-type CVN to form domain-swapped dimers is the result of a previously unidentified cis-peptide bond present in the monomeric state. The energetic analysis additionally indicates that the highest-affinity ligand for CVN characterized to date (α-Man-(1,2)-α-Man-(1,2)-α-Man) is recognized asymmetrically by the two binding sites. Finally, we are able to provide a detailed map of the role of all binding site functional groups (both backbone and side chain) to various aspects of molecular recognition: general affinity for cognate ligands, specificity for distinct oligosaccharide targets, and the asymmetric recognition of α-Man-(1,2)-α-Man-(1,2)-α-Man. Taken as a whole, these results complement past experimental characterization (both structural and thermodynamic) to provide the most complete understanding of carbohydrate recognition by CVN to date. The results also provide strong support for the application of similar approaches to the understanding of other protein-carbohydrate complexes.

Rational and computational design of stabilized variants of cyanovirin-N that retain affinity and specificity for glycan ligands.

Cyanovirin-N (CVN) is an 11 kDa pseudosymmetric cyanobacterial lectin that has been shown to inhibit infection by the human immunodeficiency virus by binding to high-mannose oligosaccharides on the surface of the viral envelope glycoprotein gp120. In this work, we describe rationally designed CVN variants that stabilize the protein fold while maintaining high affinity and selectivity for their glycan targets. Poisson-Boltzmann calculations and protein repacking algorithms were used to select stabilizing mutations in the protein core. By substituting the buried polar side chains of Ser11, Ser20, and Thr61 with aliphatic groups, we stabilized CVN by nearly 12 °C against thermal denaturation, and by 1 M GuaHCl against chemical denaturation, relative to a previously characterized stabilized mutant. Glycan microarray binding experiments confirmed that the specificity profile of carbohydrate binding is unperturbed by the mutations and is identical for all variants. In particular, the variants selectively bound glycans containing the Manα(1→2)Man linkage, which is the known minimal binding unit of CVN. We also report the slow denaturation kinetics of CVN and show that they can complicate thermodynamic analysis; in particular, the unfolding of CVN cannot be described as a fixed two-state transition. Accurate thermodynamic parameters are needed to describe the complicated free energy landscape of CVN, and we provide updated values for CVN unfolding.

Interpretation of p-cyanophenylalanine fluorescence in proteins in terms of solvent exposure and contribution of side-chain quenchers: a combined fluorescence, IR and molecular dynamics study.

The use of noncoded amino acids as spectroscopic probes of protein folding and function is growing rapidly, in large part because of advances in the methodology for their incorporation. Recently p-cyanophenylalanine has been employed as a fluorescence and IR probe, as well as a FRET probe to study protein folding, protein-membrane interactions, protein-protein interactions and amyloid formation. The probe has been shown to be exquisitely sensitive to hydrogen bonding interactions involving the cyano group, and its fluorescence quantum yield increases dramatically when it is hydrogen bonded. However, a detailed understanding of the factors which influence its fluorescence is required to be able to use this popular probe accurately. Here we demonstrate the recombinant incorporation of p-cyanophenylalanine in the N-terminal domain of the ribosomal protein L9. Native state fluorescence is very low, which suggests that the group is sequestered from solvent; however, IR measurements and molecular dynamics simulations show that the cyano group is exposed to solvent and forms hydrogen bonds to water. Analysis of mutant proteins and model peptides demonstrates that the reduced native state fluorescence is caused by the effective quenching of p-cyanophenylalanine fluorescence via FRET to tyrosine side-chains. The implications for the interpretation of p-cyanophenylalanine fluorescence measurements and FRET studies are discussed.

Optimized parameters for continuum solvation calculations with carbohydrates.

Continuum solvent models have been shown to be an efficient method for the calculation of the energetics of biomolecules in solution. However, for these methods to produce accurate results, an appropriate set of atomic radii or volumes is needed. While these have been developed for proteins and nucleic acids, the same is not true of carbohydrates. Here, a set of optimized parameters for continuum solvation calculations of carbohydrates in conjunction with the Carbohydrate Solution Force Field are presented. Explicit solvent free-energy perturbation simulations were performed on a set of hexapyranose sugars and used to fit atomic radii for Poisson-Boltzmann and generalized-Born calculations, and to fit atomic volumes for use with the analytical continuum electrostatics model. The solvation energetics computed with the optimized radii and a Poisson-Boltzmann model show remarkable agreement with explicit solvent simulation, with a root-mean-square error of 1.19 kcal/mol over a large test set of sugars in many conformations. The generalized-Born model gives slightly poorer agreement, but still correlates very strongly, with an error of 1.69 kcal/mol. The analytical continuum electrostatics model correlates well with the explicit solvent results, but gives a larger error of 4.71 kcal/mol. The remarkable agreement between the solvation free energies computed in explicit and implicit solvent provides strong motivation for the use of implicit solvent models in the simulation of carbohydrate-containing systems.

Action-at-a-distance interactions enhance protein binding affinity.

The identification of protein mutations that enhance binding affinity may be achieved by computational or experimental means, or by a combination of the two. Sources of affinity enhancement may include improvements to the net balance of binding interactions of residues forming intermolecular contacts at the binding interface, such as packing and hydrogen-bonding interactions. Here we identify noncontacting residues that make substantial contributions to binding affinity and that also provide opportunities for mutations that increase binding affinity of the TEM1 beta-lactamase (TEM1) to the beta-lactamase inhibitor protein (BLIP). A region of BLIP not on the direct TEM1-binding surface was identified for which changes in net charge result in particularly large increases in computed binding affinity. Some mutations to the region have previously been characterized, and our results are in good correspondence with this results of that study. In addition, we propose novel mutations to BLIP that were computed to improve binding significantly without contacting TEM1 directly. This class of noncontacting electrostatic interactions could have general utility in the design and tuning of binding interactions.

Design of improved protein inhibitors of HIV-1 cell entry: Optimization of electrostatic interactions at the binding interface.

Continuum electrostatic methods are a powerful tool for the analysis and design of biomolecular complexes, with methodologies that allow for the detailed analysis of the electrostatic contributions to binding affinities and procedures for computing the properties of electrostatically optimal ligands. We have applied these methods to the design of improved inhibitors of HIV-1 cell entry. HIV infection of a cell requires viral-cell membrane fusion, an event partially driven by a large-scale conformational change in the viral membrane glycoprotein gp41. This transformation involves the docking of a helix from the C-terminal region of three gp41 chains against a pre-formed trimeric-coiled coil; several protein constructs that inhibit membrane fusion act by binding to an isolated C-terminal helix and blocking the formation of the fusogenic structure. A detailed analysis of the electrostatic contributions to the binding of one such inhibitor (5-Helix) to a C-terminal helix was performed using the X-ray crystal structure of the core of the HIV-1 gp41 ectodomain as a structural model, and several residues on 5-Helix that make substantial contributions to binding, both favorable and unfavorable, were identified. An electrostatic affinity optimization methodology was applied to the side chains of 5-Helix, with the results showing that significant improvements in binding affinity are possible if the electrostatic contributions to the binding free energy are optimized. Several mutations accessible by experimental methods are suggested, with calculated improvements in binding affinity of as much as 500-fold and greater.

Escherichia coli glutaminyl-tRNA synthetase is electrostatically optimized for binding of its cognate substrates.

Natural evolution has resulted in protein molecules displaying a wide range of binding properties that include extremes of affinity and specificity. A detailed understanding of the principles underlying protein structure-function relationships, particularly with respect to binding properties, would greatly enhance molecular engineering and ligand design studies. Here, we have analyzed the interactions of an aminoacyl-tRNA synthetase for which strong evolutionary pressure has enforced high specificity for substrate binding and catalysis. Electrostatic interactions have been identified as one efficient mechanism for enhancing binding specificity; as such, the effects of charged and polar groups were the focus of this study. The binding of glutaminyl-tRNA synthetase from Escherichia coli to several ligands, including the natural substrates, was analyzed. The electrostatic complementarity of the enzyme to its ligands was assessed using measures derived from affinity optimization theory. The results were independent of the details of the calculational parameters, including the value used for the protein dielectric constant. Glutamine and ATP, two of the natural ligands, were found to be extremely complementary to their binding sites, particularly in regions seen to make electrostatic interactions in the structure. These data suggest that the optimization of electrostatic interactions has played an important role in guiding the evolution of this enzyme. The results also show that the enzyme is able to effectively select for high affinity and specificity for the same chemical moieties both in the context of smaller substrates, and in that of a larger reactive intermediate. The regions of greatest non-complementarity between the enzyme and ligands are the portions of the ligand that make few polar contacts with the binding site, as well as the sites of chemical reaction, where overly strong electrostatic binding interactions with the substrate could hinder catalysis. The results also suggest that the negative charge on the phosphorus center of glutaminyl-adenylate plays an important role in the tight binding of this intermediate, and thus that adenylate analogs that preserve the negative charge in this region may bind substantially tighter than analogs where this group is replaced with a neutral group, such as the sulfamoyl family, which can make similar hydrogen bonds but is uncharged.

X-ray structural and simulation analysis of a protein mutant: the value of a combined approach.

The effect of the mutation Arg 96 to His on the stability of bacteriophage T4 lysozyme has been previously studied by calorimetric experiments, X-ray crystallography, and free energy simulation techniques. The experimental and calculated values for the difference between the free energy of denaturation of the mutant and the wild type are in reasonable agreement. However, the two approaches led to different explanations for the loss in stability. To analyze the differences, a series of refinements based on the crystallographic data were performed, a number of aspects of the simulations were reexamined, and continuum electrostatic calculations were done to complement the latter. The results of those comparisons provide a better understanding of the origin of the free energy difference in this mutant. Furthermore, they show the importance of the combined use of simulations and crystallography for interpreting the effects of mutations on the energetics of the system.

A Statistical Framework for Hierarchical Methods in Molecular Simulation and Design.

A statistical framework for performance analysis in hierarchical methods is described, with a focus on applications in molecular design. A theory is derived from statistical principles, describing the relationships between the results of each hierarchical level by a functional correlation and an error model for how values are distributed around the correlation curve. Two key measures are then defined for evaluating a hierarchical approach-completeness and excess cost-conceptually similar to the sensitivity and specificity of dichotomous prediction methods. We demonstrate the use of this method using a simple model problem in conformational search, refining the results of an in vacuo search of glucose conformations with a continuum solvent model. Second, we show the usefulness of this approach when structural hierarchies are used to efficiently make use of large rotamer libraries with the Dead-end Elimination and A* algorithms for protein design. The framework described is applicable not only to the specific examples given but to any problem in molecular simulation or design that involves a hierarchical approach.

Energetic decomposition with the generalized-born and Poisson-Boltzmann solvent models: lessons from association of G-protein components.

Continuum electrostatic models have been shown to be powerful tools in providing insight into the energetics of biomolecular processes. While the Poisson-Boltzmann (PB) equation provides a theoretically rigorous approach to computing electrostatic free energies of solution in such a model, computational cost makes its use for large ensembles of states impractical. The generalized-Born (GB) approximation provides a much faster alternative, although with a weaker theoretical framework. While much attention has been given to how GB recapitulates PB energetics for the overall stability of a biomolecule or the affinity of a complex, little attention has been given to how the contributions of individual functional groups are captured by the two methods. Accurately capturing these individual electrostatic components is essential both for the development of a mechanistic understanding of biomolecular processes and for the design of variant sequences and structures with desired properties. Here, we present a detailed comparison of the group-wise decomposition of both PB and GB electrostatic free energies of binding, using association of various components of the heterotrimeric-G-protein complex as a model. We find that, while net binding free energies are strongly correlated in the two models, the correlations of individual group contributions are highly variable; in some cases, strong correlation is seen, while in others, there is essentially none. Structurally, the GB model seems to capture the magnitude of direct, short-range electrostatic interactions quite well but performs more poorly with moderate-range "action-at-a-distance" interactions--GB has a tendency to overestimate solvent screening over moderate distances, and to underestimate the costs of desolvating charged groups somewhat removed from the binding interface. Despite this, however, GB does seem to be quite effective as a predictor of those groups that will be computed to be most significant in a PB-based model.

Computational models explain the oligosaccharide specificity of cyanovirin-N.

The prokaryotic lectin cyanovirin-N (CV-N) is a potent inhibitor of HIV envelope-mediated cell entry, and thus is a leading candidate among a new class of potential anti-HIV microbicides. The activity of CV-N is a result of interactions with the D1 arm of high-mannose oligosaccharides on the viral glycoprotein gp120. Here, we present computationally refined models of CV-N recognition of the di- and trisaccharides that represent the terminal three sugars of the D1 arm by each CV-N binding site. These models complement existing structural data, both from NMR spectroscopy and X-ray crystallography. When used with a molecular dynamics/continuum electrostatic (MD/PBSA) approach to compute binding free energies, these models explain the relative affinity of each site for the two saccharides. This work presents the first validation of the application of continuum electrostatic models to carbohydrate-protein association. Taken as a whole, the results both provide models of CV-N sugar recognition and demonstrate the utility of these computational methods for the study of carbohydrate-binding proteins.

Rational design of new binding specificity by simultaneous mutagenesis of calmodulin and a target peptide.

Calcium-saturated calmodulin (CaM) binds and influences the activity of a varied collection of target proteins in most cells. This promiscuity underlies the role of CaM as a shared participant in calcium-dependent signal transduction pathways but imposes a handicap on popular CaM-based calcium biosensors, which display an undesired tendency to cross-react with cellular proteins. Designed CaM/target pairs that retain high affinity for one another but lack affinity for wild-type CaM and its natural interaction partners would therefore be useful as sensor components and possibly also as elements of "synthetic" cellular-signaling networks. Here, we have adopted a rational approach to creating suitably modified CaM/target complexes by using computational design methods to guide parallel site-directed mutagenesis of both binding partners. A hierarchical design procedure was applied to suggest a small number of complementary mutations on CaM and on a peptide ligand derived from skeletal-muscle light-chain kinase (M13). Experimental analysis showed that the procedure was successful in identifying CaM and M13 mutants with novel specificity for one another. Importantly, the designed complexes retained an affinity comparable to the wild-type CaM/M13 complex. These results represent a step toward the creation of CaM and M13 derivatives with specificity fully orthogonal to the wild-type proteins and show that qualitatively accurate predictions may be obtained from computational methods applied simultaneously to two proteins involved in multiple-linked binding equilibria.

Evaluation of electrostatic interactions.

Described here are several computational procedures for the analysis of electrostatic interactions in molecular complexes, all based on a continuum model of solvation. The first section describes how to compute the residual potential, a measure of how electrostatically complementary a ligand is for its receptor. The second procedure describes electrostatic component analysis, a method by which the electrostatic contribution to the binding free energy can be broken up into terms directly attributable to individual chemical groups. Finally, electrostatic affinity optimization is described. This procedure is particularly useful in determining what portions of a ligand are the most suboptimal, and thus provide the greatest opportunity for the design of improvements.