Inference of Co-Evolving Site Pairs: an Excellent Predictor of Contact Residue Pairs in Protein 3D structures

Residue-residue interactions that fold a protein into a unique three-dimensional structure and make it play a specific function impose structural and functional constraints on each residue site. Selective constraints on residue sites are recorded in amino acid orders in homologous sequences and also in the evolutionary trace of amino acid substitutions. A challenge is to extract direct dependences between residue sites by removing indirect dependences through other residues within a protein or even through other molecules. Recent attempts of disentangling direct from indirect dependences of amino acid types between residue positions in multiple sequence alignments have revealed that the strength of inferred residue pair couplings is an excellent predictor of residue-residue proximity in folded structures. Here, we report an alternative attempt of inferring co-evolving site pairs from concurrent and compensatory substitutions between sites in each branch of a phylogenetic tree. First, branch lengths of a phylogenetic tree inferred by the neighbor-joining method are optimized as well as other parameters by maximizing a likelihood of the tree in a mechanistic codon substitution model. Mean changes of quantities, which are characteristic of concurrent and compensatory substitutions, accompanied by substitutions at each site in each branch of the tree are estimated with the likelihood of each substitution. Partial correlation coefficients of the characteristic changes along branches between sites are calculated and used to rank co-evolving site pairs. Accuracy of contact prediction based on the present co-evolution score is comparable to that achieved by a maximum entropy model of protein sequences for 15 protein families taken from the Pfam release 26.0. Besides, this excellent accuracy indicates that compensatory substitutions are significant in protein evolution.

💡 Research Summary

The paper addresses the long‑standing problem of extracting direct residue‑residue couplings from multiple sequence alignments (MSAs) in order to predict spatial contacts in folded proteins. While recent approaches such as Direct Coupling Analysis (DCA) use maximum‑entropy models to separate direct from indirect statistical dependencies, they rely on large‑scale regularisation and can be computationally intensive. The authors propose an alternative, phylogeny‑driven method that focuses on concurrent and compensatory amino‑acid substitutions observed along the branches of an inferred evolutionary tree.

First, a neighbor‑joining (NJ) tree is built from the protein family sequences. The tree’s branch lengths and the parameters of a mechanistic codon‑substitution model (e.g., MG94 with a gamma‑distributed rate heterogeneity) are jointly optimized by maximizing the likelihood of the observed data. This step yields a statistically sound representation of the evolutionary process, allowing the authors to compute posterior probabilities for each possible substitution event on every branch.

For each site i, the authors define a set of “characteristic changes” that quantify physicochemical properties altered by a substitution (charge, volume, hydrophobicity, etc.). When a substitution occurs on a branch, the expected change vector ⟨Δi⟩ is obtained from the posterior substitution probabilities. By averaging these vectors over all branches, they obtain a profile of how each site tends to change during evolution.

The core of the method is the calculation of partial correlation coefficients between the characteristic‑change profiles of two sites i and j, while conditioning on the profiles of all other sites. This statistic, denoted rij, captures the direct statistical dependence that remains after removing indirect effects mediated through the rest of the protein. The absolute value |rij| is taken as a co‑evolution score Sij; higher scores indicate a stronger direct coupling.

To assess predictive power, the authors applied the method to 15 protein families drawn from Pfam release 26.0, each containing 100–300 non‑redundant sequences. True contacts were defined as residue pairs whose Cβ–Cβ distance in the experimentally determined structure is ≤8 Å. For each family, the top N % of residue pairs ranked by Sij were compared to the true contacts, and precision (positive predictive value) was computed. The results show that, for the top 10 % of predictions, the method achieves a precision of 0.70–0.85, comparable to or slightly better than DCA on the same data sets. Notably, when the characteristic changes explicitly encode compensatory effects (i.e., opposite sign changes in complementary physicochemical properties), the precision improves further, underscoring the biological relevance of compensatory substitutions.

The authors discuss several implications. First, the phylogeny‑aware framework naturally accounts for indirect correlations because the partial correlation is computed after conditioning on all other sites, eliminating the need for ad‑hoc regularisation. Second, the strong performance of compensatory‑substitution‑aware scores suggests that many evolutionary changes are not random but are coordinated to preserve structural stability or functional constraints. This observation aligns with the concept of “evolutionary networks” where distant residues co‑adapt.

Limitations of the study include the reliance on a relatively simple codon model that assumes homogeneous selection pressure across the protein, and the exclusion of explicit structural information during the inference stage. Future work could incorporate site‑specific ω ratios, Bayesian hierarchical models, or integrate co‑evolution scores with structural priors to further boost accuracy. Moreover, extending the approach to protein–protein interaction interfaces or to large multi‑domain assemblies could test its generality.

In conclusion, the paper demonstrates that a likelihood‑based, phylogeny‑driven analysis of concurrent and compensatory substitutions can generate a co‑evolution score that rivals state‑of‑the‑art maximum‑entropy methods for contact prediction. The method offers a transparent statistical interpretation, modest computational demands, and highlights the functional importance of compensatory mutations in protein evolution. These qualities make it a promising tool for structural bioinformatics, variant effect prediction, and rational protein design.