Generalized Centroid Estimators in Bioinformatics

In a number of estimation problems in bioinformatics, accuracy measures of the target problem are usually given, and it is important to design estimators that are suitable to those accuracy measures. However, there is often a discrepancy between an employed estimator and a given accuracy measure of the problem. In this study, we introduce a general class of efficient estimators for estimation problems on high-dimensional binary spaces, which representmany fundamental problems in bioinformatics. Theoretical analysis reveals that the proposed estimators generally fit with commonly-used accuracy measures (e.g. sensitivity, PPV, MCC and F-score) as well as it can be computed efficiently in many cases, and cover a wide range of problems in bioinformatics from the viewpoint of the principle of maximum expected accuracy (MEA). It is also shown that some important algorithms in bioinformatics can be interpreted in a unified manner. Not only the concept presented in this paper gives a useful framework to design MEA-based estimators but also it is highly extendable and sheds new light on many problems in bioinformatics.


Introduction
In estimation problems in bioinformatics, the space of solutions is generally large and often high-dimensional. Among them, a number of fundamental problems in bioinformatics, such as alignment of biological sequences, prediction of secondary structures of RNA sequences, prediction of biological networks, and estimation of phylogenetic trees, are classified into estimation problems whose solutions are in a high-dimensional binary space. Such problems are generally difficult to solve, and the estimates are often unreliable.
The popular solutions for these problems, such as for the secondary structure of RNA with minimum free energy, are the maximum likelihood (ML) estimators. The ML estimator maximizes the probability that the estimator is exactly correct, but that probability is generally very small. Noticing the drawbacks of the ML estimators, Carvalho and Lawrence have proposed the centroid estimator, which represents an ensemble of all the possible solutions and minimizes the expected Hamming loss of the prediction [1].
In this paper, we conduct a theoretical analysis of estimation problems in high-dimensional binary space, and present examples and solutions in bioinformatics. The theories in this paper provide a unified framework for designing superior estimators for estimation problems in bioinformatics. The estimators discussed in this paper, including the ML estimator and the centroid estimator, are formalized as maximum expected gain (MEG) estimators, which maximize the estimator-specific gain functions with respect to the given probability distribution. The objective of the estimation is not always to find the exact solution with an extremely small probability or to find the solution with the minimum Hamming loss, but rather to find the most accurate estimator. Therefore, we adopt the principle of maximum expected accuracy (MEA), which has been successfully applied to various problems in bioinformatics, such as the alignment of biological sequences [2][3][4], the secondary structure prediction of RNA [5][6][7][8] and other applications [9][10][11].
Theoretical analysis, however, shows that those MEA estimators are not always robust with respect to accuracy measures. To address this, we previously proposed the c-centroid estimator in a few specific problems [4,12]. In this paper, in order to make the ccentroid estimator easily applicable to other estimation problems, we introduce an abstract form of the c-centroid estimator, which is defined on general binary spaces and designed to fit to the commonly used accuracy measures. The c-centroid estimator is a generalization of the centroid estimator, and offers a more robust framework for estimators than the previous estimators. We extend the theory of maximum expected gain (MEG) estimators and ccentroid estimators for two advanced problems: the estimators that represent the common solutions for multiple entries, and the estimators for marginalized probability distributions.

Materials and Methods
Problem 1 (Pairwise alignment of two biological sequences) Given a pair of biological (DNA, RNA, protein) sequences x and x', predict their alignment as a point in A(x,x'), the space of all the possible alignments of x and x'.
Problem 2 (Prediction of secondary structures of RNA sequences) Given an RNA sequence x, predict its secondary structure as a point in S(x), the space of all the possible secondary structures of x.
A point in A(x,x'), can be represented as a binary vector of DxDDx'D dimensions by denoting the aligned bases across the two sequences as ''1'' and the remaining pairs of bases as ''0''. A point in S(x) can also be represented as a binary vector of DxD(DxD{1)=2 dimensions, which represent all the pairs of the base positions in x, by denoting the base pairs in the secondary structures as ''1''. In each problem, the predictive space (A(x,x') or S(x)) is a subset of a binary space (f0,1g DxDDx'D or f0,1g DxD(DxD{1)=2 ) because the combinations of aligned bases or base pairs are restricted (see ''Discrete (binary) spaces in bioinformatics'' in Appendices for more formal definitions). Therefore, Problem 1 and Problem 2 are special cases of the following more general problem: Problem 3 (Estimation problem on a binary space) Given a data set D and a predictive space Y (a set of all candidates of a prediction), which is a subset of n-dimensional binary vectors f0,1g n , that is, Y 5f0,1g n , predict a point y in the predictive space Y .
Not only Problem 1 and Problem 2 but also a number of other problems in bioinformatics are formulated as Problem 3, including the prediction of biological networks and the estimation of phylogenetic trees (Problem 4).
To discuss the stochastic character of the estimators, the following assumption is introduced.
Assumption 1 (Existence of probability distribution) In Problem 3, there exists a probability distribution p(yDD) on the predictive space Y .
For Problem 3 with Assumption 1, we have the following Bayesian maximum likelihood (ML) estimator.
Definition 1 (Bayesian ML estimator [1]) For Problem 3 with Assumption 1, the estimator y y (ML)~a rgmax which maximizes the Bayesian posterior probability p(yDD), is referred to as a Bayesian maximum likelihood (ML) estimator.
For problems classified as Problem 3, Bayesian ML estimators have dominated the field of estimators in bioinformatics for years. The classical solutions of Problem 1 and Problem 2 are regarded as Bayesian ML estimators with specific probability distributions, as seen in the following examples.
Example 1 (Pairwise alignment with maximum score) In Problem 1 with a scoring model (e.g., gap costs and a substitution matrix), the distribution p(yDD) in Assumption 1 is derived from the Miyazawa model [13] (See ''Probability distributions p (a) (hDx,x') on A(x,x')'' in Appendices), and the Bayesian ML estimator is equivalent to the alignment that has the highest similarity score.
Example 2 (RNA structure with minimum free energy) In Problem 2 with a McCaskill energy model [14], the distribution p(yDD) in Assumption 1 can be obtained with the aid of thermodynamics (See ''Probability distributions p (s) (hDx) on S(x)'' in Appendices for details), and the Bayesian ML estimator is equivalent to the secondary structure that has the minimum free energy (MFE).
When a stochastic model such as a pair hidden Markov model (pair HMM) in Problem 1 or a stochastic context-free grammar (SCFG) in Problem 2 is assumed in such problems, the distribution and the ML estimator are derived in a more direct manner.
The Bayesian ML estimator regards the solution which has the highest probability as the most likely one. To provide more general criteria for good estimators, here we define the gain function that gives the gain for the prediction, and the maximum expected gain (MEG) estimator that maximizes the expected gain.
Definition 2 (Gain function) In Problem 3, for a point h [ Y and its prediction y [ Y , a gain function is defined as G : Y |Y ?R z , G(h,y).
Definition 3 (MEG estimator) In Problem 3 with Assumption 1, the maximum expected gain (MEG) estimator is defined aŝ If the gain function is designed according to the accuracy measures of the target problem, the MEG estimator is considered as the maximum expected accuracy (MEA) estimator, which has been successfully applied in bioinformatics (e.g., [9]).Although in estimation theory a loss function that should be minimized is often used, in order to facilitate the understanding of the relationship with the MEA, in this paper, we use a gain function that should be maximized.
The MEG estimator for the gain function d(y,h) is the ML estimator. Although this means that the ML estimator maximizes the probability that the estimator is identical to the true value, there is an extensive collection of suboptimal solutions and the probability of the ML estimator is extremely small in cases where n in Problem 3 is large. Against this background, Carvalho and Lawrence proposed the centroid estimator, which takes into account the overall ensemble of solutions [1]. The centroid estimator can be defined as an MEG estimator for a pointwise gain function as follows: Definition 4 (Pointwise gain function) In Problem 3, for a point h [ Y and its prediction y~fy i g n i~1 [ Y , a gain function G(h,y) written as where F i : Y |f0,1g?R z (i~1,2, . . . ,n), is referred to as a pointwise gain function. Definition 5 (Centroid estimator [1]) In Problem 3 with Assumption 1, a centroid estimator is defined as an MEG estimator for the pointwise gain function given in Eq. (1) by defining F i (h,y i )Ĩ (h i~1 )I(y i~1 )zI(h i~0 )I(y i~0 ). Throughout this paper, I( : ) is the indicator function that takes a value of 1 or 0 depending on whether the condition constituting its argument is true or false. The centroid estimator is equivalent to the expected Hamming loss minimizer [1]. If we can maximize the pointwise gain function independently in each dimension, we can obtain the following consensus estimator, which can be easily computed.
Definition 6 (Consensus estimator [1]) In Problem 3 with Assumption 1, the consensus estimatorŷ y (c)~fŷ y (c) i g n i~1 for a pointwise gain function is defined aŝ The consensus estimator is generally not contained within the predictive space Y since the predictive space Y usually has complex constraints for each dimension (see ''Discrete (binary) spaces in bioinformatics'' in Appendices). Carvalho and Lawrence proved a sufficient condition for the centroid estimator to contain the consensus estimator (Theorem 2 in [1]). Here, we present a more general result, namely, a sufficient condition for the MEG estimator for a pointwise function to contain the consensus estimator.
Theorem 1 In Problem 3 with Assumption 1 and a pointwise gain function, let us suppose that a predictive space Y can be written as where C k is defined as for an index-set I k 5f1,2, . . . ,ng. If the pointwise gain function in Eq. (1) satisfies the condition

c-centroid estimator: generalized centroid estimator
In Problem 3, the ''1''s and the ''0''s in the binary vector of a prediction y can be interpreted as positive and negative predictions, respectively. The respective numbers of true positives (TP), true negatives (TN), false positives (FP) and false negatives (FN) for a point h and its prediction y are denoted by TP (h,y), TN (h,y), FP(h,y) and FN (h,y), respectively (See also Eqs (15)- (18)).
To design a superior MEG estimator, it is natural to use a gain function of the following form, which yields positive scores for the number of true predictions (TP and TN) and negative scores for those of false predictions (FP and FN): where a k is a positive constant (k~1,2,3,4). Note that this gain function is a pointwise gain function. This gain function is naturally compatible with commonly used accuracy measures such as sensitivity, PPV, MCC and F-score (a function of TP, TN, FP and FN; see ''Evaluation measures defined using TP, TN, FP and FN'' in Appendices for definitions). The following Definition 7 and Theorem 2 characterize the MEG estimator for this gain function.
Definition 7 (ª-centroid estimator) In Problem 3 with Assumption 1 and a fixed c §0, the c-centroid estimator is defined as the MEG estimator for the pointwise gain function given in Eq. (1) by Theorem 2 The MEG estimator for the gain function in Eq. (4) is equivalent to a c-centroid estimator with c~a 1 za 4 a 2 za 3 .
Theorem 2 (see Appendices for a formal proof) is derived from the following relations: The c-centroid estimator maximizes the expected value of TNzcTP, and includes the centroid estimator as a special case where c~1. The parameter c adjusts the balance between the gain from true negatives and that from true positives.
The expected value of the gain function of the c-centroid estimator is computed as follows (see Appendices for the derivation): where Since the second term in Eq. (6) does not depend on y, the ccentroid estimator maximizes the first term. The following theorem is obtained by assuming the additional condition described below.
Theorem 3 In Problem 3 with Assumption 1, the predictive space Y satisfies the following condition: if Then, the c-centroid estimator is equivalent to the estimator that maximizes the sum of marginalized probabilities p i that are greater than 1=(cz1) in the prediction.
The condition is necessary to obtain 0 for the i that produces negative values in the first term in Eq. (6). Problem 2, Problem 1, and many other typical problems in bioinformatics satisfy this condition. Because the pointwise gain function of the c-centroid estimator satisfies Eq. (3) in Theorem 1, we can prove the following Corollary 1.
Corollary 1 (ª-centroid estimator for 0ƒªƒ1) In Problem 3 with Assumption 1, the predictive space Y is given in the same form in Eq.
(2) of Theorem 1. Then, the c-centroid estimator for c [ ½0,1 contains its consensus estimator. Moreover, the consensus estimator is identical to the following estimator y Ã~f y Ã i g: where Here, p i is the marginalized probability of the distribution for the i-th dimension of the predictive space. In Problem 1, it is known as the alignment probability, which is defined as the probability of each pair of positions across the two sequences being aligned. In Problem 2, it is known as the base pairing probability, which is defined as the probability of each pair of positions forming a base pair in the secondary structure. These marginalized probabilities can be calculated by using dynamic programming algorithms, such as the forward-backward algorithm and the McCaskill algorithm, depending on the model of the distributions. (see ''Probability distributions on discrete spaces'' in Appendices for those distributions).
Corollary 1 does not hold for cw1, but in typical problems in bioinformatics the c-centroid estimator for cw1 can be calculated efficiently by using dynamic programming, as shown in the following examples.
Example 3 (ª-centroid estimator of pairwise alignment) In Problem 1 with Assumption 1, the c-centroid estimator maximizes the sum of the alignment probabilities which are greater than 1=(cz1) (Theorem 3), and for c [ ½0,1 it can be given as the consensus estimator calculated from Eq. (8) (Corollary 1). For cw1, the c-centroid estimator is obtained by using a dynamic programming algorithm with the same type of iterations as in the Needleman-Wunsch algorithm: where M i,k stores the optimal value of the alignment between two subsequences, x 1 Á Á Á x i and x' 1 Á Á Á x' k (see ''Secondary structure prediction of an RNA sequence (Problem 2)'' in Appendices for detailed descriptions). Example 4 (ª-centroid estimator for prediction of secondary structures) In Problem 2 with Assumption 1, the ccentroid estimator maximizes the sum of the base pairing probabilities that are greater than 1=(cz1) (Theorem 3), and for c [ ½0,1 it can be given as the consensus estimator calculated from Eq. (8) (Corollary 1). For cw1, the ccentroid estimator is obtained with the aid of a dynamic programming algorithm with the same type of iterations as in the Nussinov algorithm: where M i,j stores the best score of the sub-sequence x i x iz1 Á Á Á x j (see ''Pairwise alignment of biological sequences (Problem 1)'' in Appendices for the detail descriptions).
The c-centroid estimators are implemented in LAST [4] for Problem 1 and in CentroidFold [12,15] for Problem 2.
Problem 4 (Estimation of phylogenetic trees) Given a set of operational taxonomic units S, predict their phylogenetic trees (unrooted and multi-branched trees) as a point in T (S), the space of all the possible phylogenetic trees of S.
The phylogenetic tree in T (S) is represented as a binary vector with 2 n{1 {n{1 dimension where n is the number of units in S, based on partition of S by cutting every edge in the tree (see ''The space of phylogenetic trees: T (S)'' in Appendices for details). A sampling algorithm can be used to estimate the partitioning probabilities approximately [16].
Example 5 (ª-centroid estimator of phylogenetic estimation) In Problem 4 with Assumption 1, the c-centroid estimator maximizes the number of the partitioning probabilities which are greater than 1=(cz1) (Theorem 3), and for c [ ½0,1 it can be give as the consensus estimator calculated from Eq. (8) (Corollary 1) (see ''Estimation of phylogenetic trees (Problem 4)'' in Appendices for details).
Because the Hamming distance between two trees in T (S) is known as topological distance [17], the 1-centroid estimator minimizes the expected topological distance. In contrast to Example 3 and Example 4, it appears that no method can efficiently compute the c-centroid estimator with cw1 in Example 5. Despite the difficulties of the application to phylogenetic trees, recently, a method applying the concept of generalized centroid estimators was developed [18].

Generalized centroid estimators for representative prediction
Predictions based on probability distributions on the predictive space were discussed in the previous sections. However, there are certain even more complex problems in bioinformatics, as illustrated by the following example.
Problem 5 (Prediction of common secondary structures of RNA sequences) Given a set of RNA sequences D~fx i g,i~1, . . . K and their multiple alignment of length L and the same energy model for each RNA sequence, predict their common secondary structure as a point in S 0 (L), which is the space of all possible secondary structures of length L.
In the case of Problem 5, although the probability distribution is not implemented in the predictive space, each RNA sequence x i has a probability distribution on its secondary structure derived from the energy model. Therefore, the theories presented in the previous section cannot be applied directly to this problem. However, if we devise a new type of gain function that connects the predictive space with the parameter space of the secondary structure of each RNA sequence, we can calculate the expected gain over the distribution on the parameter spaces of RNA sequences. In order to account for this type of problem in general, we introduce Assumption 2 and Definition 8 as follows.
Assumption 2 In Problem 3 there exists a probability distribution p(hDD) on the parameter space H which might be different from the predictive space Y . It should be emphasized that the MEG estimator (Definition 3), pointwise gain function (Definition 4) and Theorem 1 can be extended to the generalized gain function.
In the case of Problem 5, for example, the parameter space is the product of the spaces of the secondary structures of each RNA sequence, and the probability distribution is the product of the distributions of secondary structures of each RNA sequence. Here, the general form of the problem of representative prediction is introduced.
Problem 6 (Representative prediction) In Problem 3 with Assumption 2, if the parameter space is represented as a product space (H~P K k~1 H (k)~Y K ) and the distribution of h [ H has the form p(hDD)~P K k~1 p (k) (h k DD), predict a point y in the predictive space Y . The generalized gain function for the representative prediction should be chosen such that the prediction reflects as much as each data entry. Therefore, it is natural to use the following generalized gain function that integrates the gain for each parameter.
Definition 9 (Homogeneous generalized gain function) In Problem 6, a homogeneous generalized gain function is defined as where G' is the gain function in Definition 2.
Definition 10 (Representative estimator) In Problem 6, given a homogeneous generalized gain function G(h,y)~P K k~1 G'(h k ,y), the MEG estimator defined aŝ is referred to as the representative estimator.

Proposition 1
The representative estimator is equivalent to an MEG estimator with averaged probability distribution on the predictive space Y : and a gain function G'. This proposition shows that a representative prediction problem with any homogeneous generalized gain function can be solved in a manner similar to Problem 3 (H~Y ) with averaged probability distribution. Therefore, the c-centroid estimator for a representative prediction satisfies Corollary 2.
Corollary 2 In Problem 6, the representative estimator where G'(h k ,y) is the gain function of the c-centroid estimator on Y , is the c-centroid estimator for the averaged probability distribution and satisfies the same properties in Theorem 2, Theorem 3, and Corollary 1.

Estimators based on marginal probabilities
In the previous section, we introduced Assumption 2, where there is a parameter space H that can be different from the predictive space Y , and we discussed the problem of representative prediction. In this section, we discuss another type of problems where H=Y . An example is presented below.
Problem 7 (Pairwise alignment using homologous sequences) Given a data set D~fx,x',hg, where x and x' are two biological sequences to be aligned and h is a sequence that is homologous to both x and x', predict a point y in the predictive space Y~A(x,x') (the space of all possible alignments of x and x').
The precise probabilistic model of this problem might include the phylogenetic tree, ancestor sequences and their alignments. Here, we assume a simpler situation where the probability distribution of all possible multiple alignments of D is given. We predict the pairwise alignment of two specific sequences according to the probability distribution of multiple alignments. Although the parameter space H, which is the space of all the possible multiple alignments, can be parametrized using the parameters of the spaces of the alignments of all pairs that can be formed from the sequences in D, H itself is not the product space of these spaces because these pairwise alignments are not independent: for x,x',h [ D, x i must be aligned to x' j if both x i and x' j are aligned to h k . This type of problems can be generalized as follows.
Problem 8 (Prediction in a subspace of the parameter space) In Problem 3 with Assumption 2, if the parameter space H is represented as H5H'|H' \ , predict a point y in the predictive space Y~H'.
For the problem of representative prediction (Problem 6), generalized gain functions on H|Y were introduced (Definition 8 and Definition 9). In contrast, in Problem 8, the values of the parameters in H' \ are not important, and a point in Y~H' is predicted. In Problem 7, for example, the optimal multiple alignment of D, the pairwise alignment of h and x, and the pairwise alignment of h and x' are irrelevant, but instead we predict the pairwise alignment of x and x'. The MEG estimator for the gain function defined on H'|Y can be written aŝ From the above MEG estimator, it might appear that Problem 8 is trivial. However, it is not a simple task to calculate the marginalized distribution in Eq. (11) in actual problems.
To reduce the computational cost, we change Problem 8 by introducing an approximated probability distribution on the product space H'|H' \ a follows.
Problem 9 (Prediction in product space) In Problem 3 with Assumption 2, if the parameter space H is represented as H~H'|H' \ and the probability distribution on H is defined as predict a point y in the predictive space Y~H'. This factorization of spaces and probability distributions creates a number of inconsistencies in the parameter space with respect to the original Problem 8. In other words, the approximated distribution yields non-zero values for a point that is not included in the original H (in Problem 8) but in H'|H' \ . To reduce these inconsistencies, a new type of gain function and a new estimator are introduced as follows.
Definition 11 (ª-type pointwise gain function) In Problem 8, a c-type pointwise gain function is defined as G(h,y) in Eq. (1) in Definition 4 having where the value d i (h') [ ½0,1 in the gain function should be designed to reduce the inconsistencies resulting from the factorization. Definition 12 (Approximated ª-type estimator) In Problem 9, with a c-type pointwise gain function with F i (h,y i ) in Eq. (13) on H|Y , an approximated c-type estimator is defined as an MEG estimator:ŷ Example 6 (PCT in pairwise alignment) We obtain the approximate estimator for Problem 7 with the following settings. The parameter space is given as H~H'|H' \ , where and the probability distribution on the parameter space H is given as The approximated c-type estimator for this c-type pointwise gain function is employed in a part of probabilistic consistency transformation (PCT) [19], which is an important step toward accurate multiple alignments. See ''Pairwise alignment using homologous sequences'' in Appendices for precise descriptions.
It is easily seen that Theorem 3 applies to the approximated ctype estimator if p i in Theorem 3 is changed as follows: Moreover, to confirm whether approximated c-type estimator contains the consensus estimator for the same gain function, it is only necessary to check if

Discussion
Properties of the c-centroid estimator In this paper, general criteria for designing estimators are given by the maximum expected gain (MEG) estimator (Definition 3). The Bayesian ML estimator is an MEG estimator with the delta function d(y,h) as the gain function, which means that only the probability for the ''perfect match'' is counted. To overcome the drawbacks of the Bayesian ML estimator, the centroid estimator [1] takes into account the overall ensemble of solutions and minimizes the expected Hamming loss. Because the Hamming loss is not the standard evaluation measures for actual problems, we have proposed an estimator of a more general type, the c-centroid estimator (Definition 7), which includes the centroid estimator as a special case, c~1. The c-centroid estimator is an MEG estimator that maximizes the expected value of TNzcTP, which generally covers all possible linear combination of the numbers of true positives (TP), true negatives (TN), false positives (FP) and false negatives (FN) (Theorem 2). Since most of the evaluation measures of the prediction accuracy are functions of these numbers [20], the c-centroid estimator is related to the principle of maximum expected accuracy (MEA). It should be noted that MEG estimators have been proposed that are similar to the c-centroid estimator for some specific problems, for example, the alignment metric accuracy (AMA) estimator [21] (see Appendices for the formal definition) for pairwise alignment (Problem 1) and the MEA-based estimator [5] (see Appendices for the formal definition) for prediction of secondary structure of RNA (Problem 2). However, these estimators display a bias with respect to the accuracy measures for the problem (see Eqs. (20) and (22)), and are therefore inappropriate from the viewpoint of the principles of MEA. Moreover, these estimators cannot be introduced in a general setting, that is, Problem 3. It has been also shown that the c-centroid estimator outperforms the MEA-based estimator [5] for various probability distributions in computational experiments [12]. (See ''Pairwise alignment of biological sequences (Problem 1)'' and ''Secondary structure prediction of an RNA sequence (Problem 2)'' in Appendices for relations between the c-centroid estimator and other estimators in Problems 1 and 2, respectively.)

How to determine the parameter in c-centroid estimator
The parameter c in c-centroid estimators adjusts sensitivity and PPV (whose relation is tradeoff). MCC or F-score is often used to obtain a balanced measure between sensitivity and PPV. In RNA secondary structure predictions, it has been confirmed that the best c (with respect to MCC) of the c-centroid estimator with CONTRAfold model was larger than that with McCaskill model [12]. It shows that the best c (with respect to a given accuracy measure) depends on not only estimation problems but also probabilistic models for predictive space. The parameter c trained by using reference structures was therefore employed as the default parameter in CentroidFold [12]. In order to select the parameter automatically (with respect to a given accuracy measure such as MCC and F-score), an approximation of maximizing expected MCC (or F-score) with the c-centroid estimator can be utilized [22].

Accuracy measures and computational efficiency
The reader might consider that it is possible to design estimators that maximize the expected MCC or F-score which balances sensitivity (SEN) and positive predictive value (PPV). However, it is much more difficult to compute such estimators in comparison with the c-centroid estimator, as described below.
The expected value of the gain function of the c-centroid estimator can be written with marginalized probabilities as in Eq. (7), which can be efficiently computed by dynamic programming in many problems in bioinformatics, for example, the forwardbackward algorithm for alignment probabilities and the McCaskill algorithm for base pairing probabilities. Under a certain condition of the predictive space, which many problems in bioinformatics satisfy, the c-centroid estimator maximizes the sum of marginalized probabilities greater than 1=(cz1) (Theorem 3). Moreover, under an additional condition of the predictive space and the pointwise gain function, which again many problems in bioinformatics satisfy, the c-centroid estimators for c [ ½0,1 can be easily calculated as the consensus estimators, which collect in the binary predictive space the components that have marginalized probabilities greater than 1=(cz1) (Corollary 1). For cw1, there often exist dynamic programming algorithms that can efficiently compute the c-centroid estimators (Examples 4 & 3), but there are certain problems, such as Problem 4, which seem to have no efficient dynamic programming algorithms.
The gain function of the estimators that maximize MCC or Fscore, and also SEN or PPV contain multiplication and/or division of TP, TN, FP and FN, while the gain function of the c-centroid estimator contains only the weighted sums of these values (i.e., TNzc : TP). Therefore, the expected gain is not written with marginalized probabilities as in Eq. (7), and it is difficult to design efficient computational algorithms for those estimators. In predicting secondary structures of RNA sequences (Problem 2), for example, it is necessary to enumerate all candidate secondary structures or sample secondary structures for an approximation in order to compute the expected MCC/F-score of a predicted secondary structure.

Probability distributions are not always defined on predictive space
After discussing the standard estimation problems on a binary space where the probability distribution is defined on the predictive space, we have proposed a new category of estimation problems where the probability distribution is defined on a parameter space that differs from the predictive space (see Assumption 2). Two types of estimators for such problems, for example, estimators for representative prediction and estimators based on marginalized distribution, have been discussed.
Prediction of the common secondary structure from an alignment of RNA sequences (Problem 5) is an example of representative prediction. The probability distribution is not implemented in the predictive space, the space of common secondary structure, but each RNA sequence has a probability distribution for its secondary structure. Because the ''correct'' reference for the common secondary structure is not known in general, direct evaluation of the estimated common secondary structure is difficult. In the popular evaluation process for this problem, the predicted common secondary structure is mapped to each RNA sequence and compared to its reference structure. Using the homogeneous generalized gain function exactly implements this evaluation process and the MEG estimator for the averaged probability distribution is equivalent to the MEG estimator for homogeneous generalized gain function. Therefore, we can use the averaged base pairing probabilities according to the alignment as the distribution for the common secondary structure (see ''Common secondary structure prediction from a multiple alignment of RNA sequences'' in Appendices for detailed discussion). The representative estimator for Problem 5 is implemented in software CentroidAlifold. Another example of representative prediction is the ''alignment of alignments'' problem, which is the fundamental element of progressive multiple alignment of biological sequences. The evaluation process using the sum of pairs score corresponds to using the homogeneous generalized gain function. (see ''Alignment between two alignments of biological sequences'' in Appendices for detailed discussion).
Estimation problems of marginalized distributions can be formalized as prediction in a subspace of the parameter space (Problem 8). If we can calculate the marginalized distribution on the predictive space from the distribution on the parameter space, all general theories apply to the predictive space and the marginalized distribution. In actual problems, such as pairwise alignment using homologous sequences (Problem 7), however, computational cost for calculation of the marginalized probability is quite high. We introduced the factorized probability distribution (Eq. (12)) for approximation, the c-type pointwise gain function (Definition 11) to reduce the inconsistency caused by the factorization, and the approximated c-type estimator (Definition 12). In Problem 7, the probability consistency transformation (PCT), which is widely used for multiple sequence alignment, is interpreted as an approximated c-type estimator. Prediction of secondary structures of RNA sequences on the basis of homologous sequences [23] (see Problem 13 in Appendices) and pairwise alignment for structured RNA sequences are further examples of this type of problems.

Application of c-centroid estimator to cluster centroid
In case probability distribution on the predictive space is multimodal, c-centroid estimators can provide unreliable solutions. For example, when there are two clusters of secondary structures in predictive spaces and those structures are exclusive, the c-centroid estimator might give a ''chimeric'' secondary structure whose free energy is quite high. To avoid this situation, Ding et al. [24] proposed a notion of the cluster centroid, which is computed by the centroid estimator with a given cluster in a predictive space. We emphasize that the extension of cluster centroid by using ccentroid estimator is straightforward and would be useful.

Conclusion
In this work, we constructed a general framework for designing estimators for estimation problems in high-dimensional discrete (binary) spaces. The theory is regarded as a generalization of the pioneering work conducted by Carvalho and Lawrence, and is closely related to the concept of MEA. Furthermore, we presented several applications of the proposed estimators (see Table 1 for summary) and the underlying theory. The concept presented in this paper is highly extendable and sheds new light on many problems in bioinformatics. In future research, we plan to investigate further applications of the c-centroid and related estimators presented in this paper.

Discrete (binary) spaces in bioinformatics
In this section, we summarize three discrete spaces that appear in this paper. These discrete spaces are often used in the definition of the predictive spaces and the parameter spaces. It should be noted that every discrete space described below is identical in form to Eq. (2).
The space of alignments of two biological sequences: A(x,x'). We define a space of the alignments of two biological (DNA, RNA and protein) sequences x and x', denoted by A(x,x'), as follows. We set I and is defined by Here I is a set of index-sets: The inclusion y [ C(I (1) i ) means that position i in the sequence x aligns with at most one position in the sequence x' in the alignment y, y [ C(I (2) j ) means that position j in the sequence x' aligns with at most one position in the sequence x and y [ C(I (3) ikjl ) means the alignment (i,k) and (j,l) is not crossing. Note that A(x,x') depends on only the length of two sequences, namely, DxD and Dx'D.
The space of secondary structures of RNA: S(x). We define a space of the secondary structures of an RNA sequence x, denoted by S(x), as follows. We set I Here I is a set of index-sets The inclusion y [ C(I (1) i ) means that position i in the sequence x belongs to at most one base-pair in a secondary structure y, and y [ C(I (2) ijkl ) means two base-pairs whose relation is pseudo-knot are not allowed in y. Note that S(x) depends on only the length of the RNA sequence x, that is, DxD.
The space of phylogenetic trees: T (S). We define a space of phylogenetic trees (unrooted and multi-branch trees) of a set of S~f1, Á Á Á ,ng, denoted by T (S), as follows. We set f DX Dw1^DX Dvn{1g, as a base index set and we define binary variables h X for X [ I (0) by Type of estimator Data Type of estimator and is defined by Note that T (S) depends on only the number of elements in S. We now give several properties of T (S) that follow directly from the definition.
Lemma 2 The topological distance [17] between two phylogenetic trees T 1 and T 2 in T (S) is where I( : ) is the indicator function.
Remark 1 If we assume the additional condition P X h X( (4n{6){2n)=2~n{3, then T (S) is a set of binary trees.

Probability distributions on discrete spaces
We use three probability distributions in this paper. Probability distributions p (a) (hDx,x') on A(x,x'). For two protein sequences x and x', a probability distribution p (a) (hDx,x') over the space A(x,x'), which is the space of pairwise alignments of x and x' defined in the previous section, is given by the following models. 1. Miyazawa model [13] and Probalign model [25]: where S(h) is the score of an alignment h under the given scoring matrix (We define S(h)~P h ij~1 s(x i ,x j ){ (penalty for gaps) where s(x i ,x j ) is a score for the correspondence of bases x i and x j ), T is the thermodynamic temperature and Z(T) is the normalization constant, which is known as a partition function.

Pair Hidden Markov Model (pair HMM) [19]
: where p(s) is the initial probability of starting in state s, a(s i ?s iz1 ) is the transition probability from s i to s iz1 and b(o i Ds i ) is the omission probability for either a single letter or aligner residue pair o i in the state s i . 3. CONTRAlign (pair CRF) model [26]: where w is a parameter vector and f (h,x,x') is a vector of features that indicates the number of times each parameter appears, V(x,x') denotes the set of all possible alignments of x and x'. We do not describe the feature vectors and refer readers to the original paper [26].
Remark 2 Strictly speaking, the alignment space in the pair hidden Markov model and the CONTRAlign model consider the patterns of gaps. In these cases, we obtain the probability space on A(x,x') by a marginalization.
Probability distributions p (s) (hDx) on S(x). For an RNA sequence x, a probability distribution p (s) (hDx) over S(x), which is the space of secondary structures of x defined in the previous section is given by the following models.
1. McCaskill model [14]: This model is based on the energy models for secondary structures of RNA sequences and is defined by where E(h,x) denotes the energy of the secondary structure that is computed using the energy parameters of Turner Lab [27], k and T are constants and Z(x) is the normalization term known as the partition function.

Stochastic Context free grammars (SCFGs) model [28]
: where p(x,s) is the joint probability of generating the parse s and is given by the product of the transition and emission probabilities of the SCFG model and V'(x) is all parses of x, V(h) is all parses for a given h.

CONTRAfold (CRFs; conditional random fields) model [5]:
This model gives us the best performance on secondary structure prediction although it is not based on the energy model.
is all parses for a given h.

Probability distributions p (t) (hDS) on T(S).
A probability distribution p (t) (hDS) on T (S) is given by probabilistic models of phylogenetic trees, for example, [29,30]. Those models give a probability distribution on binary trees and we should marginalize these distributions for multi-branch trees.
It should be noted that these measures can be written as a function of TP, TN, FP and FN. See [20] for other evaluation measures.

Schematic diagrams of representative and approximated c-type estimators
The schematic diagrams of the MEG estimator (Definition 3), the representative estimator (Definition 10) and the approximated c-type estimator (Definition 12) are shown in Figure 1, Figure 2 and Figure 3, respectively.

Applications in bioinformatics
In this section we describe several applications to bioinformatics of the general theories. Some of these applications have already been published. In those cases, we briefly explain the applications and the readers should see the original paper for further descriptions as well as the computational experiments. All of the applications in this section are summarized in Table 1.
Pairwise alignment of biological sequences (Problem 1). The pairwise alignment of biological (DNA, RNA, protein) sequences (Problem 1) is another fundamental and important problem of sequence analysis in bioinformatics (cf. [31]).
The c-centroid estimator for Problem 1 can be introduced as follows: Estimator 1 (ª-centroid estimator for Problem,:align) For Problem 1, we obtain the c-centroid estimator where the predictive space Y is equal to A(x,x') and the probability distribution on Y is taken by p (a) (hDx,x').
First, Theorem 2 and the definition of A(x,x') lead to the following property.
Property 1 (A relation of Estimator 1 with accuracy measures) The c-centroid estimator for Problem 1 is suitable for the accuracy measures: SEN, PPV, MCC and F-score with respect to the aligned-bases in the predicted alignment.
Note that accurate prediction of aligned-bases is important for the analysis of alignments, for example, in phylogenetic analysis. Therefore, the measures in above are often used in evaluations of alignments e.g. [4].
The marginalized probability p ik~p is called the aligned-base (matching) probability in this paper. The aligned-base probability matrix fp ik g i,k can be computed by the forward-backward algorithm whose time complexity is equal to O(DxDDx'D) [31]. Now, Theorem 3 leads to the following property.
Property 2 (Computation of Estimator 1) The pairwise alignment of Estimator 1 is found by maximizing the sum of aligned-base probabilities p ik (of the aligned-bases in the predicted alignment) that are larger than 1=(cz1). Therefore, it can be computed by a Needleman-Wunsch-style dynamic programming (DP) algorithm [32] after calculating the aligned-base matrix fp ik g: where M i,k stores the optimal value of the alignment between two subsequences, The time complexity of the recursion of the DP algorithm in Eq.  By using Corollary 1, we can predict the pairwise alignment of Estimator 1 with c[½0,1 without using the DP algorithm in Eq. (19).
Property 3 (Computation of Estimator 1 with 0ƒªƒ1) The pairwise alignment of the c-centroid estimator can be predicted by collecting the aligned-bases whose probabilities are larger than 1=(cz1).
The genome alignment software called LAST (http://last.cbrc. jp/) [4,33] employs the c-centroid estimator accelerated by an Xdrop algorithm, and the authors indicated that Estimator 1 reduced the false-positive aligned-bases, compared to the conventional alignment (maximum score estimator).
Relations of Estimator 1 with existing estimators are summarized as follows:   For Problem 1, Schwartz et al. [21] proposed an Alignment Metric Accuracy (AMA) estimator, which is similar to the ccentroid estimator (see also [3]). The AMA estimator is a maximum gain estimator (Definition 3) with the following gain function.
In the above equation, G f §0 is a gap factor, which is a weight for the prediction of gaps. We refer to the function G (AMA) (h,y) as the gain function of the AMA estimator. In a similar way to that described in the previous section, we obtain a relation between G (AMA) (h,y) and G (centroid) (h,y) (the gain function of the c-centroid estimator).
If we set 1=G f~c , then we obtain and C h is a value which does not depend on y. If I(h ij 1~1 )I(y ij 2~1 )~1 for j 1 =j 2 , then we obtain I(h ij 1 1)I(y ij 1~0 )~1 and I(h ij 2~0 )I(y ij 2~1 )~1, and this means that (i, j 1 ) is an aligned pair that is a false negative and (i,j 2 ) is an aligned pair that is a false positive when h is a reference alignment and y is a predicted alignment. Therefore, the terms A(h, y) (in Eq. (20)) in the gain function of AMA are not appropriate for the evaluation measures SEN, PPV, MCC and F-score for aligned bases. In summary, the c-centroid estimator is suitable for the evaluation measures: SEN, PPV and F-score with respect to the aligned-bases while the AMA estimator is suitable for the AMA.
Secondary structure prediction of an RNA sequence (Problem 2). Secondary structure prediction of an RNA sequence (Problem 2) is one of the most important problems of sequence analysis in bioinformatics. Its importance has increased due to the recent discovery of functional non-coding RNAs (ncRNAs) because the functions of ncRNAs are closely related to their secondary structures [35].
c-centroid estimator for Problem 2 can be introduced as follows: Estimator 2 (ª-centroid estimator for Problem 2) For Problem 2, we obtain the c-centroid estimator (Definition 7) where the predictive space Y is equal to S(x) and the probability distribution on Y is taken by p (s) (hDx).
The general theory of the c-centroid estimator leads to several properties. First, the following property is derived from Theorem 2 and the definition of S(x).

Property 4 (A relation of Estimator 2 with accuracy measures)
The c-centroid estimator for Problem 2 is suitable for the widely-used accuracy measures of the RNA secondary structure prediction: SEN, PPV and MCC with respect to base-pairs in the predicted secondary structure.
Because the base-pairs in a secondary structure are biologically important, SEN, PPV and MCC with respect to base-pairs are widely used in evaluations of RNA secondary structure prediction, for example, [5,12,36].
The marginalized probability p ij~p is called a base-pairing probability. The base-paring probability matrix fp ij g ivj can be computed by the Inside-Outside algorithm whose time complexity is equal to O(DxD 3 ) where DxD is the length of RNA sequence x [14,31]. Then, Theorem 3 leads to the following property.
Property 5 (Computation of Estimator 2) The secondary structure of Estimator 2 is found by maximizing the sum of the base-pairing probabilities p ij (of the base-pairs in the predicted structure) that are larger than 1=(cz1). Therefore, it can be computed by a Nussinov-style dynamic programming (DP) algorithm [37] after calculating the base-pairing probability matrix fp ij g: where M i,j stores the best score of the sub-sequence x i x iz1 Á Á Á x j . If we replace ''(cz1)p ij {1'' with ''1'' in Eq. (21), the DP algorithm is equivalent to the Nussinov algorithm [37] that maximizes the number of base-pairs in a predicted secondary structure. The time complexity of the recursion of the DP algorithm in Eq. (21) is equal to O(DxD 3 ). Hence, the total computational cost for predicting the secondary structure of the ccentroid estimator remains O(DxD 3 ), which is the same time complexity as for standard software: Mfold [38], RNAFOLD [39] and RNASTRUCTURE [40].
By using Corollary 1, we can predict the secondary structure of Estimator 2 with c [ ½0,1 without using the DP algorithm in Eq. (21).
Property 6 (Computation of Estimator 2 with 0vªƒ1) The secondary structure of the c-centroid estimator with c [ ½0,1 can be predicted by collecting the base-pairs whose probabilities are larger than 1=(cz1).
The software CENTROIDFOLD [12,15] implements Estimator 2 with various probability distributions for the secondary structures, such as the CONTRAFOLD and McCASKILL models.
Relations of Estimator 2 with other estimators are summarized as follows: For Problem 2, Do et al. [5] proposed an MEA-based estimator, which is similar to the c-centroid estimator. (The MEA-based estimator was also used in a recent paper [6].) The MEA-based estimator is defined by the maximum expected gain estimator (Definition 3) with the following gain function for h and y [ S(x).
where h Ã and y Ã are symmetric extensions of (upper triangular matrices) h and y, respectively (i.e. h Ã ij~h ij for ivj and h Ã ij~h ji for jvi; the definition of y Ã is similar.). It should be noted that, under the general estimation problem of Problem 3, the gain function of Eq. (22) cannot be introduced, and the gain function is specialized for the problem of RNA secondary structure prediction. The relation between the gain function of the c-centroid estimator (denoted by G (centroid) (h,y) and defined in Definition 7) and the one of the MEA-based estimator is where the additional term A(h,y) is positive for false predictions of base-pairs (i.e., FP and FN ) and C(h) does not depend on the prediction y (see [12] for the proof). This means the MEA- In computational experiments, the authors confirmed that the c-centroid estimator is always better than the MEA-based estimator when we used the same probability distribution of secondary structures. See [12] for details of the computational experiments.
Estimation of phylogenetic trees (Problem 4). The ccentroid estimator for Problem 4 can be introduced as follows: Estimator 3 (c-centroid estimator for Problem 4) For Problem 4, we obtain the c-centroid estimator (Definition 7) where the predictive space Y is equal to T (S) and the probability distribution on Y is taken by p (t) (hDS).
The following property is easily obtained by Theorem 2 and [17].
Property 7 (Relation of 1-centroid estimator and topological distance) The c-centroid estimator with c~1 (i.e. centroid estimator) for Problem 4 minimizes expected topological distances.
For X [ I (0) (I (0) is a set of partitions of S and is formally defined in the previous section), we call the marginalized probability p X~P h [ T (S) I(h X~1 )p (t) (hDS) partitioning probability. However, it is difficult to compute fp X g X [ I (0) as efficiently as in the prediction of secondary structures of RNA sequences, where it seems possible to compute the base-pairing probability matrix in polynomial time by using dynamic programming). Instead, a sampling algorithm can be used for estimating fp X g X [ I (0) approximately [16] for this problem. Once fp X g X [ I (0) is estimated, Theorem 3 leads to the following: Property 8 (Computaion of Estimator 3) The phylogenetic tree of Estimator 3 is found by maximizing the sum of the partitioning probabilities p X (of the partitions given by the predicted tree) that are larger than 1=(cz1).
In contrast to Estimator 1 (the c-centroid estimator for secondary structure prediction of RNA sequence) and Estimator 2 (the c-centroid estimator for pairwise alignment), it appears that there is no efficient method (such as dynamic programming algorithms) to computed Estimator 3 with cw1. Estimator 1 with c [ ½0,1, however, can be computed by using the following property, which is directly proven by Corollary 1 and the definition of the space T (S).
Alignment between two alignments of biological sequences. In this section we consider the problem of the alignment between two multiple alignments of biological sequences (Figure 4), which is often important in the multiple alignment of RNA sequences [19]. This problem is formulated as follows.
Problem 10 (Alignment between two alignments of biological sequences) The data is represented as D~fA,A'g where A and A' are alignments of biological sequences and the predictive space Y is equal to A(A,A'), that is, the space of the alignments of A and A 0 .
In the following, l(A) and n(A) denote the length of the alignment and the number of sequences in the alignment A, respectively. If both A and A 0 contain a single biological sequence (with no gap), Problem 10 is equivalent to conventional pairwise alignment of biological sequences (Problem 1). As in common secondary structure prediction, the representative estimator plays an important role in this application.
Estimator 4 (Representative estimator for Problem 10) For Problem 10, we obtain the representative estimator (Definition 10). The gain function G'(h k ,y) is the gain function of the c-centroid estimator. The parameter space H is represented as a product space HP The probability distribution on the parameter space H is given by p(hDD)P x') is given in the previous section (when x or x' contains some gaps, p (a) (hDx,x') is defined by the sequences with the gaps removed).
Corollary 2 proves the following properties of Estimator 5. Property 10 (A Relation of Estimator 4 with accuracy measures) Estimator 4 is consistent with the accuracy process for Problem 10 that is shown in Figure 5. We compare every pairwise alignment of x [ A and x' [ A' with the reference alignment. These comparisons are made using TP, TN, FP and FN with respect to the aligned-bases (e.g., using SEN, PPV and F-score).
Property 11 (Computation of Estimator 4) Estimator 4 can be given by maximizing the sum of probabilities p ik that are larger than Therefore, the pairwise alignment of Estimator 4 can be computed by the Needleman-Wunsch-type DP algorithm of Eq. (19) in which we replace p ij with Eq. (24). where p ik is defined in Eq. (24). The probability matrix fp ik g 1ƒiƒl(A),1ƒkƒl(A)' is often called an averaged aligned-base (matching) probability matrix of A and A 0 . In the iterative refinement of the ProbCons [19] algorithm, the existing multiple alignments are randomly partitioned into two groups and those two multiple alignments are re-aligned. This procedure is equivalent to Problem 10.
The estimator used in ProbCons is identical to Estimator 4 in the limit c??. Therefore, the estimator used in ProbCons is a special case of Estimator 4 and it only takes into account the SEN or SPS (sum-of-pairs score) of a predicted alignment.
Common secondary structure prediction from a multiple alignment of RNA sequences. Common secondary structure prediction from a given multiple alignment of RNA sequences plays important role in RNA research including non-coding RNA (ncRNA) [43] and viral RNAs [44], because it is useful for phylogenetic analysis of RNAs [45] and gene finding [43,[46][47][48]. In contrast to conventional secondary structure prediction of RNA sequences (Problem 2), the input of common secondary structure prediction is a multiple alignment of RNA sequences and the output is a secondary structure whose length is equal to the length of the input alignment (see Figure 6).
Problem 11 (Common secondary structure prediction) The data is represented as D~fAg where A is a multiple alignment of RNA sequences and the predictive space Y is identical to S(A) (the space of secondary structures whose length is equal to the alignment).
The representative estimator (Definition 10) directly gives an estimator for Problem 11.
Estimator 5 (The representative estimator for Problem 11) For Problem 11, we obtain the representative estimator (Definition 10) as follows. The gain function G'(h k ,y) is the gain function of the c-centroid estimator. The parameter space is equal to H~P x [ A S(x) where S(x) is the space of secondary structures. The probability distribution on H is given by is the probability distribution of the secondary structures of x [ A after observing the alignment A.
For example, p x (h x DA) can be given by extending the p (s) (hDx), although we have also proposed more appropriate probability distribution (see [49] for the details).
Corollary 2 proves the following properties of Estimator 5.  Property 13 (A relation of Estimator 5 with accuracy measures) Estimator 5 is consistent with an evaluation process for common secondary structure prediction: First, we map the predicted common secondary structure into secondary structures in the multiple alignment, and then the mapped structures are compared with the reference secondary structures based on TP, TN, FP and FN of the base-pairs using, for example, SEN, PPV and MCC (Figure 7).
Much research into common secondary structure prediction employs the evaluation process in Figure 7 (e.g., [50]).
Property 14 (Computation of Estimator 5) The common secondary structure of Estimator 5 is given by maximizing the sum of the averaged base-pairing probabilities p ij where Therefore, the common secondary structure of the estimator can be computed using the dynamic programming algorithm in Eq. (10) if we replace p ij with p ij .
Also, we can predict the secondary structure of Estimator 5 without conducting Nussinov-style DP: Property 15 (Computation of Estimator 5 with 0ƒªƒ1) The secondary structure of Estimator 5 with c [ ½0,1 can be predicted by collecting the base-pairs whose averaged base-paring probabilities are larger than 1=(cz1).
It should be noted that the tools of common secondary structure prediction, RNAALIFOLD [50], PETFOLD [8] and MCCASKILL-MEA [7] are also considered as a representative estimators (Definition 10). In [49], the authors systematically discuss those points. See [49] for details.
Pairwise alignment using homologous sequences. As in the previous application to RNA secondary structure prediction using homologous sequences, if we obtain a set of homologous sequences H for the target sequences x and x' (see Figure 8), we would have more accurate estimator for the pairwise alignment of x and x' than Estimator 1. The problem is formulated as follows.
Problem 12 (Pairwise alignment using homologous sequences) The data is represented as D~fx,x',Hg where x and x' are two biological sequences that we would like to align, and H is a set of homologous sequences for x and x'. The predictive space Y is given by Y~A(x,x') which is the space of the pairwise alignments of two sequences x and x'.
The difference between Problem 1 and this problem is that we can use other biological sequences (that seem to be homologous to x and x') besides the two sequences x and x' which are being aligned.
We can introduce the probability distribution (denoted by p (a) (hDx,x',h)) on the space of multiple alignments of three sequences x, x' and h (denoted by A(x,x',h) and whose definition is similar to that of A(x,x')) by a model such as the triplet HMM (which is similar to the pair HMM). Then, we obtain a probability distribution on the space of pairwise alignments of x and x' (i.e., A(x,x')) by marginalizing p (a) (hDx,x',h) into the space A(x,x'): where W is the projection from A(x,x',h) into A(x,x'). Moreover, by averaging these probability distributions over h [ H, we obtain the following probability distribution on A(x,x'): where DHD is the number of sequences in H.
The c-centroid estimator with the distribution in Eq. (27) directly gives an estimator for Problem 12. However, to compute the aligned-base-pairs (matching) probabilities p ik with respect to this distribution demands a lot of computational time, so we employ the approximated c-type estimator (Definition 12) of this ccentroid estimator as follows.
Estimator 6 (Approximated ª-type estimator for Problem 12) We obtain the approximated c-type estimator (Definition 12) for Problem 12 with the following settings. The parameter space is given by H~H'|H' \ where and the probability distribution on the parameter space H' is defined by Figure 7. An evaluation process for common secondary structure prediction (Problem 11). The comparison between each secondary structure and the reference secondary structure is done using TP, TN, FP and FN with respect to the base-pairs. doi:10.1371/journal.pone.0016450.g007 where DhD is the length of the sequence h. Property 16 (Computation of Estimator 6) The alignment of Estimator 6 is equal to the alignment that maximizes the sum of p ik larger than 1=(cz1) where Therefore, the recursive equation of the dynamic program to calculate the alignment of Estimator 6 is given by replacing p ik in Eq. (19) with Eq. (30).
Moreover, by using Theorem 1, we have the following proposition, which enables us to compute the proposed estimator for c [ ½0,1 without using (Needleman-Wunsch-type) dynamic programming.
Property 17 (Computation of Estimator 6 for 0ƒªƒ1) The pairwise alignment of Estimator 6 with c [ ½0,1 can be predicted by collecting the aligned-bases whose probability p ik in (30) is larger than 1=(cz1).
It should be noted that fp ik g 1ƒiƒDxD,1ƒkƒDx'D is identical to the probability consistency transformation (PCT) of x and x' [19]. In ProbCons [19], the pairwise alignment is predicted by the Estimator 6 with sufficiently large c. Therefore, the estimator for Problem 12 used in the ProbCons algorithm is a special case of Estimator 6.
RNA secondary structure prediction using homologous sequences. If we obtain a set of homologous RNA sequences for the target RNA sequence, we might have a more accurate estimator [23] for secondary structure prediction than the ccentroid estimator (Estimator 2). This problem is formulated as follows and was considered in [23] for the first time (See Figure 9). Problem 13 (RNA secondary structure prediction using homologous sequences) The data D is represented as D~fx,Hg where x is the target RNA sequence for which we would like to make secondary structure predictions and H is the set of its homologous sequences. The predictive space Y is identical to S(x), the space of the secondary structures of an RNA sequence x.
The difference between this problem and Problem 2 is that we are able to employ homologous sequence information for predicting the secondary structure of the target RNA sequence. In this problem, it is natural that we assume the target sequence x and each homologous sequence h [ H share common secondary structures. The common secondary structure is naturally modeled by a structural alignment (that considers not only the alignment between bases but also the alignment between base-pairs), and the probability distribution (denoted by p (sa) (hDx,x')) on the space of the structural alignments of two RNA sequences x and x' (denoted by SA(x,x')) is given by the Sankoff model [51]. By marginalizing the distribution p (sa) into the space of secondary structures S(x) of the target sequence x, we obtain more reliable distribution p(hDx) on S(x): where W is the projection from SA(x,h) into S(x). Moreover, by averaging these probability distributions on S(x), we obtain the following probability distribution of secondary structures of the target sequence.
where DHD is the number of sequences in H. The c-centroid estimator with the probability distribution in Eq. (32) gives a reasonable estimator for Problem 13, because Eq. (32) considers consensus secondary structures between x and h [ H. However, the calculation of the c-estimator requires huge computational cost because it requires O(nL 6 ) for computing the base-paring probability matrix fp ik g where p ik~P h [ S(x) I(h ij~1 )p(hDx) with the distribution of Eq. (32). Therefore, we employ the approximated c-type estimator (Definition 12) of the c-centroid estimator, which is equivalent to the estimator proposed in [23].
where C k is defined as (proof) It is sufficient to show that the consensus estimatorŷ y (c) is contained in the predictive space Y because G G(ŷ y)ƒ G G(ŷ y (c) ) for allŷ y in the MEG estimators, where G(y) :~E hDD ½G(h,y)~ð G(h,y)p(hDD)dh: If we assume thatŷ y (c) is not contained in the predictive space, Y that is, y y (c) = [Y , then there exists a k 0 such thatŷ y (c) = [C k0 . Becauseŷ y (c) is a binary vector, there exist indexes i,j [ I k0 such that i=j,ŷ y (c) i~1 andŷ y (c) j~1 . By the definition ofŷ y (c) , we obtain In order to prove the last inequality, we use Eq. (??). This leads to a contradiction and the theorem is proved.
Remark 3 It should be noted that the above theorem holds for an arbitrary parameter space including continuous-valued spaces.
Proof of Theorem 2. (proof) Because I(y i~1 )zI(y i~0 )~1 for arbitrary i, we obtain, using the definitions given in equations (15), (16), (17) and (18) where p i~p (h i~1 DD)~P h [ H I(h i~1 )p(hDD) is the marginalized probability. Therefore, we should always predict y i~0 whenever p i v1=(cz1), because the assumption of Theorem 3 ensures that the prediction y i~0 never violate the condition of the predictive space Y . Theorem 3 follows by using those facts.