BACKGROUND
in the field of gene clustering, gene expression data has been widely used assuming that genes that have similar expression patterns should have similar cellular functions and are likely to be involved in the same cellular processes  <cit> .
however, this assumption might be too simplistic considering the complexity of real biological systems.
it has become more and more acknowledged that different data sources offer information from different perspectives, and their combinations might make the prediction more accurate.
there are many types of biological data available besides gene expression data, such as protein-dna binding data, protein-protein interaction data, evolutionary conservation data, gene ontology information, et cetera.
however, different data types have different characteristics, and thus how to integrate multiple heterogeneous data types into a single framework and make the results more accurate has become one of the most challenging problems.
in this study, we developed a clustering algorithm that can cluster genes based on beta distributed and gaussian distributed data, which are represented by protein-dna binding probabilities  and gene expression data, respectively, in a real case study.
other possible data sources that can be naturally modeled with beta distributions include e.g.
correlations  <cit>  and pair-wise and multiple sequence similarities  <cit> , and other possible gaussian distributed data sources include various other microarray-based measurements.
many unsupervised methods have been developed and widely used in gene clustering.
they can be roughly classified into three categories, which are heuristic, iterative relocation and model-based methods  <cit> .
the first two approaches suffer from solving some basic practical issues such as 'how to define the number of clusters' and 'how to handle outliers', which can be easily handled by model-based methods.
for the first issue, the problem can be recasted as the model selection problem; and for the second question, the outliers can be handled by adding one or more components which represent a different distribution for them  <cit> .
moreover, model-based clustering methods outweigh approaches within the other two categories in their statistical nature  <cit> .
so in this study, we choose model-based clustering as the framework for the unsupervised data fusion.
expectation maximization  algorithm is often used to solve the problem of maximum likelihood estimation with incomplete data, and thus is commonly adopted in model-based clustering.
although em algorithm for gaussian distribution is well-known, less information is available about that for other distributions, not mentioning combinations of different distributions.
in this study, beta distributed data and gaussian distributed data are integrated into one combined mixture model.
we have developed three types of em algorithms, the standard em , an approximated em  and a hybrid em  algorithm for bgmm, whose comparisons were done using simulated data.
emh was used for bgmm in the simulations and real case studies.
performance tests with bgmm and its component models  were done both with simulated and real data, and the results show that our joint mixture model can yield more accurate results.
these results also demonstrate the idea that the more data that are integrated the more comprehensive the result will be.
two commonly used model selection criteria are likelihood-based methods and approximation-based methods, of which approximation-based methods are widely preferred due to their simplicity and less computational cost  <cit> .
these methods include penalized likelihood, closed-form approximations to the bayesian solution, and monte carlo sampling of the bayesian solution, among which the first two methods are most prevalent.
four well-known model selection criteria, bayesian information criterion , integrated classification likelihood-bic , akaike information criterion , and modified aic  were tested for bgmm and its two extreme models in this study.
icl is reported to work well for bmm  <cit> , and aic as well as bic are commonly used as the criterion in gmm  <cit> .
our simulation results in this study suggests that aic or icl is preferred in bgmm depending on which em algorithm is employed.
the following sections are organized as 'methods', 'results and discussion', and 'conclusions'.
in section 'methods', we introduced bgmm together with all its three types of em algorithms, and described the formulation of the four tested model selection criteria.
in section 'results and discussion', we first compared the three types of em algorithms in bgmm , and then compared the performance of bgmm with bmm and gmm.
in section 'conclusions', we first summarized the main work of this study, discussed the possible extension and limitations of the current work, and in the end briefly mentioned the related future work.
methods
in this section, bgmm and all its three types of em algorithms are first introduced, and then the approximation based model selection criteria which are compared in this study are described in detail.
mixture model based clustering
in model-based clustering methods, each observation xj, where j =  <dig> ..., n and n is the number of genes, is drawn from a finite mixture distribution with the prior probability πi, component-specific distribution  and its parameters θi.
the formula is given as  <cit>
where θ = { from  for simplicity.
bgmm
in bgmm, we define θ = t, π = t,  and ], where p <dig> and p <dig> each represents the dimension of the observations in bmm and gmm, respectively.
we also denotes y and z as the observations of beta distributed and gaussian distributed data, respectively, function f of y and f of z as the density function of beta and gaussian distribution, respectively, and x = t. y and z can be used to denote different data sources in different contexts, which for example denote tf binding probabilities and gene expression data in our bioinformatics application.
bgmm is built from bmm and gmm with the assumption that, for each component i, the beta distributed and gaussian distributed data are independent.
in the bmm part, each component is assumed to be the product of p <dig> independent beta distributions, whose probability density function is defined as
where  and .
likewise, each component is assumed to follow a gaussian distribution in the gmm part, whose probability density function of each component for each gene is defined as
where , ,  and .
we assume the commonly used gaussian model where the covariance matrix is a diagonal matrix.
this approximation is useful especially for high-dimensional data since it significantly reduces the number of parameters that need to be estimated from data.
it is worth noting that the above mixture model construction implicitly assumes that the two data sources share the same clustering structure, which is a reasonable assumption for the general problem of clustering gene expression and tf binding data .
however, this assumption does not necessarily hold in all other clustering problems, in which case our method is not applicable .
em algorithm is applied to estimate the parameters θ iteratively.
we have developed three types of em algorithms for bgmm, the standard em , an approximated em  and a hybrid em , which are described in detail in the following sections.
em algorithms
the standard em algorithm
in the standard em algorithm, the data log-likelihood  can be written as
given x = {xj: j =  <dig> ..., n}, whose direct maximization, however, is difficult.
in order to make the maximization of equation  <dig> tractable, the problem is casted in the framework of incomplete data.
since we assume that the beta distributed and gaussian distributed data are independent, lc can be factored as
if we define cj ∈ { <dig> ..., g} as the clustering membership of xj, then the complete data log-likelihood can be written as
where χ is the indicator function of whether xj is from the ith component or not.
in the em algorithm, e step computes the expectation of the complete data log-likelihood
where θ represents the parameters estimated in the mth iteration .
by computing the expectation, equation  <dig> becomes
where
according to bayes' rule.
note that  is the estimated posterior probability of xj coming from component i at iteration m, and we can assign each xj to its component based on .
equations  <dig> and  <dig> show that our assumption of the beta distributed and gaussian distributed data being independent carries over to the expected log-likelihood as well.
in the em algorithm of bgmm, αiu's and βiu's, which are the parameters of the bmm part, are estimated using newton-raphson method.
let θ1i = , then the parameters are updated by
with the constraint θ1i ≥  <dig>  where  is the hessian matrix evaluated at  and  is the lagrangian function of  .
the parameters of the gmm part, μiv's and , in bgmm can be estimated by the standard em algorithm of gmm with diagonal covariance matrix, which works by iterating over
and π's are updated by
where  is calculated from equation  <dig> .
note that {u =  <dig> ..., p1} and {v =  <dig> ..., p2}.
from the above equations, it is easy to see that the standard em for bgmm will reduce to the standard em for bmm when p <dig> goes to  <dig> and shrink to the standard em for gmm when p <dig> =  <dig>
approximated and hybrid em algorithms
we also developed an approximated em algorithm for bgmm, whose main difference compared with the standard one is that it maximizes equation  <dig> instead of equation  <dig>
in e step, τji's are first calculated with the current parameters, according to which xj's are clustered to their corresponding clusters using cj = i <dig> where i <dig> = arg maxi τji.
then in m step, the new parameters are estimated so as to maximize equation  <dig>  given the hard clusters obtained in e step.
given that the beta and gaussian distributed data are assumed to be independent, ml parameter estimates for beta and gaussian parts can be computed separately, which corresponds to the basic ml estimation using standard techniques.
in the approximated em, the new  and  are estimated with a numerical optimization method, 'betafit', which is implemented in matlab, and the new  and  are calculated by
respectively, where  is composed of all the genes in cluster i estimated from e step,  refers to the  of cluster i, and .
update of πi's and calculation of τji's remain the same with the standard em algorithm.
in the end, we developed one type of hybrid em algorithm, whose αiu's and βiu's are maximized by the approximated em, μiv's,  and πi's are updated by the standard em.
the approximate em for clustering is analogous to the viterbi training for hidden markov models .
viterbi training has been proposed as an alternative to the standard em in the cases where the standard em becomes computationally too expensive.
although there are no convergence guarantees in general, the viterbi training has been found useful due to its efficiency and, in particular, when one seeks to decode the state  via viterbi algorithm.
the same considerations apply for the clustering problem as well, where the approximate em optimizes the hard clustering and parameters iteratively.
moreover, because parameter estimates remain fixed for a given hard clustering, the optimization is a discrete process and, therefore, convergence is achieved exactly.
the hybrid method shares  the benefits from both the standard em and the approximate em.
in order to avoid the possible local maxima, we run the algorithm  multiple times with different initial values.
the parameters αiu's and βiu's for each dimension of the beta distribution u  are initialized by method-of-moments so that their means are randomly distributed within the range of y1u,..., ynu and variances are equal for all clusters , μiv's and  are obtained from the randomly initialized fuzzy c-means clustering results, and πi's are initialized with the uniform probability 1/g.
in this study, for each data set, we run each em algorithm  <dig> times with different initial values, and for all the tested models, we set the convergence threshold  and maximum number of iterations to  <dig>  and  <dig> respectively.
all the simulations have reached their convergence according to the statistics stored during the simulations.
model selection
four well-known approximation-based model selection criteria, bic  <cit> , icl  <cit> , aic  <cit> , and aic <dig>  <cit>  are compared in bgmm and its extreme models, according to which the best-performing criterion for each model is chosen.
calculations for the above criteria are defined as
where d is the number of free parameters, and m  is the total amount of the data .
note that  is the estimated entropy of the fuzzy classification matrix cji =   <cit> .
the number of free parameters d are distinct in different models.
in bmm, we have p1g free αiu's, p1g free βiu's, and g -  <dig> free πi's (), so db = 2p1g + g -  <dig>  in gmm, as we have p <dig> free σv's, p2g free μiv's, and also g -  <dig> free πi's, thus dg = p <dig> + p2g + g -  <dig>  in the joint model, the number of free parameters is the summation of those in its extreme models minus one set of free πi's, therefore we have dbg = 2p1g + p <dig> + p2g + g -  <dig>
RESULTS
in this section, we first compared the performance of bgmm with different em algorithms by artificial data, according to which one em was chosen for later simulations.
then we tested the integration idea  by comparing bgmm with its two extreme cases.
performance test of bgmm with artificial data
to evaluate the overall performance of a clustering method, we developed one scoring system to evaluate the clustering accuracy when dealing with artificial data.
it searches the best matching between the cluster labels of the results  and the ground truth clustering among all of their possible associating ways.
the score for the best match is denoted as 'e score', and is defined as
in this scoring system, tj denotes the ground truth clustering membership of data j; r stands for all possible associating ways between the estimated and the true clusters, where ri is the label of data belonging to component i predicted by the clustering algorithm, and r is chosen from labels  <dig>   <dig> ...,max{, g} .
also note that e represents the individual score of each gene, e is the average score of all the genes for each repetition, 'e score' of each repetition is the one corresponding to the optimal q, and the final 'e score' of each data set is the median of the  <dig> 'e score's.
it is worth noticing that the estimated and assumed number of clusters,  and g, vary with the model selection criteria, and thus cause different 'e scores', rendering this scoring system not only records the accuracy of the results but also reflects the influence of model selection criterion.
performance test of different em algorithms in bgmm
we first compared the performance of ems, ema and emh in bgmm.
for simplicity, we denote bgmm that employs ems, ema or emh as bgmms, bgmma or bgmmh, correspondingly.
the artificial data set for the performance test was designed according to our model, whose parameters are listed in table  <dig>  the data set was divided into high quality  and low quality  data, namely 'gb' , 'bb' , 'gg'  and 'bg'  respectively.
we also designed two kinds of 'bg's, 'bgm' and 'bgv', which were hard to be clustered compared to 'gg' with respect to close means and large variances, respectively.
the data set was designed to have three underlying clusters,  <dig> genes  and four features .
the simulation was repeated  <dig> times with randomly generated data sets, and the comparison results of the clustering accuracy were depicted in figure  <dig>
'gb', 'bb' and 'gg' stand for good, bad beta distributed data and good gaussian distributed data re-spectively; 'bgm' and 'bgv' represent bad gaussian distributed data which are hard to be clustered with respect to close means and large variances respectively.
in order to choose the best model selection criterion  for each type of bgmm, we summed up the number of hits of the correct number of clusters for each tested case.
the summation results for bic, icl, aic and aic <dig> are  <dig>   <dig>   <dig> and  <dig>  respectively, in bgmms,  <dig>   <dig>   <dig> and  <dig>  respectively, in bgmma, and  <dig>   <dig>   <dig> and  <dig>  respectively, in bgmmh.
therefore, icl is upheld by bgmms, and aic is embraced by both bgmma and bgmmh in this simulation.
we evaluated the clustering accuracy of different types of bgmm with each best model selection criterion.
simulation results show that, although different algorithms perform slightly different for different cases , the overall prediction accuracy of the three methods are similar as shown in figure  <dig>  we also compared the running time of the three methods under the same background framework, where no significant difference among them was detected.
one important application of the proposed algorithm is to cluster genes based on protein-dna binding probabilities and gene expression data, which are assumed to be of beta and gaussian distributions in bgmm.
this parametric assumption is supported by our good clustering results and additional distributional assessments.
in some cases, however, our parametric assumptions might be violated due to various reasons, especially for expression data.
for example, different platforms used to measure transcriptome might affect the distribution of expression data.
although this problem can be solved by extending the current algorithm to other parametric distributions quite easily, it is important to know how sensitive bgmm is to the violation of the parametric assumptions and how robust the algorithm is in dealing with noisy variables.
to address this, we run three additional simulations with the three em algorithms, where gene expression and protein-dna binding data are simulated from laplace and kumaraswamy distributions, respectively.
laplace and kumaraswamy distributions are used to replace gaussian and beta distributions separately in simulation  <dig> and  <dig>  and both distributions are replaced in simulation  <dig>  note that laplace and kumaraswamy distributions have the same support as gaussian and beta distributions, respectively.
means and variances used in laplace distribution are the same with those of gaussian distribution , and α 's and β 's used in kumaraswamy distribution are also the ones used in beta distribution .
as shown in figure  <dig>  all three em algorithms work similarly, and are not excessively sensitive to the parametric assumptions used in this study.
parameters of the beta distributed data and the symbols are the same with what has been shown in table  <dig>
based on the above test, the three em algorithms perform equally well.
therefore, we simply used bgmmh for the performance tests and referred to it as 'bgmm' for simplicity in the following text.
performance test of bgmm with its component models
simulations shown in this section were dedicated to test how well bgmm could integrate different data sources.
we compared the performance of bgmm  with its two extreme models, bmm with ema  and gmm with ems , for this purpose.
a slightly different data set was used, where the gaussian distributed data was designed to be less distinguishable than what has been shown in the previous section.
the parameters of the redesigned data are shown in table  <dig>  where all the rest information including the dimensions of the data  and the repetitions  remain the same.
we used the same method as what we did in the previous section to select the best criterion for bgmm, bmm and gmm.
ordered by bic, icl, aic and aic <dig>  the summations of the hits are  <dig>   <dig>   <dig>   <dig>  respectively, in bgmm,  <dig>   <dig>   <dig>   <dig>  respectively, in bmm, and  <dig>   <dig>   <dig>   <dig>  respectively, in gmm, according to which aic was chosen as the best criterion for all the three tested models in this simulation.
the comparison results of bgmm with its extreme models are shown in figure  <dig>  for expression data whose variances are not too large, the joint model can improve the clustering accuracy regardless of the quality of the data compared with either of its extreme models .
however, when gaussian distributed data has too much overlap among the clusters, bgmm does not necessarily show its superiority  when the variances are too large as shown in the case of bgv.
it is indicated that bgmm is sensitive to the variances of gaussian distributed data since bgv is designed to have similar noise level as that of bgm.
these results demonstrate that the em algorithm of bgmm has the power of reinforcing each extreme model with information from the other one, but does not necessarily outweigh both of them if the gaussian distributed data contains too much noise with respect to large variances .
performance test of bgmm with real data
we applied our methods to mouse protein-dna binding data and gene expression data.
the binding data is modeled as beta distribution, which are the binding probabilities output from a method called 'probtf'  <cit> .
probtf uses genome sequences and transcription factor sequence specificities to compute the protein-dna binding probabilities.
this method answers the question of whether the whole gene promoter has one or more binding sites for a tf.
since it processes each promoter as a whole, the computational predictions provide insights into the functional role of a tf in the regulatory program of a target gene.
the rationale for this is the fact that the higher binding probability anywhere on the promoter  implies higher probability of a regulatory relationship.
further, the method is able to make use of practically any genome-level information, such as evolutionary conservation, nucleosome positioning, chip-chip, and other prior knowledge .
the protein-dna binding data contains the probabilities of  <dig> tfs binding to  <dig> genes, calculated with mouse-specific position weight matrices from the transfac database .
the gene expression data is modeled as gaussian distribution, which is composed of  <dig> genes measured from  <dig> conditions  <cit> .
there are  <dig> genes measured in both data sets.
we removed the genes whose gene expression profiles have low absolute values  with matlab function 'genelowvalfilter', and then choose genes that have annotations available for sure with the functional classification tool of david database .
in the end, we obtained  <dig> genes for the following studies.
to see how well our data satisfy the parametric assumptions, we did the following test.
for protein-dna binding data, we grouped all the binding probabilities  into two beta-distributed clusters  and drew their pdfs, each representing the binding and unbinding cases, respectively.
figure  <dig> shows that the genome-wide binding data can remarkably well be approximated with two beta-distributed components.
similarly as shown in figure  <dig>  expression data can be fitted into a gaussian distribution.
this agrees with previous studies where gene expression data from a microarray platform is commonly assumed to be normally distributed.
although the above preliminary test does not correspond to our clustering method exactly, it demonstrates that our parametric assumptions are indeed reasonably good.
the bgmm clustering method effectively increases the number of clusters to which the data is split and further improves the fit to the data.
the binding data corresponding to two sets of tfs were chosen to cluster the genes together with its corresponding expression data by bgmm, bmm and gmm.
the clustering results were then compared and evaluated by gene ontology .
the first set of tfs was randomly chosen with respect to their biological significance , while the second set was carefully selected by our model .
there are three subsets of 'setrand', each of which was chosen based on certain criterion.
we arbitrarily choose three thresholds to be compared with the median of the binding probabilities of a certain tf, and tfs that exceed this threshold are used for clustering .
the thresholds for 'setrand1' to 'setrand3' are  <dig> ,  <dig>  and  <dig> , respectively, and the number of tfs selected are  <dig>   <dig> and  <dig>  correspondingly.
'setreal' was selected by bmm.
we first clustered the genes based on two sets of tfs by bmm, which were bach <dig> and bach <dig> combined with mafk, respectively.
this is because that the two bach proteins are both reported to interact with mafk protein.
then we compared the genes whose cluster has the lowest enrichment score from each clustering result, and the common set which contains  <dig> genes was chosen.
we further clustered all the  <dig> tfs based on the  <dig> genes by bmm, and focused on the cluster that contains bach <dig>  bach <dig> and mafk.
this cluster turns out to be composed of all the tfs that belong to the families fos, jun, maf and nf-e <dig> among our tested tfs, which are all ap- <dig> components of the leucine zipper factors class.
there are  <dig> tfs  in this cluster, all of which were chosen to form 'setreal'.
maf family proteins  heterodizes with cnc-related bzip factors which include nf-e <dig> family proteins   <cit> ; while fos family  form hetero  or homo  dimers with jun family  proteins  <cit> .
these dimers bind to dna at certain motif that contains ap- <dig> binding sites  <cit> .
the result that our bmm can cluster tfs which have similar binding profiles into one single cluster demonstrates the applicability of bmm, one extreme case of bgmm.
go was employed in this study to validate the clustering results.
in order to find the most significant annotated terms by looking at the probabilities that the terms are counted by chance, we used the hypergeometric probability distribution to calculate the p-values of gene enrichment score  for each cluster by each model with each model selection criterion .
we compared the means and medians of those p-values across all the groups clustered by each model, whose results are shown in table  <dig>  it is worth mentioning that the clustering result is obtained by running the algorithm  <dig> times and taking the one whose expected complete data log-likelihood is the maximum, and each p-value shown in table  <dig> is the mean or median of the p-values of all the ontology groups  corresponding to the best clustering result .
from this table, it is clear that, no matter whether the tfs were randomly selected or not, both means and medians of the p-values of bgmm are lower than those of either bmm or gmm, regardless of which aspect  was considered and which model selection criterion was used.
these results indicate that our bgmm can cluster the genes in a more reasonable way with respect to their biological functions, localizations and processes involved.
it is also seen from table  <dig> that, there are two cases where the four model selection criteria have different prediction results, one is in bmm of the case setrand <dig> where the results chosen by aic yields the smallest p-values, and the other is in gmm of the case setreal where aic selects the best model in terms of the smallest p-value, both of which accord well with our simulation results.
moreover, the choice of tfs whose binding probabilities are used in clustering does obviously affect the results and, therefore, tfs should be carefully chosen based on biological knowledge of a specific problem.
in this study, although binding data of randomly  chosen tfs  also give lower p-values, the obtained clusters might not provide best insight into our biological problem.
we therefore carefully studied the results obtained from setreal, which are discussed below.
statistics shown in this table are the group average of the p-values of gene enrichment score .
tfs of 'setrand1' to 'setrand3' were randomly chosen according to certain thresholds , while tfs of 'setreal' were carefully selected by bmm.
note: 'f', 'c' and 'p' stand for the three aspects of gene ontology; 'all' means all aspects are included.
'4' represents that all the four criteria indicate the same clustering result.
statistics shown here are all rounded to  <dig> decimals.
there are eight clusters obtained from setreal by bgmm, among which three groups have p-values below  <dig>  if all the aspects were taken into account .
the three clusters were named 'clu1' to 'clu3' and ordered from the highest average expression level to the lowest.
the expression patterns  and the medians of the genes  within one cluster are shown in figure  <dig>  six toll-like receptor  agonists which are cpg, pam2csk <dig>  pam3csk <dig>  lps, poly i:c and r <dig> were used as the treatments, and four gene knock-out mutants and different time points were included to increase the diversity of the tlr-stimulated gene expression data set and the number of measurements  <cit> .
the first four tlr agonists are bacterial-associated, while poly i:c is viral-associated and r <dig> is anti-viral stimuli.
they were used here to stimulate tlr-stimulated macrophages, which represent various pathogen-associated molecular patterns.
among the genes that have been deleted, adapters myd <dig> and ticam <dig>  could provide a structural platform for the recruitment of kinases and downstream effector molecules, were reported crucial for signaling by most type i il- <dig> receptor/tlr family members  <cit> .
however, bjökbacka et al.
reported that the majority of the host response to lps is regulated independently of myd <dig>  and genes appearing to be ticam1-dependent can be classified as both myd88-independent and myd88-dependent  <cit> .
CONCLUSIONS
this paper presents a novel beta-gaussian mixture model, bgmm, for gene clustering from beta distributed and gaussian distributed data.
we developed three types of em algorithms for bgmm in this study, whose overall performance are similar according to our simulations.
we simply chose emh as the core of bgmm for further performance test, which was done by comparing bgmm with its two component models, bmm and gmm, with both artificial and real data.
results from artificial data indicate that our joint model works best if the variances of the gaussian distributed data were not too large, and go validation of the real case studies show that the joint model yields more comprehensive results no matter what model selection criterion is used and whether the data is carefully chosen or not.
for the carefully selected real data, we started from limited known tfs  and ended up with all the tfs  within the tested scope that have the same common features, which demonstrates the usability of one extreme case of bgmm .
after clustering the genes with the  <dig> tfs, we obtained three distinguished gene groups which might be involved in the myd88-dependent tlr-3/ <dig> signaling cascades.
these results not only tested the performance of the joint model, but also demonstrated its usability in real cases and in some possible applications.
the main contribution of this paper is that it has proposed a framework for multiple data integration through mixture modeling that has not been addressed by anyone else before.
the proposed bgmm is designed to integrate beta distributed and gaussian distributed data.
however, the way how those data are incorporated is not limited to the data types that we have used in this study.
in principle, data of other parametric distribution can be easily integrated by combining its particular em algorithm into this framework .
therefore, the framework proposed in this paper is applicable to many other problems and not limited to the particular problem considered here.
one of the basic assumptions in this paper is that the ground truth clustering for gaussian and beta distributed data are the same.
this is because transcriptional regulation is largely controlled by the tfs that bind to the gene promoters, thus the expression profiles of genes whose regulatory regions are bound by the same/similar factors are expected to be similar.
although the above statement is generally true, it might be violated due to post-transcriptional modifications etc., in which case the method may not be directly applicable.
however, if post-transcriptional or other phenomena become a real problem, it can be compensated by integrating more information sources, such as protein-protein interactions, into the proposed clustering framework.
on the other hand, if the two data sources do not share the same clustering structure, then an alternative modeling strategy would be needed, such as a hierarchical bayes model that would model a true clustering structure but allow individual structures for both data types.
another issue that is worth mentioning is how different data pre-processing and microarray platforms can affect the distribution of gene expression data and, thereby, clustering results.
fortunately, as discussed in section, the alternative distributions we tried on bgmm had a very small effect on the clustering results, suggesting that bgmm is considerably robust to small fluctuations in the distributional assumptions.
more importantly, one major advantage of bgmm is its flexibility of easily being extended to other parametric distributions.
that is, if in a particular problem data come from different distributions, then one can relatively straightforwardly develop a similar model-based approach as proposed here to model the problem at hand in a precise way by fitting data to those specific distributions.
we employed the diagonal covariance matrix model in the em algorithm of gmm to reduce the number of parameters to be estimated so that it can be easily applied to large dimensional real data.
in particular, the number of parameters in diagonal matrix is p <dig> which is remarkably smaller than that of the full covariance matrix .
diagonal covariance matrix automatically assumes no correlations among gaussian distributed data.
so if we want to preserve the correlation information among time series data by the proposed framework without introducing too many new parameters, it is possible, e.g., to develop similar estimation algorithms for a covariance model where off-diagonal constant correlations are assumed or use a more general covariance matrix  <cit> .
in the future, we could improve the proposed model so that it can account for the correlations among gene expressions.
we could also integrate more data sources into this framework and apply it to more real problems.
in this aspect, we could either combine other data sources into the framework as a component model, or convert them into prior information which can be used to stratify the model  <cit> .
authors' contributions
xfd and hl designed the study and developed the methods.
xfd implemented the algorithms, did the performance tests, and wrote the manuscript.
te derived the standard em algorithm for bmm.
te and hl derived the standard em algorithm for bgmm.
xfd, hl, te and oy-h prepared the manuscript.
all authors have read and approved the final manuscript.
appendix
derivation of α's, β's, and π's
define
recall
where , , , ,  and .
derivation of α's and β's
define the parameter vector
thus the new estimate  is obtained as follows
where  is the hessian matrix evaluated at .
note that Ψ and Ψ' represents the digamma and trigamma functions respectively, which are the first and second logarithmic derivatives of the gamma function.
derivation of π's

 in the field of gene clustering, gene expression data has been widely used assuming that genes that have similar expression patterns should have similar cellular functions and are likely to be involved in the same cellular processes. however, this assumption might be too simplistic considering the complexity of real biological systems. it has become more and more acknowledged that different data sources offer information from different perspectives, and their combinations might make the prediction more accurate. there are many types of biological data available besides gene expression data, such as protein-dna binding data, protein-protein interaction data, evolutionary conservation data, gene ontology information, et cetera. however, different data types have different characteristics, and thus how to integrate multiple heterogeneous data types into a single framework and make the results more accurate has become one of the most challenging problems. in this study, we developed a clustering algorithm that can cluster genes based on beta distributed and gaussian distributed data, which are represented by protein-dna binding probabilities and gene expression data, respectively, in a real case study. other possible data sources that can be naturally modeled with beta distributions include e.g. correlations and pair-wise and multiple sequence similarities, and other possible gaussian distributed data sources include various other microarray-based measurements. many unsupervised methods have been developed and widely used in gene clustering. they can be roughly classified into three categories, which are heuristic, iterative relocation and model-based methods. the first two approaches suffer from solving some basic practical issues such as'how to define the number of clusters' and'how to handle outliers', which can be easily handled by model-based methods. for the first issue, the problem can be recasted as the model selection problem; and for the second question, the outliers can be handled by adding one or more components which represent a different distribution for them. moreover, model-based clustering methods outweigh approaches within the other two categories in their statistical nature. so in this study, we choose model-based clustering as the framework for the unsupervised data fusion. expectation maximization algorithm is often used to solve the problem of maximum likelihood estimation with incomplete data, and thus is commonly adopted in model-based clustering. although em algorithm for gaussian distribution is well-known, less information is available about that for other distributions, not mentioning combinations of different distributions. in this study, beta distributed data and gaussian distributed data are integrated into one combined mixture model. we have developed three types of em algorithms, the standard em, an approximated em and a hybrid em algorithm for bgmm, whose comparisons were done using simulated data. emh was used for bgmm in the simulations and real case studies. performance tests with bgmm and its component models were done both with simulated and real data, and the results show that our joint mixture model can yield more accurate results. these results also demonstrate the idea that the more data that are integrated the more comprehensive the result will be. two commonly used model selection criteria are likelihood-based methods and approximation-based methods, of which approximation-based methods are widely preferred due to their simplicity and less computational cost. these methods include penalized likelihood, closed-form approximations to the bayesian solution, and monte carlo sampling of the bayesian solution, among which the first two methods are most prevalent. four well-known model selection criteria, bayesian information criterion, integrated classification likelihood-bic, akaike information criterion, and modified aic were tested for bgmm and its two extreme models in this study. icl is reported to work well for bmm, and aic as well as bic are commonly used as the criterion in gmm. our simulation results in this study suggests that aic or icl is preferred in bgmm depending on which em algorithm is employed. the following sections are organized as'methods','results and discussion', and'conclusions'. in section'methods', we introduced bgmm together with all its three types of em algorithms, and described the formulation of the four tested model selection criteria. in section'results and discussion', we first compared the three types of em algorithms in bgmm, and then compared the performance of bgmm with bmm and gmm. in section'conclusions', we first summarized the main work of this study, discussed the possible extension and limitations of the current work, and in the end briefly mentioned the related future work. methods in this section, bgmm and all its three types of em algorithms are first introduced, and then the approximation based model selection criteria which are compared in this study are described in detail. mixture model based clustering in model-based clustering methods, each observation xj, where j =..., n and n is the number of genes, is drawn from a finite mixture distribution with the prior probability πi, component-specific distribution and its parameters θi. the formula is given as where θ = { from for simplicity. bgmm in bgmm, we define θ = t, π = t, and ], where p and p each represents the dimension of the observations in bmm and gmm, respectively. we also denotes y and z as the observations of beta distributed and gaussian distributed data, respectively, function f of y and f of z as the density function of beta and gaussian distribution, respectively, and x = t. y and z can be used to denote different data sources in different contexts, which for example denote tf binding probabilities and gene expression data in our bioinformatics application. bgmm is built from bmm and gmm with the assumption that, for each component i, the beta distributed and gaussian distributed data are independent. in the bmm part, each component is assumed to be the product of p independent beta distributions, whose probability density function is defined as where and. likewise, each component is assumed to follow a gaussian distribution in the gmm part, whose probability density function of each component for each gene is defined as where,, and. we assume the commonly used gaussian model where the covariance matrix is a diagonal matrix. this approximation is useful especially for high-dimensional data since it significantly reduces the number of parameters that need to be estimated from data. it is worth noting that the above mixture model construction implicitly assumes that the two data sources share the same clustering structure, which is a reasonable assumption for the general problem of clustering gene expression and tf binding data. however, this assumption does not necessarily hold in all other clustering problems, in which case our method is not applicable. em algorithm is applied to estimate the parameters θ iteratively. we have developed three types of em algorithms for bgmm, the standard em, an approximated em and a hybrid em, which are described in detail in the following sections. em algorithms the standard em algorithm in the standard em algorithm, the data log-likelihood can be written as given x = {xj: j =..., n}, whose direct maximization, however, is difficult. in order to make the maximization of equation tractable, the problem is casted in the framework of incomplete data. since we assume that the beta distributed and gaussian distributed data are independent, lc can be factored as if we define cj ∈ {..., g} as the clustering membership of xj, then the complete data log-likelihood can be written as where χ is the indicator function of whether xj is from the ith component or not. in the em algorithm, e step computes the expectation of the complete data log-likelihood where θ represents the parameters estimated in the mth iteration. by computing the expectation, equation becomes where according to bayes' rule. note that is the estimated posterior probability of xj coming from component i at iteration m, and we can assign each xj to its component based on. equations and show that our assumption of the beta distributed and gaussian distributed data being independent carries over to the expected log-likelihood as well. in the em algorithm of bgmm, αiu's and βiu's, which are the parameters of the bmm part, are estimated using newton-raphson method. let θ1i =, then the parameters are updated by with the constraint θ1i ≥ where is the hessian matrix evaluated at and is the lagrangian function of. the parameters of the gmm part, μiv's and, in bgmm can be estimated by the standard em algorithm of gmm with diagonal covariance matrix, which works by iterating over and π's are updated by where is calculated from equation. note that {u =..., p1} and {v =..., p2}. from the above equations, it is easy to see that the standard em for bgmm will reduce to the standard em for bmm when p goes to and shrink to the standard em for gmm when p = approximated and hybrid em algorithms we also developed an approximated em algorithm for bgmm, whose main difference compared with the standard one is that it maximizes equation instead of equation in e step, τji's are first calculated with the current parameters, according to which xj's are clustered to their corresponding clusters using cj = i where i = arg maxi τji. then in m step, the new parameters are estimated so as to maximize equation given the hard clusters obtained in e step. given that the beta and gaussian distributed data are assumed to be independent, ml parameter estimates for beta and gaussian parts can be computed separately, which corresponds to the basic ml estimation using standard techniques. in the approximated em, the new and are estimated with a numerical optimization method,'betafit', which is implemented in matlab, and the new and are calculated by respectively, where is composed of all the genes in cluster i estimated from e step, refers to the of cluster i, and. update of πi's and calculation of τji's remain the same with the standard em algorithm. in the end, we developed one type of hybrid em algorithm, whose αiu's and βiu's are maximized by the approximated em, μiv's, and πi's are updated by the standard em. the approximate em for clustering is analogous to the viterbi training for hidden markov models. viterbi training has been proposed as an alternative to the standard em in the cases where the standard em becomes computationally too expensive. although there are no convergence guarantees in general, the viterbi training has been found useful due to its efficiency and, in particular, when one seeks to decode the state via viterbi algorithm. the same considerations apply for the clustering problem as well, where the approximate em optimizes the hard clustering and parameters iteratively. moreover, because parameter estimates remain fixed for a given hard clustering, the optimization is a discrete process and, therefore, convergence is achieved exactly. the hybrid method shares the benefits from both the standard em and the approximate em. in order to avoid the possible local maxima, we run the algorithm multiple times with different initial values. the parameters αiu's and βiu's for each dimension of the beta distribution u are initialized by method-of-moments so that their means are randomly distributed within the range of y1u,..., ynu and variances are equal for all clusters, μiv's and are obtained from the randomly initialized fuzzy c-means clustering results, and πi's are initialized with the uniform probability 1/g. in this study, for each data set, we run each em algorithm times with different initial values, and for all the tested models, we set the convergence threshold and maximum number of iterations to and respectively. all the simulations have reached their convergence according to the statistics stored during the simulations. model selection four well-known approximation-based model selection criteria, bic, icl, aic, and aic are compared in bgmm and its extreme models, according to which the best-performing criterion for each model is chosen. calculations for the above criteria are defined as where d is the number of free parameters, and m is the total amount of the data. note that is the estimated entropy of the fuzzy classification matrix cji =. the number of free parameters d are distinct in different models. in bmm, we have p1g free αiu's, p1g free βiu's, and g - free πi's (), so db = 2p1g + g - in gmm, as we have p free σv's, p2g free μiv's, and also g - free πi's, thus dg = p + p2g + g - in the joint model, the number of free parameters is the summation of those in its extreme models minus one set of free πi's, therefore we have dbg = 2p1g + p + p2g + g -  in this section, we first compared the performance of bgmm with different em algorithms by artificial data, according to which one em was chosen for later simulations. then we tested the integration idea by comparing bgmm with its two extreme cases. performance test of bgmm with artificial data to evaluate the overall performance of a clustering method, we developed one scoring system to evaluate the clustering accuracy when dealing with artificial data. it searches the best matching between the cluster labels of the results and the ground truth clustering among all of their possible associating ways. the score for the best match is denoted as'e score', and is defined as in this scoring system, tj denotes the ground truth clustering membership of data j; r stands for all possible associating ways between the estimated and the true clusters, where ri is the label of data belonging to component i predicted by the clustering algorithm, and r is chosen from labels...,max{, g}. also note that e represents the individual score of each gene, e is the average score of all the genes for each repetition,'e score' of each repetition is the one corresponding to the optimal q, and the final'e score' of each data set is the median of the'e score's. it is worth noticing that the estimated and assumed number of clusters, and g, vary with the model selection criteria, and thus cause different'e scores', rendering this scoring system not only records the accuracy of the results but also reflects the influence of model selection criterion. performance test of different em algorithms in bgmm we first compared the performance of ems, ema and emh in bgmm. for simplicity, we denote bgmm that employs ems, ema or emh as bgmms, bgmma or bgmmh, correspondingly. the artificial data set for the performance test was designed according to our model, whose parameters are listed in table the data set was divided into high quality and low quality data, namely'gb','bb','gg' and'bg' respectively. we also designed two kinds of'bg's,'bgm' and'bgv', which were hard to be clustered compared to'gg' with respect to close means and large variances, respectively. the data set was designed to have three underlying clusters, genes and four features. the simulation was repeated times with randomly generated data sets, and the comparison results of the clustering accuracy were depicted in figure'gb','bb' and'gg' stand for good, bad beta distributed data and good gaussian distributed data re-spectively;'bgm' and'bgv' represent bad gaussian distributed data which are hard to be clustered with respect to close means and large variances respectively. in order to choose the best model selection criterion for each type of bgmm, we summed up the number of hits of the correct number of clusters for each tested case. the summation results for bic, icl, aic and aic are and respectively, in bgmms, and respectively, in bgmma, and and respectively, in bgmmh. therefore, icl is upheld by bgmms, and aic is embraced by both bgmma and bgmmh in this simulation. we evaluated the clustering accuracy of different types of bgmm with each best model selection criterion. simulation results show that, although different algorithms perform slightly different for different cases, the overall prediction accuracy of the three methods are similar as shown in figure we also compared the running time of the three methods under the same background framework, where no significant difference among them was detected. one important application of the proposed algorithm is to cluster genes based on protein-dna binding probabilities and gene expression data, which are assumed to be of beta and gaussian distributions in bgmm. this parametric assumption is supported by our good clustering results and additional distributional assessments. in some cases, however, our parametric assumptions might be violated due to various reasons, especially for expression data. for example, different platforms used to measure transcriptome might affect the distribution of expression data. although this problem can be solved by extending the current algorithm to other parametric distributions quite easily, it is important to know how sensitive bgmm is to the violation of the parametric assumptions and how robust the algorithm is in dealing with noisy variables. to address this, we run three additional simulations with the three em algorithms, where gene expression and protein-dna binding data are simulated from laplace and kumaraswamy distributions, respectively. laplace and kumaraswamy distributions are used to replace gaussian and beta distributions separately in simulation and and both distributions are replaced in simulation note that laplace and kumaraswamy distributions have the same support as gaussian and beta distributions, respectively. means and variances used in laplace distribution are the same with those of gaussian distribution, and α's and β's used in kumaraswamy distribution are also the ones used in beta distribution. as shown in figure all three em algorithms work similarly, and are not excessively sensitive to the parametric assumptions used in this study. parameters of the beta distributed data and the symbols are the same with what has been shown in table based on the above test, the three em algorithms perform equally well. therefore, we simply used bgmmh for the performance tests and referred to it as'bgmm' for simplicity in the following text. performance test of bgmm with its component models simulations shown in this section were dedicated to test how well bgmm could integrate different data sources. we compared the performance of bgmm with its two extreme models, bmm with ema and gmm with ems, for this purpose. a slightly different data set was used, where the gaussian distributed data was designed to be less distinguishable than what has been shown in the previous section. the parameters of the redesigned data are shown in table where all the rest information including the dimensions of the data and the repetitions remain the same. we used the same method as what we did in the previous section to select the best criterion for bgmm, bmm and gmm. ordered by bic, icl, aic and aic the summations of the hits are respectively, in bgmm, respectively, in bmm, and respectively, in gmm, according to which aic was chosen as the best criterion for all the three tested models in this simulation. the comparison results of bgmm with its extreme models are shown in figure for expression data whose variances are not too large, the joint model can improve the clustering accuracy regardless of the quality of the data compared with either of its extreme models. however, when gaussian distributed data has too much overlap among the clusters, bgmm does not necessarily show its superiority when the variances are too large as shown in the case of bgv. it is indicated that bgmm is sensitive to the variances of gaussian distributed data since bgv is designed to have similar noise level as that of bgm. these results demonstrate that the em algorithm of bgmm has the power of reinforcing each extreme model with information from the other one, but does not necessarily outweigh both of them if the gaussian distributed data contains too much noise with respect to large variances. performance test of bgmm with real data we applied our methods to mouse protein-dna binding data and gene expression data. the binding data is modeled as beta distribution, which are the binding probabilities output from a method called'probtf'. probtf uses genome sequences and transcription factor sequence specificities to compute the protein-dna binding probabilities. this method answers the question of whether the whole gene promoter has one or more binding sites for a tf. since it processes each promoter as a whole, the computational predictions provide insights into the functional role of a tf in the regulatory program of a target gene. the rationale for this is the fact that the higher binding probability anywhere on the promoter implies higher probability of a regulatory relationship. further, the method is able to make use of practically any genome-level information, such as evolutionary conservation, nucleosome positioning, chip-chip, and other prior knowledge. the protein-dna binding data contains the probabilities of tfs binding to genes, calculated with mouse-specific position weight matrices from the transfac database. the gene expression data is modeled as gaussian distribution, which is composed of genes measured from conditions. there are genes measured in both data sets. we removed the genes whose gene expression profiles have low absolute values with matlab function'genelowvalfilter', and then choose genes that have annotations available for sure with the functional classification tool of david database. in the end, we obtained genes for the following studies. to see how well our data satisfy the parametric assumptions, we did the following test. for protein-dna binding data, we grouped all the binding probabilities into two beta-distributed clusters and drew their pdfs, each representing the binding and unbinding cases, respectively. figure shows that the genome-wide binding data can remarkably well be approximated with two beta-distributed components. similarly as shown in figure expression data can be fitted into a gaussian distribution. this agrees with previous studies where gene expression data from a microarray platform is commonly assumed to be normally distributed. although the above preliminary test does not correspond to our clustering method exactly, it demonstrates that our parametric assumptions are indeed reasonably good. the bgmm clustering method effectively increases the number of clusters to which the data is split and further improves the fit to the data. the binding data corresponding to two sets of tfs were chosen to cluster the genes together with its corresponding expression data by bgmm, bmm and gmm. the clustering results were then compared and evaluated by gene ontology. the first set of tfs was randomly chosen with respect to their biological significance, while the second set was carefully selected by our model. there are three subsets of'setrand', each of which was chosen based on certain criterion. we arbitrarily choose three thresholds to be compared with the median of the binding probabilities of a certain tf, and tfs that exceed this threshold are used for clustering. the thresholds for'setrand1' to'setrand3' are, and, respectively, and the number of tfs selected are and correspondingly.'setreal' was selected by bmm. we first clustered the genes based on two sets of tfs by bmm, which were bach and bach combined with mafk, respectively. this is because that the two bach proteins are both reported to interact with mafk protein. then we compared the genes whose cluster has the lowest enrichment score from each clustering result, and the common set which contains genes was chosen. we further clustered all the tfs based on the genes by bmm, and focused on the cluster that contains bach bach and mafk. this cluster turns out to be composed of all the tfs that belong to the families fos, jun, maf and nf-e among our tested tfs, which are all ap- components of the leucine zipper factors class. there are tfs in this cluster, all of which were chosen to form'setreal'. maf family proteins heterodizes with cnc-related bzip factors which include nf-e family proteins; while fos family form hetero or homo dimers with jun family proteins. these dimers bind to dna at certain motif that contains ap- binding sites. the result that our bmm can cluster tfs which have similar binding profiles into one single cluster demonstrates the applicability of bmm, one extreme case of bgmm. go was employed in this study to validate the clustering results. in order to find the most significant annotated terms by looking at the probabilities that the terms are counted by chance, we used the hypergeometric probability distribution to calculate the p-values of gene enrichment score for each cluster by each model with each model selection criterion. we compared the means and medians of those p-values across all the groups clustered by each model, whose results are shown in table it is worth mentioning that the clustering result is obtained by running the algorithm times and taking the one whose expected complete data log-likelihood is the maximum, and each p-value shown in table is the mean or median of the p-values of all the ontology groups corresponding to the best clustering result. from this table, it is clear that, no matter whether the tfs were randomly selected or not, both means and medians of the p-values of bgmm are lower than those of either bmm or gmm, regardless of which aspect was considered and which model selection criterion was used. these results indicate that our bgmm can cluster the genes in a more reasonable way with respect to their biological functions, localizations and processes involved. it is also seen from table that, there are two cases where the four model selection criteria have different prediction results, one is in bmm of the case setrand where the results chosen by aic yields the smallest p-values, and the other is in gmm of the case setreal where aic selects the best model in terms of the smallest p-value, both of which accord well with our simulation results. moreover, the choice of tfs whose binding probabilities are used in clustering does obviously affect the results and, therefore, tfs should be carefully chosen based on biological knowledge of a specific problem. in this study, although binding data of randomly chosen tfs also give lower p-values, the obtained clusters might not provide best insight into our biological problem. we therefore carefully studied the results obtained from setreal, which are discussed below. statistics shown in this table are the group average of the p-values of gene enrichment score. tfs of'setrand1' to'setrand3' were randomly chosen according to certain thresholds, while tfs of'setreal' were carefully selected by bmm. note:'f','c' and'p' stand for the three aspects of gene ontology;'all' means all aspects are included.'4' represents that all the four criteria indicate the same clustering result. statistics shown here are all rounded to decimals. there are eight clusters obtained from setreal by bgmm, among which three groups have p-values below if all the aspects were taken into account. the three clusters were named'clu1' to'clu3' and ordered from the highest average expression level to the lowest. the expression patterns and the medians of the genes within one cluster are shown in figure six toll-like receptor agonists which are cpg, pam2csk pam3csk lps, poly i:c and r were used as the treatments, and four gene knock-out mutants and different time points were included to increase the diversity of the tlr-stimulated gene expression data set and the number of measurements. the first four tlr agonists are bacterial-associated, while poly i:c is viral-associated and r is anti-viral stimuli. they were used here to stimulate tlr-stimulated macrophages, which represent various pathogen-associated molecular patterns. among the genes that have been deleted, adapters myd and ticam could provide a structural platform for the recruitment of kinases and downstream effector molecules, were reported crucial for signaling by most type i il- receptor/tlr family members. however, bjökbacka et al. reported that the majority of the host response to lps is regulated independently of myd and genes appearing to be ticam1-dependent can be classified as both myd88-independent and myd88-dependent.  this paper presents a novel beta-gaussian mixture model, bgmm, for gene clustering from beta distributed and gaussian distributed data. we developed three types of em algorithms for bgmm in this study, whose overall performance are similar according to our simulations. we simply chose emh as the core of bgmm for further performance test, which was done by comparing bgmm with its two component models, bmm and gmm, with both artificial and real data. results from artificial data indicate that our joint model works best if the variances of the gaussian distributed data were not too large, and go validation of the real case studies show that the joint model yields more comprehensive results no matter what model selection criterion is used and whether the data is carefully chosen or not. for the carefully selected real data, we started from limited known tfs and ended up with all the tfs within the tested scope that have the same common features, which demonstrates the usability of one extreme case of bgmm. after clustering the genes with the tfs, we obtained three distinguished gene groups which might be involved in the myd88-dependent tlr-3/ signaling cascades. these results not only tested the performance of the joint model, but also demonstrated its usability in real cases and in some possible applications. the main contribution of this paper is that it has proposed a framework for multiple data integration through mixture modeling that has not been addressed by anyone else before. the proposed bgmm is designed to integrate beta distributed and gaussian distributed data. however, the way how those data are incorporated is not limited to the data types that we have used in this study. in principle, data of other parametric distribution can be easily integrated by combining its particular em algorithm into this framework. therefore, the framework proposed in this paper is applicable to many other problems and not limited to the particular problem considered here. one of the basic assumptions in this paper is that the ground truth clustering for gaussian and beta distributed data are the same. this is because transcriptional regulation is largely controlled by the tfs that bind to the gene promoters, thus the expression profiles of genes whose regulatory regions are bound by the same/similar factors are expected to be similar. although the above statement is generally true, it might be violated due to post-transcriptional modifications etc., in which case the method may not be directly applicable. however, if post-transcriptional or other phenomena become a real problem, it can be compensated by integrating more information sources, such as protein-protein interactions, into the proposed clustering framework. on the other hand, if the two data sources do not share the same clustering structure, then an alternative modeling strategy would be needed, such as a hierarchical bayes model that would model a true clustering structure but allow individual structures for both data types. another issue that is worth mentioning is how different data pre-processing and microarray platforms can affect the distribution of gene expression data and, thereby, clustering results. fortunately, as discussed in section, the alternative distributions we tried on bgmm had a very small effect on the clustering results, suggesting that bgmm is considerably robust to small fluctuations in the distributional assumptions. more importantly, one major advantage of bgmm is its flexibility of easily being extended to other parametric distributions. that is, if in a particular problem data come from different distributions, then one can relatively straightforwardly develop a similar model-based approach as proposed here to model the problem at hand in a precise way by fitting data to those specific distributions. we employed the diagonal covariance matrix model in the em algorithm of gmm to reduce the number of parameters to be estimated so that it can be easily applied to large dimensional real data. in particular, the number of parameters in diagonal matrix is p which is remarkably smaller than that of the full covariance matrix. diagonal covariance matrix automatically assumes no correlations among gaussian distributed data. so if we want to preserve the correlation information among time series data by the proposed framework without introducing too many new parameters, it is possible, e.g., to develop similar estimation algorithms for a covariance model where off-diagonal constant correlations are assumed or use a more general covariance matrix. in the future, we could improve the proposed model so that it can account for the correlations among gene expressions. we could also integrate more data sources into this framework and apply it to more real problems. in this aspect, we could either combine other data sources into the framework as a component model, or convert them into prior information which can be used to stratify the model. authors' contributions xfd and hl designed the study and developed the methods. xfd implemented the algorithms, did the performance tests, and wrote the manuscript. te derived the standard em algorithm for bmm. te and hl derived the standard em algorithm for bgmm. xfd, hl, te and oy-h prepared the manuscript. all authors have read and approved the final manuscript. appendix derivation of α's, β's, and π's define recall where,,,, and. derivation of α's and β's define the parameter vector thus the new estimate is obtained as follows where is the hessian matrix evaluated at. note that Ψ and Ψ' represents the digamma and trigamma functions respectively, which are the first and second logarithmic derivatives of the gamma function. derivation of π's 
we developed three types of estimation algorithms for bgmm, the standard expectation maximization algorithm, an approximated em and a hybrid em, and propose to tackle the model selection problem by well-known model selection criteria, for which we test the akaike information criterion, a modified aic, the bayesian information criterion, and the integrated classification likelihood-bic. a number of different approaches have been proposed for that purpose, out of which different mixture models provide a principled probabilistic framework. this paper presents a novel beta-gaussian mixture model for clustering genes based on gaussian distributed and beta distributed data. the proposed bgmm can be viewed as a natural extension of the beta mixture model and the gaussian mixture model. the proposed bgmm method differs from other mixture model based methods in its integration of two different data types into a single and unified probabilistic modeling framework, which provides a more efficient use of multiple data sources than methods that analyze different data sources separately. moreover, bgmm provides an exceedingly flexible modeling framework since many data sources can be modeled as gaussian or beta distributed random variables, and it can also be extended to integrate data that have other parametric distributions as well, which adds even more flexibility to this model-based clustering framework. applications with real mouse gene expression data and protein-dna binding probabilities also demonstrate that bgmm can yield more biologically reasonable results compared with either of its two extreme cases. performance tests with simulated data show that combining two different data sources into a single mixture joint model greatly improves the clustering accuracy compared with either of its two extreme cases, gmm or bmm. cluster analysis is increasingly often supplemented with multiple data sources nowadays, and these heterogeneous information sources should be made as efficient use of as possible.
 cluster analysis has become a standard computational method for gene function discovery as well as for more general explanatory data analysis. a number of different approaches have been proposed for that purpose, out of which different mixture models provide a principled probabilistic framework. cluster analysis is increasingly often supplemented with multiple data sources nowadays, and these heterogeneous information sources should be made as efficient use of as possible.  this paper presents a novel beta-gaussian mixture model for clustering genes based on gaussian distributed and beta distributed data. the proposed bgmm can be viewed as a natural extension of the beta mixture model and the gaussian mixture model. the proposed bgmm method differs from other mixture model based methods in its integration of two different data types into a single and unified probabilistic modeling framework, which provides a more efficient use of multiple data sources than methods that analyze different data sources separately. moreover, bgmm provides an exceedingly flexible modeling framework since many data sources can be modeled as gaussian or beta distributed random variables, and it can also be extended to integrate data that have other parametric distributions as well, which adds even more flexibility to this model-based clustering framework. we developed three types of estimation algorithms for bgmm, the standard expectation maximization algorithm, an approximated em and a hybrid em, and propose to tackle the model selection problem by well-known model selection criteria, for which we test the akaike information criterion, a modified aic, the bayesian information criterion, and the integrated classification likelihood-bic.  performance tests with simulated data show that combining two different data sources into a single mixture joint model greatly improves the clustering accuracy compared with either of its two extreme cases, gmm or bmm. applications with real mouse gene expression data and protein-dna binding probabilities also demonstrate that bgmm can yield more biologically reasonable results compared with either of its two extreme cases. one of our applications has found three groups of genes that are likely to be involved in myd88-dependent toll-like receptor 3/ signaling cascades, which might be useful to better understand the tlr-3/ signal transduction. 