1. Biophysics and Structural Biology
  2. Computational and Systems Biology
Download icon

Codon-level information improves predictions of inter-residue contacts in proteins by correlated mutation analysis

  1. Etai Jacob
  2. Ron Unger Is a corresponding author
  3. Amnon Horovitz Is a corresponding author
  1. Bar-Ilan University, Israel
  2. Weizmann Institute of Science, Israel
Research Article
  • Cited 4
  • Views 2,255
  • Annotations
Cite as: eLife 2015;4:e08932 doi: 10.7554/eLife.08932

Abstract

Methods for analysing correlated mutations in proteins are becoming an increasingly powerful tool for predicting contacts within and between proteins. Nevertheless, limitations remain due to the requirement for large multiple sequence alignments (MSA) and the fact that, in general, only the relatively small number of top-ranking predictions are reliable. To date, methods for analysing correlated mutations have relied exclusively on amino acid MSAs as inputs. Here, we describe a new approach for analysing correlated mutations that is based on combined analysis of amino acid and codon MSAs. We show that a direct contact is more likely to be present when the correlation between the positions is strong at the amino acid level but weak at the codon level. The performance of different methods for analysing correlated mutations in predicting contacts is shown to be enhanced significantly when amino acid and codon data are combined.

https://doi.org/10.7554/eLife.08932.001

eLife digest

Genes contain instructions to make proteins from building blocks called amino acids. The instructions are encoded in units called codons that each specify a single amino acid in the chain. A small mutation in a particular codon can change the amino acid found at the corresponding position in the protein. Some amino acids interact with other amino acids in the chain, thereby enabling the protein to adopt the three-dimensional shape it needs to work properly. Therefore, a mutation that affects one of these amino acids may have a large impact on the ability of the protein to work.

A mutation at one position in the protein may, however, have little effect if it is accompanied by a ‘compensatory’ mutation at another position. Such compensatory mutations are more likely to occur when the two positions in the protein are close to each other. To identify such mutations, the amino acid sequences of similar proteins from different organisms are aligned and compared.

A computational method called ‘correlated mutation analysis’ searches for pairs of positions in the alignment that display co-variation, i.e. where particular mutations at one position tend to be accompanied by certain mutations at the second position. These pairs are then ranked according to the strength of their correlation and those with the highest ranking are predicted to be in close contact. Such predictions are, however, far from perfect and can give false results.

Jacob et al. developed and tested a new technique of correlated mutation analysis by examining codon sequences as well as amino acid sequences. The rationale behind the technique relies on the fact that several different codons can encode the same amino acid, so that a mutation in a codon does not always change the amino acid it encodes. Therefore, a strong correlation at the amino acid level can be accompanied by a weak correlation at the codon level. In such cases the positions are more likely to be in contact than in cases where there is a strong correlation also at the codon level since the correlation can then be due to constraints at the DNA or RNA level.

Jacob et al. tested their approach using different methods for analyzing correlated mutations that were proposed in previous studies. This showed that the predictions obtained using both amino acid and codon data are significantly more accurate than those obtained by comparing amino acid sequences only. Future work will test whether combining amino acid and codon data can also be used to predict interactions between different proteins.

https://doi.org/10.7554/eLife.08932.002

Introduction

The effects of mutations that disrupt protein structure and/or function at one site are often suppressed by mutations that occur at other sites either in the same protein or in other proteins. Such compensatory mutations can occur at positions that are distant from each other in space, thus, reflecting long-range interactions in proteins (Horovitz et al., 1994; Lee et al., 2008). It has often been assumed, however, that most compensatory mutations occur at positions that are close in space, thus motivating the development of computational methods for identifying co-evolving positions as distance constraints in protein structure prediction (Göbel et al., 1994). These methods, which rely on multiple sequence alignments (MSA) of homologous proteins as inputs, will become increasingly more useful in the coming years owing to the explosive growth in sequence data. The output of methods for correlated mutation analysis (CMA) is a rank order of the pairs of columns in the alignment according to the statistical and/or physical signficance attached to the correlation observed for each pair. The various methods for CMA that have been developed in the past 15 years differ in the measures they employ for attaching significance to the correlations (Livesay et al., 2012; de Juan et al., 2013; Mao et al., 2015). Early measures include, for example, mutual information (MI) from information theory (Gloor et al., 2005) and observed-minus-expected-squared (OMES) in the chi-square test (Kass and Horovitz, 2002).

Statistically significant correlations in MSAs that do not reflect interactions between residues in contact, that is, false positives, can stem from (i) various indirect physical interactions and (ii) common ancestry. The extent of false positives due to the latter source is manifested in the large number of correlations between positions in non-interacting proteins that can be observed when the sequences of non-interacting proteins from the same organism are concatanated and subjected to CMA (Noivirt et al., 2005). Several methods for removing false positives owing to common ancestry were developed (Pollock et al., 1999; Wollenberg and Atchley, 2000; Noivirt et al., 2005; Dunn et al., 2008) but their success in contact prediction using CMA remained limited. False positives due to the former source, that is, indirect physical interactions, can occur when, for example, correlations corrersponding to positions i and j that are in contact and positions j and k that are in contact lead to a correlation for positions i and k that are not in contact. Methods that remove such transitive correlations have been developed in recent years and include, for example, Direct Coupling Analysis (DCA or DI for Direct Information) (Weigt et al., 2009; Morcos et al., 2011), Protein Sparse Inverse COVariance (PSICOV) (Jones et al., 2012) and Gremlin's pseudolikelihood method (Kamisetty et al., 2013). These methods have been found to be very successful in identifying contacting residues (Marks et al., 2012) and they outperform earlier methods (Mao et al., 2015). Nevertheless, their accuracy, which is ∼80% for the correlations in the top 0.1% (ranked by their scores), drops to ∼50% for the top 1% (Mao et al., 2015). Given that the number of contacts in a protein with N residues is ∼N (Faure et al., 2008), it follows that for proteins with, for example, 100 residues (i.e. with 4560 potential contacts between residues separated by at least 5 residues in the sequence) only about 25% of the contacts (i.e. 23 of the top 1% 46 predictions) will be identified by these CMA methods. In addition, these methods require large MSAs comprising thousands of sequences in order to perform well and such sequence data are not always available. Consequently, it is clear that much can be gained from further improvements in methods of CMA. Here, we describe the development and application of a new method for CMA that uses both amino acid and codon MSAs as inputs instead of relying exclusively on amino acid MSAs as done before. We show that contact prediction is improved in a meaningful manner when amino acid and codon information are combined.

Results and discussion

The key premise underlying the method introduced here is that a correlation at the amino acid level between two positions is more likely to reflect a direct interaction if the correlation at the codon level for these positions is weak (Figure 1). In other words, it is assumed that cases of strong correlations at both the amino acid and codon levels for a pair of positions are less likely to reflect selection to conserve protein contacts and more likely to reflect selection to conserve interactions involving DNA or RNA and/or common ancestry. Given this rationale in mind, we decided to test whether contact identification improves when all the pairs of positions are ranked using a score that increases with (i) increasing strength of the correlation at the amino acid level and (ii) decreasing strength of the correlation at the codon level. Such a score, Si, is given, for example, by: Si = Siα(AA)/Si(C), where Si(AA) and Si(C) are the scores generated by method i (e.g., MI) for the amino acid and codon alignments, respectively, and the value of the power α is determined empirically depending on the method (see below).

Example of a pairwise correlation in a multiple amino acid sequence alignment and two possible corresponding codon alignments.

A correlation at the amino acid level between two positions i and j may (top left) or may not (top right) be accompanied by a correlation at the codon level. The premise of the method introduced here is that a correlation at the amino acid level between two positions is more likely to reflect a direct interaction if the correlation at the codon level for these positions is weak (top right).

https://doi.org/10.7554/eLife.08932.003

The approach outlined above was tested for the OMES (Kass and Horovitz, 2002), MI (Gloor et al., 2005), MIp (Dunn et al., 2008) and DCA (Morcos et al., 2011) methods using 114 MSAs each comprising at least 2000 sequences of length between 200 and 500 residues. In the case of the PSICOV method (Jones et al., 2012), only 86 MSAs out of the 114 MSAs were used since the others didn't pass this method's threshold for amino acid sequence diversity. Each MSA also included at least one sequence with a known crystal structure at a resolution <3 Å in which at least 80% of all the residues are resolved. The mean accuracy of contact identification was plotted as a function of the top ranked number of predicted pairwise contacts (Figure 2—figure supplement 1) or as a function of the top ranked fraction of protein length, L, number of predicted pairwise contacts (Figure 2). Residues were considered as being in contact if the distance between their Cβ atoms is ≤8 Å following the definition used in CASP experiments (Ezkurdia et al., 2009) and other studies (Kamisetty et al., 2013; Skwark et al., 2014) (see also Figure 2—figure supplement 2). The results show that the PSICOV and DCA methods outperform the OMES, MI and MIp methods (Figure 2) as established before (Mao et al., 2015). They also show that combining amino acid and codon data leads to an improvement in the predictions by OMES, MI, DCA and PSICOV. In the case of MIp, however, no improvement is observed despite the fact that this method performs worse than DCA and PSICOV. In MIp, a term called average product correction (APC) is subtracted from the MI score for each pair of positions in order to reduce false positives. Removing this correction from PSICOV where it also exists and including the codon data yielded the best method (Figure 2). Hence, we can conclude that there is an overlap between the background noise reduced upon including the APC term and codon data and that including the latter can be more advantageous as we observe for PSICOV.

Figure 2 with 2 supplements see all
Plots of the mean accuracy of contact identification by various methods of correlated mutation analysis as a function of the top ranked fraction of protein length, L, number of predicted pairwise contacts.

The mean accuracies of contact identification by the OMES, MI, MIp, DCA and PSICOV methods are shown either with or without incorporating codon data. Residues were defined as being in contact if the distance between their Cβ atoms is ≤8 Å. PSICOV* indicates that it was carried out without the APC.

https://doi.org/10.7554/eLife.08932.004

The extent of improvement increases with increasing values of the power α until a maximum is reached (Figure 3A) at a value of α that depends on the method used and different values of αmax were, therefore, chosen accordingly. Cross-validation by dividing the MSA data into training and test sets showed that the values of αmax are stable, that is, they do not vary depending on the set of MSAs (Figure 3—figure supplement 1). Given these values of α max, the significance of the extent of improvement was assessed by comparing for each MSA the accuracy of the contact predictions using the different methods with and without incorporating codon data. Significance levels were determined using two non-parametric tests: (i) the Wilcoxon signed-rank test, which takes into account both the number of MSAs for which the accuracy of the contact predictions increases upon incorporating codon data (e.g., 81 in the case of DCA) and the magnitude of the improvement; and (ii) the sign test, which only considers the number of MSAs with improved accuracy. The extent of improvement achieved by incorporating codon data was found to be highly significant as indicated by the p-values obtained using both tests (Figure 3B).

Figure 3 with 1 supplement see all
The effect of the relative weights of amino acid and codon information on contact prediction improvement and its statistical significance.

(A) The median of the extent of improvement in contact prediction for 114 MSAs (86 in the case of PSICOV) is plotted as a function of the value of the power α which determines the relative weights of the amino acid and codon correlations in the score, Si (Si = Siα(AA)/Si(C), where Si(AA) and Si(C) are the respective amino acid and codon scores generated by method i). The extent of improvement was determined by calculating the difference in the areas under the curves (AUC) of prediction accuracy vs number of predictions for each method i with and without incorporation of the codon data normalized by the area under the curve generated without codon data. The analysis was done for domains of length between 200 and 500 residues and at least 2000 coding sequences in their MSA. The value of α which maximizes the median improvement was used for predictions. Maximal respective improvements of 3.9% and 4.2% were found for DCA and MI when α is 2.5, 17.6% for OMES when α is 1.7 and 1.13% for PSICOV when α is 11.2. (B) Stacked bar plots showing the number of MSAs for which including codon data improved the contact predictions using the different methods (orange) and the number of those for which it was otherwise (green). The statistical significance of the improvement achieved by incorporating codon data is indicated by the top and bottom p-values obtained using the Wilcoxon signed-rank and sign tests, respectively.

https://doi.org/10.7554/eLife.08932.007

The improvement in the predictions upon combining amino acid and codon data, when residues are defined as being in contact if the distance between their Cβ atoms is ≤8 Å, led us to examine whether direct contacts are identified better using this contact definition compared with an ‘All’ definition used by others (Morcos et al., 2011) according to which a contact exists if at least one inter-atomic distance between the residues is ≤8 Å. A non-redundant set of 2481 proteins with an available crystal structure at a resolution better than 1.6 Å was compiled and all the residue pairs in each structure that are in contact according to these two definitions were identified. We then determined for each protein what is the fraction of the residue pairs in contact according to these two definitions that are actually in direct physical contact, that is, with a distance <3.5 Å between two of their respective heavy atoms. It should be noted that atoms can interact with each other even if the distance between them is larger than 3.5 Å via, for example, weak electrostatic interactions but pairs of atoms which are closer than 3.5 Å can always be considered as interacting. It may be seen that, on average, pairs in direct contact constitute only about 10% of the pairs in contact according to the ‘All’ definition and 30% of the pairs in contact according to the Cβ-based definition (Figure 2—figure supplement 2). The better success of DCA in identifying contacts according to the Cβ-based definition when amino acid and codon data are combined is, therefore, an important result since more pairs that are in true physical contact are identified in this way.

Our finding that the Cβ-based definition of contacts is better than the ‘All’ definition but still poor (only 30% of the pairs defined as being in contact are in physical contact) prompted us to test the performance of our method for additional contact definitions. The mean of the extent of improvement in contact prediction for 114 domains (or 86 in the case of PSICOV) was, therefore, determined as a function of the distance that must exist between at least two Cβ atoms in different residues in order for them to be defined as being in contact. It may be seen that, in the cases of PSICOV, OMES and DCA, the maximum improvements in contact prediction upon combining amino acid and codon data are when these distances are about 5.5, 7 and 5.5 Å, respectively, and that, in the cases of DCA and OMES, the improvement decreases dramatically when this distance is >∼10 Å (Figure 4). In the case of MI, the extent of improvement upon combining amino acid and codon data is found to be relatively insensitive to the distance used to define a contact and is maximal when it is ∼4.5 Å (Figure 4). These data, therefore, show again that the improvement in contact prediction upon combining amino acid and codon data is greatest when the distance used for contact definition does not lead to many pairs being defined in contact when in fact they are not in direct physical contact.

Improvement in contact prediction as a function of the distance used to define a physical contact.

The mean of the extent of improvement in contact prediction for 114 domains (or 86 in the case of PSICOV) is plotted as a function of the distance that must exist between two Cβ atoms in different residues in order for them to be defined as being in contact. The extent of improvement was determined by calculating the difference in the areas under the curves of prediction accuracy vs number of predictions by OMES, MI, DCA and PSICOV with and without incorporation of the codon data normalized by the area under the curve generated without codon data. The analysis was done for domains of length between 200 and 500 residues and at least 2000 coding sequences in their MSA. The contact predictions were made for the seven sequences with available crystal structures that have the highest resolution and that in all cases is <3 Å.

https://doi.org/10.7554/eLife.08932.009

The added value in combining amino acid and codon data can be illustrated for contact prediction by DCA in the case of Kex1Δp, a prohormone-processing carboxypeptidase from Saccharomyces cerevisiae that lacks the acidic domain and membrane-spanning portion of Kex1p. The crystal structure of Kex1Δp was solved at a resolution of 2.4 Å (Shilton et al., 1997) and its MSA consists of 1877 sequences. The predictions by DCA with or without incorporating codon data are shown in the respective top and bottom halves of the Kex1Δp contact map (Figure 5A). A comparison of the predictions by the two approaches shows that those made with incorporation of codon data are more long-range (in sequence) and more spread throughout the protein structure than those made without incorporation of codon data. Examples for such long-range contacts between different secondary structure elements in Kex1Δp that are predicted only when also the codon data is used include the interactions between Thr148 with Phe185, Ala186 with Leu208 and Leu190 with Leu368 (Figure 5B). This and other examples (Figure 5—figure supplement 1) show that incorporation of codon data can yield predictions of contacts between residues that are distant in sequence and are, thus, of more value for structure prediction.

Figure 5 with 1 supplement see all
Added value of combining amino acid and codon data in contact prediction by DCA illustrated for Kex1Δp, a prohormone-processing carboxypeptidase from Saccharomyces cerevisiae.

(A) Contact map of the structure of Kex1Δp (PDB ID: 1AC5) in which all the contacts are shown as gray rectangles. Residues were defined as being in contact if the distance between their Cβ atoms (Cα for glycine) is ≤8 Å. The top 100 predictions of contacts made with or without incorporating codon data are highlighted above (in red) and below (in green) the diagonal, respectively, and those predicted by both methods by black circles. (B) The crystal structure of Kex1Δp with predicted contacts highlighted. Only true predicted contacts that were not predicted by the original method are highlighted. Each contacting pair has a different color. The contacts were predicted using an MSA with 1877 coding sequences with a length of 415 codons. The magnified region shows some long-range contacts between different secondary structure elements that are predicted only when also the codon data is used.

https://doi.org/10.7554/eLife.08932.010

Conclusions

The input for methods for analysing correlated mutations has exclusively been multiple amino acid sequence alignments. Here, we have shown that improved contact prediction can be achieved by analysing both amino acid and codon MSAs together. The premise of our approach is that direct contacts are more likely if the correlation at the amino acid level is high but at the codon level is low. The score we propose, which reflects this expectation, can be used in conjunction with different methods of CMA but other possible scores should be examined in future work. Importantly, we find cases where contacts between residues that are distant in sequence and, thus, of greatest value for structure prediction are predicted only by using the combined method. Future work should test other potential applications of combined analysis of amino acid and codon MSAs such as predicting protein–protein interactions and, more generally, in feature selection in machine learning.

Materials and methods

Collection of sequences

Protein sequence datasets were collected from Pfam version 27.0 (Finn et al., 2014) based on representative proteomes (Chen et al., 2011) at 75% co-membership threshold (RP75). Protein coding sequences (CDS) of the collected proteins from Pfam were retrieved based on Uniprot cross reference annotations (for Refseq, Ensembl, EMBL and Ensembelgenomes databases in that order of priority) using the EMBL-EBI's WSDbfetch services (McWilliam et al., 2009) and Ensembl REST API (Beta version) (Yates et al., 2015). All collected CDSs were aligned in accordance to the Pfam HMM based MSAs using tranalign tool from the EMBOSS package (Rice et al., 2000). Pfam domain families with more than 2000 successfully retrieved coding sequences were used for further analysis (total of 551 MSAs). Only families with a known crystal structure at a resolution of 3 Å or better (more than 95% of the families have at least three such structures) and with an overlap of at least 80% of the domain sequence to the ATOM sequence in the solved structure were included in the analysis (total of 460 MSAs). Our analysis was also restricted for proteins with more than 200 residues that have a large number of potential contacts for prediction (114 MSAs). PDB structures were assigned to Pfam families in accordance to the mapping in the files downloaded from http://www.rcsb.org/pdb/rest/hmmer?file=hmmer_pdb_all.txt and ftp://ftp.ebi.ac.uk/pub/databases/msd/sifts/text/pdb_chain_uniprot.lst. PDB structures were retrieved and their coordinates were extracted using the bio3D R package (Grant et al., 2006). Pairwise sequence alignments for mapping were performed using Biostrings (Pages H., Aboyoun P., Gentleman R. and DebRoy S. Biostrings: String objects representing biological sequences, and matching algorithms. R package version 2.34.1).

Evaluation of prediction accuracy

The evaluation was based on the all structures with the highest resolution (at least 3 Å) but, in cases where families have more than 30 known structures with unique sequences, only 30 with the best resolution were used (in cases of structures with the same resolution we arbitrarily chose one). The average accuracy of contact predictions for all the crystal structures of each domain family was then calculated so that domain families with many crystal structures would not be over-represented. Two definitions for a contact between two amino acids were employed: a distance of less than 8 Å between Cβ atoms and a distance of less than 8 Å between any two heavy atoms. Only pairs of residues that are separated by at least five amino acids in the protein sequence were considered. Accuracy was calculated as the proportion of true contacts from the N pairs with the highest score in that set. We evaluated the improvement of our method using the difference in the area under the curve (AUC) of the accuracy vs number of predicted pairs of our method relative to the results of the original OMES, MI, MIp, PSICOV and DCA methods. AUC was calculated using the auc function in MESS package in R with the default parameters.

Determination of the number of pairs in physical contact using different contact definitions

A non-redundant set of 2481 PDB entries with a percentage identity cutoff of 20%, resolution better than 1.6 Å and an R-factor cutoff of 0.25 was downloaded from the pre-compiled CullPDB lists (Wang and Dunbrack, 2003) at http://dunbrack.fccc.edu/Guoli/pisces_download.php on February 25, 2015. Two residues were defined to be in a physical contact if they have at least one pair of atoms with a distance ≤3.5 Å. The number of true physical contacts, that is, those with a distance ≤3.5 Å between two of their respective heavy atoms, was determined for each protein in the set and divided by the number of residue pairs defined to be in contact if at least one inter-atomic distance between them is ≤8 Å (designated ‘All’) or if the distance between their Cβ atoms is ≤8 Å. Only pairs of residues that are separated by at least five amino acids along the protein sequence were considered.

Methods for analysing correlated mutations

The score for a pair of positions i and j, S (i,j), for the OMES (Observed Minus Expected Squared) method was calculated, as follows (Kass and Horovitz, 2002; Fodor and Aldrich, 2004):

SOMES(i,j)=ab(OBSaibjEXPaibj)2

where OBSaibj and EXPaibj are the respective observed and expected number of sequences in the MSA with residue type a at position i and residue type b at position j. The score for the mutual information, MI, method was calculated as follows (Gloor et al., 2005):

SMI(i,j)=a=121b=121f(i,a;j,b)logf(i,a;j,b)f(i,a)f(j,b)

where f(i,a) and f(j,b) denote the respective frequencies of occurrence of residue type a at position i and residue type b at position j and f(i,a; j,b) denotes the joint probability of occurrence of residue type a at position i and type b at position j. In the case of the MIp method (Dunn et al., 2008), an average product correction (APC) term is subtracted from the MI score for each pair of positions. The APC term, which is a measure of the background MI shared by positions i and j, is given by:

APC(i,j)=MI(i,x¯)MI(j,x¯)MI¯

where terms in the nominator are the respective average MI values of positions i and j with all other positions in the alignment and the term in the denominator is the average background MI of all the positions in the alignment. The MIp score is given by:

SMIp(i,j)=SMI(i,j)APC(i,j)

The Direct Coupling Analysis (DCA) method (Morcos et al., 2011) was implemented in R for amino acid and codon MSAs based on a Matlab source code provided by Weigt et al. (http://dca.rice.edu/portal/dca/download). The PSICOV code was downloaded from http://bioinfadmin.cs.ucl.ac.uk/downloads/PSICOV/ and used for the predictions based on amino acid MSAs with the default parameters for faster options as recommended by the authors (-p -r 0.001 and with the -l option in order to avoid using the APC term). The PSICOV code was modified in order to carry out the same analysis for codon MSAs and a python script was implemented to perform the whole analysis as done for the other methods using Pfam MSA files in Stockholm format and fasta MSA files as inputs. PSICOV was used here either with the APC for amino acid MSAs or without the APC for the predictions based on both amino acid and codon MSAs.

Available software

The R and Python source codes for the contact prediction by all methods, C source code modifications to PSICOV V2.1b3, R source code for structure-domain sequence mapping and python scripts for generating codon MSAs are available at https://etaijacob.github.io/. Details on the relevant R packages that will be available on CRAN will also be provided at: https://etaijacob.github.io/.

References

  1. 1
    Protein Data Bank
    1. HM Berman
    2. J Westbrook
    3. Z Feng
    4. G Gilliland
    5. TN Bhat
    6. H Weissig
    7. IN Shindyalov
    8. PE Bourne
    (2015)
    Protein Data Bank, RCSB Protein Data Bank, http://www.rcsb.org/.
  2. 2
  3. 3
  4. 4
  5. 5
  6. 6
  7. 7
  8. 8
  9. 9
  10. 10
  11. 11
  12. 12
  13. 13
  14. 14
  15. 15
  16. 16
  17. 17
  18. 18
  19. 19
  20. 20
  21. 21
  22. 22
  23. 23
  24. 24
  25. 25
  26. 26
  27. 27
  28. 28
  29. 29
  30. 30
  31. 31
  32. 32
  33. 33
  34. 34

Decision letter

  1. Michael Levitt
    Reviewing Editor; Stanford University, United States

eLife posts the editorial decision letter and author response on a selection of the published articles (subject to the approval of the authors). An edited version of the letter sent to the authors after peer review is shown, indicating the substantive concerns or comments; minor concerns are not usually shown. Reviewers have the opportunity to discuss the decision before the letter is sent (see review process). Similarly, the author response typically shows only responses to the major concerns raised by the reviewers.

Thank you for submitting your work entitled “Codon-level information improves predictions of inter-residue contacts in proteins by correlated mutation analysis” for peer review at eLife. Your submission has been favorably evaluated by Aviv Regev (Senior Editor) and two reviewers, one of whom is a member of our Board of Reviewing editors. The Reviewing Editor has drafted this decision to help you prepare a revised submission.

The paper presents a potentially important advance in the area of deducing spatial contacts between amino acids in globular proteins, by proposing a simple yet elegant modification, which combines the analysis of amino acid and codon MSAs. While the results are potentially interesting, there are some important statistical analyses which are still lacking and need to addressed in the revision. In particular (1) the authors should incorporate statistical testing to assess the significance of the differences in accuracy with and without incorporating codon information. One would want to know if the improvement is robust, and it could be also that there are specific sub-classes where it is more/less beneficial. The scope of alternative methods to compare to should also be expanded. (2) The authors should show compelling analyses that address the question of over-fitting vs. generalization. This can be accomplished by cross-validation (the correct way to choose a set of parameters), and then with an additional test set (for seeing how the parameters chosen by cross validation perform on new data).

Essential revisions:

1) In Figure 1, there are differences shown for accuracy with and without incorporating codon information, but it is unclear if these differences are significant. Figure 1 shows the mean accuracy across the alignments, but it would be straightforward to also show or report the variance or standard deviation across the alignments. Then an explicit P-value could be put on a null hypothesis that incorporation of codon information increases or decreases accuracy. In general, the apparent result in Figure 1A that incorporation of codon information into SDCA either helps or hurts depending on how contacts are defined reduces confidence that incorporation of codon information is robustly improving the underlying algorithms.

2) In Figure 1B, the incorporation of codon information appears to help OMES and MI, but a more appropriate comparison might be to other algorithms that attempt to correct for phylogenetic artifacts (see for example http://bioinformatics.oxfordjournals.org/content/24/3/333.short; Mutual information without the influence of phylogeny or entropy dramatically improves residue contact prediction).

3) The data in Figure 2 should be on one plot. The figure needs to clearly emphasize that one set of lines are without the codon data and the other set are with it. This makes a huge difference and it the important conclusion of the work. The effect of how close distances are determined (Cβ or All) is minor and confusing.

4) In tuning a free parameter, such as is being performed in Figure 3, it is not surprising that a value for the free parameter alpha can be found that improves the performance. But does this just reflect over-fitting? It would perhaps be more informative to split the data into a training set and a test set, and evaluate the improvement on the test set which was not used to estimate alpha. By repeating this process on permutations of randomly chosen test sets and training sets, it might be possible to determine how stable estimates of alpha are or, alternatively, to demonstrate that the optimal value of alpha is highly variable from training set to training set.

https://doi.org/10.7554/eLife.08932.012

Author response

Essential revisions:

1) In Figure 1, there are differences shown for accuracy with and without incorporating codon information, but it is unclear if these differences are significant. Figure 1 shows the mean accuracy across the alignments, but it would be straightforward to also show or report the variance or standard deviation across the alignments. Then an explicit P-value could be put on a null hypothesis that incorporation of codon information increases or decreases accuracy. In general, the apparent result in Figure 1A that incorporation of codon information into SDCA either helps or hurts depending on how contacts are defined reduces confidence that incorporation of codon information is robustly improving the underlying algorithms.

We agree with the reviewers that an assessment of the significance of our results should have been provided. We have, therefore, added a new panel to Figure 3 that shows stacked bar plots of the number of MSAs for which including codon data improved the contact predictions using the different methods and the number of MSAs for which it was otherwise. The improvement achieved by incorporating codon data is shown to be significant as indicated by P-values obtained using the Wilcoxon signed-rank and sign tests.

2) In Figure 1B, the incorporation of codon information appears to help OMES and MI, but a more appropriate comparison might be to other algorithms that attempt to correct for phylogenetic artifacts (see for example http://bioinformatics.oxfordjournals.org/content/24/3/333.short; Mutual information without the influence of phylogeny or entropy dramatically improves residue contact prediction).

In the original manuscript, we showed that incorporation of codon information improves predictions by the OMES, MI and DCA methods. Following the comment of the reviewers, we tested whether incorporation of codon information improves predictions also by the MIp method and were surprised to see no significant improvement. In MIp, the background MI shared by positions i and j (termed APC) is subtracted from their MI score. It appears, therefore, that there is an overlap between the improvements generated by including the APC term and codon data. The PSICOV method, which appears to outperform all other methods, also includes the APC term. We, therefore, tested whether including codon data instead of the APC term improves the PSICOV method and found that there is a significant improvement. The effects of including codon data on the predictions by all 5 methods (OMES, MI, MIp, DCA and PSICOV) are shown in the revised Figure 2. Although incorporation of codon data does not improve the MIp method, it does lead to an improvement in the two state-of-the-art methods, DCA and PSICOV, which outperform all the older methods including MIp.

3) The data in Figure 2 should be on one plot. The figure needs to clearly emphasize that one set of lines are without the codon data and the other set are with it. This makes a huge difference and it the important conclusion of the work. The effect of how close distances are determined (Cβ or All) is minor and confusing.

Figure 2 was altered as suggested by the reviewers so that the comparison of all the methods with and without incorporating codon data is shown in the same figure using only the Cβ definition which is better at identifying contacts and is the most widely used. The comparison between Cβ and All has been removed as we agree that it is minor and confusing.

4) In tuning a free parameter, such as is being performed in Figure 3, it is not surprising that a value for the free parameter alpha can be found that improves the performance. But does this just reflect over-fitting? It would perhaps be more informative to split the data into a training set and a test set, and evaluate the improvement on the test set which was not used to estimate alpha. By repeating this process on permutations of randomly chosen test sets and training sets, it might be possible to determine how stable estimates of alpha are or, alternatively, to demonstrate that the optimal value of alpha is highly variable from training set to training set.

We followed the reviewers’ suggestion and split the MSAs into 1000 different training and test tests. The alpha was determined for each learning set and then used for predictions carried out for its corresponding test set. The results, which are presented in a new supplementary Figure (Figure 3–figure supplement 1), show that the estimates of alpha are indeed stable and that the prediction accuracies for the learning and test sets are similar.

https://doi.org/10.7554/eLife.08932.013

Article and author information

Author details

  1. Etai Jacob

    1. The Mina and Everard Goodman Faculty of Life Sciences, Bar-Ilan University, Ramat-Gan, Israel
    2. Department of Structural Biology, Weizmann Institute of Science, Rehovot, Israel
    Contribution
    EJ, Conception and design, Acquisition of data, Analysis and interpretation of data, Drafting or revising the article
    Competing interests
    The authors declare that no competing interests exist.
  2. Ron Unger

    The Mina and Everard Goodman Faculty of Life Sciences, Bar-Ilan University, Ramat-Gan, Israel
    Contribution
    RU, Conception and design, Analysis and interpretation of data, Drafting or revising the article
    For correspondence
    ron@biomodel.os.biu.ac.il
    Competing interests
    The authors declare that no competing interests exist.
  3. Amnon Horovitz

    Department of Structural Biology, Weizmann Institute of Science, Rehovot, Israel
    Contribution
    AH, Conception and design, Analysis and interpretation of data, Drafting or revising the article
    For correspondence
    Amnon.Horovitz@weizmann.ac.il
    Competing interests
    The authors declare that no competing interests exist.

Funding

Israel Science Foundation (ISF) (772/13)

  • Ron Unger

Israel Science Foundation (ISF) (158/12)

  • Amnon Horovitz

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

Acknowledgements

This work was supported by grants 158/12 (to A. H.) and 772/13 (to R. U.) of the Israel Science Foundation. A.H. is an incumbent of the Carl and Dorothy Bennett Professorial Chair in Biochemistry.

Reviewing Editor

  1. Michael Levitt, Reviewing Editor, Stanford University, United States

Publication history

  1. Received: May 22, 2015
  2. Accepted: September 13, 2015
  3. Accepted Manuscript published: September 15, 2015 (version 1)
  4. Accepted Manuscript updated: September 25, 2015 (version 2)
  5. Version of Record published: October 13, 2015 (version 3)

Copyright

© 2015, Jacob et al.

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.

Metrics

  • 2,255
    Page views
  • 502
    Downloads
  • 4
    Citations

Article citation count generated by polling the highest count across the following sources: Crossref, PubMed Central, Scopus.

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Download citations (links to download the citations from this article in formats compatible with various reference manager tools)

Open citations (links to open the citations from this article in various online reference manager services)

Further reading