S. Rackovsky

A united-residue force field for off-lattice protein-structure simulations. I. Functional forms and parameters of long-range side-chain interaction potentials from protein crystal data

A two‐stage procedure for the determination of a united‐residue potential designed for protein simulations is outlined. In the first stage, the long‐range and local‐interaction energy terms of the total energy of a polypeptide chain are determined by analyzing protein‐crystal data and averaging the all‐atom energy surfaces. In the second stage (described in the accompanying article), the relative weights of the energy terms are optimized so as to locate the native structures of selected test proteins as the lowest energy structures. The goal of the work in the present study is to parameterize physically reasonable functional forms of the potentials of mean force for side‐chain interactions. The potentials are of both radial and anisotropic type. Radial potentials include the Lennard‐Jones and the shifted Lennard‐Jones potential (with the shift parameter independent of orientation). To treat the angular dependence of side‐chain interactions, three functional forms of the potential that were designed previously to describe anisotropic systems are evaluated: Berne‐Pechukas (dilated Lennard‐Jones); Gay‐Berne (shifted Lennard‐Jones with orientation‐dependent shift parameters); and Gay‐Berne‐Vorobjev (the same as the preceding one, but with one more set of variable parameters). These functional forms were used to parameterize, within a short‐distance range, the potentials of mean force for side‐chain pair interactions that are related by the Boltzmann principle to the pair correlation functions determined from protein‐crystal data. Parameter determination was formulated as a generalized nonlinear least‐squares problem with the target function being the weighted sum of squares of the differences between calculated and “experimental” (i.e., estimated from protein‐crystal data) angular, radial‐angular, and radial pair correlation functions, as well as contact free energies. A set of 195 high‐resolution nonhomologous structures from the Protein Data Bank was used to calculate the “experimental” values. The contact free energies were scaled by the slope of the correlation line between side‐chain hydrophobicities, calculated from the contact free energies, and those determined by Fauchere and Pliška from the partition coefficients of amino acids between water and n‐octanol. The methylene group served to define the reference contact free energy corresponding to that between the glycine methylene groups of backbone residues. Statistical analysis of the goodness of fit revealed that the Gay‐Berne‐Vorobjev anisotropic potential fits best to the experimental radial and angular correlation functions and contact free energies and therefore represents the free‐energy surface of side‐chain‐side‐chain interactions most accurately. Thus, its choice for simulations of protein structure is probably the most appropriate. However, the use of simpler functional forms is recommended, if the speed of computations is an issue.

A united-residue force field for off-lattice protein-structure simulations. II. Parameterization of short-range interactions and determination of weights of energy terms by Z-score optimization

Continuing our work on the determination of an off‐lattice united‐residue force field for protein‐structure simulations, we determined and parameterized appropriate functional forms for the local‐interaction terms, corresponding to the rotation about the virtual bonds (Utor), the bending of virtual‐bond angles (Ub), and the energy of different rotameric states of side chains (Urot). These terms were determined by applying the Boltzmann principle to the distributions of virtual‐bond torsional and virtual‐bond angles and side‐chain rotameric states, respectively, calculated from a data base of 195 high‐resolution nonhomologous proteins. The complete energy function was constructed by combining the individual energy terms with appropriate weights. The weights were determined by optimizing the so‐called Z‐score value (which is the normalized difference between the energy of the native structure and the mean energy of non‐native structures) of the histidine‐containing phosphocarrier protein from Streptococcus faecalis (1PTF; an 88‐residue α + β protein). To accomplish this, a database of Cα patterns was created using high‐resolution nonhomologous protein structures from the Protein Data Bank, and the distributions of energy components of 1PTF were obtained by threading its sequence through ∼500 randomly chosen Cα‐patterns from the X‐ray structures in the PDB, followed by energy minimization, with the energy function incorporating initially guessed weights. The resulting minimized energies were used to optimize the Z‐score value of 1PTF as a function of the weights of the various energy terms, and the new weights were used to generate new energy‐component distributions. The process was iterated, until the weights used to generate the distributions and the optimized weights were self‐consistent. The potential function with the weights of the various energy terms obtained by optimizing the Z‐score value for 1PTF was found to locate the native structures of other test proteins (within an average RMS deviation of 3 Å): calcium‐binding protein (4ICB), ubiquitin (1UBQ), α‐spectrin (1SHG), major cold‐shock protein (1MJC), and cytochrome b5 (3B5C) (which included α and β structures) as distinctively lowest in energy in similar threading experiments.

United-residue force field for off-lattice protein-structure simulations: III. Origin of backbone hydrogen-bonding cooperativity in united-residue potentials

Based on the dipole model of peptide groups developed in our earlier work [Liwo et al., Prot. Sci., 2, 1697 (1993)], a cumulant expansion of the average free energy of the system of freely rotating peptide‐group dipoles tethered to a fixed α‐carbon trace is derived. A graphical approach is presented to find all nonvanishing terms in the cumulants. In particular, analytical expressions for three‐ and four‐body (correlation) terms in the averaged interaction potential of united peptide groups are derived. These expressions are similar to the cooperative forces in hydrogen bonding introduced by Koliński and Skolnick [J. Chem. Phys., 97, 9412 (1992)]. The cooperativity arises here naturally from the higher order terms in the power‐series expansion (in the inverse of the temperature) for the average energy. Test calculations have shown that addition of the derived four‐body term to the statistical united‐residue potential of our earlier work [Liwo et al., J. Comput. Chem., 18, 849, 874 (1997)] greatly improves its performance in folding poly‐l‐alanine into an α‐helix. © 1998 John Wiley & Sons, Inc. J Comput Chem 19: 259–276, 1998

On the properties and sequence context of structurally ambivalent fragments in proteins.

The goal of this work is to characterize structurally ambivalent fragments in proteins. We have searched the Protein Data Bank and identified all structurally ambivalent peptides (SAPs) of length five or greater that exist in two different backbone conformations. The SAPs were classified in five distinct categories based on their structure. We propose a novel index that provides a quantitative measure of conformational variability of a sequence fragment. It measures the context‐dependent width of the distribution of (ϕ,ξ) dihedral angles associated with each amino acid type. This index was used to analyze the local structural propensity of both SAPs and the sequence fragments contiguous to them. We also analyzed type‐specific amino acid composition, solvent accessibility, and overall structural properties of SAPs and their sequence context. We show that each type of SAP has an unusual, type‐specific amino acid composition and, as a result, simultaneous intrinsic preferences for two distinct types of backbone conformation. All types of SAPs have lower sequence complexity than average. Fragments that adopt helical conformation in one protein and sheet conformation in another have the lowest sequence complexity and are sampled from a relatively limited repertoire of possible residue combinations. A statistically significant difference between two distinct conformations of the same SAP is observed not only in the overall structural properties of proteins harboring the SAP but also in the properties of its flanking regions and in the pattern of solvent accessibility. These results have implications for protein design and structure prediction.

Improvement of the Treatment of Loop Structures in the UNRES Force Field by Inclusion of Coupling between Backbone- and Side-Chain-Local Conformational States

The UNited RESidue (UNRES) coarse-grained model of polypeptide chains, developed in our laboratory, enables us to carry out millisecond-scale molecular-dynamics simulations of large proteins effectively. It performs well in ab initio predictions of protein structure, as demonstrated in the last Community Wide Experiment on the Critical Assessment of Techniques for Protein Structure Prediction (CASP10). However, the resolution of the simulated structure is too coarse, especially in loop regions, which results from insufficient specificity of the model of local interactions. To improve the representation of local interactions, in this work we introduced new side-chain-backbone correlation potentials, derived from a statistical analysis of loop regions of 4585 proteins. To obtain sufficient statistics, we reduced the set of amino-acid-residue types to five groups, derived in our earlier work on structurally optimized reduced alphabets, based on a statistical analysis of the properties of amino-acid structures. The new correlation potentials are expressed as one-dimensional Fourier series in the virtual-bond-dihedral angles involving side-chain centroids. The weight of these new terms was determined by a trial-and-error method, in which Multiplexed Replica Exchange Molecular Dynamics (MREMD) simulations were run on selected test proteins. The best average root-mean-square deviations (RMSDs) of the calculated structures from the experimental structures below the folding-transition temperatures were obtained with the weight of the new side-chain-backbone correlation potentials equal to 0.57. The resulting conformational ensembles were analyzed in detail by using the Weighted Histogram Analysis Method (WHAM) and Wards minimum-variance clustering. This analysis showed that the RMSDs from the experimental structures dropped by 0.5 Å on average, compared to simulations without the new terms, and the deviation of individual residues in the loop region of the computed structures from their counterparts in the experimental structures (after optimum superposition of the calculated and experimental structure) decreased by up to 8 Å. Consequently, the new terms improve the representation of local structure.

Improvement of statistical potentials and threading score functions using information maximization

We show that statistical potentials and threading score functions, derived from finite data sets, are informatic functions, and that their performance depends on the manner in which data are classified and compressed. The choice of sequence and structural parameters affects estimates of the conditional probabilities P(C|S), the quantification of the effect of sequence S on conformation C, and determines the amount of information extracted from the data set, as measured by information gain. The mathematical link between information gain and mean conformational energy, established in this work using the local backbone potential as model, demonstrates that manipulation of descriptive parameters also alters the “energy” values assigned to native conformation and to decoy structures in the test pool, and consequently, the performance of such statistical potential functions in fold recognition exercises. We show that sequence and structural partitions that maximize information gain also minimize the mean energy of the ensemble of native conformations. Moreover, we establish an informatic basis for the placement of the native score within an energy spectrum given by the decoy pool in a threading exercise. We discover that, among all informatic quantities, information gain is the best predictor of threading success, even better than the standard Z‐score. Consequently, the choices of sequence and structural descriptors, extent of compression, and levels of discretization that maximize information gain must also produce the best potential functions. Strategies to optimize these parameters with respect to information extraction are therefore relevant to building better statistical potentials. Last, we demonstrate that the backbone torsion potential, defined by the trimer sequence, can be an effective tool in greatly reducing the set of possible conformations from a vast decoy pool. Proteins 2006.

Information and discrimination in pairwise contact potentials

We examine the information‐theoretic characteristics of statistical potentials that describe pairwise long‐range contacts between amino acid residues in proteins. In our work, we seek to map out an efficient information‐based strategy to detect and optimally utilize the structural information latent in empirical data, to make contact potentials, and other statistically derived folding potentials, more effective tools in protein structure prediction. Foremost, we establish fundamental connections between basic information‐theoretic quantities (including the ubiquitous Z‐score) and contact “energies” or scores used routinely in protein structure prediction, and demonstrate that the informatic quantity that mediates fold discrimination is the total divergence. We find that pairwise contacts between residues bear a moderate amount of fold information, and if optimized, can assist in the discrimination of native conformations from large ensembles of native‐like decoys. Using an extensive battery of threading tests, we demonstrate that parameters that affect the information content of contact potentials (e.g., choice of atoms to define residue location and the cut‐off distance between pairs) have a significant influence in their performance in fold recognition. We conclude that potentials that have been optimized for mutual information and that have high number of score events per sequence–structure alignment are superior in identifying the correct fold. We derive the quantity “information product” that embodies these two critical factors. We demonstrate that the information product, which does not require explicit threading to compute, is as effective as the Z‐score, which requires expensive decoy threading to evaluate. This new objective function may be able to speed up the multidimensional parameter search for better statistical potentials. Lastly, by demonstrating the functional equivalence of quasi‐chemically approximated “energies” to fundamental informatic quantities, we make statistical potentials less dependent on theoretically tenuous biophysical formalisms and more amenable to direct bioinformatic optimization. Proteins 2008.

Computational study of packing a collagen‐like molecule: Quasi‐hexagonal vs “Smith” collagen microfibril model

The lateral packing of a collagen‐like molecule, CH3CO‐(Gly‐L‐Pro‐L‐Pro)4‐NHCH3, has been examined by energy minimization with the ECEPP/3 force field. Two current packing models, the Smith collagen microfibril twisted equilateral pentagonal model and the quasi‐hexagonal packing model, have been extensively investigated. In treating the Smith microfibril model, energy minimization was carried out on various conformations including those with the symmetry of equivalent packing, i.e., in which the triple helices were arranged equivalently with respect to each other. Both models are based on the experimental observation of the characteristic axial periodicity, D = 67 nm, of light and dark bands, indicating that, if any superstructure exists, it should consist of five triple helices. The quasi‐hexagonal packing structure is found to be energetically more favorable than the Smith microfibril model by as much as 31.2 kcal/mol of five triple helices. This is because the quasi‐hexagonal packing geometry provides more nonbonded interaction possibilities between triple helices than does the Smith microfibril geometry. Our results are consistent with recent x‐ray studies with synthetic collagen‐like molecules and rat tail tendon, in which the data were interpreted as being consistent with either a quasi‐hexagonal or a square‐triangular structure.

Discriminative ability with respect to amino acid types: Assessing the performance of knowledge-based potentials without threading

We present a novel method designed to analyze the discriminative ability of knowledge‐based potentials with respect to the 20 residue types. The method is based on the preference of amino acids for specific types of protein environment, and uses a virtual mutagenesis experiment to estimate how much information a given potential can provide about environments of each amino acid type. This allows one to test and optimize the performance of real potentials at the level of individual amino acids, using actual data on residue environments from a dataset of known protein structures. We have applied our method to long‐range and medium‐range pairwise distance‐dependent potentials. The results of our study indicate that these potentials are only able to discriminate between a very limited number of residue types, and that discriminative ability is extremely sensitive to the choice of parameters used to construct the potentials, and even to the size of the training dataset. We also show that different types of pairwise distance potentials are dominated by different types of interactions. These dominant interactions strongly depend on the type of approximation used to define residue position. For each potential, our methodology is able to identify a potential‐specific amino acid distance matrix and a reduced amino acid alphabet of any specified size, which may have implications for sequence alignment and multibody models. Proteins 2002;49:266–284.

Property‐based sequence representations do not adequately encode local protein folding information

We examine the informatic characteristics of amino acid representations based on physical properties. We demonstrate that sequences rewritten using contracted alphabets based on physical properties do not encode local folding information well. The best four‐character alphabet can only encode ∼57% of the maximum possible amount of structural information. This result suggests that property‐based representations that operate on a local length scale are not likely to be useful in homology searches and fold‐recognition exercises. Proteins 2007.