Abstract
Large classes of problems can be formulated as inverse problems, where the goal is to find parameter values that best explain some observed measures. The relationship between parameters and observations is typically highly non-linear, with relatively high dimensional observations and correlated multidimensional parameters. To deal with these constraints via inverse regression strategies, we consider the Gaussian Local Linear Mapping (GLLiM) model, a special instance of mixture of expert models. We propose a general scheme to design a Bayesian nonparametric GLLiM model to avoid any commitment to an arbitrary number of experts. A tractable estimation algorithm is designed using variational Bayesian expectation-maximisation. We establish posterior consistency for the number of mixture components after the merge-truncate-merge algorithm post-processing. Illustrations on simulated data show good results in terms of recovering the true number of experts and the regression function.
Keywords: B
ayesian nonparametrics; mixture of experts; inverse problems; Gaussian locally-linear mapping models; linear cluster-weighted models; variational inference; clustering; regression; model selection.
1 Introduction
Inverse problems.
Many problems can be formulated as inverse problems, where the goal is to find parameter values that best explain an observed phenomenon. Typical constraints in practice are that the relationships between parameters and observations are highly non-linear, with relatively high dimensional observations and correlated multivariate parameters. To handle these constraints, we consider probabilistic mixtures of locally linear models, namely the Gaussian locally linear mapping (GLLiM) approach of Deleforge et al. (2015), which includes affine instances of the mixture of experts (MoE) model (Xu et al., 1995) and the classical inverse linear regression model (Hoadley, 1970; Li, 1991), as special cases. We propose to address inverse problems in a Bayesian framework, making use of the availability of simulations from forward models of interest. The GLLiM approach has been used in many applications, e.g., in medical imaging (Boux et al., 2021), planetary science (Kugler et al., 2022; Nguyen et al., 2024a), head-pose estimation problem in computer vision (Lathuilière et al., 2017) and quantitative trait prediction from biological data (Blein-Nicolas et al., 2024).
Mixture of experts models (MoE)
are generalizations of neural network architectures proposed by Jacobs et al. (1991). Further, these flexible models also generalize the classical mixture models (MMs) and mixture of regression models (McLachlan and Peel, 2000). Their flexibility comes from the fact that they allow the mixture weights (or the gating functions) to depend on the explanatory variables, together with the component densities (or the experts). In regression, MoE models with Gaussian experts and softmax or normalized Gaussian gating functions (as in GLLiM) are the most popular choices. These models are powerful tools for modelling complex nonlinear relationships between outputs (responses) and inputs (predictors) that arise from different subpopulations. The popularity of these conditional mixture density models is largely due to their universal approximation properties (Norets, 2010; Nguyen et al., 2016, 2019, 2021a) as well as their good convergence rate, see, e.g., for mixture of regression in Ho et al. (2022) and for MoE in Jiang and Tanner (1999); Nguyen et al. (2024c, 2023a, b). It is worth noting that these results improve the approximation capabilities and convergence rates of unconditional MMs, as discussed in Genovese and Wasserman (2000); Rakhlin et al. (2005); Nguyen (2013); Shen et al. (2013); Ho and Nguyen (2016a, b); Nguyen et al. (2020, 2023b). At a high level, universal approximation theorems state that given a large enough number of components, MMs and MoE models can approximate a large class of unconditional and conditional probability density functions (cPDF), respectively, to any degree of accuracy. See, e.g., Yuksel et al. (2012); Nguyen and Chamroukhi (2018); Do et al. (2023); Nguyen (2021); Chen et al. (2022), for further detailed reviews of practical and theoretical aspects of MoE models in statistics and in diverse domains, e.g., natural language processing and computer vision.
Model selection in MoE models.
Although universal approximation allows us to conclude that, given a sufficient number of components, a finite MoE can approximate any other cPDF to an arbitrary degree of accuracy, it is not clear how to choose a large enough number of components for realistic problems. This motivates a careful study of interesting and important model selection problems for MMs and MoE models, which have attracted much attention in statistics and machine learning over the last 50 years, see, e.g., Celeux et al. (2019); Gormley and Frühwirth-Schnatter (2019) for a recent comprehensive review.
When selecting the best data-driven number of components for MoE models, there are several approaches to controlling and accounting for model complexity. Typically, model selection in MoE models is performed using an information criterion, such as the Akaike information criterion (AIC; Akaike, 1974; Frühwirth-Schnatter et al., 2018), the Bayesian information criterion (BIC; Schwarz, 1978; Berrettini et al., 2024, BIC-GLLiM; Forbes et al., 2022a, b), the BIC-like approximation of integrated classification likelihood (ICL-BIC; Biernacki et al., 2000; Frühwirth-Schnatter et al., 2012), extended BIC for mixture Gaussian graphical models (eBIC; Nguyen and Li, 2024; Foygel and Drton, 2010), or Sin and White information criterion (SWIC; Sin and White, 1996; Westerhout et al., 2024) for stationary sequences that are not independently and identically distributed. However, an important limitation of these criteria is that they are only asymptotically valid. This means that there are no finite sample guarantees when using AIC, BIC, ICL-BIC or SWIC to choose between different levels of complexity. Therefore, their use in small sample settings is ad hoc.
To overcome such difficulties, and to partially support the so-called slope heuristic approach (Birgé and Massart, 2007, also see Arlot, 2019 for a recent review), Nguyen et al. (2021b, 2022c, 2022b, 2023d, 2023c) recently established non-asymptotic risk bounds in the form of weak oracle inequalities, provided that lower bounds on the penalties hold, in high-dimensional regression scenarios for a variety of MoE models, including GLLiM. Another approach, as proposed by Nguyen et al. (2022a), is based on the closed testing principle, which leads to a sequential testing procedure that enables confidence statements regarding the order of a finite mixture model. These methods facilitate an optimal data-driven selection of the number of components in finite-sample settings. Additionally, Bayesian model selection, which employs marginal likelihoods as a criterion, provides another strategy for determining the number of components in MoE models (Frühwirth-Schnatter, 2019; Zens, 2019). However, all previous approaches require that a range of models with different values be trained and compared, which can be a computational bottleneck in a high-dimensional framework.
Recently, Kock et al. (2022) proposed computationally efficient variational inference approaches for architecture selection in high-dimensional deep Gaussian mixture models using overfitted mixtures (see, e.g., Rousseau and Mengersen, 2011; Forbes et al., 2019), where unnecessary components are dropped in the estimation. The same idea is applied in the context of overfitting sparse Bayesian factor analysis in Frühwirth-Schnatter et al. (2024), where the authors employed the Markov Chain Monte Carlo (MCMC) procedure, followed by post-processing of posterior draws, to facilitate practical Bayesian inference. However, we are interested here in the more general context of the Bayesian nonparametric (BNP) approach, (see, e.g., Hjort et al., 2010; Ghosal and Van der Vaart, 2017), where it is not necessary to know an upper bound of the true number of components as in Bayesian overfitted MM. This is one motivation for the BNP priors that are considered here for GLLiM.
Dirichlet process mixture models (DP-MMs) and Pitman-Yor process mixture models (PYP-MMs) are among the most popular BNP models, particularly suitable for density estimation and probabilistic clustering. However, the posterior of the DP-MMs or PYP-MMs are inconsistent for the number of components if the true number of components is finite and the concentration parameter is fixed (see, e.g., Miller and Harrison, 2014, and Alamichel et al., 2024 for a review). This is because a BNP prior such as DP or PYP places zero probability on mixing measures with a finite number of supporting atoms. An interesting recent result in Ascolani et al. (2022) is that consistency for the number of components can be achieved if a prior is placed on the concentration parameter of the DP-MM, under some assumptions on this prior.
It appears that BNP-GLLiM tends to produce many small extraneous components around the true clusters in our numerical experiment Section 6. This makes it difficult to use them to infer the true number of components when this becomes a quantity of interest (Maceachern and Müller, 1998; Green and Richardson, 2001). This encourages the use of a novel, simple post-processing algorithm in the spirit of the Merge-Truncate-Merge (MTM) introduced by Guha et al. (2021). This post-processing procedure consistently estimates the number of components for any general Bayesian prior, even without knowing its exact structure, as long as the posterior for that prior contracts to the true mixing distribution at a known rate.
Contributions.
To address the challenges of highly nonlinear inverse problems with relatively high dimensional observations and correlated parameters, we propose a novel BNP-GLLiM model and inference procedure, which is computationally efficient and avoids any commitment to an arbitrary number of components. Although BNP mixture models (BNP-MM), which are special cases of the BNP-GLLiM model, have been extensively studied in the literature (Escobar and West, 1995; Maceachern and Müller, 1998; Neal, 2000; Arbel et al., 2021; Li et al., 2022; Durand et al., 2022), the extension to MoE models in an inverse regression framework has not been covered. In addition, we establish theoretical properties such as posterior consistency for recovering the true number of components in BNP-GLLiM using the post-processing merge-truncate-merge (MTM) algorithm. Finally, our illustrations on simulated data show good results in terms of recovering the true number of components and mean regression functions. It is worth emphasising that, for the first time, we provide evidence that MTM consistency holds not only for the MMs results of Guha et al. (2021), but also in the more general context of MoE models for cPDFs.
Notations.
Throughout this paper, is abbreviated as for , where denotes the positive natural numbers. The notation refers to a definition. It is used to simplify the notation or expression. For a parametric model , refers to its dimension, i.e., the total number of parameters to be estimated. Furthermore, denotes vertical vector concatenation. Throughout, we use the following colour rule for observations and parameters: observations are represented in green; latent, random or unknown parameters in red; and (fixed) hyperparameters in blue.
Outline.
The paper is organized as follows. In Section 2, we first discuss how to construct the BNP-GLLiM model. A VBEM algorithm and the corresponding ELBO are described in Section 3, and predictive cPDFs in Section 4. Next, Section 5 shows how we can integrate the Merge-Truncate-Merge (MTM) post-processing procedure and prove consistency for the MTM output. This is useful to perform regression, clustering and model selection, simultaneously. We experimentally evaluate our new results on simulated datasets in Section 6. Some perspectives are provided in Section 7. We recall the standard Bayesian nonparametric priors and variational Bayesian expectation-maximisation principle in Appendices A and B, respectively. All details of VBEM for the BNP-GLLiM model, evidence lower-bound, and technical proofs not included in the main paper are relegated to Appendices C, D and E, respectively. Appendix F presents a more general model with a hyperprior on the gating parameters, called BNP-GLLiM2.
2 BNP-GLLiM model
In Sections 2.1 and 2.2, we present the advantage of adopting a GLLiM model, which uses an inverse regression approach to estimate nonlinear high-to-low dimensional mappings. Such a strategy allows to greatly reduce the number of required parameters. To avoid any commitment to an arbitrary number of components, we then construct the BNP-GLLiM model in Section 2.3.
2.1 Inverse regression framework
We are interested in estimating the law of a low-dimensional random variable conditionally on a high-dimensional , where typically . We follow an inverse regression framework as in e.g., Li (1991); Deleforge et al. (2015). Therefore, in training, the low-dimensional variable plays the role of the regressor, while the response is a function of , possibly corrupted by noise through inverse cPDF ), where is an inverse parameter. The low dimension of the regressor allows to drastically reduce the number of parameters to be estimated. In addition, the forward parameter and cPDF are tractable after estimating the inverse parameter . Therefore, this density can be used to predict the low-dimensional response of a high-dimensional test point . This inverse-then-forward regression strategy justifies the unconventional notation: for the high-dimensional input and for the low-dimensional response. Here and subsequently, we refer to the low-dimensional data sample as , the high-dimensional data sample as . We denote the observed values as , which are independently and identically distributed (i.i.d.) samples from the random variables .
2.2 Nonlinear high-to-low dimensional mapping via GLLiM model
The GLLiM models, as originally introduced in Deleforge et al. (2015), are used to capture the non-linear relationship between the response and the set of covariates from a high-to-low heterogeneous data. Specifically, Deleforge et al. (2015) overcame the difficulty of high-to-low regression by tackling the problem the other way round, i.e., low-to-high. This means that the roles of input and response variables are swapped so that the low-dimensional variable becomes the regressor as in Section 2.1. GLLiM then relies on a piecewise linear model in the following way. The high-dimensional response is approximated by the local affine mappings :
Here, is an indicator function and is a latent variable that captures a cluster relationship, such that if comes from cluster . Matrices and vectors define cluster-specific affine transformations. In addition, are error terms that capture both the reconstruction error due to the local affine approximations as well as the observation noise in .
Following the usual assumption that is a zero-mean Gaussian variable with a covariance matrix , it follows that
| (1) |
where is the vector of model parameters, is the Gaussian cPDF of dimension . In order to enforce the affine transformations to be local, is defined as a mixture of Gaussian components as follows:
where , , and belongs to a probability simplex, defined as . Then, via the conditional property of Gaussian variables and hierarchical decomposition, given any , we obtain the following inverse conditional density:
| (2) |
Using the inverse regression framework, in (1), the roles of input and response variables should be reversed so that becomes the covariate and plays the role of the multivariate response. Therefore, based on a similar previous hierarchical one in (2), its corresponding forward conditional density from to is defined by
| (3) |
A useful feature of GLLiM models is described in the following Lemma 2.1, established for multivariate Gaussian and Student components in Deleforge et al. (2015); Perthame et al. (2018) and which can be straightforwardly extended to Gaussian scale mixtures and elliptical distributions (Nguyen et al., 2022c; Ingrassia et al., 2012).
Lemma 2.1.
Remark 1 (GLLiM models are computationally efficient).
Without assuming anything about the structure of the parameters, the dimension of GLLiM is It is worth noting that can be very large compared to the sample size (see, e.g., Deleforge et al. 2015 for real data sets) whenever is large and . Furthermore, under the assumption that the transformations are affine, it is more realistic to make the assumption on the residual covariance matrices of the error vectors rather than on (cf., Section 3 Deleforge et al., 2015). This justifies using the inverse regression trick from Deleforge et al. (2015), drastically reducing the number of parameters to be estimated. For instance, can be modelled with equal isotropic Gaussian noise, so we have , with some positive . The number of parameters to be estimated is then . For example, it is if , , . However, if a high-to-low regression is estimated directly instead, the size of the parameter vector will be , which is in the previous example.
A notable recent illustration of the GLLiM good features is described in Blein-Nicolas et al. (2024), tackling the challenge of nonlinear quantitative trait prediction using biological data. In this work, the authors focused on predicting a small set of continuous quantitative traits () from a large set of biomarkers (). Their inverse regression approach is not only computationally efficient, but also has the advantage of preserving all covariates while assuming a specific structure for the covariance matrix, which aids dimensionality reduction and improves prediction accuracy. This method operates by estimating an inverse model through several low-dimensional regressions and then inverting these estimators to solve the initial high-dimensional regression problem using Lemma 2.1.
2.3 Construction of BNP-GLLiM model
We propose the following hierarchical representation of BNP-GLLiM model to generate a data point within our BNP-GLLiM model:
-
1.
BNP prior: .
(4) (5) Define , , where , and finally define .
-
2.
BNP-GLLiM model: for each , .
(6) (7)
3 Variational inference for BNP-GLLiM
For a brief summary of variational Bayesian expectation maximisation (VBEM) and its notation, see Appendix B. The task of conditional density estimation and clustering using the BNP-GLLiM model is mainly to estimate the unknown labels from the observed data , whose joint distribution is determined by a set of BNP prior parameters , namely the stick-breaking construction of Pitman–Yor process (PYP) (Pitman and Yor, 1997) in Appendix A; and by additional hyperparameters . Then the desired joint distribution is given by:
| (8) |
In most variational approximations, the posterior for the stick-breaking variables is approximated in a factorized form (mean-field approximation). Following the same approach, by factorizing the latent variables and the parameters, we choose the following variational distribution: .
In particular, the intractable posterior on is approximated as that factorizes so as to handle intractability, namely
| (9) |
Then the infinite state space for each is dealt with by choosing a truncation of the state space to a maximum label , see, e.g., Blei and Jordan (2006); Wang et al. (2011). In practice, this consists of assuming that the variational distributions for , satisfy for and that the variational distribution on also factorizes as with the additional condition that . Thus the truncated variational posterior of parameters is given by
In practice, for tractability reasons, we have to restrict to and to Dirac distributions, , which is equivalent to treating the as fixed unknown hyperparameters, as illustrated graphically in Figure 1. These forms of and lead to our three VB-E steps and four VB-M steps, summarized below and with more detail in Appendix C. Set the initial value of to . Then repeat the following steps, iteratively. The iteration index is omitted in the update formulas for simplicity. Note that a more complex version with a normal-inverse-Wishart (NIW) distribution on the gating parameters , referred to as BNP-GLLiM2, is presented in Appendix F.
3.1 VB-E steps
VB-E- step.
The VB-E- step corresponds to a variational approximation in the exponential family case and results in a posterior from the same family as the prior. More precisely, to achieve this, we use (5), (6), (8), and are only interested in the functional dependence of (22) on the variable . Thus, any terms that do not depend on can be included in the additive normalization constant. Then, given for that corresponds to the weight of the cluster , see more details in Section C.1, it holds that
VB-E- step.
The variational posterior is more complex, but has a simple form in the DP case (). Specifically, we have to compute
| (10) |
where is the digamma function defined by .
When then is a gamma distribution with . Otherwise (PYP case), is only identified up to a normalizing constant but the required and can be computed by importance sampling, see Section C.2 for more details.
We next consider the derivation of the update equation for the factor .
VB-E- step.
By using the mean-field approximation (9) and the truncation (see Section C.3), for all and for all , this step consists in computing
| (11) | ||||
Note that in the above formula, symbols are the hyperparameters Specifically defined in the following Sections 3.2 and 3.2. It is important to note that (11) provides assignment probabilities rather than intermediate commitments to hard assignments of . However, the hard assignments can be postponed to the end if desired to obtain a segmentation by the following MAP estimation .
3.2 VB-M steps
The maximisation step consists of updating the hyperparameters , where , by maximizing the free energy as follows:
The VB-M-step can therefore be divided into 4 independent sub-steps, as listed below. From the conditional independence of and given , the solution for the VB-M- (in the DP case) step is straightforward. Only the M- (in the PYP case) and M- steps are more involved.
VB-M- step.
This step is straightforward in the DP case (). It can be expressed easily using the fact that both the prior and the variational posterior are Gamma distributions, and using the cross-entropy properties,
where is given in (10). We can also solve this step numerically using importance sampling in the more general case of PYP (). For more details, see Appendix A.7 in Lü et al. (2020).
VB-M- step.
This step is divided into solving optimisation problems:
We can then update the Gaussian gating parameters as follows:
The technical details will be left to the Section C.4.
VB-M- step.
Using the same idea, this step is divided into sub-steps, which include the following optimisation problems
Given the following quantities:
We can update the parameters for the Gaussian experts as follows (cf. Section C.5):
3.3 Evidence lower-bound (ELBO)
Evaluating the ELBO in (20) allows us to not only monitor the bound during the re-estimation to test for convergence but also to check both the mathematical expressions for the solutions and their software implementation. Indeed, the value of this bound (20) at each step of the iterative re-estimation procedure should not decrease (Svensén and Bishop, 2005), in particular, see recent results for Bayesian nonparametric mixture models in Appendix A of Durand et al. (2022). Recall that . Here, in order to keep the notation uncluttered, we will sometimes omit the subscripts on the expectation operators because each expectation is taken with respect to all of the random variables in its argument, and the hat superscript on the hyperparameters of distribution. If and there are enough training data, the ELBO can be evaluated via the fact that the integral reduces to a point evaluation at the posterior mean of each parameter, see, e.g., Yuan and Neubauer (2008); Luo and Sun (2017); Wu and Ma (2018); Nguyen and Bonilla (2014). When , we can analytically compute the ELBO in the BNP-GLLiM via Proposition 3.1 which is proved in Section E.2.
Proposition 3.1.
If , the ELBO in the BNP-GLLiM is decomposed as follows:
| (12) |
where all the terms have a closed-form expression.
The closed-form expressions for the terms of the right-hand side of Equation 12 are provided in Appendix D. Note that if the free energy is computed at the end of each VBEM iteration, as in Section 3.2, we have .
4 Predictive conditional density
The most popular uses of BNP-GLLiM with discrete random probability measures, such as the one displayed in (7), relate to conditional density estimation and data clustering. Specifically, we are interested in the predicted conditional density for a new value of the observed variables. Note that there will be a corresponding latent variable associated with these observations. If , we can use the previous remark in Section 3.3, where the integral reduces to a point evaluation at the posterior mean of each parameter. When , we can analytically approximate such densities via several following theorems. In Theorems 4.1, 4.2 and 4.3, the notation “” means that we approximate the desired densities of the BNP-GLLiM by a mixture of Gaussians using factorized variational approximation posteriors and a truncation of .
4.1 Joint density
We first compute the joint density via Theorem 4.1, which is proved in Section E.3.
Theorem 4.1.
With , we have
| (13) | ||||
| (14) | ||||
4.2 Inverse conditional density
We then show how to approximate the inverse conditional density . This predictive density in BNP-GLLiM is approximated by a GLLiM via Theorem 4.2 with the proof in Section E.4.
Theorem 4.2.
We approximate the inverse conditional density and its conditional expectation for prediction by:
| (15) | ||||
4.3 Forward conditional density
Given the inverse conditional density , we approximate the following forward conditional density via Theorem 4.3, whose proof is provided in Section E.5.
Theorem 4.3.
We approximate the forward conditional density and its corresponding conditional expectation and variance for prediction and uncertainty estimation by
5 Bayesian nonparametric model selection
Notations. A coupling between and is a joint distribution on , which is expressed as a matrix with marginal probabilities and , for any and . We use to denote the space of all such couplings. Regarding the space of mixing measures, let and respectively denote the space of all mixing measures with exactly and at most support points, all in some parameter space . Additionally, with to denote the set of all discrete measures with finite supports on . Moreover, denotes the space of all discrete measures (including those with countably infinite supports) on . Finally, stands for the space of all probability measures on .
5.1 Posterior contraction rate in Bayesian infinite mixtures
Problem setup.
We first recall the GMM where we have i.i.d. samples coming from a true but unknown distribution with given PDF
| (16) |
where is a true but unknown mixing distribution with exactly number of support points, where is also unknown. Furthermore, is a chosen parameter space, to which we believe that the true parameters belong. In a well-specified setting, all support points of reside in , but this may not be the case in a misspecified setting. In this section, we assume that the GMM is well-specified, i.e., the data are i.i.d. samples from the mixture density , where mixing measure has support atoms in compact parameter space .
A Bayesian modeller places a prior distribution on a suitable subspace of . Then, the posterior distribution over is given by:
Here, the GMM is defined in (16) with unknown number of support points. We are interested in the posterior contraction behaviour of toward , in addition to recovering the true number of components .
We next recall the notion of Wasserstein distance for mixing measures that prove useful in the next sections.
Wasserstein distance for MM.
It is useful to analyze the identifiability and convergence of parameter estimation in mixture models using the notion of Wasserstein distance, as in Nguyen (2013); Ho and Nguyen (2016b). This distance can be defined as the optimal cost of moving masses in the transformation from one probability measure to another (Villani, 2003, 2009).
Definition 5.1.
Suppose is equipped with a metric . The Wasserstein distance between two discrete measures and is given by
where couplings are defined at the beginning of this section. See Delon and Desolneux (2020) for more details.
It should be emphasized that if a sequence of probability measures converges to under the metric at a rate for some , then there exists a subsequence of such that the set of atoms of converges to the atoms of , up to a permutation of the atoms, at the same rate .
Posterior contraction rate in infinite mixtures.
With a similar idea as in Guha et al. (2021), our starting point is the availability of a mixing measure sample that is drawn from the posterior distribution , where are i.i.d. samples of the mixing density . Under certain conditions on the kernel density, it can be established that for some Wasserstein metric ,
| (17) |
for all constant , while is a vanishing rate. Thus can be assumed to be (slightly) slower than the actual rate of posterior contraction of the mixture measure. We can also write that is a rate such that, under the posterior distribution , . See Nguyen (2013); Gao and Vaart (2016); Ho and Nguyen (2016b) for concrete examples of posterior contraction rates in infinite and (overfitted) finite mixtures.
5.2 Merge-Truncate-Merge (MTM) algorithm for BNP-GLLiM
Link between GLLiM and joint GMM.
We start by noting that a GLLiM model on , see (2), with unconstrained parameters , is equivalent to a GMM on the joint variable with unrestricted parameters, via Lemma 5.2, which is briefly proved in Section E.1.
Lemma 5.2.
A GLLiM model on with unconstrained parameters , defined in (2), is equivalent to a GMM on the joint variable with unconstrained parameters , i.e.,
The parameter can be expressed as a function of by:
| (18) |
Here, we have defined
Note that the symmetry of the covariance matrix implies that and are symmetric, while . The parameter vector can be expressed as a function of by:
| (19) |
Merge-Truncate-Merge (MTM) algorithm
is a post-processing procedure applied to a posterior sample of the mixing measure in BNP-MM, essential for achieving posterior contraction rates under the Wasserstein metric (Guha et al., 2021). We propose in Algorithm 1 an algorithm for BNP joint GMM which follows the same steps as the original MTM algorithm for BNP-MM from Guha et al. (2021). This algorithm involves two main stages: the Merge procedure and the Truncate-Merge procedure. In the Merge stage, atoms are reordered by simple random sampling without replacement, ensuring random permutation. Sequentially, atoms are merged based on a distance threshold , updating weights and removing merged atoms from the set, resulting in a new measure with reordered weights. In the Truncate-Merge stage, atoms are divided into two sets based on a weight threshold . For each atom in the significant weight set, if another atom within the threshold distance exists, it is moved to the negligible weight set. Atoms in the negligible set are then merged with the nearest significant atom, resulting in the final measure and the number of its supporting atoms . See Algorithm 1 for a pseudo-code description.
As a consequence of our MTM algorithm, we obtain the theoretical guarantee of Theorem 5.3 for the outcome of Algorithm 1.
Theorem 5.3 (MTM consistency for BNP joint GMM).
Let be a posterior sample from the posterior distribution of any Bayesian procedure, namely, according to which the upper bound (17) holds for all . Let and be the outcome of Algorithm 1 applied to , for an arbitrary constant . Then the following hold
-
(a)
in -probability.
-
(b)
For all , in -probability.
Proof of Theorem 5.3.
Lemma 5.2 implies that BNP joint GMM and BNP-GLLiM are considered equivalent with respect to the number of components in the model selection problem. Therefore, using Theorem 3.2 from Guha et al. (2021) and Lemma 5.2, it follows that the result of the MTM Algorithm 1 for BNP joint GMM is a consistent estimate of both the number of components and the mixing measure. The latter also admits the upper bound of the posterior contraction rate , which leads to the desired Theorem 5.3. ∎
Remark 2.
Regarding the above theorem, we provide the following comments on posterior consistency for the number of components in BNP-GLLiM after the MTM algorithm post-processing.
-
(i)
As a complementary result to Guha et al. (2021), the aim of this paper is to study the practical viability of MTM Algorithm 1 and Theorem 5.3 in the context of high-to-low dimensional inverse regression via BNP-GLLiM model. In order to do this, we first need to specify the metric in , e.g., Here, denotes either the -norm elements in or the entrywise -norm for matrices in .
-
(ii)
In practice, one may not have a mixing measure sampled from the posterior , but rather a sample of itself. In particular, to deal with large data sets, we need to use VBEM. Therefore, in BNP-GLLiM, instead, we only obtain a sample from the variational posterior . Here, and are defined in (14) and (13), respectively. However, as long as is sufficiently close to in the sense that , we can still apply the MTM algorithm to , instead. This requires an extension of the above Theorem 5.3 to cover this scenario and verify this approximation condition, which we leave for future work.
6 Numerical experiments
The code to reproduce our simulation study is publicly available111https://github.com/Trung-TinNGUYEN/BNP-GLLiM and all simulations below were performed in Python 3.9.13 on a standard Unix machine. For proof-of-concept numerical experiments, we consider only the simple data generating mechanism with and and demonstrate that BNP-GLLiM performs well in model selection, clustering and cPDF estimation with the MTM procedure. Real world examples are postponed to future work.
6.1 Data generation
We illustrate our theoretical results on simulated datasets in a more general setting for the BNP approach compared to those considered by Chamroukhi et al. (2010); Montuelle and Le Pennec (2014); Nguyen et al. (2022c). Specifically, we consider the following true inverse cPDF from GLLiM model as follows:
Here , , and , where
Figure 2a shows typical realisations of the true inverse cPDF from GLLiM, representing a -shape simulation with three clusters without labels.
6.2 Model selection, clustering and regression tasks via MTM-BNP-GLLiM
Our goal is to evaluate the inverse and forward cPDF, as well as the conditional means, to investigate the empirical performance of our MTM-BNP-GLLiM in the previous simulation. In Figure 2 it is clear that with the help of MTM Algorithm 1, MTM-BNP-GLLiM can simultaneously perform regression, clustering and model selection well. Without the MTM procedure, BNP-GLLiM performs poorly in model selection, clustering and cPDF estimation, except for conditional expectations as shown in Figure 3.
Next, we illustrate the performance of the MTM algorithm when applied to the variational posterior from BNP-GLLiM. Specifically, the samples in our trials are drawn from , wherere and are defined in (14) and (13), respectively. We know that for some constant , which depends on the covariance matrix , the location parameters and the weights , the contraction rate of mixing measures under the location Gaussian DP-MM is with respect to the -norm. Similar to Guha et al. (2021), our first attempt to choose to satisfy (17) is . In fact, we can choose any , as long as , in order for to satisfy (17).
Since we only work with finite sample , it is not expected that the posterior probability for is close to and the input to Algorithm 1 should be chosen so that . Furthermore, based on Equation (26) from Guha et al. (2021), with a useful lower bound on the posterior mass the mode, for any , , we hope to identify via the posterior mode with a reasonable estimate. To guarantee consistently using the posterior mode safety, we have to choose , with satisfying
Therefore, we can choose
In particular, it is unrealistic to obtain the exact computation of the upper bound and the lower bound . However, a reasonable estimate may be possible by considering a large range of , and show that there is a range where we can robustly identify the true number of components via the posterior mode. Guha et al. (2021) also used the same setting in their experiments. Figure 4 indicates that leads to a quite good posterior mode in our experiments.
Although we do not have a theoretical result for the convergence rate of the variational posterior of BNP-GLLiM to the true data generating process, Figure 4 seems to suggest that MTM-BNP-GLLiM gives a comparable good result to the location Gaussian DP-MM in the simulation studies in Guha et al. (2021).
7 Perspectives
To address the dimensionality issue when is less than , we could incorporate sparsity penalty terms as described in the GLLiM context in Chamroukhi et al. (2019). Alternatively, we could impose additional structural restrictions, such as block-diagonal covariance matrices as in Blein-Nicolas et al. (2024)), which extends the GLLiM method to account for hidden module-structured regulatory networks of predictors through block-diagonal structures. In particular, Blein-Nicolas et al. (2024) applied their method to the prediction of drought-related traits from protein abundance in maize genotypes. We leave for future research the intriguing but challenging questions of how to integrate BNP priors into the GLLiM model and how to establish posterior convergence theory in high-dimensional settings where .
As indicated in Remark 2, there is a crucial need to formally establish general conditions on the prior, the likelihood and the variational class to characterise the convergence rate of the variational posterior of BNP-GLLiM to the true data generating process. Using the similar “prior mass and testing” conditions as in Ghosal et al. (2000), we believe that an interesting but challenging extension of the work on variational posterior unconditional distributions for MMs (Zhang and Gao, 2020) and on adaptive Bayesian estimation for MMs and MoE models but for true posterior distribution (Kruijer et al., 2010; Shen et al., 2013; Norets and Pati, 2017) can help shed some light and answer this important question. Furthermore, it is important to establish an extensional convergence property of our VBEM algorithm for BNP-GLLiM. This property is only known for GMM from Titterington and Wang (2006). A potential improvement of the VBEM algorithm developed for BNP-GLLiM can be achieved by combining it with MCMC, taking advantage of both inference approaches as in Ruiz and Titsias (2019). Finally, as mentioned in Section 6.2, the selection of a good data-driven tuning parameter as the same idea from the slope heuristic of Birgé and Massart (2007) is crucial for the success of the MTM procedure for any BNP model. We leave these interesting but challenging questions for future research.
Acknowledgments
We would like to express our sincere gratitude to the Reviewers and Editors for their valuable feedback, which has greatly enhanced the quality of this paper. All authors acknowledge funding from the Australian Research Council grant DP230100905, and from Inria Project WOMBAT.
References
- A new look at the statistical model identification. IEEE Transactions on Automatic Control 19 (6), pp. 716–723. External Links: Document Cited by: §1. ↩ a
- Bayesian mixture models (in)consistency for the number of clusters. Scandinavian Journal of Statistics. Cited by: §1. ↩ a
- BNPdensity: Bayesian nonparametric mixture modelling in R. Australian & New Zealand Journal of Statistics 63 (3), pp. 542–564. External Links: Link, Document Cited by: §1. ↩ a
- Minimal penalties and the slope heuristics: a survey. Journal de la Société Française de Statistique 160 (3), pp. 1–106. Cited by: §1. ↩ a
- Clustering consistency with Dirichlet process mixtures. Biometrika 110 (2), pp. 551–558. External Links: ISSN 1464-3510, Link, Document Cited by: §1. ↩ a
- The variational Bayesian EM algorithm for incomplete data: with application to scoring graphical model structures. Bayesian statistics 7 (453-464), pp. 210. Cited by: Appendix B. ↩ a
- Variational algorithms for approximate Bayesian inference. Ph.D. Thesis, University of London, University College London (United Kingdom). Cited by: Appendix B, Appendix B. ↩ a b
- Identifying Brexit voting patterns in the British house of commons: an analysis based on Bayesian mixture models with flexible concomitant covariate effects. Journal of the Royal Statistical Society Series C: Applied Statistics 73 (3), pp. 621–638. External Links: ISSN 0035-9254, Link, Document Cited by: §1. ↩ a
- Assessing a mixture model for clustering with the integrated completed likelihood. IEEE Transactions on Pattern Analysis and Machine Intelligence 22 (7), pp. 719–725. External Links: Document Cited by: §1. ↩ a
- Minimal penalties for Gaussian model selection. Probability Theory and Related Fields 138 (1), pp. 33–73. External Links: Link, Document Cited by: §1, §7. ↩ a b
- Pattern Recognition and Machine Learning (Information Science and Statistics). Springer-Verlag, Berlin, Heidelberg. External Links: ISBN 0-387-31073-8 Cited by: Appendix B, §E.3, §E.5, §F.1.2, §F.4, §F.6. ↩ a b c d e f
- Variational inference for Dirichlet process mixtures. Bayesian analysis 1 (1), pp. 121–143. Cited by: Appendix A, §E.3, §E.3, §F.4, §3. ↩ a b c d e
- Nonlinear network-based quantitative trait prediction from biological data. Journal of the Royal Statistical Society Series C: Applied Statistics, pp. qlae012. Cited by: §1, §7, Remark 1. ↩ a b c d
- Bayesian inverse regression for vascular magnetic resonance fingerprinting. IEEE Transactions on Medical Imaging 40 (7), pp. 1827–1837. Cited by: §1. ↩ a
- Model selection for mixture models–perspectives and strategies. In Handbook of mixture analysis, pp. 117–154. Cited by: §1. ↩ a
- Fast joint detection-estimation of evoked brain activity in event-related fMRI using a variational approach. IEEE transactions on Medical Imaging 32 (5), pp. 821–837. Cited by: Appendix B. ↩ a
- Regularized Estimation and Feature Selection in Mixtures of Gaussian-Gated Experts Models. In Statistics and Data Science, H. Nguyen (Ed.), Singapore, pp. 42–56. External Links: ISBN 978-981-15-1960-4 Cited by: §7. ↩ a
- A hidden process regression model for functional data description. Application to curve discrimination. Neurocomputing 73 (7-9), pp. 1210–1221. Cited by: §6.1. ↩ a
- Towards Understanding the Mixture-of-Experts Layer in Deep Learning. In Advances in Neural Information Processing Systems, A. H. Oh, A. Agarwal, D. Belgrave, and K. Cho (Eds.), External Links: Link Cited by: §1. ↩ a
- Are Gibbs-type priors the most natural generalization of the Dirichlet process?. IEEE transactions on pattern analysis and machine intelligence 37 (2), pp. 212–229. Cited by: Appendix A. ↩ a
- High-dimensional regression with gaussian mixtures and partially-latent response variables. Statistics and Computing 25 (5), pp. 893–911. External Links: ISSN 1573-1375, Link, Document Cited by: §E.1, §1, §2.1, §2.2, §2.2, Remark 1. ↩ a b c d e f g h i
- A Wasserstein-type distance in the space of Gaussian mixture models. SIAM Journal on Imaging Sciences 13 (2), pp. 936–970. Cited by: Definition 5.1. ↩ a
- Maximum Likelihood from Incomplete Data Via the EM Algorithm. Journal of the Royal Statistical Society: Series B (Methodological) 39 (1), pp. 1–22. External Links: ISSN 0035-9246, Link, Document Cited by: Appendix B. ↩ a
- HyperRouter: Towards Efficient Training and Inference of Sparse Mixture of Experts. In Proceedings of the 2023 Conference on Empirical Methods in Natural Language Processing, Singapore. Cited by: §1. ↩ a
- Bayesian non-parametric spatial prior for traffic crash risk mapping: A case study of Victoria, Australia. Australian & New Zealand Journal of Statistics 64 (2), pp. 171–204. External Links: Link, Document Cited by: Appendix A, §1, §3.3. ↩ a b c
- Bayesian Density Estimation and Inference Using Mixtures. Journal of the American Statistical Association 90 (430), pp. 577–588. External Links: Link, Document Cited by: §1. ↩ a
- On the stick-breaking representation for homogeneous NRMIs. Bayesian Analysis 11 (3), pp. 697–724. Cited by: Appendix A. ↩ a
- A Bayesian Analysis of Some Nonparametric Problems. The Annals of Statistics 1 (2), pp. 209 – 230. External Links: Link, Document Cited by: Appendix A. ↩ a
- Component elimination strategies to fit mixtures of multiple scale distributions. In RSSDS 2019 - Research School on Statistics and Data Science, Communications in Computer and Information Science, Vol. 1150, Melbourne, Australia, pp. 81–95. External Links: Link, Document Cited by: §1. ↩ a
- Mixture of expert posterior surrogates for approximate Bayesian computation. In JDS 2022 - 53èmes Journées de Statistique de la Société Française de Statistique (SFdS), Lyon, France. Cited by: §1. ↩ a
- Summary statistics and discrepancy measures for approximate Bayesian computation via surrogate posteriors. Statistics and Computing 32 (5), pp. 85. External Links: ISSN 1573-1375, Link, Document Cited by: §1. ↩ a
- Extended Bayesian Information Criteria for Gaussian Graphical Models. In Advances in Neural Information Processing Systems, J. Lafferty, C. Williams, J. Shawe-Taylor, R. Zemel, and A. Culotta (Eds.), Vol. 23. External Links: Link Cited by: §1. ↩ a
- Sparse Bayesian Factor Analysis When the Number of Factors Is Unknown. Bayesian Analysis, pp. 1 – 48. Note: Publisher: International Society for Bayesian Analysis External Links: Link, Document Cited by: §1. ↩ a
- Labor market entry and earnings dynamics: Bayesian inference using mixtures-of-experts Markov chain clustering. Journal of Applied Econometrics 27 (7), pp. 1116–1137. External Links: Link, Document Cited by: §1. ↩ a
- Analysing plant closure effects using time-varying mixture-of-experts Markov chain clustering. The Annals of Applied Statistics 12 (3), pp. 1796 – 1830. Note: Publisher: Institute of Mathematical Statistics External Links: Link, Document Cited by: §1. ↩ a
- Keeping the balance—Bridge sampling for marginal likelihood estimation in finite mixture, mixture of experts and Markov mixture models. Brazilian Journal of Probability and Statistics 33 (4), pp. 706 – 733. External Links: Link, Document Cited by: §1. ↩ a
- Posterior contraction rates for deconvolution of Dirichlet-Laplace mixtures. Electronic Journal of Statistics 10 (1), pp. 608 – 627. External Links: Link, Document Cited by: §5.1. ↩ a
- Rates of convergence for the Gaussian mixture sieve. The Annals of Statistics 28 (4), pp. 1105 – 1127. External Links: Link, Document Cited by: §1. ↩ a
- Convergence rates of posterior distributions. The Annals of Statistics 28 (2), pp. 500 – 531. External Links: Link, Document Cited by: §7. ↩ a
- Fundamentals of nonparametric Bayesian inference. Vol. 44, Cambridge University Press. External Links: ISBN 0-521-87826-8 Cited by: Appendix A, §E.3, §1. ↩ a b c
- Mixture of experts models. In Handbook of mixture analysis, pp. 271–307. Cited by: §1. ↩ a
- Modelling Heterogeneity With and Without the Dirichlet Process. Scandinavian Journal of Statistics 28 (2), pp. 355–375. External Links: Link, Document Cited by: §1. ↩ a
- On posterior contraction of parameters and interpretability in Bayesian mixture modeling. Bernoulli 27 (4), pp. 2159 – 2188. External Links: Link, Document Cited by: §1, §1, item i, §5.1, §5.2, §5.2, §6.2, §6.2, §6.2. ↩ a b c d e f g h i j k
- Bayesian Nonparametrics. Cambridge Series in Statistical and Probabilistic Mathematics, Cambridge University Press. External Links: ISBN 978-0-521-51346-3, Link, Document Cited by: §1. ↩ a
- Convergence rates of parameter estimation for some weakly identifiable finite mixtures. The Annals of Statistics 44 (6), pp. 2726 – 2755. External Links: Link, Document Cited by: §1. ↩ a
- On strong identifiability and convergence rates of parameter estimation in finite mixtures. Electronic Journal of Statistics 10 (1), pp. 271–307. Cited by: §1, §5.1, §5.1. ↩ a b c
- Convergence Rates for Gaussian Mixtures of Experts. Journal of Machine Learning Research 23 (323), pp. 1–81. External Links: Link Cited by: §1. ↩ a
- A Bayesian Look at Inverse Linear Regression. Journal of the American Statistical Association 65 (329), pp. 356–369. External Links: Link, Document Cited by: §1. ↩ a
- Local Statistical Modeling via a Cluster-Weighted Approach with Elliptical Distributions. Journal of Classification 29 (3), pp. 363–401. External Links: ISSN 1432-1343, Link, Document Cited by: §2.2. ↩ a
- Gibbs Sampling Methods for Stick-Breaking Priors. Journal of the American Statistical Association 96 (453), pp. 161–173. External Links: ISSN 0162-1459, Link, Document Cited by: Appendix A. ↩ a
- Adaptive mixtures of local experts. Neural computation 3 (1), pp. 79–87. Cited by: §1. ↩ a
- Hierarchical mixtures-of-experts for exponential family regression models: approximation and maximum likelihood estimation. Annals of Statistics, pp. 987–1011. Cited by: §1. ↩ a
- Variational inference and sparsity in high-dimensional deep Gaussian mixture models. Statistics and Computing 32 (5), pp. 70. External Links: ISSN 1573-1375, Link, Document Cited by: §1. ↩ a
- Adaptive Bayesian density estimation with location-scale mixtures. Electronic Journal of Statistics 4 (none), pp. 1225 – 1257. External Links: Link, Document Cited by: §7. ↩ a
- Fast Bayesian Inversion for high dimensional inverse problems. Statistics and Computing 32 (2). Cited by: §1. ↩ a
- Deep mixture of linear inverse regressions applied to head-pose estimation. In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, pp. 4817–4825. Cited by: §1. ↩ a
- Sliced Inverse Regression for Dimension Reduction. Journal of the American Statistical Association 86 (414), pp. 316–327. External Links: ISSN 01621459, Link, Document Cited by: §1, §2.1. ↩ a b
- Deep Dirichlet process mixture models. In Proceedings of the Thirty-Eighth Conference on Uncertainty in Artificial Intelligence, J. Cussens and K. Zhang (Eds.), Proceedings of Machine Learning Research, Vol. 180, pp. 1138–1147. External Links: Link Cited by: §1. ↩ a
- Bayesian nonparametric priors for hidden Markov random fields. Statistics and Computing 30 (4), pp. 1015–1035. External Links: ISSN 1573-1375, Link, Document Cited by: Appendix A, Appendix A, §3.2. ↩ a b c
- Variational Mixtures of Gaussian Processes for Classification. In Proceedings of the Twenty-Sixth International Joint Conference on Artificial Intelligence, IJCAI-17, pp. 4603–4609. External Links: Link, Document Cited by: §3.3. ↩ a
- Estimating Mixture of Dirichlet Process Models. Journal of Computational and Graphical Statistics 7 (2), pp. 223–238. External Links: Link, Document Cited by: §1, §1. ↩ a b
- Finite Mixture Models. John Wiley & Sons. Cited by: §1. ↩ a
- The EM algorithm and extensions. Wiley, New York, USA. Cited by: Appendix B. ↩ a
- Inconsistency of Pitman-Yor Process Mixtures for the Number of Components. Journal of Machine Learning Research 15 (96), pp. 3333–3370. External Links: Link Cited by: §1. ↩ a
- Mixture of Gaussian regressions model with logistic weights, a penalized maximum likelihood approach. Electronic Journal of Statistics 8 (1), pp. 1661–1695. Cited by: §6.1. ↩ a
- A View of the EM Algorithm that Justifies Incremental, Sparse, and other Variants. In Learning in Graphical Models, M. I. Jordan (Ed.), pp. 355–368. External Links: ISBN 978-94-011-5014-9, Link, Document Cited by: Appendix B. ↩ a
- Markov Chain Sampling Methods for Dirichlet Process Mixture Models. Journal of Computational and Graphical Statistics 9 (2), pp. 249–265. External Links: Link, Document Cited by: §1. ↩ a
- Unveiling the Characteristics of the Lunar Surface by Massive Inversion of a Photometric Model. In LPI Contributions, LPI Contributions, Vol. 3040, pp. 1998. Cited by: §1. ↩ a
- Joint learning of Gaussian graphical models in heterogeneous dependencies of high-dimensional transcriptomic data. In The 16th Asian Conference on Machine Learning (Conference Track), External Links: Link Cited by: §1. ↩ a
- Approximation results regarding the multiple-output Gaussian gated mixture of linear experts model. Neurocomputing 366, pp. 208–214. External Links: ISSN 0925-2312, Link, Document Cited by: §1. ↩ a
- Practical and theoretical aspects of mixture-of-experts modeling: An overview. Wiley Interdisciplinary Reviews: Data Mining and Knowledge Discovery 8 (4), pp. e1246. Cited by: §1. ↩ a
- A universal approximation theorem for mixture-of-experts models. Neural computation 28 (12), pp. 2585–2593. Cited by: §1. ↩ a
- Order selection with confidence for finite mixture models. Journal of the Korean Statistical Society. External Links: ISSN 2005-2863, Link, Document Cited by: §1. ↩ a
- Approximations of conditional probability density functions in Lebesgue spaces via mixture of experts models. Journal of Statistical Distributions and Applications 8 (1), pp. 13. External Links: ISSN 2195-5832, Link, Document Cited by: §1. ↩ a
- A General Theory for Softmax Gating Multinomial Logistic Mixture of Experts. In Proceedings of The 41st International Conference on Machine Learning, Cited by: §1. ↩ a
- Demystifying Softmax Gating in Gaussian Mixture of Experts. In Advances in Neural Information Processing Systems, Cited by: §1. ↩ a
- Towards Convergence Rates for Parameter Estimation in Gaussian-gated Mixture of Experts. In Proceedings of The 27th International Conference on Artificial Intelligence and Statistics, S. Dasgupta, S. Mandt, and Y. Li (Eds.), Proceedings of Machine Learning Research, Vol. 238, pp. 2683–2691. External Links: Link Cited by: §1. ↩ a
- Fast Allocation of Gaussian Process Experts. In Proceedings of the 31st International Conference on Machine Learning, E. P. Xing and T. Jebara (Eds.), Proceedings of Machine Learning Research, Vol. 32, Bejing, China, pp. 145–153. External Links: Link Cited by: §3.3. ↩ a
- Approximation of probability density functions via location-scale finite mixtures in Lebesgue spaces. Communications in Statistics - Theory and Methods 52 (14), pp. 5048–5059. External Links: Link, Document Cited by: §1. ↩ a
- Non-asymptotic model selection in block-diagonal mixture of polynomial experts models. Preprint. arXiv:2104.08959. External Links: Link Cited by: §1. ↩ a
- Model selection by penalization in mixture of experts models with a non-asymptotic approach. In JDS 2022 - 53èmes Journées de Statistique de la Société Française de Statistique (SFdS), Lyon, France. Cited by: §1. ↩ a
- A non-asymptotic theory for model selection in high-dimensional mixture of experts via joint rank and variable selection. In AJCAI Australasian Joint Conference on Artificial Intelligence 2023, External Links: Link Cited by: §1. ↩ a
- Non-asymptotic oracle inequalities for the Lasso in high-dimensional mixture of experts. arXiv:2009.10622. External Links: Link Cited by: §1. ↩ a
- A non-asymptotic approach for model selection via penalization in high-dimensional mixture of experts models. Electronic Journal of Statistics 16 (2), pp. 4742 – 4822. External Links: Link, Document Cited by: §1, §2.2, §6.1. ↩ a b c
- Approximation by finite mixtures of continuous density functions that vanish at infinity. Cogent Mathematics & Statistics 7 (1), pp. 1750861. External Links: ISSN null, Link, Document Cited by: §1. ↩ a
- Model Selection and Approximation in High-dimensional Mixtures of Experts Models: from Theory to Practice. PhD Thesis, Normandie Université, (en). External Links: Link Cited by: §1. ↩ a
- Convergence of latent mixing measures in finite and infinite mixture models. The Annals of Statistics 41 (1), pp. 370–400. External Links: Link, Document Cited by: §1, §5.1, §5.1. ↩ a b c
- Adaptive Bayesian estimation of conditional densities. Econometric Theory 33 (4), pp. 980–1012. External Links: Document Cited by: §7. ↩ a
- Approximation of conditional densities by smooth mixtures of regressions. The Annals of Statistics 38 (3), pp. 1733 – 1766. External Links: Link, Document Cited by: §1. ↩ a
- Inverse regression approach to robust nonlinear high-to-low dimensional mapping. Journal of Multivariate Analysis 163 (C), pp. 1–14. Cited by: §2.2. ↩ a
- The two-parameter Poisson-Dirichlet distribution derived from a stable subordinator. The Annals of Probability 25 (2), pp. 855–900. External Links: Document Cited by: Appendix A, Appendix A, §3. ↩ a b c
- Risk bounds for mixture density estimation. ESAIM: PS 9, pp. 220–229. External Links: Link, Document Cited by: §1. ↩ a
- Asymptotic behaviour of the posterior distribution in overfitted mixture models. Journal of the Royal Statistical Society: Series B (Statistical Methodology) 73 (5), pp. 689–710. External Links: Link, Document Cited by: §1. ↩ a
- A Contrastive Divergence for Combining Variational Inference and MCMC. In Proceedings of the 36th International Conference on Machine Learning, K. Chaudhuri and R. Salakhutdinov (Eds.), Proceedings of Machine Learning Research, Vol. 97, pp. 5537–5545. External Links: Link Cited by: §7. ↩ a
- Estimating the dimension of a model. The Annals of Statistics 6 (2), pp. 461–464. Cited by: §1. ↩ a
- A constructive definition of Dirichlet priors. Statistica Sinica 4, pp. 639–650. Cited by: Appendix A. ↩ a
- Adaptive Bayesian multivariate density estimation with Dirichlet mixtures. Biometrika 100 (3), pp. 623–640. Cited by: §1, §7. ↩ a b
- Information criteria for selecting possibly misspecified parametric models. Journal of Econometrics 71 (1), pp. 207–225. External Links: ISSN 0304-4076, Link, Document Cited by: §1. ↩ a
- Robust Bayesian mixture modelling. Neurocomputing 64, pp. 235–252. External Links: ISSN 0925-2312, Link, Document Cited by: §3.3. ↩ a
- Convergence properties of a general algorithm for calculating variational Bayesian estimates for a normal mixture model. Bayesian Analysis 1 (3), pp. 625 – 650. External Links: Link, Document Cited by: §7. ↩ a
- Topics in optimal transportation. Graduate Studies in Mathematics, Vol. 58, American Mathematical Society. External Links: ISBN 0-8218-3312-X Cited by: §5.1. ↩ a
- Optimal transport: old and new. Vol. 338, Springer. Cited by: §5.1. ↩ a
- Online Variational Inference for the Hierarchical Dirichlet Process. In Proceedings of the Fourteenth International Conference on Artificial Intelligence and Statistics, G. Gordon, D. Dunson, and M. Dudík (Eds.), Proceedings of Machine Learning Research, Vol. 15, Fort Lauderdale, FL, USA, pp. 752–760. External Links: Link Cited by: §3. ↩ a
- On the Asymptotic Distribution of the Minimum Empirical Risk. In Forty-first International Conference on Machine Learning, External Links: Link Cited by: §1. ↩ a
- A Two-Layer Mixture Model of Gaussian Process Functional Regressions and Its MCMC EM Algorithm. IEEE Transactions on Neural Networks and Learning Systems 29 (10), pp. 4894–4904. External Links: Document Cited by: §3.3. ↩ a
- An Alternative Model for Mixtures of Experts. In Advances in Neural Information Processing Systems, G. Tesauro, D. Touretzky, and T. Leen (Eds.), Vol. 7. External Links: Link Cited by: §1. ↩ a
- Variational Mixture of Gaussian Process Experts. In Advances in Neural Information Processing Systems, D. Koller, D. Schuurmans, Y. Bengio, and L. Bottou (Eds.), Vol. 21. External Links: Link Cited by: §3.3. ↩ a
- Twenty Years of Mixture of Experts. IEEE Transactions on Neural Networks and Learning Systems 23 (8), pp. 1177–1193. External Links: ISSN 2162-2388 VO - 23, Document Cited by: §1. ↩ a
- Bayesian shrinkage in mixture-of-experts models: identifying robust determinants of class membership. Advances in Data Analysis and Classification 13 (4), pp. 1019–1051. External Links: ISSN 1862-5355, Link, Document Cited by: §1. ↩ a
- Convergence rates of variational posterior distributions. The Annals of Statistics 48 (4), pp. 2180 – 2207. External Links: Link, Document Cited by: §7. ↩ a
Supplementary Materials for
“Bayesian nonparametric mixture of experts for inverse problems”
In this supplementary material, we first recall the standard Bayesian nonparametric priors and variational Bayesian expectation-maximisation principle in Appendices A and B, respectively. Then, all specifications of the VBEM for the BNP-GLLiM model, the evidence lower bound and the technical proofs that are not included in the main paper are placed in Appendices C, D and E. . Finally, Appendix F proposes a more general model with a hyper prior on the gating parameters, referred to as BNP-GLLiM2.
Appendix A Bayesian nonparametric priors
Stick-breaking construction of Dirichlet process.
Note that the Dirichlet process (DP) (Ferguson, 1973) is a central BNP prior and is the infinite-dimensional generalization of the Dirichlet distribution. Therefore, for the sake of completeness, let us first recall the definition of the DP. A DP on the space is defined as a random process characterized by a concentration parameter and a base distribution , denoted by , such that for any finite partition of , the random vector is Dirichlet distributed:
We use the stick-breaking construction of the DP (SBDP), due to Sethuraman (1994):
Pitman–Yor process.
As a generalized version of the Dirichlet process, in the Pitman–Yor process (PYP) (Pitman and Yor, 1997), the ’s are independent () but not identically distributed. Specifically,
Here is a discount parameter and is a concentration parameter . The PYP is a two-parameter generalisation of the DP that allows one to control the tail behaviour when modelling data with either exponential or power-law tails (Ishwaran and James, 2001; Pitman and Yor, 1997). The PYP reduces to a DP when . More general stick-breaking representations are possible, e.g., Gibbs-type priors (De Blasi et al., 2015; Ghosal and Van der Vaart, 2017) or homogeneous normalised random measures with independent increments (Favaro et al., 2016). The PYP has a power-law behaviour for the number of clusters. This can make it more suitable for a number of applications. In other words, the number of clusters grows as for PYP, while growing more slowly as for DP.
Since the hyperparameters and can have a significant effect on the growth of the number of clusters with data sample size, it is possible to specify priors for them. For the DP case obtained with , it is suggested in Blei and Jordan (2006) to use a gamma prior, , where the hyperparameters and can be estimated or fixed. A natural question is whether one can also find a tractable prior for the discount parameter . Following the work of Lü et al. (2020), we use the following prior that satisfies the constraints and ,
where is a shifted gamma distribution and is a distribution depending on some parameter which is not specified at the moment but which can typically be assumed to be a uniform distribution on the interval . Such a shifted gamma distribution is the distribution of a variable , where is considered fixed and follows a gamma distribution . The PDF of this shifted gamma distribution is obtained from the standard gamma distribution as .
Hierarchical representation of BNP-MM.
Appendix B Variational Bayesian expectation-maximisation principle
The clustering task consists mainly of estimating the unknown labels from the observed data , whose joint distribution is determined by a set of parameters denoted by and often by additional hyperparameters .
The expectation-maximisation (EM) algorithm (Dempster et al., 1977; McLachlan and Krishnan, 1997) is a generative technique for maximum likelihood estimation (MLE) in the presence of unobserved latent variables or missing data. An EM iteration consists of two steps usually referred to as the E-step in which the expectation of the so-called complete log-likelihood is computed and the M-step in which this expectation is maximized over . An equivalent way to define EM is the following. As discussed in Neal and Hinton (1998), EM can be viewed as an alternating maximisation procedure of a function defined, for any probability distribution over labels , by
where is the entropy of and is the expectation with respect to . The function depends on the observations , which are fixed throughout and are therefore omitted from the notation.
When prior knowledge on the parameters is available, an alternative approach consists of replacing the MLE by a maximum a posteriori (MAP) estimation of using the prior knowledge encoded in a distribution . More precisely, the MLE of is then replaced by a point estimation . In this paper, instead of considering only point estimation of , we carry out a fully Bayesian approach. That is, we integrate out as follows
This integration requires the computation of the density , which is usually not available in closed-form. As an alternative to costly simulation-based methods (e.g., Markov chain Monte Carlo (MCMC)), an EM-like procedure using variational approximation can provide approximations of the marginal posterior distributions and . This approach is referred to as VBEM for variational Bayesian EM, as introduced by Beal and Ghahramani (2003).
To deal with the BNP-GLLiM model, we need to use the VBEM with hyperparameter optimisation of Beal (2003, Figure 2.5 and Algorithm 5.3). Let and denote the distributions over and , respectively, which will serve as approximations to the true posteriors. Specifically, in the Bayesian setting, the intractable posterior is approximated by the variational posterior .
Similar to standard EM, VBEM maximizes the following evidence lower bound (often abbreviated ELBO, and sometimes called the variational lower bound or negative variational free energy), defined for arbitrary and distributions by
| (20) | ||||
alternatively over and . Here, KL stands for Kullback-Leibler divergence. It is worth noting that adding a prior on is formally equivalent to considering as missing variables, while the hyperparameters play the role of the parameters of interest in MLE.
The alternate maximisation over leads to the VBEM algorithm, which can be decomposed into three steps. It is easy to show, using the KL divergence properties, that the maximisation over and leads to the following E-steps, see, e.g., Chaari et al. (2013, Appendix A), Beal (2003, Theorem 2.1) and Bishop (2006, Section 10.1.1), which is essentially coordinate ascent in the function space of variational distributions. Furthermore, the following update rules for E-steps converge to a local maximum of . At the th iteration, using current values and , we get the following updating,
| (21) | |||
| (22) | |||
In practice, we can decide which parameters to treat as genuine parameters or as hyperparameters , depending on whether some prior knowledge is available for only a subset of the parameters, or whether the model has hyperparameters for which no prior information is available. Furthermore, for complex models, and may need to be further restricted to simpler forms, such as factorised forms, to ensure tractable VB-E steps. This is illustrated in the next Section F.1 for the BNP-GLLiM inference.
Appendix C Details of VBEM for BNP-GLLiM model
C.1 VB-E- step from Section 3.1
To achieve results from Section 3.1, we make use of (5), (6), (40), and are only interested in the functional dependence of the right-hand side of (22) on the variable . Thus, any terms that do not depend on can be absorbed into the additive normalization constant, giving
Here,
Furthermore, we used the fact that
Finally, we have for , let correspond to the weight of cluster , then
C.2 VB-E- step from Section 3.1
In the PY case, to achieve results from Section 3.1, we make use of (5), (4), (40), (41), and are only interested in the functional dependence of the right-hand side of (22) to the variables . Thus, any terms that do not depend on can be included in the additive normalization constant. This results in being proportional to
where we used the fact that . Except in the DP-GLLiM case, i.e., , the normalizing constant, , for is not tractable. However, to perform VBEM in Section F.1, we do not need the full distribution, but only the means and . One solution, therefore, is to use importance sampling or MCMC to compute these expectations by means of Monte Carlo sums. Via the prior on given in (4), it holds that
| (23) |
Here, given that is the digamma function defined by , we have
We propose to use the following important distribution where is the uniform distribution on , denoted as . Then we obtain an expression for the importance weights,
The importance sampling scheme then consists of the following steps
-
•
For , first simulate independently from and then simulate conditionally with the -shifted gamma . This later simulation is easily obtained by simulating a standard and then subtracting from the result.
-
•
Compute the importance weights .
-
•
Approximate the means
Note that this complication is due to the PY. In the DP-GLLiM case, by substituting in (23), the step is much simpler, as it reduces to computing the approximate posterior expectation of , namely,
C.3 VB-E- step from Section 3.1
In some situations, it is useful to use a -of- binary vector to represent the latent variable for each observation . To be more precise, we introduce a -dimensional binary random variable , with a -of- representation in which a particular element is equal to , i.e., , and all other elements are equal to . The values of thus satisfy and . If there is no confusion, we also denote as the latent matrix . It is worth to mentioning that when using a -of- representation of , we can also write down marginal the conditional distributions of and , corresponding to (7), in the form
| (24) | ||||
| (25) | ||||
| (26) | ||||
The observations and are therefore i.i.d. and generated from the same GMM (25) and infinite GLLiM (26), respectively. Similarly, (6) can be written down in the form
By using the decomposition (8), the representation (24), (26) and absorbing any terms that are independent on into the additive normalization constant, we obtain
Here, we used the fact that
By taking exponential of both sides and taking into account the normalized constant, it holds that
Note also that if and only if the latent matrix reduces to a sparse matrix which has only one position different from , namely . This leads to the following simplified notation:
C.3.1 Updating
By using matrix derivatives, the derivative of the log likelihood with respect to is given by
Finally, setting to zero yields
C.4 VB-M- step from Section 3.2
This step divides into sub-steps that involve the following optimisations
By definition, we have
We aim to solve the following optimisation
Similarly with Sections C.5.1 and C.3.1, we obtain the following update:
C.5 VB-M- step from Section 3.2
By definition, we have
We aim to solve the following optimisation
C.5.1 Updating
The derivative of the log likelihood with respect to is given by
Setting this derivative to zero, we obtain the solution for VB-M- step given by
| (27) |
C.5.2 Updating
The derivative of the log likelihood with respect to is given by
Then, we set this derivative w.r.t. equal to zero, giving
Here, the last equality is obtained by firstly define the following quantities,
Then, we used the fact that
and
Here, we also use the equalities
and for each , it holds that
Appendix D Details on the ELBO
In this section, we provide the closed-form expressions for the ELBO stated in Proposition 3.1. Let us recal that when , the ELBO in the BNP-GLLiM is derived as follows:
The terms of the right-hand side of the above equation have the following closed-form expressions:
| (28) | ||||
| (29) | ||||
| (30) | ||||
| (31) | ||||
| (32) | ||||
| (33) | ||||
Appendix E Technical proofs
E.1 Proof of Lemma 5.2
We first want to prove (5.2). Using the partition of a joint Gaussian with , we obtain
| (34) |
Recall that
| (35) |
By identifying the parameters of (E.1) and (E.1), it holds that
E.2 Proof of Proposition 3.1
Using the sum and product rules for both discrete and continuous variables, the ELBO in BNP-GLLiM (20) is given by
| (36) |
Next, we evaluate the various terms in the ELBO (36).
Proof of (28)
Via the mean field approximation and the truncation, we have the following computations:
where
Proof of (29)
Similarly to the previous proof, we obtain
where
Proof of (30)
Via calculation, it follows the expressions of the following quantities,
| (37) |
Via (37), it holds that
Proof of (31)
Given a chosen truncated value , it holds that
Here, we have
where we have defined
Next, for the sake of simplicity, for , we use a uniform prior so that parameter does not have to be taken into account. Then it holds that
When , the normalizing constant for is not tractable. Nevertheless, to compute the ELBO, we do not need the full distribution but only the means , , , and . One solution is therefore to use importance sampling or MCMC to compute these expectations via Monte Carlo sums.
When , using integration by parts, it holds that and hence . Furthermore, the posterior is again a gamma distribution with and . Therefore, we have the following tractable formulas:
Proof of (32)
Due to the mean-field approximation (9) and truncation, this step is analytically computed as follows:
Proof of (33)
We have
Note that these terms involving expectations of the logs of the distributions simply represent the negative entropies of those distributions.
Since is not tractable, when , we cannot calculate analytically . Furthermore, it is also difficult to approximate it using MCMC or importance sampling.
When , the posterior is again a gamma distribution with
Since we had , its differential entropy is given by
E.3 Proof of Theorem 4.1
Recall that . Then,
| (38) |
Note that in (E.3), is in fact the (unknown) true posterior distribution of the parameters given a sample . Because the integrations w.r.t. true posterior distribution are intractable, we approximate the predictive conditional density by replacing the true posterior distribution with its truncated variational posterior of parameters given by
Recall that the infinite state space for each is dealt with by choosing a truncation of the state space to a maximum label (Blei and Jordan, 2006). In practice, this consists of assuming that the variational distributions for , satisfy for and that the variational distribution on also factorizes as with the additional condition that . In particular, here we choose , where is estimated from some suitable procedures.
For simplicity, we consider the case when . Then we have
Here, by defining , we used the fact that
Indeed, we made use of the following result for the joint Gaussian, see, e.g., Bishop (2006, Eq. (2.115), page 93). Given a marginal Gaussian distribution for and a conditional Gaussian distribution for given in the form
then the joint distribution of is given by
| (39) |
In our situation, the desired result is obtained via using , .
Furthermore, we also used the fact that
Next, we aim to prove that
Indeed, recall that we have defined
and to deal with the infinite state space for each , we considered a truncation of the state space to a maximum label (Blei and Jordan, 2006). In practice, this consists of assuming that the variational distributions for , satisfy for and that the variational distribution on also factorizes as with the additional condition that . Based on the proof from Ghosal and Van der Vaart (2017, Lemma 3.4), it holds that a necessary and sufficient condition to guarantee that these ’s sum to almost surely, i.e.,
is that the expectation tends to as tends to . In particular, if are i.i.d., e.g., when , it suffices that . Then
E.4 Proof of Theorem 4.2
From the product rule of probability, we see that this conditional distribution can be evaluated from the joint and marginal distributions. Furthermore, by integrating out and , the predictive conditional density is then given by
Next, with a similar step as in the proof of Theorem 4.1, we also obtain
Therefore, we obtain
which is a mixture of Gaussian experts since we have
belongs to a dimensional probability simplex.
E.5 Proof of Theorem 4.3
To deal with high-dimensional regression data, namely high-to-low regression, given the inverse conditional density , we want to compute the following forward conditional density
Then, we have to compute or approximate . Using Theorem 4.1, we obtain
Indeed, we made use of the following results for marginal and conditional Gaussians, see, e.g., Bishop (2006, Eq. (2.115), page 93). Given a marginal Gaussian distribution for and a conditional Gaussian distribution for given in the form
then the marginal distribution of and the conditional distribution of given are given by
In our situation, the desired result is obtained via using , .
Finally, we obtain
where
Here, we used the fact that , namely,
with
When required, it is straightforward to approximate the expectation and covariance matrix of as follows:
where we used the following definitions
Appendix F BNP-GLLiM2: a model with an hyperprior on the gating parameters
F.1 VBEM for BNP-GLLiM2
A more general BNP-GLLiM model, referred to as BNP-GLLiM2 can be considered by specifying a prior on the gating parameters as a normal-inverse-Wishart (NIW) distribution parameterized by with a PDF
The assumptions on the other parameters are not changed, so that hyperparameters and parameters are now as follows:
BNP-GLLiM2 can be represented graphically as in Figure 5.
The joint distribution of the observed data and all latent variables can be expressed hierarchically as
| (40) |
Following the same idea as in Section 3, we only consider the truncated variational posterior of parameters as follows
| (41) |
These forms of and lead to our four VB-E steps and three VB-M steps, summarized below with details in Appendix C. Set the initial value of to . Then, repeat iteratively the following steps. The iteration index is omitted in the update formulas for simplicity.
VB-E steps
Note that the VB-E-, VB-E- steps are the same as in Section 3. We only highlight the modified steps as follows.
We first consider the derivation of the update equation for the factor .
F.1.1 VB-E- step
By using the mean-field approximation (9) and the truncation, see Section C.3 for more details, for all and for , this step consists of computing
| (42) |
Here, given represents the neighbors of , we define by
| (43) | ||||
Note that in the above formula, symbols and are the hyperparameters Specifically defined in the following Sections F.1.2 and 3.2.
Proof of (43).
With respect to the VBEM for BNP-GLLiM2 model from Section F.1.1, it is almost similar to the previous step in Section C.3, except that we have to take into account the randomness of and . Namely, we have
| (44) |
F.1.2 VB-E- step
This step is divided into parts where the computation is similar to that in standard Bayesian GMM with a choice of conjugate prior. Hence, for each , the variational posterior is a Normal-inverse-Wishart density defined as
| (45) |
Here, the hyperparameters are updated as follows (see, e.g., Bishop (2006, Section 10.2.1)):
| (46) | ||||
VB-M steps
The maximisation step consists of updating the hyperparameters , where , by maximizing the free energy, if they are not set heuristically:
| (47) |
The VB-M-step can therefore be divided into four independent sub-steps as listed below. From the conditional independence of and given , the solutions for the VB-M- (in the DP case) and VB-M- steps are straightforward. Only the M- step (in the PYP case) and are more involved.
Note that the VB-M-, VB-M- steps are the same as in Section 3. We only highlight the modified step below.
F.1.3 VB-M- step
F.1.4 ELBO for BNP-GLLiM2
Proposition F.1.
When , the ELBO in BNP-GLLiM2 is determined analytically as follows:
| (48) |
Here, we have the following update formulas:
| (49) | |||
| (50) | |||
| (51) | |||
| (52) |
| (53) | |||
| (54) | |||
F.2 Predictive conditional density for BNP-GLLiM2
F.2.1 Joint density
We first show how to compute the joint density via Theorem F.2, which is proved in Section F.4.
Theorem F.2.
We approximate the joint density of BNP-GLLiM2 by a mixture of product between Gaussian and Student’s t-distributions as follows:
| (55) |
Here, the positive semidefinite shape matrices of Student’s t-distributions are given by
| (56) |
F.2.2 Inverse conditional density
We then show how to approximate the inverse conditional density . This predictive density in BNP-GLLiM2 is approximated by a MoE via Theorem F.3 with the proof in Section F.5.
Theorem F.3.
We approximate the inverse conditional density of BNP-GLLiM2 by a MoE as follows:
Here, the gating posteriors are defined as
Furthermore, for any , it holds that
The prediction task is carried out via the following approximation
F.2.3 Forward conditional density
To deal with high-dimensional regression data, namely high-to-low regression, given the inverse conditional density , we want to approximate the following forward conditional density via Theorem F.4, whose proof is provided in Section F.6.
Theorem F.4.
It holds that
which is a mixture of Gaussian experts, where, for all ,
Here, , are chosen via discretizing -space, , into a grid, e.g., uniform. Note that for simplicity, we evaluate the integrand as a Riemann integral with a truncated value and a number of point for approximating the integration but we can use any scheme to approximate such -dimensional integration.
F.3 Proof of Proposition F.1
Using the sum and product rules for both discrete and continuous variables, the ELBO in BNP-GLLiM (20) is given by
| (57) | ||||
Proof of (50)
| (58) |
Lemma F.5.
We can compute the expectations w.r.t. the variational distributions of the parameters as follows:
| (59) |
Furthermore, for each , it holds that
| (60) |
Proof of (52)
Given a chosen truncated value , it holds that
Note that and are calculated in the same way as in Proposition 3.1.
Finally, we have to compute the remaining term
where
Proof of (54)
We have
Note that these terms involving expectations of the logs of the distributions simply represent the negative entropies of those distributions. In particular, the first two terms are calculated in the same way as in Proposition 3.1.
Similarly, we obtain
F.4 Proof of Theorem F.2
Recall that we defined , . Then,
| (61) |
Note that in (F.4), is in fact the (unknown) true posterior distribution of the parameters given a sample . Because the integrations w.r.t. true posterior distribution are intractable, we approximate the predictive conditional density by replacing the true posterior distribution with its truncated variational posterior of parameters given by
Recall that the infinite state space for each is dealt with by choosing a truncation of the state space to a maximum label (Blei and Jordan, 2006). In practice, this consists of assuming that the variational distributions for , satisfy for and that the variational distribution on also factorizes as with the additional condition that . In particular, here we choose , where is estimated from some suitable procedures.
Lemma F.6.
For each , it holds that
Proof of Lemma F.6
By definition, we obtain
| (62) |
When the size of the data set is large, i.e., , this predictive distribution (F.4) becomes a mixture of Gaussians with component means and precisions . In particular, we made use of the following results for marginal and conditional Gaussians, see, e.g., Bishop (2006, Eq. (2.115), page 93). Given a marginal Gaussian distribution for and a conditional Gaussian distribution for given in the form
then the marginal distribution of and the conditional distribution of given are given by
where
In our situation, via using , we obtain
Notice that the Wishart distribution is a conjugate prior for the Gaussian distribution with known mean and unknown precision. Therefore, it holds that the product of
is again a Wishart distribution without normalized. This can be verified by focusing on the dependency on . More precisely, by using the trace trick of quadratic form, , we obtain
Here, , and
Via the normalization constant we have
Here,
and is the squared Mahalanobis distance defined by
Then, the last equality holds since we have
F.5 Proof of Theorem F.3
From the product rule of probability, we see that this conditional distribution can be evaluated from the joint and marginal distributions. Furthermore, by integrating out and , the predictive conditional density is then given by
Next, with a similar step as in the proof of Theorem F.2, we also obtain
Therefore
which is a mixture of Gaussian experts since we have
belongs to a dimensional probability simplex.
F.6 Proof of Theorem F.4
To deal with high-dimensional regression data, namely high-to-low regression, given the inverse conditional density , we want to compute the following forward conditional density
Then, we have to compute or numerically approximate . Using Theorem F.2 and definition of Student’s t-distribution, we obtain
Then, by definition of Student’s t-distribution, it holds that
Furthermore, we used the fact that
Indeed, we made use of the following results for marginal and conditional Gaussians, see, e.g., Bishop (2006, Eq. (2.115), page 93). Given a marginal Gaussian distribution for and a conditional Gaussian distribution for given in the form
then the marginal distribution of and the conditional distribution of given are given by
In our situation, the desired result is obtained via using , .
Therefore, we obtain
where, for all ,
Here, we used the fact that , namely,
where
The last approximation is deduced by using the fact that one simplistic strategy for evaluating integration would be to discretize -space (-dimensional) into a uniform grid and to evaluate the integrand as a Riemann integral with a truncated value and a number of point for approximating the integration.