Linear programming optimization and a double statistical filter for protein threading protocols.

The design of scoring functions (or potentials) for threading, differentiating native-like from non-native structures with a limited computational cost, is an active field of research. We revisit two widely used families of threading potentials: the pairwise and profile models. To design optimal scoring functions we use linear programming (LP). The LP protocol makes it possible to measure the difficulty of a particular training set in conjunction with a specific form of the scoring function. Gapless threading demonstrates that pair potentials have larger prediction capacity compared with profile energies. However, alignments with gaps are easier to compute with profile potentials. We therefore search and propose a new profile model with comparable prediction capacity to contact potentials. A protocol to determine optimal energy parameters for gaps, using LP, is also presented. A statistical test, based on a combination of local and global Z-scores, is employed to filter out false-positives. Extensive tests of the new protocol are presented. The new model provides an efficient alternative for threading with pair energies, maintaining comparable accuracy. The code, databases, and a prediction server are available at http://www.tc.cornell.edu/CBIO/loopp.

Enzymatic circularization of a malto-octaose linear chain studied by stochastic reaction path calculations on cyclodextrin glycosyltransferase.

Cyclodextrin glycosyltransferase (CGTase) is an enzyme belonging to the alpha-amylase family that forms cyclodextrins (circularly linked oligosaccharides) from starch. X-ray work has indicated that this cyclization reaction of CGTase involves a 23-A movement of the nonreducing end of a linear malto-oligosaccharide from a remote binding position into the enzyme acceptor site. We have studied the dynamics of this sugar chain circularization through reaction path calculations. We used the new method of the stochastic path, which is based on path integral theory, to compute an approximate molecular dynamics trajectory of the large (75-kDa) CGTase from Bacillus circulans strain 251 on a millisecond time scale. The result was checked for consistency with site-directed mutagenesis data. The combined data show how aromatic residues and a hydrophobic cavity at the surface of CGTase actively catalyze the sugar chain movement. Therefore, by using approximate trajectories, reaction path calculations can give a unique insight into the dynamics of complex enzyme reactions.

Distance-dependent, pair potential for protein folding: results from linear optimization.

The results of an optimization of a folding potential are reported. The complete energy function is modeled as a sum of pairwise interactions with a flexible functional form. The relevant distance between two amino acids (2 - 9 A) is divided into 13 intervals, and the energy of each interval is optimized independently. We show, in accord with a previous publication (Tobi et al., Proteins 2000;40:71-85) that it is impossible to find a pair potential with the above flexible form that recognizes all native folds. Nevertheless, a potential that rates correctly a subset of the decoy structures was constructed and optimized. The resulting potential is compared with a distance-dependent statistical potential of Bahar and Jernigan. It is further tested against decoy structures that were created in the Levitt's group. On average, the new potential places native shapes lower in energy and provides higher Z scores than other potentials.

fw2.2: a quantitative trait locus key to the evolution of tomato fruit size.

Domestication of many plants has correlated with dramatic increases in fruit size. In tomato, one quantitative trait locus (QTL), fw2.2, was responsible for a large step in this process. When transformed into large-fruited cultivars, a cosmid derived from the fw2.2 region of a small-fruited wild species reduced fruit size by the predicted amount and had the gene action expected for fw2.2. The cause of the QTL effect is a single gene, ORFX, that is expressed early in floral development, controls carpel cell number, and has a sequence suggesting structural similarity to the human oncogene c-H-ras p21. Alterations in fruit size, imparted by fw2.2 alleles, are most likely due to changes in regulation rather than in the sequence and structure of the encoded protein.

Chemical development of latent fingerprints: computational design of ninhydrin analogues.

The design of chemical compounds for development of latent fingerprints is explored computationally. Our main findings are: (a) We show why past attempts to improve the widely used ninhydrin gave relatively small improvements (referring to color only). The optical transition is connected with a "transition core" and therefore is influenced little by substitution on the aromatic rings. (b) We propose new analogues of ninhydrin with a significant potential such as thiono derivatives.

On the design and analysis of protein folding potentials.

Pairwise interaction models to recognize native folds are designed and analyzed. Different sets of parameters are considered but the focus was on 20 x 20 contact matrices. Simultaneous solution of inequalities and minimization of the variance of the energy find matrices that recognize exactly the native folds of 572 sequences and structures from the protein data bank (PDB). The set includes many homologous pairs, which present a difficult recognition problem. Significant recognition ability is recovered with a small number of parameters (e.g., the H/P model). However, full recognition requires a complete set of amino acids. In addition to structures from the PDB, a folding program (MONSSTER) was used to generate decoy structures for 75 proteins. It is impossible to recognize all the native structures of the extended set by contact potentials. We therefore searched for a new functional form. An energy function U, which is based on a sum of general pairwise interactions limited to a resolution of 1 angstrom, is considered. This set was infeasible too. We therefore conjecture that it is not possible to find a folding potential, resolved to 1 angstrom, which is a sum of pair interactions.

Unit-vector RMS (URMS) as a tool to analyze molecular dynamics trajectories.

The Unit-vector RMS (URMS) is a new technique to compare protein chains and to detect similarities of chain segments. It is limited to comparison of C(alpha) chains. However, it has a number of unique features that include exceptionally weak dependence on the length of the chain and efficient detection of substructure similarities. Two molecular dynamics simulations of proteins in the neighborhood of their native states are used to test the performance of the URMS. The first simulation is of a solvated myoglobin and the second is of the protein MHC. In accord with previous studies the secondary structure elements (helices or sheets) are found to be moving relatively rigidly among flexible loops. In addition to these tests, folding trajectories of C peptides are analyzed, revealing a folding nucleus of seven amino acids.

Coupling the folding of homologous proteins.

The empirical observation that homologous proteins fold to similar structures is used to enhance the capabilities of an ab initio algorithm to predict protein conformations. A penalty function that forces homologous proteins to look alike is added to the potential and is employed in the coupled energy optimization of several homologous proteins. Significant improvement in the quality of the computed structures (as compared with the computational folding of a single protein) is demonstrated and discussed.

Computer simulations of carbon monoxide photodissociation in myoglobin: structural interpretation of the B states.

The early diffusion processes of a photodissociated ligand (carbon monoxide) in sperm whale myoglobin and its Phe29 mutant are studied computationally. An explicit solvent model is employed in which the protein is embedded in a box of at least 2300 water molecules. Electrostatic interactions are accounted for by using the particle mesh Ewald. Two hundred seventy molecular dynamics trajectories are computed for 10 ps. Different models of solvation and the ligand, and their influence on the diffusion are examined. The two B states of the CO are identified as "docking" sites in the heme pocket. The sites have a similar angle with respect to the heme normal, but differ in the orientation in the plane. The computational detection of the B states is stable under a reasonable variation of simulation conditions. However, in some trajectories only one of the states is observed. It is therefore necessary to use extensive simulation data to probe these states. Comparison to diffraction experiments and spectroscopy is performed. The shape of the experimental infrared spectra is computed. The overall linewidth is in an agreement with experiment. The contributions to the linewidth (van der Waals and electrostatic interactions) are discussed.

Knowledge-based structure prediction of MHC class I bound peptides: a study of 23 complexes.

BACKGROUND: The binding of T-cell antigenic peptides to MHC molecules is a prerequisite for their immunogenicity. The ability to identify binding peptides based on the protein sequence is of great importance to the rational design of peptide vaccines. As the requirements for peptide binding cannot be fully explained by the peptide sequence per se, structural considerations should be taken into account and are expected to improve predictive algorithms. The first step in such an algorithm requires accurate and fast modeling of the peptide structure in the MHC-binding groove. RESULTS: We have used 23 solved peptide-MHC class I complexes as a source of structural information in the development of a modeling algorithm. The peptide backbones and MHC structures were used as the templates for prediction. Sidechain conformations were built based on a rotamer library, using the 'dead end elimination' approach. A simple energy function selects the favorable combination of rotamers for a given sequence. It further selects the correct backbone structure from a limited library. The influence of different parameters on the prediction quality was assessed. With a specific rotamer library that incorporates information from the peptide sidechains in the solved complexes, the algorithm correctly identifies 85% (92%) of all (buried) sidechains and selects the correct backbones. Under cross-validation, 70% (78%) of all (buried) residues are correctly predicted and most of all backbones. The interaction between peptide sidechains has a negligible effect on the prediction quality. CONCLUSIONS: The structure of the peptide sidechains follows from the interactions with the MHC and the peptide backbone, as the prediction is hardly influenced by sidechain interactions. The proposed methodology was able to select the correct backbone from a limited set. The impairment in performance under cross-validation suggests that, currently, the specific rotamer library is not satisfactorily representative. The predictions might improve with an increase in the data.

Kinetics of peptide folding: computer simulations of SYPFDV and peptide variants in water.

The folding of Ser-Tyr-Pro-Phe-Asp-Val (SYPFDV), and sequence variants of this peptide (SYPYD and SYPFD) are studied computationally in an explicit water environment. An atomically detailed model of the peptide is embedded in a sphere of TIP3P water molecules and its optimal structure is computed by simulated annealing. At distances from the peptide that are beyond a few solvation shells, a continuum solvent model is employed. The simulations are performed using a mean field approach that enhances the efficiency of sampling peptide conformations. The computations predict a small number of conformations as plausible folded structures. All have a type VI turn conformation for the peptide backbone, similar to that found using NMR. However, some of the structures differ from the experimentally proposed ones in the packing of the proline ring with the aromatic residues. The second most populated structure has, in addition to a correctly folded backbone, the same hydrophobic packing as the conformation measured by NMR. Our simulations suggest a kinetic mechanism that consists of three separate stages. The time-scales associated with these stages are distinct and depend differently on temperature. Electrostatic interactions play an initial role in guiding the peptide chain to a roughly correct structure as measured by the end-to-end distance. At the same time or later the backbone torsions rearrange due to local tendency of the proline ring to form a turn: this step depends on solvation forces and is helped by loose hydrophobic interactions. In the final step, hydrophobic residues pack against each other. We also show the existence of an off the pathway intermediate, suggesting that even in the folding of a small peptide "misfolded" structures can form. The simulations clearly show that parallel folding paths are involved. Our findings suggest that the process of peptide folding shares many of the features expected for the significantly larger protein molecules.

Simultaneous and coupled energy optimization of homologous proteins: a new tool for structure prediction.

BACKGROUND: Homology-based modeling and global optimization of energy are two complementary approaches to prediction of protein structures. A combination of the two approaches is proposed in which a novel component is added to the energy and forces similarity between homologous proteins. RESULTS: The combination was tested for two families: pancreatic hormones and homeodomains. The simulated lowest-energy structure of the pancreatic hormones is a reasonable approximation to the native fold. The lowest-energy structure of the homeodomains has 80% of the native contacts, but the helices are not packed correctly. The fourth lowest energy structure of the homeodomains has the correct helix packing (RMS 5.4 A and 82% of the correct contacts). Optimizations of a single protein of the family yield considerably worse structures. CONCLUSIONS: Use of coupled homologous proteins in the search for the native fold is more successful than the folding of a single protein in the family.

Novel methods for molecular dynamics simulations.

In the past year, significant progress was made in the development of molecular dynamics methods for the liquid phase and for biological macromolecules. Specifically, faster algorithms to pursue molecular dynamics simulations were introduced and advances were made in the design of new optimization algorithms guided by molecular dynamics protocols. A technique to calculate the quantum spectra of protein vibrations was introduced.

Anharmonic wave functions of proteins: quantum self-consistent field calculations of BPTI.

The harmonic approximation for the potential energy of proteins is known to be inadequate for the calculation of many protein properties. To study the effect of anharmonic terms on protein vibrations, the anharmonic wave functions for the ground state and low-lying excited states of the bovine pancreatic trypsin inhibitor (BPTI) were calculated. The results suggest that anharmonic treatments are essential for protein vibrational spectroscopy. The calculation uses the vibrational self-consistent field approximation, which includes anharmonicity and interaction among modes in a mean-field sense. Properties obtained include the quantum coordinate fluctuations, zero-point energies, and the vibrational absorption spectrum.

Computer determination of peptide conformations in water: different roads to structure.

Fragments of proteins (short peptides) that "fold" suggest a mechanism of how complete conformational search in protein folding is avoided. We used a computational method to determine structures of two foldable peptides in explicit water: RVEW and CSVTC. The optimization starts from random structures and no experimental constraints are used. In agreement with NMR data, the simulations find a hydrophobic pair (Val/Trp) in REVW. The structure of CSVTC is induced by a surface water that bridges two amide hydrogens, a drive to structure hypothesized by Ben-Naim [Ben-Naim, A. (1990) J. Chem. Phys. 93, 8196-8210] that is largely ignored in studies of folding. Tendency to structure in short peptide chains suggests a mechanism for the formation of short-range nucleation sites in protein folding.

Sodium in gramicidin: an example of a permion.

The reaction path and free energy profile of Na+ were computed in the interior of the channel protein gramicidin, with the program MOIL. Gramicidin was represented in atomic detail, but surrounding water and lipid molecules were not included. Thus, only short range interactions were investigated. The permeation path of the ion was an irregular spiral, far from a straight line. Permeation cannot be described by motions of a single Na+ ion. The minimal energy path includes significant motion of water and channel atoms as well as motion of the permeating ion. We think of permeation as motion of a permion, a quasi-particle that includes the many body character of the permeation process, comparable with quasi-particles like holes, phonons, and electrons of solid-state physics. Na+ is accompanied by a plug of water molecules, and motions of water, Na+, and the atoms of gramicidin are highly correlated. The permion moves like a linear polymer made of waters and ion linked and moving coherently along a zigzag line, following the reptation mechanism of polymer transport. The effective mass, free energy, and memory kernel (of the integral describing time-dependent friction) of short range interactions were calculated. The effective mass of the permion (properly normalized) is much less than Na+. Friction varies substantially along the path. The free energy profile has two deep minima and several maxima. In certain regions, the dominant motions along the reaction path are those of the channel protein, not the permeating ion: there, ion waits while the other atoms move. At these waiting sites, the permion's motion along the reaction path is a displacement of the atoms of gramicidin that prepare the way for the Na+ ion.

Nitric oxide recombination to double mutants of myoglobin: role of ligand diffusion in a fluctuating heme pocket.

Picosecond recombination of nitric oxide to the double mutants of myoglobin, His64Gly-Val68Ala and His64Gly.Val68Ile, at E7 and E11, has been studied experimentally and by computation. It is shown that distal residues have a profound effect on NO recombination. Recombination in the mutants may be explained in terms of fluctuating free volume and structure of the heme pocket. The double mutants provide insight into the effects of free volume and steric hindrance on rates of ligand rebinding following photolysis. Water molecules of the first solvation shell replace surface residues deleted by mutation and can block apparent holes in the protein structure. Thus, water molecules extend the time required for ligands to escape significantly to a nanosecond time scale, which is much longer than would be expected for an open heme pocket. Both nearly exponential (G64A68) and markedly nonexponential (native and G64I68) kinetics are observed, a result at variance with expectation from the model of Petrich et al. [Petrich, J.W., Lambry, J.C., Kuczera, K., Karplus, M., Poyart, C., & Martin, J.L. (1991) Biochemistry 30, 3975-3987], which attributes nonexponential kinetics to proximal effects.

Molecular dynamics simulation of NO recombination to myoglobin mutants.

Molecular dynamics simulations on two coupled electronic surfaces are employed to investigate the geminate recombination of nitric oxide to mutants of sperm whale myoglobin. A model for the ground and the excited states is constructed based on experimental data. The crossing between the surfaces is treated using the Landau-Zener formula. The reaction probability and the recombination curves are calculated directly by histogramming the results of an ensemble of trajectories. The experimental trend is reproduced in which the picosecond recombination rate of different mutants increases in the order Phe29 > Leu29 > Val29 > Ala29. Furthermore, in accord with the experiment on significantly longer time scales an opposite trend is obtained, in which the recombination rate for Ala29 is larger than for Phe29. These results are explained by constrained diffusion of the ligand in the heme pocket. The average and the transient volume of the heme pocket is modified by the 29 mutants.

Ligand binding and conformation change in the dimeric hemoglobin of the clam Scapharca inaequivalvis.

The reaction with carbon monoxide of the cooperative dimeric hemoglobin from Scapharca inaequivalvis has been examined by flash photolysis. In the nanosecond time range, geminate rebinding of 5% of dissociated CO occurs with a rate constant of 1.4 x 10(7) s-1. There is a change in absorbance of deoxyhemoglobin following photolysis at a rate of 1.2 x 10(6) s-1, consistent with a shift in the position of the Soret band to longer wavelengths. The amplitude of the change is proportional to the population of deoxydimer. In much of the Soret region this change is greater than the absorbance excursion associated with geminate recombination. There is at least one other slower change associated with the singly liganded species. Geminate rebinding of NO has components of 50, 8, and 0.035 ns-1, accounting for 75%, 25%, and less than 1% of the total reaction observed after a 35-ps photolysis flash. Simulation of diffusion of NO by molecular dynamics shows the ligands moving from the heme pocket to a subsidiary space between the edge of the heme and the surface of the protein.