A Bayesian Sampler for Optimization of Protein Domain Hierarchies

Abstract

The process of identifying and modeling functionally divergent subgroups for a specific protein domain class and arranging these subgroups hierarchically has, thus far, largely been done via manual curation. How to accomplish this automatically and optimally is an unsolved statistical and algorithmic problem that is addressed here via Markov chain Monte Carlo sampling. Taking as input a (typically very large) multiple-sequence alignment, the sampler creates and optimizes a hierarchy by adding and deleting leaf nodes, by moving nodes and subtrees up and down the hierarchy, by inserting or deleting internal nodes, and by redefining the sequences and conserved patterns associated with each node. All such operations are based on a probability distribution that models the conserved and divergent patterns defining each subgroup. When we view these patterns as sequence determinants of protein function, each node or subtree in such a hierarchy corresponds to a subgroup of sequences with similar biological properties. The sampler can be applied either de novo or to an existing hierarchy. When applied to 60 protein domains from multiple starting points in this way, it converged on similar solutions with nearly identical log-likelihood ratio scores, suggesting that it typically finds the optimal peak in the posterior probability distribution. Similarities and differences between independently generated, nearly optimal hierarchies for a given domain help distinguish robust from statistically uncertain features. Thus, a future application of the sampler is to provide confidence measures for various features of a domain hierarchy.

Evaluating, Comparing, and Interpreting Protein Domain Hierarchies

“…The MAP, which corresponds to the maximum LLR (denoted here as LLR MAP ), is the value to which the sampler (with simulated annealing) tends to converge, though this is not guaranteed. As illustrated by Neuwald (2014), however, obtaining the same maximum LLR in multiple sampling runs using distinct random seeds suggests that a nearly optimal hierarchy is typically found. Of course, when comparing such MAP-associated LLRs for distinct hierarchies, it is essential that the same input alignment be used throughout.…”

“…We define as corresponding those pairs of nodes or subtrees, one from each hierarchy, with ‡ 90% overlap in their sequence sets; as roughly corresponding such pairs with > 50% but < 90% overlap; and as noncorresponding the remaining nodes. We represent sequence sets using a previously described data structure (Neuwald and Green, 1994) facilitating efficient bitwise set operations.…”

“…Moreover, our recent article in this journal (Neuwald, 2014) describes a procedure for obtaining an optimal or nearly optimal domain hierarchy via Bayesian Markov chain Monte Carlo (MCMC) sampling (Liu, 2008). This approach takes as input a multiple alignment of, in principle, every available sequence for a given protein domain.…”

“…The statistical model for protein domain hierarchies is described in a previous article (Neuwald, 2014) and thus is merely summarized here. Such a model is illustrated in Figure 1, as follows: Figure 1A shows a protein domain hierarchy consisting of N = 12 nodes and of M = 12 functionally divergent subgroups, each of which corresponds to a subtree.…”

“…The inner product of two vectors is denoted by AE Á‚Á ae. The third through seventh terms in Equation 1 correspond to the logarithm of the product of the prior probabilities, as defined in Neuwald (2014). …”