PROFcon: novel prediction of long-range contacts

Marco Punta ^1,2,* & Burkhard Rost ^1,2,3

This article is published in (Bioinformatics, issue, 2005 and pages) © copyright Oxford University Press (2005). OUP is the only authorised source. All copying of this article including placing on another website requires the written permission of the copyright owner.

Abstract

Motivation: Despite the continuing advance in the experimental determination of protein structures, the gap betweenthe number of known protein sequences and structures continues to increase. Prediction methods can bridge this sequence-structure gap only partially. Better predictions of non-local contacts between residues could improve comparative modeling, fold recognition, and could assist the experimental structure determination.

Results:  Here, we introduced PROFcon, a novel contact prediction method that combines information from alignments, from predictions of secondary structure and solvent accessibility, from the region between two residues, and from the average properties of the entire protein. In contrast to some other methods, PROFcon predicted short and long proteins at similar levels of accuracy. As expected, PROFcon was clearly less accurate when tested on sparse evolutionary profiles, i.e. on families with few homologues. Prediction accuracy was highest for proteins belonging to the SCOP alpha/beta class. PROFcon compared favorably with state-of-the-art prediction methods at the CASP6 meeting. While the performance may still be perceived as low, our method clearly pushed the mark higher. Furthermore, predictions are already accurate enough to seed predictions of global features of protein structure.

Availability:   http://www.predictprotein.org/submit_profcon.html

Contact:  punta@cubic.bioc.columbia.edu

Key words: protein structure prediction, inter-residue contacts, evolutionary information, neural networks.

Abbreviations used

Introduction

Protein three-dimensional (3D) structure is one key to understanding biological function. Structures unravel details needed to engineer residue mutations and to design protein-specific ligands. Thus, the knowledge of protein structure can impact medical and clinical research. Structural genomics is a large-scale effort that has the experimental determination of most unknown folds as one aim [13, 14, 15, 16, 17]. Structural genomics consortia are progressing rapidly and already contribute every third (J Liu & B Rost, unpublished) of the sequence-unique structures added to the Protein Data Bank (PDB [9]) of 3D structures. Nevertheless, the sequence-structure gap, i.e. the difference between the number of proteins with known sequence and those with known structure is increasing much more rapidly. Methods that predict aspects of protein structure continue to be a crucial means of obtaining structural information that helps in the unraveling of protein function [18, 19, 20, 21].

Despite significant advances over the last years, computational biology can still not reliably generate biologically meaningful 3D models for proteins that have no detectable homology to proteins of known structure [1]. In the absence of a reliable solution to the protein structure prediction problem, developers have addressed simplified problems, such as the prediction of protein secondary structure and solvent accessibility; such methods have evolved into successful, automatic tools that continue to significantly impact experimental and computational biology [10, 22]. One problem with these methods is that they only predict local features of 3D structure; these local features are further simplified by projecting 3D information onto a 1D representation of protein structure (e.g. strings of secondary structure) that neither captures the 2D information of contacts between residues, nor local "irregularities" such as bends of helices. In contrast, two-dimensional (2D) maps of distances between residues reduce the dimensionality of the problem in a way that, in principle, allows the reconstruction of the full 3D structure [23, 24, 25]. NMR spectroscopy exploits this fact in determining structures from distance constraints. Most of the available 2D prediction methods simplify the description from real-valued distances to binary-valued contacts (below an assigned distance threshold; here we defined residues as in contact when their C-beta atoms were closer than 8Å, i.e. 0.8nm, Methods). Contact prediction methods extract information from correlated mutations [26, 27, 28], use neural networks with [29] or without correlated mutations [30], Hidden Markov Models [31, 32], Support Vector Machines [33] and genetic programming [34]. The prediction of 2D maps continues to be a very difficult problem. Consequently, performance is rather limited. Nevertheless, automatic 2D predictions have been used successfully for the prediction of protein structure [28, 35, 36]. Furthermore, no matter how inaccurate 2D predictions, they are still better than constraints from the best de novo 3D prediction methods [4].

Here, we introduced PROFcon, a new method for predicting inter-residue contacts through a simple neural network. For the network input we mixed different sources of information most of which had been used separately in some way before (Methods). We considered information from two "windows" around two residues i and j for which the probability of a spatial contact was predicted. Each sequence position k within the two windows (k ∈ {i-n, …, i+n, j-n, …, j+n}) was characterized by evolutionary substitution profiles from multiple sequence alignments, conservation weights, predicted secondary structure, and predicted solvent accessibility. Additionally, we used the complexity of residues i and j and classified the pair ij into one of seven classes based on their physico-chemical properties. Information from the sequence segment that connects a pair ij has been shown to correlate with the probability of contact formation [37]. The segment length - usually referred to as sequence separation - has been used successfully to predict contacts [38, 33]. However, characterizing the connecting segment in more detail is likely to further improve predictions for residues that are not too far apart [37]. Therefore, we added a third window describing the region at the center of a segment between i and j (analogous to the windows around i and j). Finally, we introduced global features such as the overall compositions of amino acids, predicted secondary structure composition, and protein length. We evaluated the sustained level of performance of our method on a dataset of experimentally determined X-ray structures from the PDB [8, 9]. We benchmarked proteins of different lengths, different number of aligned homologous sequences, and different structural classes assigned according to SCOP [39]. PROFcon performed favorably compared to other state-of-the-art contact prediction methods in the recent CASP6 (Dec. 2004) assessment of blind predictions (O Grana & A Valencia, Madrid, manuscript in preparation). The method is available as an Internet prediction server at http://www.predictprotein.org/submit_profcon.html ; all datasets used for this work are at http://www.rostlab.org/results/2004/profcon/ .

System and Methods

Data sets and cross-validation. All proteins used for the development of our methods were taken from the PDB [9], i.e. have known structure. To avoid biasing methods according to the accidental composition of proteins in the PDB, we extracted a subset of proteins that are not clearly related in sequence. The EVA server evaluates structure prediction methods [40] and maintains a continuously updated subset of sequence-unique PDB chains (no pair of proteins in this set has HSSP-value above 0 [6, 7]). In particular, we used the December 2003 EVA release, a set of 3201 protein chains of known structure. We removed all non-X-ray structures, all membrane and coiled-coil proteins and proteins with physical chain breaks [37]. We divided the proteins into three data sets: (1) for training, we selected structures with high resolution (≤2.0Å), (2) for cross-training, i.e. the optimization of all free neural network parameters such as “stop training” structures with low resolution in the interval (2.5-3.0 Å), and (3) for testing structures with medium resolution (2.0-2.5 Å). Due to CPU limitations, we had to reduce the test set further by excluding all proteins longer than 400 residues. Training, cross-training, and test set contained 748, 466 and 633 proteins, respectively.

Definition of contact.  More for the sake of enabling the direct comparison with other methods than for any other reason, we choose the standard threshold for considering a pair of residues to be in contact [41, 42, 43, 26, 44, 45, 46], namely a maximal distance of 8Å between their C-beta atoms (C-alpha for glycines).

Neural network architecture. We trained standard feed-forward neural networks with back-propagation and momentum term [47]. We addressed the extremely unequal distribution of true (contact) and false (non-contact) samples by balanced training [47]. Symmetry between the contact probabilities for the prediction between ij and ji was enforced through a simple post-processing by averaging over both raw output values [30]. In total, we used 738 input, 100 hidden, and 2 output units (contact, non-contact). The input features corresponded to three different ways of describing each pair of residues ( Table 1S in Supplementary Materials), we used: (1) information from the local environment of both residues, (2) information from the segment connecting i and j, and (3) global information from the entire protein.

(1) Local information from immediate residue environment.  For each residue pair ij in a protein, the network incorporates information from all residues in two windows of size 9 centered around i and j (corresponding to the intervals {i-4;i+4} and {j-4;j+4}). Each residue position within the two windows was characterized by 29 input units: 20 for the evolutionary profile (i.e. frequency of occurrence of the 20 amino acid types at that position, as obtained from multiple sequence alignments [48, 49])), one additional unit served as a spacer accounting for the N- and C-terminal residues [50, 51], 4 units coded for the predicted secondary structure (three for helix/strand/other and one for the reliability of the secondary structure prediction at that residue), 3 units for the predicted solvent accessibility (two units for buried/exposed and one value for prediction reliability) and, finally, 1 for the conservation weight [48]. Alignments were obtained through PSI-BLAST [52] using our standard protocol of three automatic iterations [49] and then filtering the aligned sequences at 80% sequence identity, i.e. any two sequences with >80% percentage pairwise sequence identity were removed at the end. We used PROFphd [10, 11, 12] to predict secondary structure and solvent accessibility. Note that we trained and tested on predicted rather than observed values for 1D structure to account for the fact that secondary structure predictions are more correlated to each other than they are to observed secondary structure [53]. As a consequence, using predictions both in training and testing can result in a more coherent input to the networks and hence can help to ease the classification task (note that similar advantages hold for the prediction of secondary structure itself [47, 48, 12]). As the two local windows together accounted for 18 residue positions, we needed a total of 522 input units for their description (18*29). We also introduced additional features to better characterize the central residues i and j, namely a coarse-grained bio-physical classification [54] (7 input units: hydrophobic-hydrophobic, polar-polar, charged-polar, opposite charges, same charges, aromatic-aromatic, other) and we specified whether or not i and j were in low-complexity regions (according to SEG program [55], 2 input units).

(2) Local information from connecting segment. Since the centre of segments that connect two residues i and j has been shown to be most informative for the contact formation between these two residues [37], we introduced another window of five consecutive residues that spanned the interval {int(|i-j|/2)-2; int(|i-j|/2)+2}. Each residue in this window was characterized by the same information as used for the windows around i and j, i.e. we used 29 input units for each residue. Also, it has been shown that the probability for contact formation decreases as sequence separation increases [56, 43, 44, 57]. Therefore, we also had to encode the length of the segment that connects i and j; for this, we used 11 input units that corresponded to sequence separations of 6, 7, 8, 9, 10-14, 15-19, 20-24, 25-29, 30-39, 40-49, and >49 (values chosen by intuition instead of by optimization). Finally, we added features that described the entire segment, namely its amino acid composition (20 units), its secondary structure composition (3 units) and the fraction of SEG-low-complexity [55] residues in that segment (1 unit). Overall, we used 180 input units to describe the connecting segment.

(3) Global information. The use of global information can help the network to ‘decide’ whether or not two residues are in contact. For example, knowing that the protein is very short should increase the probability of having a contact between two cysteine residues (disulfide bridges); knowing that the protein is longer than the average domain length (~100 residues) should decrease the probability of long-range contacts (fewer inter-domains contacts). Here, we used only very coarse-grained global features, namely 20+3 units to describe the composition of amino acids and secondary structure of the entire protein, and 4 units to describe the protein length (intervals 1-61, 61-120, 121-240, >241; again, values were not optimized but chosen identically to our PHD methods [47, 58, 48]).

Measuring performance. Many measures have been introduced to evaluate the performance of 2D predictions. We applied criteria that were basically identical to those used at CASP/CAFASP [2, 29, 3, 4]. Here, we only briefly sketched these scores that are described in detail elsewhere [26, 28, 29, 4]. Accuracy (also referred to as "specificity") was defined by:

where NC_ok is the number of correctly predicted contacts, i.e. the true positives (TP), NC_prd is the number of predicted contacts, which corresponds to the sum of true positives (TP) and false positives (FP). We also followed the tradition to evaluate performance on a number of predicted contacts that is proportional to a fraction of the protein length L. The rationale behind this choice is that the overall number of contacts in a protein is linearly correlated to the protein length. The advantage is that accuracy estimates are related to a quantity that can be evaluated from sequence alone (hence, it is known a priori). However, in isolation accuracy alone does not suffice, instead we need to contrast it with the coverage (also referred to as "sensitivity"), defined as:

where NC_ok is the number of correctly predicted contacts, i.e. the true positives (TP), NC_obs is the number of observed contacts, which corresponds to the sum of true positives (TP) and false negatives (FN). The following example illustrates the importance of considering coverage. If we fix the number of predictions to a fraction of the length of a protein (e.g. L/2), and if we assume a linear relation between the number of contacts (NC_obs) and the protein length (L), i.e. NC_obs=alpha+beta*L, we can write the coverage as:

where NC_prd is the number of predicted contacts. A simple regression on our test set (633 proteins with L≤400 and considering only sequence separations s≥6) estimated the two free parameters to be: alpha~-220 and beta~5 (valid only for L≥45, see denominator in eqn. 3). For a protein of 100 residues, we would obtain Cov@0.18*Acc, while a protein with 400 residues would yield Cov~0.11*Acc. In other words, the same level of accuracy corresponds to different coverage for proteins of different length. Therefore, even if we predict a number of contacts proportional to the length of the protein L, we need to report Cov along with Acc in order to capture the performance of a method.

For each score we reported the average over all proteins in the test set and the associated estimates for standard deviations of the averages obtained from bootstrapping [59]. In Table 1, Table 2, Table 3, and Table 4 we reported performance for two particular values in the minimal sequence separation (s≥6 and s≥24), i.e. two different definitions of ‘non-local’ contacts, and for a number L/2 of predictions. In Fig. 1 we reported accuracies for different sequence separation values, and in Fig. 1 S the accuracy for different numbers of predicted contacts ( Fig. 1S , panels A and B) and for different values in the prediction reliability, i.e. the normalized network output (Fig. 1 "δ-evaluation" [35] of the performance of the method ( Table 2S in Supplementary Material). In this case, a contact between two residues i and j is considered as correctly predicted if at least one inter-contact is observed between residues in the interval {i-d, i+d} and residues in the interval {j-d, j+d}, i.e. between residues in two windows of size 2d+1 around i and j. This gives an idea of how far off misplaced predictions are. Note that we never reported scores for data set used for optimizing any free parameter; instead all scores given estimate the performance for proteins of unknown structure.

SCOP classes.  A coarse-grained classification groups proteins according to their structural class [60]. Here, we used the SCOP [61] classification (release 1.65). Of the 633 proteins in our test set, we could assign 522 to one of the four major classes (131 all-alpha, 103 all-beta, 119 alpha+beta, 169 alpha/beta). The remaining proteins were either in other structural classes (peptides, small and multi domain proteins), or were not assigned by SCOP, yet.

Contact density. We defined contact density, i.e. the density of non-local contacts in a protein as:

where L was the length of the protein, NC_obs(s) was the number of experimentally observed contacts (C-beta ≤8 Angstrom) at sequence separation s or greater, and N(s) was the number of residue pairs in that protein separated by s or more sequence positions (counting only the upper diagonal of the symmetric matrix). We chose s=6 for all reports of contact density for simplification. D depends on at least two factors, namely the protein length and structural class of a protein. For instance, in our test set, D(6)=0.035 for proteins with L<150 (199 proteins) and D(6)=0.0155 for proteins with 250≤L≤400 (226 proteins). In other words, longer proteins have lower contact density as is well known. Contact densities for different structural classes were reported below (Results).

Biophysical features. We analyzed the performance separately for different biophysical contexts. Secondary structure was taken from DSSP [62], with DSSP states G, I and H treated as helix (H), B and E as strand (E) and all other DSSP states as "other" (L). Accessibility was also taken from DSSP. We considered nine amino acids as hydrophobic (alanine, leucine, isoleucine, valine, tryptophane, phenylalanine, proline, cysteine, and methionine).

Results and Discussion

Connecting segment very informative for contact formation.  First, we confirmed [37] that the information from the segment connecting two residues i and j improves the prediction of contacts (Fig. 1 Table 1). For sequence separations >20, accuracy was always lower than 20% (Fig. 1), while it reached almost 40% (s=6 and L/20, Fig. 1) for "less global" residues. This difference could partially be explained by the background probability (more contacts at shorter separations as illustrated by the curve for random in Fig.1). On the other hand, it may be easier for the neural networks to learn what determines the formation of more local contacts (more samples, stronger sequence signatures, e.g. for beta turns and beta hairpins). In order to take this strong dependency of performance on sequence separation into account, we always reported values for performance for two different values of minimal sequence separation (s≥6 and s≥24). Note that we provided more details about the explicit dependency of performance in the supplementary material (in particular Fig. 1S ).

^a  Two different methods are compared: networks that use only information from the local sequence environment of two contacting residues (local only) and our system that also uses information from the connecting segment (PROFcon).
^b  Minimal sequence separation: n (=6 or 24) means that only contacts between pairs i, j minimally n residues apart are considered, i.e. |i-j| ≥ n.
^c  .Body { color: #000000; font-family: 'Helvetica', 'Helvetica'; font-size: 12.00pt; font-style: normal; font-variant: normal; font-weight: normal; letter-spacing: 0; line-height: 1.17; margin-bottom: 0.000000pt; margin-left: 0.00pt; margin-right: 0.00pt; margin-top: 0.000000pt; padding-bottom: 12.000000pt; padding-top: 0.000000pt; text-align: left; text-decoration: none; text-indent: 0.00pt; text-transform: none; vertical-align: 0.000000em; } .Bullet { list-style: disc; margin-bottom: 0.000000pt; margin-left: 0.000000pt; margin-top: 0.000000pt; padding-left: 7pt; text-indent: 6.000000pt; } Number of proteins in respective test set.
^d  Acc: Accuracy (aka Specificity, eqn. 1) and Cov: coverage (aka Selectivity, eqn. 2). Note that we report numbers valid for the L/2 strongest predictions.
^e  Standard error from bootstrapping.

Fig. 1

Fig. 1: PROFcon accuracy vs. sequence separation. Accuracy (aka specificity, Eqn. 1) for three particular values of the number of predicted contacts that were considered, as common practice, these values were chosen depending on protein length L (L/20, L/10, L/5; black symbols) and the random baseline (gray circles); the random accuracy is defined as Rnd(s)=NCobs(s)/(L-s), where NCobs(s) is the number of contacts separated by a sequence separation of exactly s residues; (L-s) is the total number of residue pairs in a protein of length L, which are separated by s residues. Trivially, if fewer contacts are predicted (L/20 < L/10 < L/5), accuracy is higher. However, this increased accuracy comes at the expense of a low coverage. Note that all thresholds shown implied that the number of predicted contacts was smaller than the number of observed contacts.

Evolutionary profiles were crucial for performance.  Due to the large size of our datasets and our limited resources we could not systematically test the relevance of all the input features that we used through a leave-one-out type of test. One aspect that we investigated by training separate networks was the contribution of non-local input information, in general: networks using only local features were less accurate than our final system (Table 1). Detailed analyses of our results and preliminary work on smaller datasets suggested which of the remaining input features were most important for performance. Evolutionary information was clearly most relevant, as had been noted before [57]. Even the simplest measure for the information in a multiple sequence alignment, namely the number of proteins aligned, clearly correlated with performance, e.g. the accuracy dropped from 37% for alignments with over 200 proteins to 23% for alignments with fewer than 15 proteins at sequence separations ≥6 and from 24% to 13% at separations ≥24 (Column Acc in Table 2). Note that, although accuracy correlated dramatically with the number of aligned sequences, differences in coverage (Column Cov in Table 2) were not statistically significant. This was probably related to the different protein composition of the reported subsets, in terms of length and structural class.

^a  Number of proteins in multiple sequence alignment used to extract evolutionary profiles.
^b  As in Table 1 ; note all values compiled for the first L/2 predictions (Methods).

Contact density dependent on type of protein.  It is well known that the contact density decreases with increasing protein length. Thus, contact predictions are more difficult for longer proteins [29, 30]. We observed that the contact density (eqn. 4) also depends on the structural class: all-beta proteins had the highest density (D(6)=0.040), while all-alpha proteins had the lowest contact density (D(6)=0.022). Since it is more difficult to predict low-density than high-density contacts, most existing methods for contact prediction strongly depend on protein length [38, 30] and structural class [33, 34].

Similar accuracy but better performance for short proteins. PROFcon reached surprisingly similar levels of accuracy for proteins of very different length (Table 3 column Acc). However, when also considering the coverage/selectivity of our predictions (Table 3 column Cov), we noted the length-dependency of our method: short proteins (especially of length < 100) had, by far, the highest coverage (1.5-2 times higher than for long proteins).

^a Protein length.
^b  As in Table 1; note all values compiled for the first L/2 predictions (Methods).

All-alpha worst and alpha/beta best. Our test set was large enough to distinguish between the four major structural classes in SCOP [39], namely all-alpha, all-beta, alpha/beta, and alpha+beta (at least 100 proteins in each class). We found that PROFcon performed clearly worst on all-alpha, especially for shorter sequence separations (Table 4 ); we verified that this effect was true at levels of identical coverage (data not shown). Performance was largely similar for the other three classes with the exception of very long-range contacts (s≥24) that were predicted best in alpha/beta proteins (Table 4). A similar trend has previously been reported on a much smaller data set [34]. This trend might originate from particular strand-helix-strand modules that are abundant in alpha/beta proteins, namely those with two flanking strands that contact each other. This type of structural motif may be easy to predict, especially for a system that relies on information about the connecting segments.

^a As in Table 1; note all values compiled for the first L/2 predictions (Methods).

50% of predicted contacts within two residues of an observed contact. "δ-analysis" (Methods) shows that many predicted contacts fall very close to observed contacts (Table 2S in the Supplementary Material). For example, for sequence separation s≥6, 50% of all predicted contacts (L/2 predictions) are within 2 residues of an experimental contact ( Table 2S ).

Correct for core, hydrophobic, and regular secondary structure. At least half of the contacts correctly predicted by PROFcon were between residues in identical regular secondary structures (helix-helix and strand-strand, Table 3S in the Supplementary Material); this was independent of the structural class and of sequence separation ( Table 3S ). Although most correctly predicted contacts were between regular secondary structures, all residues contacting between helices were, on average, predicted the least accurately (slightly worse than mixed). Strand-strand contacts were by far the most accurately predicted (>40% for s≥6 and > 20% for s≥24, in all classes, Table 3S ). In alpha/beta proteins long-range strand-strand contacts (s≥24) were predicted at levels of accuracy as high as 42% to be compared to 20% and 24% in all-beta and alpha+beta, respectively. This indicated that PROFcon captured a strong signal from long-range preferences determining the formation of sheets in alpha/beta proteins that may be related to strand-helix-strand modules. In analogy to contacts between regular secondary structures, contacts between hydrophobic residues also constitute most of the correctly predicted contacts ( Table 4S in the Supplementary Material). Unlike for secondary structure, the contribution of hydrophobic pairs is slightly increasing for higher sequence separations. As expected, predicted contacts are on the average more buried (core residues) and less distant in sequence than observed contacts (averages between 11 and 22 AngstromŲ for correctly predicted contacts, Table 4S ; L/2 predictions).

CASP6 and comparisons to other methods.  Comparing the performance of PROFcon to that of other contact prediction methods is not an easy task; different groups use different datasets and, as shown (Table 2, Table 3, Table 4), the dataset composition (number of sequences in alignments, protein length, structural class composition) significantly alters average scores. Furthermore, different groups use different scores, often even different definitions for what is considered a long-range inter-residue contact. The only reasonable comparison of the performance of methods is based on the same scores and the same data set. Such a data set must be sufficiently large to contain a representative sequence-unique subset of proteins [47, 58, 63, 64, 22]. Furthermore, the set should not overlap with any of the proteins used for development. While cross-validation provides some clues, it does not suffice for rigorous comparisons. At the moment, there is no set available that meets all conditions for a comprehensive, meaningful comparison of our method to others. The best approximation might be the data from CASP6 (Dec. 2004), with the caveat that this set was far too small to draw definite conclusions. PROFcon appeared to be one of the top three contact prediction methods at CASP6, as judged by the assessor (A Valencia, CASP6 web site; note that the limitation of the data set did not allow any distinction in performance between the top three methods).

Contact predictions capture relevant information and are useful! More than many other structure prediction methods, contact predictions continue to suffer from the fact that performance appears to be so low. Are 2D predictions of any use? Our predictions clearly capture important globular information as demonstrated by a method that succeeds in predicting folding rates exclusively based on PROFcon predictions [65]. Ortiz, Skolnick and colleagues have shown that even more noisy contact predictions provide important constraints for the prediction of 3D structures [35], and used by experts, contact predictions have also been shown to aid in fold recognition [28, 66]. Two particular examples from CASP6 may help elucidating what exactly is captured by 2D predictions.

Example 1: CASP6 target T0230 (Fig. 2 A) is a small protein (~100 residues) characterized by the following topology: two helices, two interacting strands (labeled A and B in Fig. 2 A), one helix (labeled 1), another strand (labeled C, which interacts with helix A) and one last helix. At CASP6, T0230 was classified as fold recognition analogous target (FR/A), i.e. a protein for which a structure was known in PDB, however, the similarity between the template and the known structure could not be identified through sequence homology. PROFcon strongly predicted the interaction between the two anti-parallel strands A and B that are separated by a short loop (Fig. 2). It also correctly identified the sparse cluster of interactions between helix 1 and strand C. However, it wrongly predicted the main interaction of strand C. In fact, while predicting only a few interactions between parallel strands C and B (that are in contact in the structure), it suggested a strong contact between C and A that is not observed.

Fig. 2

Fig. 2: Predicted and observed contacts for two CASP6 targets. Examples of predictions for two targets from CASP6: T230 (A) and T216_2 (B, second domain of target T216). Right side: simple sketch of 3D structure using VMD [68], helices shown as red cylinders, strands as yellow arrows; other residues are marked in blue. Left side: contact maps; the upper triangles contain the best 2*L predictions of PROFcon (red dots); the lower triangles contain the experimental contact maps (black dots). All contacts between residues closer than six residues (sequence separation <6) are removed, and the lines parallel to the diagonals indicate sequence separations of 24. The blue boxes highlight specific clusters of interactions between selected secondary structure elements in the observed structure; gray boxes highlight over-predicted (predicted, not observed) clusters of contacts between selected secondary structure elements. The labels on the boxes are taken from the 3D sketches on the right. For example, on the top panel (A) the blue box label a-b means that the highlighted area corresponds to contacts occurring between strands a and b of target T0230.

Example 2: The substantially longer (~210 residues) domain labeled T0216_2 in CASP6 (Fig. 2 B) is part of a two-domain protein. CASP6 assessors classified this target as a novel fold (NF), i.e. a domain for which no other domain from PDB had a structural similarity. We focused on some of the strands because this domain has a rather complicated architecture. Specifically, we considered two groups of strands labeled A, B, C and 1, 2, 3, 4 in Fig. 2 B. PROFcon correctly identified contacts between pairs of strands A-B, 2-3 and 3-4, while it completely missed the interactions between A-C and 1-2. In both examples (Fig. 2), incorrect contact predictions between regular secondary structure segments often occurred far from the main diagonal of the map (i.e. at large sequence separations, see lines marking s=24). The overall accuracy for 2*L predictions and s≥6 was close to average for both T0230 (26%) and T0216_2 (21%).

These two examples seem to suggest that even when the overall prediction accuracy is rather low, PROFcon still correctly identifies contacts between regular secondary structures that are separated by less than 20-30 residues. This ability could be exploited by integrating the predictions into fold recognition and/or de novo prediction methods.

Unique combination of information makes the difference. Our approach to contact prediction was by no means "radically new". We simply combined sources of information slightly differently from what other groups did. In doing so, we ended up with a simple neural network the size of which is not unusual for the predictions of 1D structure [10, 22], but is slightly larger than what has previously been used to predict 2D structure. As a result we needed a very large training set with almost 400,000 positive (contact) samples. Preliminary tests demonstrated that we needed at least these many samples to fully profit from all the input features that we combined. The downside was the increase in computational complexity: the development required several terabytes and many CPU years. These constraints also impacted our ability to separately test - and possibly optimize - the various input features that we considered. One outstanding feature of PROFcon is its consistent performance across a wide range of protein lengths and for both less and more long-range contacts. This consistency was clearly related to the introduction of information from the segment connecting two contacting residues (Table 1).

Future. Major improvements from here may have to introduce ways of post-processing raw predictions. Our method did not exploit any of the constraints imposed on a contact map of an entire protein, e.g. that residues can form only a limited number of contacts. In the past, an intricate post-processing has been proposed for beta-strand pairing [67]. Although the concepts embedded in HMMSTR also address this task [32], no solution exists that comprehensively refines inter-residue contact predictions.

Conclusions

Better predictions of 2D information captured by inter-residue contacts from sequence could help predicting important aspects of protein structures. However, the difficulty of the task and the perceived lack of acceptable performance have so far hampered progress. We presented a new method that exploits evolutionary information in form of multiple sequence alignments and other sequence information relevant for predicting contacts through simple neural networks. While none of our ideas revolutionized the field, the particular combination of information chosen made a significant difference in sustained prediction performance; the major novelty was the particular way in which we successfully used information other than from the sequence environment of the two contacting residues (Table 1). Our method, PROFcon, was particularly successful in its consistent predictions of contacts across a wide range of protein lengths, as well as for residues closer in sequence (separated by at least 6 residues) as for residues very far apart in sequence (separated by at least 24 residues). Nevertheless, PROFcon performed better for short proteins (Table 3), for proteins for which sequence alignment methods detected many homologues (Table 2), and for proteins with beta-strands (Table 4); it was particularly successful for alpha/beta proteins (Table 4). Overall, visual inspections for individual contact maps suggest that the predictions contain more useful information than might be expected from the low levels of accuracy and coverage.

Acknowledgements

Thanks to Jinfeng Liu and Megan Restuccia (Columbia) for computer assistance; to the EVA team, in particular, to Dariusz Przybylski, Ingrid Koh, and Volker Eyrich (all Columbia), Osvaldo Grana and Alfonso Valencia (CNB Madrid). Many thanks for insightful discussions to Yanay Ofran (Columbia), Søren Brunak (CBS Copenhagen), Piero Fariselli and Rita Casadio (both Bologna Univ.), Alfonso Valencia (CNB Madrid), Reinhard Schneider (LION), and Chris Sander (Sloan Kettering, NYC). This work was supported by the grant R01-GM64633-01 from the NIH. Last, not least, thanks to Amos Bairoch (SIB, Geneva), Rolf Apweiler (EBI, Hinxton), Phil Bourne (San Diego Univ.), and their crews for maintaining excellent databases and to all experimentalists who enabled this analysis by making their data publicly available.

References