Better prediction of protein contact number using a support vector regression analysis of amino acid sequence
- Zheng Yuan^{1}Email author
DOI: 10.1186/1471-2105-6-248
© Yuan; licensee BioMed Central Ltd. 2005
Received: 04 July 2005
Accepted: 13 October 2005
Published: 13 October 2005
Abstract
Background
Protein tertiary structure can be partly characterized via each amino acid's contact number measuring how residues are spatially arranged. The contact number of a residue in a folded protein is a measure of its exposure to the local environment, and is defined as the number of C_{ β }atoms in other residues within a sphere around the C_{ β }atom of the residue of interest. Contact number is partly conserved between protein folds and thus is useful for protein fold and structure prediction. In turn, each residue's contact number can be partially predicted from primary amino acid sequence, assisting tertiary fold analysis from sequence data. In this study, we provide a more accurate contact number prediction method from protein primary sequence.
Results
We predict contact number from protein sequence using a novel support vector regression algorithm. Using protein local sequences with multiple sequence alignments (PSI-BLAST profiles), we demonstrate a correlation coefficient between predicted and observed contact numbers of 0.70, which outperforms previously achieved accuracies. Including additional information about sequence weight and amino acid composition further improves prediction accuracies significantly with the correlation coefficient reaching 0.73. If residues are classified as being either "contacted" or "non-contacted", the prediction accuracies are all greater than 77%, regardless of the choice of classification thresholds.
Conclusion
The successful application of support vector regression to the prediction of protein contact number reported here, together with previous applications of this approach to the prediction of protein accessible surface area and B-factor profile, suggests that a support vector regression approach may be very useful for determining the structure-function relation between primary protein sequence and higher order consecutive protein structural and functional properties.
Background
Prediction of protein three-dimensional structure from primary sequence is the central problem in structural bioinformatics. One approach is to use known structur∈homolog proteins as templates to determine the tertiary structures of novel proteins of unknown structure. Approaches include comparative modelling, threading and fold recognition methods. One protein structural feature is of particular interest here, namely, residue contact number (CN) which can be used to enhance protein fold recognition [1]. This measure has also been regarded as the conserved solvent exposure descriptor of similar folds without a common evolutionary origin [2]. Furthermore, contact number may be used to determine the energy function allowing molecular dynamics simulations of protein structures [3]. Here, we seek to use protein contact number to assist with the tertiary fold prediction of novel proteins for which an accurate functional relationship between a protein's primary sequence and its residues' contact numbers must be determined. To fulfil the task, we use a new method, the so-called support vector regression, to approximate the sequence-contact number relationship. We demonstrate that, as a result, we achieve more accurate predicted contact numbers than have been achieved to date.
The contact number, or coordination number, of a given residue of a folded protein is defined as the number of C_{ β }(or C_{ α }) atoms in other residues within a sphere around the C_{ β }(or C_{ α }) atom of that given residue. Previous approaches to the prediction of protein contact number fall into two categories: classification and regression. In the classification approach, residue contact numbers were first classified into two populations allowing a subsequent use of machine learning methods such as recurrent neural networks to perform predictions [4, 5]. Unfortunately, decomposing contact numbers into two states via an arbitrary threshold oversimplifies the problem and much original information is lost. In contrast, the regression approach provides a direct and more accurate way to determine a functional relationship matching contact numbers and protein sequence and thus to provide more accurate contact number predictions. A recent study of Kinjo et al. [3], followed this approach but used a simple linear regression scheme to determine the functional relationship. They reported that the predicted and observe contact numbers had a correlation coefficient (CC) of 0.627. However, most functions in nature are non-linear and cannot be accurately approximated by linear formulas. Under the reasonable expectation that the sequence-contact number is indeed nonlinear, we use a more complicated machine learning method to determine the relationship and expect thereby to obtain more accurate predictions. In particular, we adopt a support vector regression (SVR) algorithm fully capable of determining a non-linear sequence-contact number relationship.
In our former work, we studied the dependence of protein accessible surface area (ASA) [6, 7] and B-factor [8] on primary sequence. These works established that ASAs can be predicted and match observed values with a correlation coefficient of 0.69, while B-factors can be predicted and match observed values with a correlation coefficient of 0.53. These approaches established that multiple sequence inputs outperform single sequence inputs significantly. The importance of using multiple sequences was also observed in prior predictions of contact numbers [3]. Thus, in this present work, we focus on multiple sequence inputs. For completeness, we examine a range of different definitions of contact number ("consecutive" and "discrete"), and also examine whether including further information such as protein molecular weight and amino acid composition allows improved predictions. As a result, we are able to make predictions which match observed values with a correlation coefficient of 0.73, a significant improvement on earlier studies.
Results
Contact numbers according to different r_{ d }values
The mean ($\overline{N}$) and standard deviation (SD) of contact numbers according to different radius (r_{ d }) cutoffs. All results are expressed as ($\overline{N}$, SD).
r_{ d }= 8 Å | r_{ d }= 10 Å | r_{ d }= 12 Å | r_{ d }= 14 Å | |
---|---|---|---|---|
Discrete | 6.14, 3.29 | 12.90, 6.19 | 23.50, 10.46 | 35.39,15.39 |
Consecutive | 6.27, 3.25 | 13.07, 6.14 | 23.56, 10.41 | 35.53, 15.39 |
Estimating the sequence-contact number relationship
Correlation coefficient (CC) and root mean square error (RMSE) for different contact number predictions. All results are expressed as mean ± standard deviation.
r_{ d }= 8 Å | r_{ d }= 10 Å | r_{ d }= 12 Å | r_{ d }= 14 Å | ||
---|---|---|---|---|---|
Discrete | CC | 0.64 ± 0.01 | 0.66 ± 0.01 | 0.69 ± 0.01 | 0.69 ± 0.01 |
RMSE | 0.77 ± 0.01 | 0.75 ± 0.01 | 0.72 ± 0.01 | 0.72 ± 0.02 | |
Consecutive | CC | 0.66 ± 0.01 | 0.67 ± 0.01 | 0.70 ± 0.01 | 0.70 ± 0.01 |
RMSE | 0.75 ± 0.01 | 0.74 ± 0.01 | 0.72 ± 0.01 | 0.72 ± 0.02 |
For all radius cutoffs, predictions using consecutive contact numbers are slightly better than predictions using their discrete counterparts. However, when the larger thresholds (e.g. 12 Å and 14 Å) are used, the accuracy difference decreases to insignificance. Previous work has shown that CNs with a radius cutoff of 12 Å or 14 Å are more useful for protein fold recognition [1]. Likewise, in this work, we also find that these cutoffs give better predictions. But, compared with the discrete contact numbers, the consecutive contact numbers give only a very slight improvement in predictions. The best accuracies are for consecutive contact numbers with thresholds of 12 Å and 14 Å, in which case the correlation coefficient between predicted and observed values can reach 0.70. If we convert the normalized contact numbers to their original absolute ones, the RMSE of 0.72 is equal to an actual error of 7.5 for a threshold of 12 Å, and an actual error of 11.1 for a threshold of 14 Å.
Sequence weight, as a feature, can improve prediction accuracy significantly
Prediction accuracy can be improved by taking account of protein size as measured by its weight. Given a sequence, the weight is the sum of individual weights of consisting amino acids. We use discrete contact numbers (radius cutoff = 12 Å) and divide all proteins into three groups with equal number of proteins, according to their weights. For the three groups, the weights (mean ± standard deviation) are 10611 ± 1637, 17421 ± 2666 and 40073 ± 5952 Daltons. Their average correlation coefficient between predicted and observed contact numbers are 0.61, 0.67 and 0.72, while their average RMSEs are 7.49, 7.09 and 7.35, respectively. These results suggest that smaller molecules are worst predicted.
Correlation coefficients (CCs) and root mean square errors (RMSEs) for individual proteins in different weight groups. The results are given as (mean, median).
LS | LS+W | LS+AA | LS+W+AA | ||
---|---|---|---|---|---|
W≤3485 | CC | 0.61, 0.65 | 0.64, 0.67 | 0.62, 0.66 | 0.64, 0.67 |
RMSE | 7.49, 6.95 | 6.68, 6.41 | 7.24, 6.82 | 6.71, 6.45 | |
13485<W≤22750 | CC | 0.67, 0.70 | 0.68, 0.71 | 0.68, 0.71 | 0.69, 0.72 |
RMSE | 7.09, 6.79 | 6.76, 6.54 | 7.05, 6.79 | 6.76, 6.60 | |
W>22750 | CC | 0.72, 0.73 | 0.72, 0.74 | 0.72, 0.73 | 0.73, 0.74 |
RMSE | 7.35, 7.07 | 7.12, 6.95 | 7.24, 6.94 | 7.10, 6.90 | |
All | CC | 0.67, 0.70 | 0.68, 0.71 | 0.68, 0.71 | 0.68, 0.71 |
RMSE | 7.31, 6.94 | 6.86, 6.66 | 7.18, 6.86 | 6.86, 6.66 |
For all cases, it was determined that amino acid composition information can improve prediction performance. However, sequence weight can give yet more significant improvements. For example, in the group of small molecules, data about amino acid composition can increase the CC mean to 0.62 and decrease the RMSE mean to 7.24, while sequence weight data can increase the CC mean to 0.64 and decrease the RMSE mean to 6.68. When information about both sequence weight and amino acid composition is used together, we find still further improvement compared with using each data individually, although this may not be reflected by all measures. Particularly, the difference between "LS+W" and "LS+W+AA" is very minor. However, all the results clearly show that sequence weight is more important than amino acid composition in the prediction of contact numbers.
Correlation coefficients (CCs) and root mean square errors (RMSEs) are calculated for all residues when performing support vector regression algorithm using different input information.
CC | RMSE | |
---|---|---|
LS | 0.69 ± 0.01 | 0.72 ± 0.01 |
LS+W | 0.733 ± 0.005 | 0.68 ± 0.01 |
LS+AA | 0.708 ± 0.009 | 0.71 ± 0.01 |
LS+W+AA | 0.734 ± 0.006 | 0.68 ± 0.01 |
Examining performance by formulating regression as a two-class problem
The least prediction accuracy (%) for two-class problems according to different contact number definitions. Only local sequence information is used.
r_{ d }= 8 Å | r_{ d }= 10 Å | r_{ d }= 12 Å | r_{ d }= 14 Å | |
---|---|---|---|---|
Discrete | 74.1 | 74.5 | 75.8 | 75.2 |
Consecutive | 75.6 | 75.7 | 76.1 | 75.6 |
Discussion
Protein structural properties such as secondary structure, solvent accessibility and contact number provide valuable information for prediction of protein tertiary structures. How to improve the prediction accuracy of these parameters is still a challenging problem. Following Rost and Sander's pioneering work [11] on how to find a conserved and useful prediction index, Hamelryck [2] examined the conservation of nine solvent-exposure measures and found that contact number is the most conserved (correlation coefficient 0.72). His study suggested that CN is more suitable for fold recognition than other descriptors such as ASA. However, difficulties in accurately expressing the prediction problem (for example, it was previously framed as a two class problem using an arbitrary threshold) limited its further application. Recent work on contact number [3] formulating the problem for a regression analysis has enhanced studies in this area. From our work here, we confirm the utility of a regression analysis, and more specifically, establish that allowing for non-linearity via support vector regression allows a more accurate determination of the sequence-contact number relation which further illuminates relationships between protein structural and functional properties and their primary sequence and other features.
Conclusion
We provide a new method for the prediction of protein contact number. Using protein local sequence information generated by multiple sequence alignments, the correlation coefficient between predicted and observed contact numbers can reach 0.70, with normalized root mean square error less than 0.72. The addition of information about sequence weight and amino acid composition as input features can increase the correlation coefficient to 0.734 and decrease the root mean square error to 0.68. This improvement is mainly attributed to the information about sequence weight while the information about amino acid composition only contributes slightly. Moreover, more than half of the proteins are predicted with correlation coefficients greater than 0.71. The prediction accuracies in the two-class problems, regardless of the cutoff thresholds, are greater than 77.0%. The successful application of SVR approach in this study suggests that it can more accurately describe the relationship between protein contact numbers and primary sequence.
Methods
Residue contact number
We take two definitions of contact number in this study, namely, that of "discrete" and "consecutive" contact number. The "discrete" contact number, N_{ d }, is defined by the number of C_{ β }atoms on other residues located within a sphere of radius r_{ d }centred on the C_{ β }atom of the residue of interest. The discrete contact number for i-th residue in a sequence with M residues is given by
${N}_{d}^{i}={\displaystyle \sum _{j:|j-i|>2}^{M}\sigma ({r}_{i,j})\{\begin{array}{c}\sigma ({r}_{i,j})=1\text{if}{r}_{i,j}<{r}_{d}\\ \sigma ({r}_{i,j})=0\text{if}{r}_{i,j}\ge {r}_{d}\end{array},\left(1\right)}$
where r_{i,j}is the distance between the C_{ β }atoms of the i th and j th residues which are understood to be separated in sequence by at least two amino acids. Note that ${N}_{d}^{i}$ is a discrete integer. By replacing the step function σ(r_{i,j}) with a sigmoid function, ${N}_{d}^{i}$ becomes a real number. This procedure was previously adopted by Kinjo et al. [3] to smooth the discrete contact numbers. A particular sigmoid function is given by
σ(r_{i,j}) = 1/{1 + exp [3(r_{i,j}- r_{ d })]}. (2)
We have tried four values of r_{ d }(8 Å, 10 Å, 12 Å and 14 Å) with discrete and consecutive definitions and thus have 8 combinations all of which will be used in our SVR approach.
Normalization of contact number
The distributions of contact numbers can be approximated by normal distributions, as shown in Fig. 1. With respect to a certain r_{ d }, we calculate the mean ($\overline{N}$) and standard deviation (SD). So, the normalized contact number N_{ norm }is determined by the following formula:
${N}_{norm}=\frac{N-\overline{N}}{SD}.\left(3\right)$
At the first step, we predict the normalized contact number because 1) it is easy to handle the data, and 2) it is easy to compare the results for different r_{ d }thresholds. At the second step, we recover the absolute contact numbers from their predicted normalized values using this equation.
Sequence coding
We predict contact number from protein local sequence. For a given residue, the local sequence contains its N-terminal and C-terminal seven nearest-neighbour residues. Thus, the local sequence makes a window of fifteen amino acids. We code each residue in the window using the PSI-BLAST position-specific scoring matrix [12]. The matrices are obtained by querying the input sequence using PSI-BLAST against the NCBI non-redundant protein sequence database with three rounds, masking coil-coiled and low-complexity regions [13]. The elements in the row of the matrix reflect the probabilities for 20 amino acids occurring at this position. All the elements are divided by 10 for normalization and thus each residue is represented by a 20-dimesional vector. Since the residues in coil-coil and low-complexity regions do not have meaningful scores, we encode the residue with an orthogonal scheme. In the 20-dimensional vector coding a given residue, only the entry representing this type of amino acid is assigned as 0.5 with all other entries set as zeros. To consider the terminal residues, we expend the 20-dimensional vector to being 21-dimensional for all residues. When the last entry is set as 0.5 and other entries have zeros, it represents a blank residue added to the N-terminal or the C-terminal to make a local sequence of 15-residue length. For all other residues, the 21-st entries are set to zero. In summery, a residue is coded by a 315-dimensional vector.
Support vector regression
To find the function between protein local sequence and normalized contact number, we use ∈-insensitive support vector regression (∈-SVR) [14, 15]. The expected function can be formulated as
f(X_{ i }) = 〈W, Φ(X_{ i })〉 + b, (4)
where W is the weight and b is the bias. Φ(X_{ i }) is a non-linear function mapping a data point from the input space to the feature space, so consequently, SVR is able to perform non-linear regression. The goal of the regression is to find the optimal W and b using some optimisation criteria. In ε-SVR, errors greater than ε are penalized, where two positive variables ξ and ξ* are used to measure the deviation of samples outside the ε-insensitive tube. The optimisation problem can be expressed as
$Minimize\frac{1}{2}{\Vert W\Vert}^{2}+C{\displaystyle \sum _{i=1}^{M}({\xi}_{i}+{\xi}_{i}^{*})},$
$subjectto\{\begin{array}{c}f({X}_{i})-{y}_{i}\le \epsilon +{\xi}_{i}\\ {y}_{i}-f({X}_{i})\le \epsilon +{\xi}_{i}^{*}\\ {\xi}_{i},{\xi}_{i}^{*}\ge 0\text{for}i=1,\dots ,M\end{array}\left(5\right)$
where C is the regularization constant that determines the trad∈off between the norm and the error penalty.
The solution of the above problem was given by the authors of ∈SVR [14, 15] as follows,
$f(X)={\displaystyle \sum _{i=1}^{M}({\alpha}_{i}-{\alpha}_{i}^{*})\u3008\Phi ({X}_{i}),\Phi (X)\u3009+b,\left(6\right)}$
where α_{ i }and ${\alpha}_{i}^{*}$ are Lagrange multipliers. We can replace 〈Φ(X_{ i }), Φ(X)〉, the inner product of Φ(X_{ i }) and Φ(X), by a kernel function K(X_{ i }, X), if K(X_{ i }, X) = 〈Φ(X_{ i }), Φ(X)〉. The radial basis function are used in our study, as given by
K(X_{ i }, X) = exp(-γ||X_{ i }- X||^{2}), (7)
where γ is a parameter to be tuned by the user.
We constantly set ε as 0.01, γ as 0.01 and C as 5.0, because this set of parameters yielded the best performance in our previous work [6, 8]. A number of software packages can be used to find the solution such as SVMlight [16].
Dataset preparation and prediction evaluation
To test our approach, we selected 945 unique protein chains, which were previously used for prediction of protein ASA, and were prepared by PDB-REPRDB [17]. The structures solved by X-ray crystallography were with resolution less than 2.0 Å and with an R-factor less than 0.2. All chains are at least 60 amino acids or longer, and the pair-wise identity is less than 25%. The protein names can be found in the additional file 1 (supplementary material).
The proteins are randomly divided into three groups with each group having 315 chains. Each group is in turn used for training with the remaining two groups used for testing. Therefore, each group is tested twice by the two functions derived from the other groups, and as a result we have six groups of examination results.
Pearson's correlation coefficients and root mean square errors are calculated with respects to all residues and individual proteins. In addition, the absolute errors are calculated for the residues with different contact numbers. In order to compare with previous classification methods, we use different thresholds to classify contact numbers as "contacted" or "non-contacted" and compute the overall accuracy. The accuracy is defined as the ratio between the number of correctly predicted residues and the total number.
Declarations
Acknowledgements
The author thanks Dr. Rohan Teasdale for kind assistance and Dr. Michael Gagen for critical reading of this article and some helpful suggestions. The work was supported by funds from the Australian Research Council (ARC) and The University of Queensland Early Research Grant. The computer simulations were performed at the High Performance Computing Facility at The University of Queensland.
Authors’ Affiliations
References
- Karchin R, Cline M, Karplus K: Evaluation of local structure alphabets based on residue burial. Proteins 2004, 55: 508–518. 10.1002/prot.20008View ArticlePubMed
- Hamelryck T: An amino acid has two sides: A new 2D measure provides a different view of solvent exposure. Proteins 2005, 59: 38–48. 10.1002/prot.20379View ArticlePubMed
- Kinjo AR, Horimoto K, Nishikawa K: Predicting absolute contact numbers of native protein structure from amino acid sequence. Proteins 2005, 58: 158–165. 10.1002/prot.20300View ArticlePubMed
- Pollastri G, Baldi P, Fariselli P, Casadio R: Prediction of coordination number and relative solvent accessibility in proteins. Proteins 2002, 47: 142–153. 10.1002/prot.10069View ArticlePubMed
- Pollastri G, Baldi P, Fariselli P, Casadio R: Improved prediction of the number of residue contacts in proteins by recurrent neural networks. Bioinformatics 2001, 17 Suppl 1: S234–42.View ArticlePubMed
- Yuan Z, Huang B: Prediction of protein accessible surface areas by support vector regression. Proteins 2004, 57: 558–564. 10.1002/prot.20234View ArticlePubMed
- Yuan Z, Bailey TL: Prediction of protein solvent profile using SVR. Proceedings of the 26th Annual International Conference of the IEEE EMBS 2004, 2889–2892.
- Yuan Z, Bailey TL, Teasdale RD: Prediction of protein B-factor profiles. Proteins 2005, 58: 905–912. 10.1002/prot.20375View ArticlePubMed
- Kabsch W, Sander C: Dictionary of protein secondary structure: pattern recognition of hydrogen-bonded and geometrical features. Biopolymers 1983, 22: 2577–2637. 10.1002/bip.360221211View ArticlePubMed
- Thompson MJ, Goldstein RA: Predicting solvent accessibility: higher accuracy using Bayesian statistics and optimized residue substitution classes. Proteins 1996, 25: 38–47. Publisher Full Text 10.1002/(SICI)1097-0134(199605)25:1<38::AID-PROT4>3.3.CO;2-HView ArticlePubMed
- Rost B, Sander C: Conservation and prediction of solvent accessibility in protein families. Proteins 1994, 20: 216–226. 10.1002/prot.340200303View ArticlePubMed
- Altschul SF, Madden TL, Schaffer AA, Zhang J, Zhang Z, Miller W, Lipman DJ: Gapped BLAST and PSI-BLAST: a new generation of protein database search programs. Nucleic Acids Res 1997, 25: 3389–3402. 10.1093/nar/25.17.3389PubMed CentralView ArticlePubMed
- Jones DT: Protein secondary structure prediction based on position-specific scoring matrices. J Mol Biol 1999, 292: 195–202. 10.1006/jmbi.1999.3091View ArticlePubMed
- Vapnik V: The Nature of Statistical Learning Theory. Springer, New York, ; 2000.View Article
- Smola A, Scholkopf B: A tutorial on support vector regression. NeuroCOLT Technical Report Series, NC-TR-1998–030 1998. [http://www.neurocolt.com]
- Joachims T: Making large-Scale SVM Learning Practical. In Advances in Kernel Methods - Support Vector Learning. Edited by: Schölkopf B, Burges C and Smola A. , MIT Press; 1999.
- Noguchi T, Akiyama Y: PDB-REPRDB: a database of representative protein chains from the Protein Data Bank (PDB) in 2003. Nucleic Acids Res 2003, 31: 492–493. 10.1093/nar/gkg022PubMed CentralView ArticlePubMed
Copyright
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.