# Bayesian outcome-guided multi-view mixture models with applications in
molecular precision medicine

Paul D. W. Kirk
and
Filippo Pagani
and
Sylvia Richardson
MRC Biostatistics Unit, Cambridge, UK

###### Abstract

Clustering is commonly performed as an initial analysis step for uncovering
structure in ’omics datasets, e.g. to discover molecular subtypes of disease.
The high-throughput, high-dimensional nature of these datasets means that they
provide information on a diverse array of different biomolecular processes and
pathways. Different groups of variables (e.g. genes or proteins) will be
implicated in different biomolecular processes, and hence undertaking analyses
that are limited to identifying just a single clustering partition of the
whole dataset is therefore liable to conflate the multiple clustering
structures that may arise from these distinct processes. To address this, we
propose a multi-view Bayesian mixture model that identifies groups of
variables (“views”), each of which defines a distinct clustering structure. We
consider applications in stratified medicine, for which our principal goal is
to identify clusters of patients that define distinct, clinically actionable
disease subtypes. We adopt the semi-supervised, outcome-guided mixture
modelling approach of Bayesian profile regression that makes use of a response
variable in order to guide inference toward the clusterings that are most
relevant in a stratified medicine context. We present the model, together with
illustrative simulation examples, and examples from pan-cancer proteomics. We
demonstrate how the approach can be used to perform integrative clustering,
and consider an example in which different ’omics datasets are integrated in
the context of breast cancer subtyping.

##  1 Introduction

Clustering is ubiquoitously used in the analysis of omics data as a means to
uncover structure and patterns in these large, high-dimensional datasets (e.g.
Eisen et al., 1998; Heyer et al., 1999; Alon et al., 1999; Ben-Dor et al.,
1999; Son et al., 2005). Here we are particularly interested in molecular
precision medicine applications, in which the analysis objective is to
identify molecular subtypes of disease (e.g. Golub et al., 1999; Perou et al.,
2000; Sørlie et al., 2001; Cancer Genome Atlas Network, 2012; Kuijjer et al.,
2018). In this context, the aim is to identify clusters of patients on the
basis of a diverse range of molecular variables, measurements of which are
obtained using high-throughput ’omics technologies. One feature of these
datasets is that they are typically high-dimensional, which frequently
necessitates the use of variable screening or selection strategies (e.g.
Witten and Tibshirani, 2010; Fop and Murphy, 2018; Crook et al., 2018), or
other dimension reduction techniques (e.g. Yeung and Ruzzo, 2001; McLachlan et
al., 2002; Taschler et al., 2019). However, a commonly overlooked challenge is
that ’omics datasets often define multiple clustering structures in the
patient population, as a consequence of different subsets of variables (e.g.
those corresponding to functional groups of genes or proteins) being
implicated in a variety of different biomolecular processes. Thus, depending
on the subset of variables we consider, we can identify different patient
clusters.

A number of papers have proposed methods for identifying multiple clustering
structures (e.g. Cui et al., 2007; Niu et al., 2010; Guan et al., 2010; Li and
Shafto, 2011; Niu et al., 2014). In the literature, a set of variables that
define the same clustering structure has been termed a view, while the task of
identifying views and their associated clustering structures has been referred
to as either multi-view clustering (Cui et al., 2007) or cross-clustering (Li
and Shafto, 2011). One potential challenge faced by these approaches is how to
decide which of the identified clustering structures is the most useful or
relevant for a given task. This is particularly important in stratified
medicine and disease subtyping applications, where we seek clusters (strata)
that define groups of patients who have, for example, similar prognoses or
respond similarly to treatment (e.g. Perou et al., 2000; Sørlie et al., 2001).
In practice, to determine if a given clustering structure is “relevant”, it is
common to make use of a left out outcome variable $y$ (e.g. survival data) and
to assess whether or not different clusters are associated with different
distributions of $y$ (e.g. Curtis et al., 2012). An alternative approach is to
adopt a semi-supervised (outcome guided) approach that makes use of $y$ when
performing the clustering analysis.

Bayesian profile regression is one such semi-supervised mixture modelling
approach that makes use of an outcome/response in order to guide inference
toward relevant clustering structures (Molitor et al., 2010). Informally, a
clustering is said to be relevant (for a given response) if individuals
allocated to the same cluster tend to have similar values for the response;
or, more generally, if the responses of individuals in the same cluster can be
accurately described by a common model. More precisely, a clustering is
defined to be relevant (for a given response) if the value taken by an
individual’s response is not independent of their cluster allocation. It is
clear from this definition that the relevance of a clustering can only be
specified relative to a given response – and, in particular, that for
different responses, different clusterings might be relevant.

For example, we could use Bayesian profile regression to retrospectively
cluster patients on the basis of genetic or genomics data, using their
survival times as a response to guide the clustering toward prognostically
relevant disease subtypes. If we were to use a different response (e.g.
height), we might end up with a completely different clustering structure.
Neither one of these clustering structures would necessarily be “wrong” – they
are just relevant with respect to different responses. Crucially, if we were
to adopt an unsupervised clustering approach (which is commonly the default
analysis choice), there is no guarantee that this would identify a clustering
structure that was relevant for either response.

As we demonstrate in Section 5.1, a limitation of the Bayesian profile
regression model is that, with increasing data dimension, the influence
exerted by the (typically low dimensional) response on the inference of the
clustering structure grows weaker. Here we propose a semi-supervised multi-
view Bayesian clustering model that simultaneously addresses both this
limitation, as well as the challenge faced by existing multi-view approaches
of picking out the (most) relevant clustering structure.

##  2 Profile regression

We suppose that we have data comprising observations on a vector of clustering
variables, ${\bf x}$, and responses, ${\bf y}$. We denote the concatenated
vector of clustering variables and response by ${\bf v}=[{\bf x};{\bf y}]$,
and model the data using a mixture model with $K$ components (where $K$ could
be finite or infinite), as follows:

$\displaystyle p({\bf v}|\boldsymbol{\rho},\boldsymbol{\pi})$
$\displaystyle=\sum_{k=1}^{K}\pi_{k}f_{\bf v}({\bf v}|\boldsymbol{\rho}_{k})$
(1) $\displaystyle=\sum_{k=1}^{K}\pi_{k}f_{\bf y}({\bf
y}|\boldsymbol{\theta}_{k},{\bf x})f_{\bf x}({\bf x}|\boldsymbol{\phi}_{k}),$
(2)

where $\pi_{k}$ is the mixture weight associated with the $k$-th component,
$\boldsymbol{\rho}_{k}$ denotes the parameters associated with the $k$-th
component, and we write $\boldsymbol{\rho}$ and $\boldsymbol{\pi}$ to denote
$\\{\boldsymbol{\rho}_{k}\\}_{k=1}^{K}$ and $\\{\pi_{k}\\}_{k=1}^{K}$
respectively. We assume that the joint density $f_{\bf v}$ can be factorised
into $f_{\bf x}$ and $f_{{\bf y}}$ as shown in Equation (2), with
$\boldsymbol{\phi}_{k}$ denoting the parameters of the model for ${\bf x}$,
and $\boldsymbol{\theta}_{k}$ denoting the parameters of the model for ${\bf
y}$. We will write $\boldsymbol{\Phi}$ and $\boldsymbol{\Theta}$ to denote
$\\{\boldsymbol{\phi}_{k}\\}_{k=1}^{K}$ and
$\\{\boldsymbol{\theta}_{k}\\}_{k=1}^{K}$ respectively. In the profile
regression model, it is further assumed that $f_{\bf y}({\bf
y}|\boldsymbol{\theta}_{k},{\bf x})=f_{\bf y}({\bf
y}|\boldsymbol{\theta}_{k})$; i.e. that ${\bf y}$ is conditionally independent
of ${\bf x}$ given $\boldsymbol{\theta}_{k}$ (Molitor et al., 2010). A related
model, in which this conditional independence assumption is not made, is given
in Shahbaba and Neal (2009).

In Molitor et al. (2010), the authors allow the model for ${\bf y}$ to include
a dependence upon additional adjustment covariates, ${\bf w}$, that may be
predictive of the response but that we do not wish to contribute to the
clustering (e.g. confounders that we wish to control for, such as age or sex),
together with associated “global” (i.e. not component-specific) parameters
$\boldsymbol{\beta}$ (see Appendix for details). The general profile
regression model is then:

$\displaystyle p({\bf x},{\bf
y}|\boldsymbol{\Phi},\boldsymbol{\Theta},\boldsymbol{\pi},\boldsymbol{\beta},{\bf
w})=\sum_{k=1}^{K}\pi_{k}f_{\bf y}({\bf y}|\boldsymbol{\theta}_{k},{\bf
w},\boldsymbol{\beta})f_{\bf x}({\bf x}|\boldsymbol{\phi}_{k}).$ (3)

The original formulation of the profile regression model, which we also adopt
here, is specifically in terms of infinite mixture models using Dirichlet
process priors (Molitor et al., 2010); however, we note that the model is
equally applicable in the case of finite $K$.

###  2.1 Dirichlet process formulation

In Molitor et al. (2010), a finite approximation to the Bayesian nonparametric
case was considered, using a truncated stick breaking construction of the
Dirichlet process (Ishwaran and James, 2001) to define the prior on the
$\pi_{k}$’s in Equation (3). Subsequent Bayesian profile regression papers
(Hastie et al., 2014) and implementations (Liverani et al., 2015) also
considered stick breaking constructions, with inference performed via slice
sampling (Walker, 2007; Kalli et al., 2011). Following the derivations of Neal
(2000) and Rasmussen (2000) in the unsupervised case, here we instead consider
the Dirichlet process mixture model as a limiting case of a (finite) $K$
component mixture model, in which a symmetric Dirichlet prior with parameter
$\alpha/K$ is placed on the mixture weights, $\pi_{k}$ (see also Ishwaran and
Zarepour, 2002). This is closely related to the Pólya urn (Blackwell and
MacQueen, 1973) and Chinese restaurant process (Aldous et al., 1985)
constructions for the Dirichlet process, and permits inference via a collapsed
Gibbs sampler in which the mixture weights are marginalised. Details are
provided in the Appendix, with sampling performed as in Neal (2000) – although
we note that alternative approaches for performing inference in Bayesian
mixture models could also be employed in this context (e.g. Richardson and
Green, 1997; Green and Richardson, 2001; Jain and Neal, 2004, 2007; Walker,
2007; Kalli et al., 2011; Miller and Harrison, 2017)

To provide a very brief overview, let $D=\\{({\bf x}_{i},{\bf y}_{i},{\bf
w}_{i})\\}_{i=1}^{n}$ denote a dataset comprising $({\bf x},{\bf y},{\bf w})$
triples for $n$ individuals, such that $({\bf x}_{i},{\bf y}_{i},{\bf w}_{i})$
corresponds to the $i$-th individual. As is common for mixture models, we
introduce latent component allocation variables, $z_{i}$, where $z_{i}=k$ if
the $i$-th individual is associated with the $k$-th component, and
$p(z_{i}=k|\boldsymbol{\pi})=\pi_{k}$. We define ${\bf
z}=\\{z_{1},\ldots,z_{n}\\}$ to be the multiset of all $n$ component
allocations. The component-conditional likelihood associated with the $i$-th
individual is then:

$\displaystyle p({\bf
x}_{i},\mathbf{y}_{i}|z_{i}=k,\boldsymbol{\phi},\boldsymbol{\theta},\boldsymbol{\beta},{\bf
w}_{i})=f_{\bf y}({\bf y}_{i}|\boldsymbol{\theta}_{k},{\bf
w}_{i},\boldsymbol{\beta})f_{\bf x}({\bf x}_{i}|\boldsymbol{\phi}_{k}).$ (4)

Independent priors, $p(\boldsymbol{\theta}),p(\boldsymbol{\phi})$ and
$p(\boldsymbol{\beta})$, are taken for the component-specific and global
parameters. The collapsed Gibbs sampler then iterates between:

  * •

Step (i): updating the component allocations, $z_{i}$, given the data $D$, the
most recently sampled parameters
$\boldsymbol{\Theta},\boldsymbol{\Phi},\boldsymbol{\beta}$, and the most
recently sampled values for the other allocation variables, ${\bf z}_{-i}={\bf
z}\backslash\\{z_{i}\\}$; and

  * •

Step (ii): updating the parameters
$\boldsymbol{\Theta},\boldsymbol{\Phi},\boldsymbol{\beta}$, given the data $D$
and the most recently sampled component allocations ${\bf z}$.

In the finite $K$ case, if we take a symmetric Dirichlet prior with parameter
$\alpha/K$ for the mixture weights,
$\pi_{1},\ldots,\pi_{K}\sim\mbox{Dir}(\alpha/K)$, then the conditional
posterior probability of allocating individual $i$ to the $k$-th component –
required for Step (i) above – is given by:

$\displaystyle p(z_{i}=k|D,$ $\displaystyle{\bf
z}_{-i},\boldsymbol{\Phi},\boldsymbol{\Theta},\boldsymbol{\beta},{\bf
w}_{i},\alpha)=b\frac{n_{-i,k}}{n-1+\alpha}f_{\bf y}({\bf
y}_{i}|\boldsymbol{\theta}_{k},{\bf w}_{i},\boldsymbol{\beta})f_{\bf x}({\bf
x}_{i}|\boldsymbol{\phi}_{k}),$ (5)

where $b$ is a normalising constant that ensures that
$\sum_{k=1}^{K}p(z_{i}=k|D,{\bf z}_{-i},\ldots)=1$, and $n_{-i,k}$ is the
number of individuals currently allocated to component $k$, excluding the
$i$-th individual. That is, if we let $\mathbbm{1}_{a}$ denote the indicator
function (which is equal to 1 if $a$ is true and zero otherwise), then
$n_{-i,k}=\sum_{z_{j}\in{\bf z}_{-i}}\mathbbm{1}_{z_{j}=k}$.

The Dirichlet process (DP) mixture model may be derived by considering the
limit $K\rightarrow\infty$, in which case the conditional posterior
probability of allocating individual $i$ to an existing component (i.e. one to
which other individuals are currently allocated) is given by:

$\displaystyle p(z_{i}=k,\mbox{ where }k\in{\bf z}_{-i}|D,$ $\displaystyle{\bf
z}_{-i},\boldsymbol{\Phi},\boldsymbol{\Theta},\boldsymbol{\beta},{\bf
w}_{i},\alpha)=b\frac{n_{-i,k}}{n-1+\alpha}f_{\bf y}({\bf
y}_{i}|\boldsymbol{\theta}_{k},{\bf w}_{i},\boldsymbol{\beta})f_{\bf x}({\bf
x}_{i}|\boldsymbol{\phi}_{k}),$ (6)

where $n_{-i,k}$ is as before and $b$ is again a normalising constant (see
Equation (8) below). For the DP, we also require the conditional posterior
probability of allocating individual $i$ to a new component, which is:

$\displaystyle p(z_{i}\not\in{\bf z}_{-i}|D,$ $\displaystyle{\bf
z}_{-i},\boldsymbol{\Phi},\boldsymbol{\Theta},\boldsymbol{\beta},{\bf
w}_{i},\alpha)=b\frac{\alpha}{n-1+\alpha}\int_{\boldsymbol{\phi},\boldsymbol{\theta}}f_{\bf
y}({\bf y}_{i}|\boldsymbol{\theta},{\bf w}_{i},\boldsymbol{\beta})f_{\bf
x}({\bf
x}_{i}|\boldsymbol{\phi})p(\boldsymbol{\phi})p(\boldsymbol{\theta})d\boldsymbol{\phi}d\boldsymbol{\theta}.$
(7)

The normalising constant, $b$, in Equations (6) and (7) is chosen to ensure
that:

$p(z_{i}\not\in{\bf z}_{-i}|D,{\bf z}_{-i},\ldots)+\sum_{k\in{\bf
z}_{-i}}p(z_{i}=k,\mbox{ where }k\in{\bf z}_{-i}|D,{\bf z}_{-i},\ldots)=1.$
(8)

Given values for the component allocation variables, ${\bf z}$, it is
relatively straightforward to perform Step (ii); i.e. to sample the parameters
$\boldsymbol{\Theta},\boldsymbol{\Phi},\boldsymbol{\beta}$. See Appendix for
full details, where we also describe how to sample the DP hyperparameter,
$\alpha$, according to the method described in Escobar and West (1995).

##  3 Variable selection in profile regression

The multi-view model that we propose in Section 4 may be regarded as an
extended form of variable selection for clustering. We therefore start by
presenting a version of the profile regression model that permits variable
selection in the same manner as in the models proposed by Law et al. (2003,
2004); Tadesse et al. (2005); Papathomas et al. (2012).

We again consider the component-conditional likelihood associated with the
$i$-th individual, making the additional assumption of conditionally
independent clustering variables:

$\displaystyle
f(\mathbf{x}_{i},\mathbf{y}_{i}|z_{i}=k,\boldsymbol{\phi},\boldsymbol{\theta},\boldsymbol{\beta},{\bf
w}_{i})$ $\displaystyle=f_{\mathbf{y}}({\bf
y}_{i}|\boldsymbol{\theta}_{k},\boldsymbol{\beta},\mathbf{w}_{i})f_{\mathbf{x}}(\mathbf{x}_{i}|\boldsymbol{\phi}_{k})$
$\displaystyle=f_{\bf y}({\bf
y}_{i}|\boldsymbol{\theta}_{k},\boldsymbol{\beta},\mathbf{w}_{i})\prod_{j=1}^{P}f_{x}(x_{ij}|\boldsymbol{\phi}_{k}),\vspace{-0.5cm}$
(9)

where $P$ is the dimension of ${\bf x}$, so ${\bf
x}_{i}=[x_{i1},\ldots,x_{iP}]^{\top}$. In the following, we assume throughout
a common parametric form, $f_{x}$, for all clustering variables; however, we
could straightforwardly extend to allow different parametric forms for
different clustering variables in order to allow, for example, modelling of
mixed (continuous and categorical) data.

In order to perform variable selection, we follow the same approach as in Law
et al. (2003, 2004) and Tadesse et al. (2005), and introduce latent binary
indicator variables $\gamma_{j}$, such that $\gamma_{j}=1$ if the $j$-th
clustering variable contributes to the clustering structure, and 0 otherwise.
We model the clustering variables for which $\gamma_{j}=0$ as having density
$q_{x}(x_{j}|\boldsymbol{\phi}_{0j})$, where $\boldsymbol{\phi}_{0j}$ is a
“global” (i.e. not component-specific) parameter and $q_{x}$ need not be of
the same parametric form as $f_{x}$. We denote the collection of all
$\boldsymbol{\phi}_{0j}$ parameters by
$\boldsymbol{\phi}_{0\cdot}=\\{\boldsymbol{\phi}_{0j}\\}_{j=1}^{P}$. The
allocation-conditional likelihood is then:

$\displaystyle
f(\mathbf{x}_{i},y_{i}|z_{i}=k,\boldsymbol{\phi},\boldsymbol{\phi}_{0\cdot},\boldsymbol{\theta},{\bf
w}_{i},\boldsymbol{\beta},\boldsymbol{\gamma})$
$\displaystyle=f_{\mathbf{y}}({\bf
y}_{i}|\boldsymbol{\theta}_{k},\boldsymbol{\beta},\mathbf{w}_{i})f_{\mathbf{x}}(\mathbf{x}_{i}|\boldsymbol{\phi}_{k})$
$\displaystyle=f_{\bf y}({\bf
y}_{i}|\boldsymbol{\theta}_{k},\boldsymbol{\beta},\mathbf{w}_{i})\prod_{j=1}^{P}q_{x}(x_{ij}|\boldsymbol{\phi}_{0j})^{\mathbbm{1}_{\gamma_{j}=0}}f_{x}(x_{ij}|\boldsymbol{\phi}_{k})^{\mathbbm{1}_{\gamma_{j}=1}}.\vspace{-0.5cm}$
(10)

The introduction of the latent indicator variables in Equation (10) results in
the grouping of clustering variables into two disjoint sets. Adopting the
terminology of Cui et al. (2007), we will refer to these sets as views. We
define View 1 to be the set of clustering variables (for which $\gamma_{j}=1$)
that define a clustering structure, and View 0 to be the set (for which
$\gamma_{j}=0$) that do not define a clustering structure. We will refer to
this latter set of clustering variables as the null view. The variable
selection model may therefore be considered to be one which not only groups
individuals together (into clusters), but also groups clustering variables
together (into 2 views). Inference of the view allocation variables,
$\gamma_{j}$, can be performed via Gibbs sampling, and amounts to sampling
according to the posterior probabilities associated with each of the models
(i.e. the “no clustering structure” model, $q_{x}$, and the “clustering
structure” model, $f_{x}$) given the observed data for the $j$-th clustering
variable, $\\{x_{ij}\\}_{i=1}^{n}$. See Appendix for details.

##  4 Multi-view Bayesian profile regression

A natural extension to the variable selection model is to allow there to be
$L$ views (with $L\in\mathbb{Z}^{+}$ now allowed to be more than 2) by
allowing the indicator variables, $\gamma_{j}$, to be categorical with $L$
categories, $0,1,2,\ldots,L-1$. It is then clear that the $\gamma_{j}$
variables are view allocation variables, which serve an analogous role to the
component allocation variables, $z_{i}$, but which act to group together the
clustering variables rather than the individuals.

As in the variable selection model, the introduction of multiple views
necessitates the introduction of view-specific models. In the variable
selection case, the clustering variables in View 1 contribute to a mixture
model (which defines clusters among the individuals), while the clustering
variables in the null view contribute to a single density (corresponding to
all individuals being in a single cluster). In general, different choices are
possible for the view-specific models. Here we consider an $L$-view case in
which View 0 is a null view, while Views $1,\ldots,L-1$ are each associated
with a mixture model, such that each of these views defines a different
clustering structure among the individuals. We associate the response with
View 1 only, and refer to this view as the relevant view. Thus, the View 1
model is a (semi-supervised) Bayesian profile regression model, while the
models for Views $2,\ldots,L-1$ are all unsupervised mixture models.

###  4.1 Allocations-conditional likelihood

For each of the non-null views, we introduce component allocation variables
$z_{i\ell}$, such that $z_{i\ell}$ denotes the component that is responsible
for the $i$-th individual in the $\ell$-th view’s mixture model. We denote the
number of mixture components in the $\ell$-th view’s mixture model by
$K_{\ell}$, where – as previously – $K_{\ell}$ may be finite or infinite. We
moreover denote the parameters associated with the $k$-th component in the
$\ell$-th view by $\boldsymbol{\phi}_{\ell k}$ and define
$\boldsymbol{\phi}_{\ell\cdot}=\\{\boldsymbol{\phi}_{k\ell}\\}_{k=1}^{K_{\ell}}$
to be the full complement of component-specific parameters associated with the
$\ell$-th view. The $\boldsymbol{\theta}$ parameters of the model for ${\bf
y}$ are associated with view $\ell=1$ only.

The likelihood associated with the $i$-th individual, conditioned on these
component allocations, is then:

$\displaystyle
f(\mathbf{x}_{i},y_{i}|\\{z_{i\ell}=k_{\ell},\boldsymbol{\phi}_{\ell\cdot}\\}_{\ell=1}^{L-1},\boldsymbol{\phi}_{0\cdot},\boldsymbol{\theta},$
$\displaystyle{\bf w}_{i},\boldsymbol{\beta},\boldsymbol{\gamma})=$
$\displaystyle f_{\bf y}({\bf
y}_{i}|\boldsymbol{\theta}_{k_{1}},\boldsymbol{\beta},\mathbf{w}_{i})\prod_{j=1}^{P}q_{x}(x_{ij}|\boldsymbol{\phi}_{0j})^{\mathbbm{1}_{\gamma_{j}=0}}\prod_{\ell=1}^{L-1}f_{x}(x_{ij}|\boldsymbol{\phi}_{\ell
k_{\ell}})^{\mathbbm{1}_{\gamma_{j}=\ell}}.\vspace{-0.5cm}$ (11)

Defining $\Gamma_{\ell}=\\{j:\gamma_{j}=\ell\\}$ to be the index set of
clustering variables allocated to the $\ell$-th view, the above may
alternatively be written as follows:

$\displaystyle f(\mathbf{x}_{i},y_{i}|$
$\displaystyle\\{z_{i\ell}=k_{\ell},\boldsymbol{\phi}_{\ell\cdot}\\}_{\ell=1}^{L-1},\boldsymbol{\phi}_{0\cdot},\boldsymbol{\theta},{\bf
w}_{i},\boldsymbol{\beta},\boldsymbol{\gamma})=$
$\displaystyle\definecolor{currentcolor}{rgb}{0,0,0}\color[rgb]{0,0,0}\underbrace{\color[rgb]{0,0,0}\left(\prod_{j\in\Gamma_{0}}q_{x}(x_{ij}|\boldsymbol{\phi}_{0j})\right)}_{\color[rgb]{0,0,0}\mathclap{\text{View
0}}}\color[rgb]{0,0,0}\definecolor{currentcolor}{rgb}{0,0,0}\color[rgb]{0,0,0}\underbrace{\color[rgb]{0,0,0}\left(f_{\bf
y}({\bf
y}_{i}|\boldsymbol{\theta}_{k_{1}},\boldsymbol{\beta},\mathbf{w}_{i})\prod_{j\in\Gamma_{1}}f_{x}(x_{ij}|\boldsymbol{\phi}_{1k_{1}})\right)}_{\color[rgb]{0,0,0}\mathclap{\text{View
1}}}\color[rgb]{0,0,0}\prod_{\ell=2}^{L-1}\definecolor{currentcolor}{rgb}{0,0,0}\color[rgb]{0,0,0}\underbrace{\color[rgb]{0,0,0}\left(\prod_{j\in\Gamma_{\ell}}f_{x}(x_{ij}|\boldsymbol{\phi}_{\ell
k_{\ell}})\right)}_{\color[rgb]{0,0,0}\mathclap{\text{View
$\ell$}}}\color[rgb]{0,0,0},\vspace{-0.5cm}$ (12)

where each of the bracketed terms corresponds to a different view, as shown.
This expression makes clear that (conditioned on the allocation of clustering
variables to views) each of the views is modelled independently, as also
illustrated in Figure 1, with View 0 being modelled as a single cluster, View
1 being modelled with a profile regression model (compare to Equation (9)),
and the remaining views each modelled by their own (unsupervised) mixture
model.

Figure 1: Schematic illustration of semi-supervised multi-view model.
Conditioned on the allocation of clustering variables to views, each view is
modelled independently, with View 0 (the null view) modelled as possessing no
clustering structure, View 1 modelled with a profile regression (semi-
supervised clustering) model that links the clustering structure defined by
the clustering variables to a response, $\bf{y}$, and the remaining views each
modelled by their own unsupervised mixture model. The level of the response,
${\bf y}$, is indicated by a column shown to the left of the data matrix in
the figure above.

###  4.2 Inference

Conditioned on the view allocations, the models describing each of the views
are independent, and hence we can perform inference for the component
allocations and parameters within each view using existing approaches either
for Bayesian profile regression models (Molitor et al., 2010) in the case of
View 1, or for unsupervised Bayesian mixture models (e.g. Richardson and
Green, 1997; Neal, 2000) for Views 2, …, $L-1$. We therefore adopt a Gibbs
sampling approach, in which we iterate between sampling the view allocation
variables and performing inference for the models in each view. We provide
details of the update for the view allocation indicators below. Updates for
the within-view parameters and latent variables are performed as in Neal
(2000); see Appendix for further details.

####  4.2.1 Updating the view allocation indicators

Define $\nu_{\ell}$ to be the prior probability that a clustering variable is
allocated to the $\ell$-th view. Given the component allocations
$z_{\cdot\ell}=\\{z_{i\ell}\\}_{i=1}^{n}$ for the $\ell$-th view (for
$\ell=1,\ldots,L-1$), the posterior probability that the $j$-th clustering
variable is allocated to the $\ell$-th view is then:

$\displaystyle
p(\gamma_{j}=\ell|z_{\cdot\ell},\boldsymbol{\phi}_{\ell\cdot})=\frac{1}{Z}\nu_{\ell}\prod_{i=1}^{n}f_{x}(x_{ij}|\boldsymbol{\phi}_{\ell
z_{i\ell}})\mbox{ for $\ell=1,\ldots,L-1$},$ (13)

where $Z$ is a normalising constant that ensures that the posterior view
allocation probabilities sum to 1, and, as before, $\boldsymbol{\phi}_{\ell
z_{i\ell}}$ denotes the parameter associated with the $z_{i\ell}$-th component
in the $\ell$-th view, and $\boldsymbol{\phi}_{\ell\cdot}$ is the full
complement of component-specific parameters associated with the $\ell$-th
view.

For the null view, we have:

$\displaystyle
p(\gamma_{j}=0|\boldsymbol{\phi}_{0j})=\frac{1}{Z}\nu_{0}\prod_{i=1}^{n}q_{x}(x_{ij}|\boldsymbol{\phi}_{0j}),$
(14)

and it is now clear that

$Z=\nu_{0}\prod_{i=1}^{n}q_{x}(x_{ij}|\boldsymbol{\phi}_{0j})+\sum_{\ell=1}^{L-1}\left(\nu_{\ell}\prod_{i=1}^{n}f_{x}(x_{ij}|\boldsymbol{\phi}_{\ell
z_{i\ell}})\right).$

In the case where conjugate priors are taken for the component-specific
parameters, $\boldsymbol{\phi}_{\ell k}$, and/or the parameters of the null
view model, $\boldsymbol{\phi}_{0.}$, these parameters may be integrated out
and the likelihood functions $f_{x}$ and $q_{x}$ may be replaced with marginal
likelihoods; see Appendix.

###  4.3 Choice of $L$, the number of views

It may be noted that inference for the view allocations is somewhat analogous
to inference for the mixture component allocations within each (non-null)
view. In the latter case, individuals are allocated with a higher probability
to components in which individuals with similar clustering variable profiles
are allocated; whereas in the former case, clustering variables are allocated
with a higher probability to views in which clustering variables defining a
similar clustering structure are allocated.

Similarly, the choice of the number of views in a multi-view model, $L$, is
analogous to the choice of the number of components, $K$, in a conventional
mixture model. Moreover, in much the same way that the prior component
allocation probabilities in a mixture model may be treated as parameters,
$\pi_{k}$, to be inferred (e.g. taking a Dirichlet or Dirichlet process
prior), it is also possible to treat the prior view allocation probabilities
in the multi-view model, $\nu_{\ell}$, as parameters. By adopting a Dirichlet
(or Dirichlet process) prior, and taking $L$ to be large (or infinite), the
number of views may, in principle, be inferred automatically. However, this
comes with an associated computational expense, since each additional view
brings with it a mixture model whose parameters and latent component
allocations must be inferred. Approximate inference procedures for these
models, such as variational Bayes (previously considered in the context of
multi-view clustering by Guan et al., 2010) will be an important direction for
future research, but in the present work we focus upon small values of $L\geq
2$, for which inference via Markov chain Monte Carlo (MCMC) is feasible.

###  4.4 Initialisation

We have found that a good initialisation strategy is to start with all
variables in the relevant view, so that irrelevant and null variables are
“selected out” at subsequent iterations. We adopt this initialisation strategy
in all examples

##  5 Examples

Although until now we have deliberately kept the exposition general, in the
examples that follow we restrict our attention to Dirichlet process mixture
models for the non-null views. In Section 5.1 we present simulation examples
that allow us to illustrate how clustering approaches that do not model
multiple clustering structures can fail, even if they exploit response
information and employ variable selection. In Section 5.2 we consider an
integrative clustering example in the context of breast cancer subtype
characterisation, in which we fit a multi-view model with 3 views (including
one null view).

###  5.1 Simulation study to illustrate the limitations of methods that ignore
multiple views

We construct a simulation example to demonstrate the limitations of existing
semi-supervised clustering approaches that fail to account for multiple
clustering structures. We consider simulated datasets with $n=300$
individuals, $p=10$ categorical clustering variables (each of which has 3
categories), and a univariate binary response, $y$. The clustering variables
define 2 views, with the first $q$ clustering variables ($q\in\\{3,5,7,9\\}$)
defining a relevant clustering structure (which relates to the response) and
the remaining $p-q$ defining an irrelevant clustering structure (unrelated to
the response). The clustering variables in the relevant view (View 1) define 6
equally-sized clusters, which are related to the response according to
$P(y_{i}=1|z_{i1}=k)=\theta_{k}$, where
$\theta_{k}\in\\{0.01,0.15,0.40,0.60,0.85,0.90\\}$. The clustering variables
in the irrelevant view (View 2) also define 6 clusters, but there is no link
between these clusters and the response; i.e.
$P(y_{i}=1|z_{i2}=k)=P(y_{i}=1)=0.485$, regardless of the component allocation
in View 2. An illustration of a dataset for $q=5$ is provided in Figure 2,
with Figure 2(a) showing the clustering structure defined by the clustering
variables in View 1, and Figure 2(b) showing the clustering structure defined
by the clustering variables in View 2 (which is irrelevant for the response
shown).

(a) Rows ordered to emphasise the relevant clustering structure.

(b) Rows ordered to emphasise the irrelevant clustering structure.

Figure 2: The same simulated dataset with 2 different row orderings. The data
within each view are simulated according to Equation (15), with $w=0.8$. (a)
The rows are ordered to show the clusters defined by the clustering variables
in the relevant view (variables 1 to 5). The probability that an individual’s
response is 0/1 depends upon their membership of these clusters. (b) As in
(a), but ordering the rows to highlight the clusters defined by the clustering
variables in the irrelevant view (variables 6 to 10). As described in the main
text, there is no link between these (irrelevant) clusters and the response.

Within each view, the categorical data are simulated according to a mixture
distribution. With probability $w\in[0,1]$ we simulate $x_{ij}$ according to:

$\displaystyle x_{ij}|z_{i}=k,\gamma_{j}=\ell$
$\displaystyle\sim\textrm{Categorical}(\boldsymbol{\phi}_{\ell kj}),$
$\displaystyle\mbox{ where }\boldsymbol{\phi}_{\ell kj}=[{\phi}_{\ell
kj1},{\phi}_{\ell kj2},{\phi}_{\ell kj3}]$
$\displaystyle\sim\textrm{Dirichlet}(0.01,0.01,0.01),$ (15)

and with probability $1-w$ we simulate $x_{ij}$ according to a discrete
uniform distribution on the three categories. The parameter $w$ thereby allows
us to control cluster separability. The clusters are easily separable when $w$
is close to 1, but the clusters are increasingly noisy and hard to separate as
$w$ approaches 0.

####  5.1.1 Results

We applied an existing semi-supervised clustering approach to each dataset, as
implemented in the PReMiuM R package for Bayesian profile regression (Liverani
et al., 2015). For each dataset, we ran PReMiuM both with and without variable
selection. To run with variable selection, we specified the “varSelectType”
option in PReMiuM to be “continuous”, which performs variable selection with
latent selection weights as described in Papathomas et al. (2012). To obtain a
“gold standard” that reflects the performance that could be achieved if only
the relevant clustering variables were selected, we additionally applied
PReMiuM to each dataset after having removed all of the irrelevant clustering
variables. In practical examples, identifying the relevant clustering
variables a priori in this way will not be possible; however, here we suppose
we have access to an oracle that can provide this information. In all cases,
we used PReMiuM to perform 10,000 Gibbs sampling iterations, discarding the
first 1,000 as burn in and thinning to retain every 5-th draw.

To summarise the MCMC output, we calculated the adjusted Rand index (ARI)
between the true clustering structure in the relevant view and each retained
clustering structure (partition) sampled from the posterior. Thus, for each
distinct PReMiuM run, we obtained a distribution of ARI values, which assesses
how well the inferred clustering structure matches the clustering structure in
the relevant view. The ARI can take a value of at most 1 (indicating a perfect
match), while a value of 0 indicates that the match is no better than we would
expect by random chance.

Figure 3 illustrates typical output for 4 datasets, each having a different
value for $q$ (the number of relevant clustering variables). Corresponding
figures for other simulated datasets are provided in Supplementary Results,
and are qualitatively similar. As we might expect, as $q$ diminishes, so too
does our ability to infer the correct clustering structure. Perhaps less
intuitively, we also see that – as a consequence of failing to model the
multiple views – running PReMiuM with variable selection does not necessarily
help to improve inference of the relevant clustering structure, and can
actually be damaging. Figure 3 shows that variable selection can be useful
when the number of relevant clustering variables is large relative to the
number of irrelevant clustering variables (e.g. if 7 out of 10 are relevant).
In these cases, the irrelevant clustering variables are discarded (i.e.
strongly down-weighted), resulting in performance that is comparable with that
achieved when the irrelevant clustering variables are artificially removed
(e.g. consider $q=7$ in Figure 3). However, when the number of relevant
clustering variables is small relative to the number of irrelevant clustering
variables (e.g. if only 3 out of 10 are relevant), variable selection can
diminish performance; as, in these cases, it is the relevant clustering
variables that are discarded. Since PReMiuM models the data as possessing only
one (non-null) clustering structure, it tends to home in on a single
“dominant” clustering structure (here, the one that is defined by the majority
of clustering variables). Variable selection reinforces this by removing or
down-weighting any clustering variables that define a different clustering
structure, resulting in performance that is better if the dominant clustering
also happens to be the relevant one, and worse if not.

Figure 3: Posterior distributions of Adjusted Rand index (ARI) values obtained
when applying PReMiuM to datasets comprising a mix of relevant and irrelevant
clustering variables. PReMiuM is applied either with (blue) or without (gold)
variable selection to datasets comprising both relevant and irrelevant
variables, or to the same datasets with the irrelevant variables removed
(grey; “gold standard” performance). Solid lines: applying PReMiuM without
variable selection to each dataset. Dashed lines: applying PReMiuM with
variable selection to each dataset. Dotted lines: applying PReMiuM to each
dataset after having artificially removed the irrelevant clustering variables.

Crucially, we stress that the results seen here are not specific either to the
inference or variable selection procedures implemented in PReMiuM, but are a
consequence of failing to model multiple (non-null) views in the data, when
they exist.

###  5.2 Breast cancer subtyping

In recent years, many authors (e.g. Shen et al., 2009; Kirk et al., 2012; Lock
and Dunson, 2013; Savage et al., 2013) have considered the problem of how to
perform integrative clustering, using multiple ’omics datasets in order to
better characterise cancer subtypes at the molecular level (see also
Kristensen et al., 2014, for a reveiw). Although multi-view approaches have
not previously been applied for this purpose, they straightforwardly permit
integrative clustering. In the multi-view model, clustering variables are
allocated to the same view if they define the same clustering structure.
Whether or not these clustering variables come from the same dataset or are of
the same data type is irrelevant; all that matters is the clustering structure
that they define.

Here we apply semi-supervised multiview modelling to identify clusters among
breast cancer tumour samples on the basis of reverse phase protein array
(RPPA) and micro-RNA (miRNA) data from TCGA (The Cancer Genome Atlas, 2012). A
great deal of existing work has considered the use of mRNA expression data for
identifying cancer subtypes, including the PAM50 predictive model for
classifying breast cancer tumour samples on the basis of the expression of 50
genes (Parker et al., 2009). We consider $n=108$ breast cancer tumour samples,
66 of which have been classified on the basis of mRNA expression data as
basal-like and 42 as Luminal A. We use this classification as a binary
response, to guide the clustering on the basis of the RPPA and miRNA data from
TCGA. The RPPA data comprise measurements on 171 proteins, while the miRNA
data comprise measurements on 423 miRNAs.

Before clustering, we process the miRNA and RPPA datasets as in Lock and
Dunson (2013). To robustify against misspecification of models for continuous
data, we take the additional pre-processing step of using tertile
discretisation within each tumour sample, and treat the resulting data as
categorical. After pre-processing, the datasets are concatenated, so that our
final working dataset has $p=594$ clustering variables (corresponding to 171
proteins and 423 miRNAs).

####  5.2.1 Results

We fitted our semi-supervised multi-view model to the concatenated TCGA miRNA
and RPPA data, assuming $L=3$ views (1 relevant, 1 irrelevant, and 1 null
view). We performed 10,000 Gibbs sampling iterations, removing the first 5,000
as burn in, and then thinning to retain every 5-th draw. For each clustering
variable, we calculated Monte Carlo estimates of the probability of being
selected into each view. Bar plots of these probabilities are shown in Figure
4. The majority of clustering variables are selected with high probability
into the null view, with relatively few selected with high probability into
the other two views.

(a)

(b)

(c)

Figure 4: Probability of each clustering variable being selected into the (a)
relevant; (b) irrelevant; and (c) null views. In each plot, the clustering
variables are ordered along the $x$-axis according to increasing probability
of selection into the relevant view.

To summarise the clustering structure within each non-null view, we first
calculated the posterior similarity matrix (PSM) for each view. A PSM is an
$n\times n$ matrix whose $i,j$-entry indicates the proportion of clusterings
sampled from the posterior in which tumour sample $i$ and $j$ had the same
cluster label (Fritsch and Ickstadt, 2009). Since PSMs summarise pairwise co-
clustering probabilities, they provide a summary of the MCMC output that
avoids challenges associated with label-switching. The PSMs for the non-null
views are provided in Figure 5. It is clear from Figure 5a that, as desired,
the clustering structure in the relevant view has a strong association with
the subtype label. As shown in Figure 5b, the irrelevant view also possesses a
strong clustering structure, but one which is not associated with the subtype
label.

(a) Posterior similarity matrix: relevant view

(b) Posterior similarity matrix: irrelevant view

Figure 5: Posterior similarity matrices summarising the allocation of
individuals to clusters in the (a) relevant; and (b) irrelevant views. To aid
visualisation, hierarchical clustering has been applied to both the rows and
columns of the PSMs (as indicated by the dendrograms).

To visualise the clusters in the relevant view, we thresholded the selection
probabilities to retain only those clustering variables selected with
probability at least 0.90 of being selected into the relevant view. This left
only 53 clustering variables, of which 32 were miRNAs and 21 were proteins. A
heatmap representation of the data for these 53 clustering variables is
provided in Figure 6a, from which both the clustering structure and its
association with the response (subtype) is clear. We similarly visualised the
clusters in the irrelevant view by retaining only those clustering variables
selected with probability at least 0.90 of being selected into the irrelevant
view. This left 76 clustering variables, all of which were proteins. The
resulting heatmap representation is provided in Figure 6b. In this case, while
there is an evident clustering structure, it is not associated with tumour
subtype.

(a) Clustering structure: relevant view

(b) Clustering structure: irrelevant view

Figure 6: Visualisation of the clustering structures in (a) relevant; and (b)
irrelevant views. Rows correspond to tumour samples, and the columns in each
plot correspond to clustering variables selected with probability at least
0.90 of being selected into the (a) relevant; and (b) irrelevant views.
Hierarchical clustering has been applied to both the rows and columns of each
plot, as indicated by the dendrograms. The PAM50 subtypes for each tumour
sample are indicated by an additional column shown to the left of each plot.

##  6 Discussion

We have demonstrated that, when there are multiple clustering structures
present in data, existing (single view) clustering approaches can fail to
recover the most relevant clustering structure, even when guided by an
appropriate response (Section 5.1.1). Moreover, traditional variable selection
approaches for clustering do not necessarily improve matters, since they tend
to select variables that define the dominant clustering structure, regardless
of whether or not it is associated with a response of interest. In Section
5.2.1, we have shown that real molecular datasets can and do possess multiple
clustering structures, and that our semi-supervised multi-view model can allow
both relevant and irrelevant structures to be identified.

While multi-view approaches clearly provide advantages relative to existing
(single view) alternatives, computational considerations are a notable
challenge. Dirichlet process mixture models are already computationally
costly, and the multi-view approach proposed here introduces an additional
Dirichlet process mixture model for each additional view. While restricting
the number of views provided adequate results in the examples considered here,
this need not be the case if the “true” number of views in the data is greater
than the selected value for $L$. Both computational approaches (such as
parallelisation, as in Suchard et al., 2010) and fast approximate inference
procedures (such as variational Bayes, as in Guan et al., 2010) will be
important considerations for future work.

## Funding

PDWK and SR were supported by the MRC (MC_UU_00002/13 and MC_UU_00002/10)
respectively. FP and SR were supported by EPSRC project EP/R018561/1.

##  7 Appendix

###  7.1 Unsupervised Bayesian mixture models

We start by considering unsupervised mixture models of the following general
form:

$\displaystyle p({\bf x}|\boldsymbol{\phi},\boldsymbol{\pi})$
$\displaystyle=\sum_{k=1}^{K}\pi_{k}f_{\bf x}({\bf x}|\boldsymbol{\phi}_{k}),$
(16)

where $\boldsymbol{\pi}=[\pi_{1},\ldots,\pi_{K}]$ is the vector of mixture
weights and
$\boldsymbol{\phi}=\\{\boldsymbol{\phi}_{1},\ldots,\boldsymbol{\phi}_{K}\\}$
is the collection of all component-specific parameters. We shall initially
consider the case of finite $K$, and return later to the infinite mixture
model.

As is common for mixture models, we introduce latent component allocation
variables, $z_{i}$, where $z_{i}=k$ if the $i$-th observation ${\bf x}_{i}$ is
associated with the $k$-th component, and $p(z_{i}=k|\pi)=\pi_{k}$. Then,

$\displaystyle p({\bf x}_{i}|z_{i},\boldsymbol{\phi})=f_{\bf x}({\bf
x}_{i}|\boldsymbol{\phi}_{z_{i}}),$ (17)

and hence

$\displaystyle p({\bf x}_{i},z_{i}=k|\boldsymbol{\phi},\pi)$
$\displaystyle=f_{\bf x}({\bf
x}_{i}|\boldsymbol{\phi}_{k})p(z_{i}=k|\boldsymbol{\pi})$ (18)
$\displaystyle=f_{\bf x}({\bf x}_{i}|\boldsymbol{\phi}_{k})\pi_{k}.$ (19)

Integrating out $z_{i}$ by summing over all $K$ possible values, we obtain (as
we would hope):

$\displaystyle p({\bf x}_{i}|\boldsymbol{\phi},\pi)$
$\displaystyle=\sum_{k=1}^{K}\pi_{k}f_{\bf x}({\bf
x}_{i}|\boldsymbol{\phi}_{k}).$ (20)

Making the usual conditional independence assumptions, the full joint model
for ${\bf x}_{i},z_{i},\boldsymbol{\phi},\boldsymbol{\pi}$ is:

$\displaystyle p({\bf x}_{i},z_{i},\boldsymbol{\phi},\pi)$
$\displaystyle=f_{\bf x}({\bf
x}_{i}|\boldsymbol{\phi}_{z_{i}})p(z_{i}|\boldsymbol{\pi})p(\boldsymbol{\pi})p(\boldsymbol{\phi})$
(21) $\displaystyle=f_{\bf x}({\bf
x}_{i}|\boldsymbol{\phi}_{z_{i}})p(z_{i}|\boldsymbol{\pi})p(\boldsymbol{\pi})\prod_{k=1}^{K}p(\boldsymbol{\phi}_{k}),$
(22)

where we assume independent priors for the component-specific parameters,
$\phi_{k}$, and a symmetric Dirichlet prior for the mixture weights,
$\pi_{1},\ldots,\pi_{K}\sim\mbox{Dir}(\alpha/K)$.

For the full dataset, we have:

$\displaystyle p({\bf x}_{1},\ldots,{\bf
x}_{n},z_{1},\ldots,z_{n},\boldsymbol{\phi},\boldsymbol{\pi})$
$\displaystyle=\left(\prod_{i=1}^{n}f_{\bf x}({\bf
x}_{i}|\boldsymbol{\phi}_{z_{i}})p(z_{i}|\boldsymbol{\pi})\right)p(\boldsymbol{\pi})p(\boldsymbol{\phi})$
(23) $\displaystyle=\left(\prod_{i=1}^{n}f_{\bf x}({\bf
x}_{i}|\boldsymbol{\phi}_{z_{i}})p(z_{i}|\boldsymbol{\pi})\right)p(\boldsymbol{\pi})\prod_{k=1}^{K}p(\boldsymbol{\phi}_{k}).$
(24)

####  7.1.1 Inference via Gibbs sampling (finite $K$ case)

Given Equation (24), it is straightforward to write down the conditionals for
Gibbs sampling. For the time being, we assume finite $K$ (from which we will
later derive the infinite limit).

##### Conditional for $\boldsymbol{\phi}_{k}$

By examination of the RHS of Equation 24, we have:

$p(\boldsymbol{\phi}_{k}|{\bf x}_{1},\ldots,{\bf
x}_{n},z_{1},\ldots,z_{n},\boldsymbol{\phi}_{-k},\pi_{1},\ldots,\pi_{K})\propto
p(\boldsymbol{\phi}_{k})\prod_{i:z_{i}=k}f_{\bf x}({\bf
x}_{i}|\boldsymbol{\phi}_{z_{i}}),$ (25)

where $\boldsymbol{\phi}_{-k}$ denotes the set comprising all
$\boldsymbol{\phi}_{j}$ for which $j\neq k$. Thus the conditional for
$\boldsymbol{\phi}_{k}$ is the posterior density for $\boldsymbol{\phi}_{k}$
given all ${\bf x}_{i}$ for which $z_{i}=k$. If conjugate priors are taken,
this posterior is available analytically, otherwise samples may be drawn by,
for example, the Metropolis-Hastings algorithm. Note that if there are no
${\bf x}_{i}$ for which $z_{i}=k$ (i.e. if the $k$-th component has no
observations associated with it), then $\boldsymbol{\phi}_{k}$ is simply
sampled from the prior, $p(\boldsymbol{\phi}_{k})$.

##### Conditional for $\pi$

By examination of the RHS of Equation 24, we have:

$\displaystyle p(\pi_{1},\ldots,\pi_{K}|{\bf x}_{1},\ldots,{\bf
x}_{n},z_{1},\ldots,z_{n},\boldsymbol{\phi}_{1},\ldots,\boldsymbol{\phi}_{K})$
$\displaystyle\propto\left(\prod_{i=1}^{n}p(z_{i}|\pi_{1},\ldots,\pi_{K})\right)p(\pi_{1},\ldots,\pi_{K}).$
(26)

Hence, the conditional for $\pi$ is the posterior for $\pi$ given the values
taken by the categorical latent allocation variables, $z_{i}$, $i=1,\ldots,n$.
If we take a conjugate Dirichlet prior, this posterior is available in closed
form.

##### Conditional for $z_{i}$

By examination of the RHS of Equation 24, we have:

$\displaystyle p(z_{i}=k|{\bf x}_{1},\ldots,{\bf
x}_{n},\boldsymbol{\phi}_{1},\ldots,\boldsymbol{\phi}_{K},\pi_{1},\ldots,\pi_{K},z_{-i})$
$\displaystyle\propto p(z_{i}=k|\pi)f_{\bf x}({\bf x}_{i}|\phi_{k}),$ (27)
$\displaystyle=\pi_{k}f_{\bf x}({\bf x}_{i}|\phi_{k}),$ (28)

where $z_{-i}$ denotes the set comprising all $z_{j}$ for which $j\neq i$.
Since $\sum_{k=1}^{K}p(z_{i}=k|{\bf x}_{1},\ldots,{\bf
x}_{n},\boldsymbol{\phi},\pi,z_{-i})=1$, it follows that the conditional is:

$\displaystyle p(z_{i}=k|{\bf x}_{1},\ldots,{\bf
x}_{n},\boldsymbol{\phi},\pi,z_{-i})$ $\displaystyle=\frac{\pi_{k}f_{\bf
x}({\bf x}_{i}|\phi_{k})}{\sum_{k=1}^{K}\pi_{k}f_{\bf x}({\bf
x}_{i}|\phi_{k})},$ (29)

which may be straightforwardly evaluated for finite $K$.

##### Marginalising $\pi$

Taking a conjugate Dirichlet prior for $\pi$, an alternative strategy is to
marginalise $\boldsymbol{\pi}$ rather than to sample it. We assume a symmetirc
Dirichlet prior with concentration parameter $\alpha/K$.

Note that the $z_{i}$’s are only conditionally independent of one another
given $\boldsymbol{\pi}$, so if we marginalise $\boldsymbol{\pi}$ then we must
be careful to model the dependence of $z_{i}$ on $z_{-i}$ in our conditional
for $z_{i}$.

We have

$\displaystyle p(z_{i}=k|{\bf x}_{1},\ldots,{\bf
x}_{n},\boldsymbol{\phi},\boldsymbol{\pi},z_{-i},\alpha)$
$\displaystyle\propto p(z_{i}=k|z_{-i},\pi,\alpha)f_{\bf x}({\bf
x}_{i}|\phi_{k})\mbox{\quad\quad\quad[cf. Equation \eqref{cicond1}]}.$ (30)

To marginalise $\boldsymbol{\pi}$, we must therefore evaluate
$\int_{\pi}p(z_{i}=k|z_{-i},\boldsymbol{\pi},\alpha)p(\boldsymbol{\pi}|\alpha)d\boldsymbol{\pi}=p(z_{i}=k|z_{-i},\alpha)$,
which is the conditional prior for $z_{i}$ given the values for the other
latent allocation variables, $z_{-i}$.

We have,

$\displaystyle p(z_{i}=k|z_{-i},\alpha)$
$\displaystyle=\frac{p(z_{i}=k,z_{-i}|\alpha)}{p(z_{-i}|\alpha)}$ (31)
$\displaystyle=\frac{\int_{\boldsymbol{\pi}}p(z_{i}=k,z_{-i}|\boldsymbol{\pi},\alpha)p(\boldsymbol{\pi}|\alpha)d\boldsymbol{\pi}}{\int_{\boldsymbol{\pi}}p(z_{-i}|\boldsymbol{\pi},\alpha)p(\boldsymbol{\pi}|\alpha)d\boldsymbol{\pi}}$
(32)
$\displaystyle=\frac{\int_{\boldsymbol{\pi}}p(z_{i}=k,z_{-i}|\boldsymbol{\pi})p(\boldsymbol{\pi}|\alpha)d\boldsymbol{\pi}}{\int_{\boldsymbol{\pi}}p(z_{-i}|\boldsymbol{\pi})p(\boldsymbol{\pi}|\alpha)d\boldsymbol{\pi}},$
(33)

where in the final line we exploit the fact that the $z_{i}$’s are
conditionally independent of $\alpha$, given $\pi$.

In order to proceed, we must evaluate this fraction. To do this we require a
standard result about Dirichlet distributions, which says that moments of
random variables distributed according to a symmetric Dirichlet distribution
with parameter $\alpha/K$ can be expressed as follows:

$E\left[\prod_{k=1}^{K}\pi_{k}^{m_{k}}\right]=\frac{\Gamma(\sum_{k=1}^{K}(\alpha/K))}{\Gamma(\sum_{k=1}^{K}((\alpha/K)+m_{k}))}\times\prod_{k=1}^{K}\frac{\Gamma((\alpha/K)+m_{k})}{\Gamma(\alpha/K)},$
(34)

where the $m_{k}$’s are any natural numbers.

Moreover, we note the following two equalities:

$p(z_{i}=k,z_{-i}|\pi)=\pi_{k}^{n_{-i,k}+1}\prod_{\begin{subarray}{c}c=1,\ldots,K\\\
c\neq k\end{subarray}}\pi_{c}^{n_{-i,c}},$

and

$p(z_{-i}|\pi)=\pi_{k}^{n_{-i,k}}\prod_{\begin{subarray}{c}c=1,\ldots,K\\\
c\neq k\end{subarray}}\pi_{c}^{n_{-i,c}},$

where $n_{-i,c}$ is the number of $z_{j}$’s with $j\neq i$ for which
$z_{j}=c$. It then follows that we may use the result given in Equation (34)
in order to evaluate the numerator and denominator in the RHS of Equation
(33). After some algebra, and exploiting the property of Gamma functions that
$\Gamma(t+1)=t\Gamma(t)$, we obtain:

$\displaystyle p(z_{i}=k|z_{-i},\alpha)$
$\displaystyle=\frac{n_{-i,k}+\alpha/K}{n-1+\alpha}.$ (35)

Hence,

$\displaystyle p(z_{i}=k|{\bf x}_{1},\ldots,{\bf
x}_{n},\boldsymbol{\phi},z_{-i},\alpha)$
$\displaystyle\propto\frac{n_{-i,k}+\alpha/K}{n-1+\alpha}\times f_{\bf x}({\bf
x}_{i}|\phi_{k}).$ (36)

Moreover, since $K$ is finite, we may straightforwardly evaluate the equality:

$\displaystyle p(z_{i}=k|{\bf x}_{1},\ldots,{\bf
x}_{n},\boldsymbol{\phi},z_{-i},\alpha)$
$\displaystyle=\frac{1}{Z}\frac{n_{-i,k}+\alpha/K}{n-1+\alpha}\times f_{\bf
x}({\bf x}_{i}|\phi_{k}),$ (37)

where

$Z=\sum_{c=1}^{K}\left(\frac{n_{-i,c}+\alpha/K}{n-1+\alpha}\times f_{\bf
x}({\bf x}_{i}|\phi_{c})\right).$ (38)

##### Marginalising $\phi$

Similarly, if a conjugate prior is available for the $\phi_{k}$’s, then these
may be marginalised too. Note that (similar to the case with the $z_{i}$’s
when we marginalised $\pi$) the ${\bf x}_{i}$’s are only conditionally
independent of one another given the $\phi_{k}$’s and the $z_{i}$’s, so if we
marginalise $\phi$ then we must be careful to model the dependence of $v_{i}$
on $v_{-i}$ in our conditional for $z_{i}$.

After some algebra, it is straightforward to show that marginalising $\phi$
gives the following for the conditional for $z_{i}$:

$\displaystyle p(z_{i}=k|{\bf x}_{1},\ldots,{\bf x}_{n},z_{-i},\alpha)$
$\displaystyle=\frac{1}{Z}\frac{n_{-i,k}+\alpha/K}{n-1+\alpha}\times\int_{\phi_{k}}f_{\bf
x}({\bf x}_{i}|\phi_{k})p(\phi_{k}|{\bf x}_{-i,k})d\phi_{k},$ (39)

where ${\bf x}_{-i,k}$ denotes all observations ${\bf x}_{j}$ for which $j\neq
i$ and $z_{j}=k$, and hence $p(\phi_{k}|{\bf x}_{-i,k})$ is the posterior for
$\phi_{k}$ given all of the observations currently associated with component
$k$ (excluding ${\bf x}_{i}$). If there are no ${\bf x}_{j}$ for which $j\neq
i$ and $z_{j}=k$ (i.e. if $k$ is a component to which no other observations
have been allocated), then we say that the $k$-th component is empty and
define $p(\phi_{k}|{\bf x}_{-i,k}):=p(\phi_{k})$ to be the prior for
$\phi_{k}$.

When implementing the sampler, it is useful to observe that

$p(\phi_{k}|{\bf x}_{-i,k})=\frac{f_{\bf x}({\bf
x}_{-i,k}|\phi_{k})p(\phi_{k})}{\int_{\phi_{k}}f_{\bf x}({\bf
x}_{-i,k}|\phi_{k})p(\phi_{k})d\phi_{k}},\mbox{ if the $k$-th component is not
empty.}$

Hence, still assuming that the $k$-th component is not empty, the integral in
Equation (39) is

$\int_{\phi_{k}}f_{\bf x}({\bf x}_{i}|\phi_{c})p(\phi_{k}|{\bf
x}_{-i,k})d\phi_{k}=\frac{\int_{\phi_{k}}f_{\bf x}({\bf x}_{i},{\bf
x}_{-i,k}|\phi_{k})p(\phi_{k})d\phi_{k}}{\int_{\phi_{k}}f_{\bf x}({\bf
x}_{-i,k}|\phi_{k})p(\phi_{k})d\phi_{k}},$ (40)

which is a ratio of marginal likelihoods: one in which we include ${\bf
x}_{i}$ amongst the observations associated with component $k$, and one in
which we exclude ${\bf x}_{i}$ from the observations associated with component
$k$.

This expression aids the interpretation of the sampler: at each iteration, and
for each component, we weigh the evidence that ${\bf x}_{i}$ is associated
with component $k$ against the evidence that ${\bf x}_{i}$ is not associated
with component $k$ (given the other observations currently associated with
that component, ${\bf x}_{-i,k}$). Intuitively, this expression ensures that
we are more likely to allocate ${\bf x}_{i}$ to a component to which similar
observations have previously been allocated.

Note also that the term $\frac{n_{-i,k}+\alpha/K}{n-1+\alpha}$ in Equation
(39) represents the conditional prior probability that ${\bf x}_{i}$ should be
allocated to component $k$ represents our prior belief that we should allocate
${\bf x}_{i}$ to component $k$, given the allocation of all of the other
observations. Since $n_{-i,k}$ is in the numerator, this expresses a “rich-
get-richer” prior belief; i.e. that, a priori, we are more likely to assign
${\bf x}_{i}$ to a component that already has many observations assigned to
it, rather than to one with fewer.

##### Final note

Note that, having marginalised the $\pi_{k}$’s and $\phi_{k}$’s, we may use
Equation (39) to sample just the $z_{i}$’s, without having to sample any other
parameters. The one exception is the $\alpha$ hyperparameter, which we may
either fix or sample (using, for example, the approach described in Escobar
and West, 1995).

####  7.1.2 Conditionals for Gibbs sampling ($K$ infinite)

The infinite case can be derived from the finite case by considering the limit
$K\rightarrow\infty$. The key change is to the conditional priors for the
component allocation variables, $p(z_{i}=k|z_{-i},\alpha)$. There are two
cases to consider:

##### 1\. $k$ is the label for an existing (non-empty) component:

This is the case where $k=z_{j}$ for some $j\neq i$. By considering the limit
as $K\rightarrow\infty$ of Equation (35), we obtain the following:

$\displaystyle\mbox{If $k=z_{j}$ for some $j\neq i$,
then\qquad}p(z_{i}=k|z_{-i},\alpha)$
$\displaystyle=\frac{n_{-i,k}}{n-1+\alpha}.$ (41)

##### 2\. $k$ is a new label:

This is the case where $k\neq z_{j}$ for all $j\neq i$. Hence, the probability
we seek is $p(z_{i}\neq z_{j}\mbox{ for all }j\neq i|z_{-i},\alpha)$.
Calculation of this probability is aided by using the following identity:

$1=p(z_{i}\neq z_{j}\mbox{ for all }j\neq
i|z_{-i},\alpha)+\sum_{k\in\mathcal{C}_{-i}}p(z_{i}=k|z_{-i},\alpha),$

where $\mathcal{C}_{-i}=\\{k:z_{j}=k\mbox{ for some }j\neq i\\}$ is the set of
all current component labels (excluding $z_{i}$). Using this identity, we
obtain (after a little algebra) the following:

$p(z_{i}\neq z_{j}\mbox{ for all }j\neq
i|z_{-i},\alpha)=\frac{\alpha}{n-1+\alpha}.$ (42)

##### Conditional for $z_{i}$

Using Equations (41) and (42), it is straightforward to show (cf. Equation
(39), and making use of Equation (40)) that we have the following conditionals
for $z_{i}$:

$\displaystyle\mbox{If $k=z_{j}$ for some $j\neq i$, then\qquad}p(z_{i}=k|{\bf
x}_{1},\ldots,{\bf x}_{n},z_{-i},\alpha)$
$\displaystyle=b\frac{n_{-i,k}}{n-1+\alpha}\frac{\int_{\phi_{k}}f_{\bf x}({\bf
x}_{i},{\bf x}_{-i,k}|\phi_{k})p(\phi_{k})d\phi_{k}}{\int_{\phi_{k}}f_{\bf
x}({\bf x}_{-i,k}|\phi_{k})p(\phi_{k})d\phi_{k}},$ (43)
$\displaystyle\mbox{and \qquad}p(z_{i}\neq z_{j}\mbox{ for all }j\neq i|{\bf
x}_{1},\ldots,{\bf x}_{n},z_{-i},\alpha)$
$\displaystyle=b\frac{\alpha}{n-1+\alpha}\int_{\phi_{k}}f_{\bf x}({\bf
x}_{i}|\phi_{k})p(\phi_{k})d\phi_{k},$ (44)

where $b$ is a normalising constant that ensures that

$p(z_{i}\neq z_{j}\mbox{ for all }j\neq i|{\bf x}_{1},\ldots,{\bf
x}_{n},z_{-i},\alpha)+\sum_{k\in\mathcal{C}_{-i}}p(z_{i}=k|{\bf
x}_{1},\ldots,{\bf x}_{n},z_{-i},\alpha)=1.$

Up to changes in notation and expansion of some terms, these expressions are
identical to those given in Equation (3.7) of Neal (2000).

###  7.2 Profile regression for categorical data

Recall that our data comprises input-output (i.e. clustering variable vector-
response) pairs, $D=\\{{\bf v}_{i}\\}_{i=1}^{n}$, where ${\bf v}_{i}=[{\bf
x}_{i},y_{i}]^{\top}$, with ${\bf x}_{i}$ and $y_{i}$ denoting the $i$-th
clustering variable vector and response, respectively. We now consider the
case in which the clustering variables and response are all categorical.

####  7.2.1 Modelling the clustering variables

We now assume that each clustering variable (i.e. each element of the vector
${\bf x}_{i}$) is categorical, with the $j$-th clustering variable having
$R_{j}$ categories, which we label as $1,2,\ldots,R_{j}$. We model the data
using categorical distributions. We define $\phi_{k,j,r}$ to be the
probability that the $j$-th clustering variable takes value $r$ in the $k$-th
component, and write
$\Phi_{k,j}=[\phi_{k,j,1},\phi_{k,j,2},\ldots,\phi_{k,j,R_{j}}]$ for the
collection of probabilities associated with the $j$-th clustering variable in
the $k$-th component. We further define
$\Phi_{k}=\\{\Phi_{k,1},\Phi_{k,2},\ldots,\Phi_{k,J}\\}$ to be the collection
of all probabilities (over all $J$ clustering variables) associated with the
$k$-th component, and $\Phi=\\{\Phi_{k}\\}_{k\in\mathcal{C}}$ to be the set of
all $\Phi_{k}$’s that are associated with non-empty components (here,
$\mathcal{C}=\\{k:z_{i}=k\mbox{ for some }i\in\\{1,\ldots,n\\}\\}$).

We assume that the clustering variables are conditionally independent, given
their component allocation, so that

$\displaystyle f_{\bf x}({\bf
x}_{i}=[x_{i1},x_{i2},\ldots,x_{iJ}]|\boldsymbol{\Phi},z_{i}=k)$
$\displaystyle=\phi_{k,j,x_{i1}}\phi_{k,j,x_{i2}}\ldots\phi_{k,j,x_{iJ}}$ (45)
$\displaystyle=\prod_{j=1}^{J}\phi_{k,j,x_{j}}$ (46)

##### Conditional for $\Phi_{k,j}$

From Equation (25), the conditional that we require for Gibbs sampling is the
posterior for $\Phi_{k,j}$, given the observations associated with the $k$-th
component. For each $j$, we adopt a conjugate Dirichlet prior for
$\Phi_{k,j}$,

$\Phi_{k,j}\sim\mbox{Dirichlet}({\bf a}_{j}),$

where ${\bf a}_{j}=[a_{j,1},\ldots,a_{j,R_{j}}]$ is the vector of
concentration parameters. The posterior is then:

$\displaystyle\Phi_{k,j}|{x}_{i_{1},j},{x}_{i_{2},j},\ldots,{x}_{i_{n_{k}},j},{\bf
a}_{j}\sim\mbox{Dirichlet}({\bf
a}_{j}+[s_{k,j,1},s_{k,j,2},\ldots,s_{k,j,R_{j}}]),$ (47)

where ${\bf x}_{i_{1}},{\bf x}_{i_{2}},\ldots,{\bf x}_{i_{n_{k}}}$ are the
observations associated with component $k$, and $s_{k,j,r}$ is defined to be
the number of observations associated with component $k$ for which the $j$-th
clustering variable is in category $r$.

##### Marginalising $\Phi_{k,j}$

We may also integrate out $\Phi_{k,j}$ in order to write down the marginal
likelihood associated with
${x}_{i_{1},j},{x}_{i_{2},j},\ldots,{x}_{i_{n_{k}},j}$. Note that the marginal
likelihood is (by definition) the prior expectation of the product
$\phi_{k,j,1}^{s_{k,j,1}}\ldots\phi_{k,j,{R_{j}}}^{s_{k,j,{R_{j}}}}$. We may
therefore use the same standard result that was used to derive Equation (34)
in order to immediately write down the marginal likelihood. Still assuming
that ${\bf x}_{i_{1}},{\bf x}_{i_{2}},\ldots,{\bf x}_{i_{n_{k}}}$ are the
observations associated with component $k$, we have:

$p({x}_{i_{1},j},{x}_{i_{2},j},\ldots,{x}_{i_{n_{k}},j}|{\bf
a}_{j})=\frac{\Gamma(\sum_{r=1}^{R_{j}}a_{j,r})}{\Gamma(\sum_{r=1}^{R_{j}}(a_{j,r}+s_{k,j,r}))}\times\prod_{r=1}^{R_{j}}\frac{\Gamma(a_{j,r}+s_{k,j,r})}{\Gamma(a_{j,r})}.$
(48)

To shorten notation, define $X_{k}=\\{{\bf x}_{i_{1}},{\bf
x}_{i_{2}},\ldots,{\bf x}_{i_{n_{k}}}\\}$ to be the set of observations
associated with component $k$, and
$X_{k,j}=\\{{x}_{i_{1},j},{x}_{i_{2},j},\ldots,{x}_{i_{n_{k}},j}\\}$ to be the
set containing the $j$-th elements of the vectors in $X_{k}$. Since we assume
that the clustering variables are conditionally independent, given their
component allocation, it follows that:

$p(X_{k}|{\bf a}_{1},\ldots,{\bf a}_{J})=\prod_{j=1}^{J}p(X_{k,j}|{\bf
a}_{j}),$ (49)

where $p(X_{k,j}|{\bf
a}_{j})=p({x}_{i_{1},j},{x}_{i_{2},j},\ldots,{x}_{i_{n_{k}},j}|{\bf a}_{j})$
is as given in Equation (48).

####  7.2.2 Modelling the response

We assume that the response, $y_{i}$, is also categorical (this time with
$R_{y}$ categories), and model using a categorical distribution. Similar to
before, we define $\theta_{k,r}$ to be the probability that $y$ takes value
$r$ in the $k$-th component, and write
$\Theta_{k}=[\theta_{k,1},\theta_{k,1},\ldots,\theta_{k,{R_{y}}}]$ to be the
collection of all probabilities associated with the $k$-th component. Also
define $Y_{k}=\\{y_{i_{1}},y_{i_{2}},\ldots,y_{i_{n_{k}}}\\}$ to be the set of
$y$’s allocated to component $k$.

##### Conditional for $\theta_{k}$

We adopt a conjugate Dirichlet prior for $\theta_{k}$,

$\theta_{k}\sim\mbox{Dirichlet}({\bf a}_{y}),$

where ${\bf a}_{y}=[a_{y,1},\ldots,a_{j,R_{y}}]$ is the vector of
concentration parameters. Similar to previously, the posterior for
$\theta_{k}$ is then

$\displaystyle\theta_{k}|Y_{k},{\bf a}_{y}\sim\mbox{Dirichlet}({\bf
a}_{y}+[s_{k,y,1},s_{k,y,2},\ldots,s_{k,y,R_{y}}]),$ (50)

where $s_{k,y,r}$ is the number of observations in $Y_{k}$ that are in the
$r$-th category.

##### Marginalising $\theta_{k}$

As before, we may also calculate the marginal likelihood, $p(Y_{k}|{\bf
a}_{j})$. We have

$p(Y_{k}|{\bf
a}_{y})=\frac{\Gamma(\sum_{r=1}^{R_{y}}a_{y,r})}{\Gamma(\sum_{r=1}^{R_{y}}(a_{y,r}+s_{k,y,r}))}\times\prod_{r=1}^{R_{y}}\frac{\Gamma(a_{y,r}+s_{k,y,r})}{\Gamma(a_{y,r})}.$
(51)

####  7.2.3 Joint marginal likelihood

In order to proceed, we need an expression for the marginal likelihood
associated with $V_{k}=\\{X_{k},Y_{k}\\}$. We assume that $y$ and ${\bf x}$
are conditionally independent, given their component allocation. It follows
that

$p(V_{k}|{\bf a}_{1},\ldots,{\bf a}_{J},{\bf a}_{y})=p(Y_{k}|{\bf
a}_{y})\prod_{j=1}^{J}p(X_{k,j}|{\bf a}_{j}),$ (52)

where the expressions for $p(Y_{k}|{\bf a}_{y})$ and $p(X_{k,j}|{\bf a}_{j})$
are as given previously. Note that:

  1. 1.

Setting $V_{k}=\\{{\bf x}_{i},{\bf x}_{-i,k}\\}$, we can evaluate the marginal
likelihood in the numerator in Equation (43).

  2. 2.

Setting $V_{k}=\\{{\bf x}_{-i,k}\\}$, we can evaluate the marginal likelihood
in the denominator in Equation (43).

  3. 3.

Setting $V_{k}=\\{{\bf x}_{i}\\}$, we can evaluate the marginal likelihood in
Equation (44).

Thus, we may evaluate all of the terms required for the conditionals for
$z_{i}$, and hence (leaving aside, for the time being, the problem of sampling
$\alpha$) we have everything we need in order to perform inference for our
model.

###  7.3 Variable selection

####  7.3.1 A null model for the clustering variable data

We introduce a null model for the ${\bf x}$ data, under the assumption that
there is no clustering structure. Under the null model, we again model the
clustering variable data using categorical distributions, and – similar to
previously – define
$\Phi_{0,j}=[\phi_{0,j,1},\phi_{0,j,2},\ldots,\phi_{0,j,R_{j}}]$ and
$\Phi_{0}=\\{\Phi_{0,1},\Phi_{0,2},\ldots,\Phi_{0,J}\\}$. However, under the
null model, we assume that there is no clustering structure (or, equivalently,
that the data form a single large cluster). The likelihood associated with
${\bf x}_{i}$ therefore does not involve a component allocation variable,
$z_{i}$, and is instead:

$\displaystyle f_{\bf x}({\bf
x}_{i}=[x_{i1},x_{i2},\ldots,x_{iJ}]^{\top}|\boldsymbol{\Phi}_{0})=\prod_{j=1}^{J}\phi_{0,j,x_{ij}}$
(53)

In our null model, the observations ${\bf x}_{1},\ldots,{\bf x}_{n}$ are
assumed to be conditionally independent given $\boldsymbol{\Phi}_{0}$. It
follows that, under the null model, the likelihood associated with the full
clustering variable dataset is:

$\displaystyle f_{\bf x}({\bf x}_{1},\ldots,{\bf
x}_{n}|\boldsymbol{\Phi}_{0})=\prod_{i=1}^{n}\prod_{j=1}^{J}\phi_{0,j,x_{ij}}.$
(54)

####  7.3.2 Variable selection

We consider the possibility that only some of the clustering variables are
relevant for determining the clustering structure, while the others are
irrelevant. We introduce binary indicators $\gamma_{j}$, such that
$\gamma_{j}$ = 1 means that the $j$-th clustering variable is relevant for the
clustering structure, while $\gamma_{j}$ = 0 means that the variable is
irrelevant.

Informally, we wish to use the full mixture model for the relevant variables,
and the null model for the irrelevant variables. More precisely, we consider
the following model for ${\bf x}$:

$\displaystyle f_{\bf x}({\bf
x}_{i}=[x_{i1},x_{i2},\ldots,x_{iJ}]^{\top}|\boldsymbol{\Phi},\boldsymbol{\Phi}_{0},z_{i}=k)$
$\displaystyle=\prod_{j=1}^{J}\phi_{k,j,x_{ij}}^{\mathbb{I}(\gamma_{j}=1)}\phi_{0,j,x_{ij}}^{\mathbb{I}(\gamma_{j}=0)}.$
(55)

####  7.3.3 Inference

We consider how the introduction of the $\gamma_{j}$’s affects how we infer
the model parameters. As previously, we start by writing down the full joint
model:

$\begin{split}p({\bf x}_{1},\ldots,{\bf
x}_{n},&z_{1},\ldots,z_{n},\boldsymbol{\theta},\boldsymbol{\Phi},\boldsymbol{\Phi}_{0},\pi,\alpha,\gamma)=\\\
&p(\pi|\alpha)p(\Phi_{0})\left(\prod_{j=1}^{J}p(\gamma_{j})\right)\left(\prod_{k=1}^{K}p(\Phi_{k})p(\theta_{k})\right)\left(\prod_{i=1}^{n}p(z_{i}|\pi)\theta_{z_{i},y_{i}}\left(\prod_{j=1}^{J}\phi_{z_{i},j,x_{ij}}^{\mathbb{I}(\gamma_{j}=1)}\phi_{0,j,x_{ij}}^{\mathbb{I}(\gamma_{j}=0)}\right)\right).\end{split}$
(56)

##### Inferring $\Phi_{k,j}$

As before, we adopt independent Dirichlet priors for $\Phi_{k,j}$,

$\Phi_{k,j}\sim\mbox{Dirichlet}({\bf a}_{j}),$ (57)

where ${\bf a}_{j}=[a_{j,1},\ldots,a_{j,R_{j}}]$ is the vector of
concentration parameters.

It is clear from Equation (56) that the conditional for $\Phi_{k,j}$ is

$p(\Phi_{k,j}|{\bf x}_{1},\ldots,{\bf
x}_{n},z_{1},\ldots,z_{n},\boldsymbol{\theta},\boldsymbol{\Phi}_{-{k,j}},\boldsymbol{\Phi}_{0},\pi,\alpha,\gamma)\propto
p(\Phi_{k,j})\prod_{i=1}^{n}\prod_{j=1}^{J}\phi_{z_{i},j,x_{ij}}^{\mathbb{I}(\gamma_{j}=1)}.$
(58)

If $\gamma_{j}=1$, then – as before – the conditional for $\Phi_{k,j}$ is the
posterior for $\Phi_{k,j}$ given all observations currently associated with
component $k$, i.e.

$\displaystyle\Phi_{k,j}|\gamma_{j}=1,{x}_{i_{1},j},{x}_{i_{2},j},\ldots,{x}_{i_{n_{k}},j},{\bf
a}_{j}\sim\mbox{Dirichlet}({\bf
a}_{j}+[s_{k,j,1},s_{k,j,2},\ldots,s_{k,j,R_{j}}]),$ (59)

where the $s_{k,j,r}$’s are as previously defined

If $\gamma_{j}=0$, then the conditional for $\Phi_{k,j}$ is just the prior,
$p(\Phi_{k,j})$, as provided in Equation (57). [But note that, in practice, we
would never need to sample from the conditional in this case].

##### Inferring $\Phi_{0,j}$

As we did for $\Phi_{k,j}$, we adopt a Dirichlet prior for $\Phi_{0,j}$,

$\Phi_{0,j}\sim\mbox{Dirichlet}({\bf a}_{j}^{(0)}),$ (60)

where ${\bf a}_{j}^{(0)}=[a_{j,1}^{(0)},\ldots,a_{j,R_{j}}^{(0)}]$ is the
vector of concentration parameters. It is clear from Equation (56) that the
conditional for $\Phi_{0,j}$ is

$p(\Phi_{0,j}|{\bf x}_{1},\ldots,{\bf
x}_{n},z_{1},\ldots,z_{n},\boldsymbol{\theta},\boldsymbol{\Phi}_{-{k,j}},\boldsymbol{\Phi}_{0},\pi,\alpha,\gamma)\propto
p(\Phi_{0,j})\prod_{i=1}^{n}\prod_{j=1}^{J}\phi_{0,j,x_{ij}}^{\mathbb{I}(\gamma_{j}=1)}.$
(61)

If $\gamma_{j}=1$, then the conditional for $\Phi_{0,j}$ is just the prior,
$p(\Phi_{0,j})$, as provided in Equation (60). [But note that, in practice, we
would never need to sample from the conditional in this case].

If $\gamma_{j}=0$, then the conditional for $\Phi_{0,j}$ is the posterior for
$\Phi_{0,j}$ given all observations, i.e.

$\displaystyle\Phi_{0,j}|\gamma_{j}=1,{x}_{1,j},{x}_{2,j},\ldots,{x}_{n,j},{\bf
a}_{j}^{(0)}\sim\mbox{Dirichlet}({\bf
a}_{j}^{(0)}+[s_{j,1},s_{j,2},\ldots,s_{j,R_{j}}]),$ (62)

where $s_{j,r}$ is defined to be the number of observations for which the
$j$-th variable is in category $r$.

####  7.3.4 Conditional for $z_{i}$

Let us temporarily return to the case of finite $K$, and recall Equation (37),
which provides the conditional for $z_{i}$ (assuming that the $\theta_{k}$’s
are sampled, rather than integrated out). Introducing the $\gamma_{j}$’s (to
allow for variable selection), and exploiting the conditional independence of
the variables given $z_{i}$, Equation (37) becomes:

$\displaystyle p(z_{i}=k|{\bf x}_{1},\ldots,{\bf
x}_{n},\boldsymbol{\theta},z_{-i},\alpha,\gamma)$
$\displaystyle=\frac{1}{Z}\frac{n_{-i,k}+\alpha/K}{n-1+\alpha}f_{y}({y}_{i}|\theta_{k})\left(\prod_{j:{\gamma_{j}}=1}f_{\bf
x}(x_{ij}|\Phi_{k,j})\right)\left(\prod_{j:{\gamma_{j}}=0}f_{\bf
x}(x_{ij}|\Phi_{0,j})\right)$ (63)
$\displaystyle=\frac{1}{Z}\frac{n_{-i,k}+\alpha/K}{n-1+\alpha}\theta_{k,y_{i}}\left(\prod_{j:{\gamma_{j}}=1}\phi_{k,j,x_{ij}}\right)\left(\prod_{j:{\gamma_{j}}=0}\phi_{0,j,x_{ij}}\right).$
(64)

##### Marginalising $\theta_{k}$

After some algebra, it is straightforward to show that marginalising the
component-specific parameters, $\theta_{k}$, gives the following for the
conditional for $z_{i}$:

$\displaystyle p(z_{i}=k|{\bf x}_{1},\ldots,{\bf x}_{n},z_{-i},\alpha)$
$\displaystyle=\frac{1}{Z}\frac{n_{-i,k}+\alpha/K}{n-1+\alpha}\times\int_{\theta_{k}}f_{y}(y_{i}|\theta_{k})p(\theta_{k}|\\{{y}_{m}\\}_{m\in\mathcal{M}_{-i,k}})d\theta_{k}$
$\displaystyle\times\prod_{j:{\gamma_{j}}=1}\left(\int_{\Phi_{k,j}}f_{{\bf
x}}(x_{ij}|\Phi_{k,j})p(\Phi_{k,j}|\\{{x}_{mj}\\}_{m\in\mathcal{M}_{-i,k}})d\Phi_{k,j}\right)$
$\displaystyle\times\prod_{j:{\gamma_{j}}=0}\left(\int_{\Phi_{0,j}}f_{{\bf
x}}(x_{ij}|\Phi_{0,j})p(\Phi_{0,j}|\\{{x}_{m}\\}_{m\in\mathcal{M}_{-i}})d\Phi_{0,j}\right),$
(65)

where

  * •

$\mathcal{M}_{-i}=\\{m:m\in\\{1,\ldots,n\\}\mbox{ and }m\neq i\\}$; and

  * •

$\mathcal{M}_{-i,k}=\\{m:m\in\\{1,\ldots,n\\}\mbox{ and }m\neq i\mbox{ and
}z_{m}=k\\}$.

Thus, as previously, $p(\theta_{k}|\\{{y}_{m}\\}_{m\in\mathcal{M}_{-i,k}})$ is
the posterior for $\theta_{k}$ given all observations (excluding $y_{i}$)
currently associated with the $k$-th component (and similarly for
$\Phi_{k,j}$). Also as before, if the $k$-th component is empty we define
$p(\theta_{k}|\\{{y}_{m}\\}_{m\in\mathcal{M}_{-i,k}}):=p(\theta_{k})$ to be
the prior for $\theta_{k}$ (and similarly for $\Phi_{k,j}$).

Note that the expression involving $\Phi_{0,j}$ is slightly different: here,
we have the posterior for $\Phi_{0,j}$ given all observations except the
$i$-th. Since this expression does not involve $k$, we may absorb this term
into the normalising constant, $Z$, so that:

$\displaystyle p(z_{i}=k|{\bf x}_{1},\ldots,{\bf x}_{n},z_{-i},\alpha)$
$\displaystyle=\frac{1}{Z^{\prime}}\frac{n_{-i,k}+\alpha/K}{n-1+\alpha}\times\int_{\theta_{k}}f_{y}(y_{i}|\theta_{k})p(\theta_{k}|\\{{y}_{m}\\}_{m\in\mathcal{M}_{-i,k}})d\theta_{k}$
$\displaystyle\times\prod_{j:{\gamma_{j}}=1}\left(\int_{\Phi_{k,j}}f_{{\bf
x}}(x_{ij}|\Phi_{k,j})p(\Phi_{k,j}|\\{{x}_{mj}\\}_{m\in\mathcal{M}_{-i,k}})d\Phi_{k,j}\right).$
(66)

##### Conditional for $z_{i}$ ($K$ infinite)

Proceeding as previously, we consider the limit as $K\rightarrow\infty$. In
this case, the conditionals for $z_{i}$ are as follows:

If $k=z_{j}$ for some $j\neq i$, then

$\displaystyle p(z_{i}=k|{\bf x}_{1},\ldots,{\bf x}_{n},z_{-i},\alpha)$
$\displaystyle=b\frac{n_{-i,k}}{n-1+\alpha}\times\frac{\int_{\theta_{k}}f_{y}(y_{i},\\{{y}_{m}\\}_{m\in\mathcal{M}_{-i,k}}|\theta_{k})p(\theta_{k})d\theta_{k}}{\int_{\theta_{k}}f_{y}(\\{{y}_{m}\\}_{m\in\mathcal{M}_{-i,k}}|\theta_{k})p(\theta_{k})d\theta_{k}}$
$\displaystyle\times\frac{\prod_{j:{\gamma_{j}}=1}\left(\int_{\Phi_{k,j}}f_{{\bf
x}}(x_{ij},\\{{x}_{mj}\\}_{m\in\mathcal{M}_{-i,k}}|\Phi_{k,j})p(\Phi_{k,j})d\Phi_{k,j}\right)}{\prod_{j:{\gamma_{j}}=1}\left(\int_{\Phi_{k,j}}f_{{\bf
x}}(\\{{x}_{mj}\\}_{m\in\mathcal{M}_{-i,k}}|\Phi_{k,j})p(\Phi_{k,j})d\Phi_{k,j}\right)}.$
(67)

Moreover,

$\displaystyle p(z_{i}\neq z_{j}\mbox{ for all }j\neq i|{\bf
x}_{1},\ldots,{\bf x}_{n},z_{-i},\alpha)$
$\displaystyle=b\frac{\alpha}{n-1+\alpha}\times\int_{\theta_{k}}f_{y}(y_{i}|\theta_{k})p(\theta_{k})d\theta_{k}$
$\displaystyle\times\prod_{j:{\gamma_{j}}=1}\left(\int_{\Phi_{k,j}}f_{{\bf
x}}(x_{ij}|\Phi_{k,j})p(\Phi_{k,j})d\Phi_{k,j}\right),$ (68)

where $b$ is a normalising constant that ensures that

$p(z_{i}\neq z_{j}\mbox{ for all }j\neq i|{\bf x}_{1},\ldots,{\bf
x}_{n},z_{-i},\alpha)+\sum_{k\in\mathcal{C}_{-i}}p(z_{i}=k|{\bf
x}_{1},\ldots,{\bf x}_{n},z_{-i},\alpha)=1.$

Note that these expressions are essentially the same as those provided in
Equations (43) and (44), except that the irrelevant variables do not
contribute to the conditionals (as we would hope). The marginal likelihoods
may be evaluated as previously described (see Sections 7.2.1 and 7.2.2).

####  7.3.5 Sampling $\gamma_{j}$

By examination of Equation (56), the conditional for $\gamma_{j}$ is
immediately

$p(\gamma_{j}|{\bf x}_{1},\ldots,{\bf
x}_{n},z_{1},\ldots,z_{n},\boldsymbol{\theta},\boldsymbol{\Phi},\boldsymbol{\Phi}_{0},\pi,\alpha,\gamma_{-j})\propto
p(\gamma_{j})\left(\prod_{i=1}^{n}\phi_{z_{i},j,x_{ij}}^{\mathbb{I}(\gamma_{j}=1)}\phi_{0,j,x_{ij}}^{\mathbb{I}(\gamma_{j}=0)}\right).$
(69)

Since $\gamma_{j}\in\\{0,1\\}$, we have:

$p(\gamma_{j}=1|{\bf x}_{1},\ldots,{\bf
x}_{n},z_{1},\ldots,z_{n},\boldsymbol{\theta},\boldsymbol{\Phi},\boldsymbol{\Phi}_{0},\pi,\alpha,\gamma_{-j})=\frac{1}{z}p(\gamma_{j}=1)\prod_{i=1}^{n}\phi_{z_{i},j,x_{ij}},$
(70)

and

$p(\gamma_{j}=0|{\bf x}_{1},\ldots,{\bf
x}_{n},z_{1},\ldots,z_{n},\boldsymbol{\theta},\boldsymbol{\Phi},\boldsymbol{\Phi}_{0},\pi,\alpha,\gamma_{-j})=\frac{1}{z}p(\gamma_{j}=0)\prod_{i=1}^{n}\phi_{0,j,x_{ij}},$
(71)

where

$z=p(\gamma_{j}=1)\prod_{i=1}^{n}\phi_{z_{i},j,x_{ij}}+p(\gamma_{j}=0)\prod_{i=1}^{n}\phi_{0,j,x_{ij}}.$

##### Marginalising $\Phi_{k,j}$ and $\Phi_{0,j}$

As previously, we may also marginalise the $\Phi_{k,j}$ and $\Phi_{0,j}$
parameters. Define $X_{k}=\\{{\bf x}_{i_{1}},{\bf x}_{i_{2}},\ldots,{\bf
x}_{i_{n_{k}}}\\}$ to be the set of observations associated with component
$k$, and $X_{k,j}=\\{{x}_{i_{1},j},{x}_{i_{2},j},\ldots,{x}_{i_{n_{k}},j}\\}$
to be the set containing the $j$-th elements of the vectors in $X_{k}$ (as
before). Also define $\mathcal{C}=\\{k:z_{i}=k\mbox{ for some }i\\}$ to be the
set of all current component labels, and
$X_{\cdot,j}=\cup_{k\in\mathcal{C}}\\{X_{k,j}\\}$ to be the set containing the
$j$-th elements of all observations.

We then obtain

$p(\gamma_{j}=1|{\bf x}_{1},\ldots,{\bf
x}_{n},z_{1},\ldots,z_{n},\boldsymbol{\theta},\pi,\alpha,\gamma_{-j})=\frac{1}{z^{\prime}}p(\gamma_{j}=1)\prod_{k\in\mathcal{C}}p(X_{k,j}|{\bf
a}_{j}),$ (72)

and

$p(\gamma_{j}=0|{\bf x}_{1},\ldots,{\bf
x}_{n},z_{1},\ldots,z_{n},\boldsymbol{\theta},\pi,\alpha,\gamma_{-j})=\frac{1}{z^{\prime}}p(\gamma_{j}=0)p(X_{\cdot,j}|{\bf
a}_{j}^{(0)}),$ (73)

where $z^{\prime}$ is a normalising constant that ensures that the
probabilities sum to 1. The marginal likelihoods may be evaluated as
previously described (see Section 7.2.1).

###  7.4 Multi-view clustering

The method for variable selection described in the previous section may be
considered to be a special case of multi-view clustering, in which we allow
the data to possess multiple clustering structures (depending upon which group
of variables we consider to be “relevant”). In the variable selection case, we
model the data as possessing 2 clustering structures: one non-trivial
clustering structure (to which the “relevant” variables contribute), and one
trivial clustering structure (comprising just a single cluster, to which the
“irrelevant” variables contribute). This may be straightforwardly extended to
$L$ clustering structures, by allowing $\gamma_{j}$ to be a categorical
variable whose value indicates the clustering structure to which the $j$-th
variable contributes.

####  7.4.1 Notational conventions

We assume that $\gamma_{j}\in\\{0,1,2,\ldots,V-1\\}$ is a categorical variable
with $L$ categories. As previously, we assume that if $\gamma_{j}=0$ then the
$j$-th variable is completely irrelevant, and does not contribute to any
clustering structure (except the trivial structure comprising one cluster).

$\displaystyle f({\bf x}_{i},y_{i}|\cdots)$
$\displaystyle=\left(\prod_{j=1}^{J}\left(\prod_{v=1}^{V-1}f_{\bf
x}(x_{ij}|\theta^{(v)}_{z_{i}^{(v)}})^{\mathbb{I}(\gamma_{j}=v)}\right)f_{\bf
x}(x_{ij}|\theta_{0})^{\mathbb{I}(\gamma_{j}=0)}\right)f_{\bf
y}({y}_{i}|\phi_{z_{i}^{(1)}},{\bf w}_{i},\beta)$

$\displaystyle p(\gamma_{j}=0|\cdots)$ $\displaystyle\propto
p(\gamma_{j}=0)\prod_{i=1}^{n}f_{\bf x}(x_{ij}|\theta_{0})$ $\displaystyle
p(\gamma_{j}=v|\cdots)$ $\displaystyle\propto
p(\gamma_{j}=v)\prod_{i=1}^{n}f_{\bf
x}(x_{ij}|\theta_{z_{i}^{(v)}}^{(v)})\mbox{\qquad\qquad for $v=1,\ldots,V-1$}$

For $\ell=1,\ldots,L-1$, we assume that there are non-overlapping variable
sets $\mathcal{S}_{\ell}$, such that each $\mathcal{S}_{\ell}$ possesses its
own clustering structure. Within each $\mathcal{S}_{\ell}$, we model the data
using a mixture model. To this end, we introduce latent component allocation
variables, $z_{i}^{(\ell)}$, which indicate the component in the model for
$\mathcal{S}_{\ell}$ that is responsible for generating the $i$-th
observation. Similarly, we define $\Phi_{k}^{(\ell)}$ to be the parameters
associated with the $k$-th component in the $\ell$-th mixture model.

We moreover define $\gamma_{j}=\ell$ if the $j$-th variable is in
$\mathcal{S}_{\ell}$. Note that, a priori, we do not know which variables
belong to which $\mathcal{S}_{\ell}$, so we must perform inference for the
$\gamma_{j}$’s.

Finally, we introduce another categorical variable,
$\nu\in\\{1,2,\ldots,L-1\\}$, such that the clustering structure present in
$\mathcal{S}_{\nu}$ is the one that is relevant for profile regression, i.e.
the likelihood for $y_{i}$,

$f_{\bf
y}(y_{i}|\boldsymbol{\theta},z_{i}^{(1)},\ldots,z_{i}^{(L-1)},\nu)=f_{\bf
y}(y|\boldsymbol{\theta}_{z_{i}^{(\ell)}}),$

depends only on the $z_{i}^{(\ell)}$ component allocation variable.

Our conditional model for ${\bf x}_{i}=[{\bf x}_{i},y_{i}]^{\top}$ (given the
$z_{i}^{(\ell)}$’s, $\gamma_{j}$’s, and $\nu$) is then

$\begin{split}p({\bf
x}_{i}|z_{i}^{(1)},\ldots,z_{i}^{(L-1)},&\boldsymbol{\theta},\Phi_{0},\Phi^{(1)},\ldots,\Phi^{(L-1)},\gamma,\nu)\\\
&=f_{\bf
y}(y|\boldsymbol{\theta}_{z_{i}^{(\ell)}})\prod_{j=1}^{J}\left(f(x_{ij}|\phi_{0,j,x_{ij}})^{\mathbb{I}(\gamma_{j}=0)}\prod_{\ell=1}^{L-1}\left(f(x_{ij}|\phi^{(\ell)}_{z_{i}^{(\ell)},j,x_{ij}})^{\mathbb{I}(\gamma_{j}=\ell)}\right)\right).\end{split}$
(74)

####  7.4.2 Conditionals for Gibbs sampling

We may follow the same method of argument presented in Section 3 in order to
derive the required conditionals. For brevity, we provide the conditionals
obtained after integrating out the component-specific parameters.

##### Conditionals for $z_{i}^{(\ell)}$: $\nu=\ell$ case

In this case, we are dealing with the clustering structure that is relevant
for profile regression, so must include the contribution of the responses
variable, $y$. We have (cf. Equations (67) and (68)): If $k=z_{j}^{(\ell)}$
for some $j\neq i$, then

$\displaystyle p(z_{i}^{(\ell)}=k|{\bf x}_{1},\ldots,{\bf
x}_{n},z_{-i}^{(\ell)},\alpha^{(\ell)},\nu=\ell)$
$\displaystyle=b^{(\ell)}\frac{n_{-i,k}^{(\ell)}}{n-1+\alpha^{(\ell)}}\times\frac{\int_{\theta_{k}}f_{y}(y_{i},\\{{y}_{m}\\}_{m\in\mathcal{M}^{(\ell)}_{-i,k}}|\theta_{k})p(\theta_{k})d\theta_{k}}{\int_{\theta_{k}}f_{y}(\\{{y}_{m}\\}_{m\in\mathcal{M}^{(\ell)}_{-i,k}}|\theta_{k})p(\theta_{k})d\theta_{k}}$
$\displaystyle\times\frac{\prod_{j:{\gamma_{j}}=1}\left(\int_{\Phi^{(\ell)}_{k,j}}f_{{\bf
x}}(x_{ij},\\{{x}_{mj}\\}_{m\in\mathcal{M}^{(\ell)}_{-i,k}}|\Phi^{(\ell)}_{k,j})p(\Phi^{(\ell)}_{k,j})d\Phi^{(\ell)}_{k,j}\right)}{\prod_{j:{\gamma_{j}}=1}\left(\int_{\Phi^{(\ell)}_{k,j}}f_{{\bf
x}}(\\{{x}_{mj}\\}_{m\in\mathcal{M}^{(\ell)}_{-i,k}}|\Phi^{(\ell)}_{k,j})p(\Phi^{(\ell)}_{k,j})d\Phi^{(\ell)}_{k,j}\right)},$
(75)

where $n_{-i,k}^{(\ell)}$ is the number of observations, ${\bf x}_{m}$, for
which $z_{m}^{(\ell)}=k$ (excluding the $i$-th observation), and
$\mathcal{M}_{-i,k}^{(\ell)}=\\{m:m\in\\{1,\ldots,n\\}\mbox{ and }m\neq
i\mbox{ and }z_{m}^{(\ell)}=k\\}$.

Moreover,

$\displaystyle p(z_{i}^{(\ell)}\neq z_{j}^{(\ell)}\mbox{ for all }j\neq i|{\bf
x}_{1},\ldots,{\bf x}_{n},c^{(\ell)}_{-i},\alpha^{(\ell)},\nu=\ell)$
$\displaystyle=b^{(\ell)}\frac{\alpha^{(\ell)}}{n-1+\alpha^{(\ell)}}\times\int_{\theta_{k}}f_{y}(y_{i}|\theta_{k})p(\theta_{k})d\theta_{k}$
$\displaystyle\times\prod_{j:{\gamma_{j}}=1}\left(\int_{\Phi^{(\ell)}_{k,j}}f_{{\bf
x}}(x_{ij}|\Phi^{(\ell)}_{k,j})p(\Phi^{(\ell)}_{k,j})d\Phi^{(\ell)}_{k,j}\right),$
(76)

where $b^{(\ell)}$ is a constant that ensures that the probabilities sum to 1,
and $\alpha^{(\ell)}$ is the concentration parameter of the Dirichlet process
prior on the mixture weights for the $\nu$-th mixture model.

##### Conditionals for $z_{i}^{(\ell)}$: $\nu\neq\ell$ case

In this case, the response variable, $y$, does not contribute to the
conditionals. We have (cf. Equations (67) and (68)): If $k=z_{j}^{(\ell)}$ for
some $j\neq i$, then

$\displaystyle p(z_{i}^{(\ell)}=k|{\bf x}_{1},\ldots,{\bf
x}_{n},z_{-i}^{(\ell)},\alpha^{(\ell)},\nu\neq\ell)$
$\displaystyle=b^{(\ell)}\frac{n_{-i,k}^{(\ell)}}{n-1+\alpha^{(\ell)}}\frac{\prod_{j:{\gamma_{j}}=\ell}\left(\int_{\Phi^{(\ell)}_{k,j}}f_{{\bf
x}}(x_{ij},\\{{x}_{mj}\\}_{m\in\mathcal{M}_{-i,k}}|\Phi^{(\ell)}_{k,j})p(\Phi^{(\ell)}_{k,j})d\Phi^{(\ell)}_{k,j}\right)}{\prod_{j:{\gamma_{j}}=\ell}\left(\int_{\Phi^{(\ell)}_{k,j}}f_{{\bf
x}}(\\{{x}_{mj}\\}_{m\in\mathcal{M}_{-i,k}}|\Phi^{(\ell)}_{k,j})p(\Phi^{(\ell)}_{k,j})d\Phi^{(\ell)}_{k,j}\right)}.$
(77)

where $n_{-i,k}^{(\ell)}$ is the number of observations, ${\bf x}_{m}$, for
which $z_{m}^{(\ell)}=k$ (excluding the $i$-th observation). Moreover,

$\displaystyle p(z_{i}^{(\ell)}\neq z_{j}^{(\ell)}\mbox{ for all }j\neq i|{\bf
x}_{1},\ldots,{\bf x}_{n},c^{(\ell)}_{-i},\alpha^{(\ell)},\nu\neq\ell)$
$\displaystyle=b^{(\ell)}\frac{\alpha^{(\ell)}}{n-1+\alpha^{(\ell)}}\prod_{j:{\gamma_{j}}=\ell}\left(\int_{\Phi^{(\ell)}_{k,j}}f_{{\bf
x}}(x_{ij}|\Phi^{(\ell)}_{k,j})p(\Phi^{(\ell)}_{k,j})d\Phi^{(\ell)}_{k,j}\right),$
(78)

where $b^{(\ell)}$ is a constant that ensures that the probabilities sum to 1,
and $\alpha^{(\ell)}$ is the concentration parameter of the Dirichlet process
prior on the mixture weights for the $\ell$-the mixture model.

####  7.4.3 Conditional for $\gamma_{j}$

Define $X_{k}^{(\ell)}=\\{{\bf x}_{i_{1}},{\bf x}_{i_{2}},\ldots,{\bf
x}_{i_{n_{k}^{(\ell)}}}\\}$ to be the set of observations associated with
component $k$ of the $\ell$-th mixture model, and
$X_{k,j}^{(\ell)}=\\{{x}_{i_{1},j},{x}_{i_{2},j},\ldots,{x}_{i_{n_{k}^{(\ell)}},j}\\}$
to be the set containing the $j$-th elements of the vectors in
$X_{k}^{(\ell)}$. Also define
$\mathcal{C}^{(\ell)}=\\{k:z_{i}^{(\ell)}=k\mbox{ for some }i\\}$ to be the
set of all current component labels for the $\ell$-th mixture model.

We have

$p(\gamma_{j}=0|\ldots)=\frac{1}{z^{\prime}}p(\gamma_{j}=0)p(X_{\cdot,j}|{\bf
a}_{j}^{(0)}),$ (79)

and, for all other $\ell\neq 0$,

$p(\gamma_{j}=\ell|\ldots)=\frac{1}{z^{\prime}}p(\gamma_{j}=\ell)\prod_{k\in\mathcal{C}^{(\ell)}}p(X_{k,j}^{(\ell)}|{\bf
a}_{j}^{(\ell)}),$ (80)

where $z^{\prime}$ is a normalising constant that ensures that the
probabilities sum to 1. The marginal likelihoods may be evaluated as
previously described (see Section 7.2.1).

####  7.4.4 Conditional for $\nu$

The conditional for $\nu$ is rather similar to the conditional for
$\gamma_{j}$. Define
$Y_{k}^{(\ell)}=\\{{y}_{i_{1}},{y}_{i_{2}},\ldots,{y}_{i_{n_{k}^{(\ell)}}}\\}$.
Then, for $\ell\in\\{1,\ldots,L-1\\}$, we have:

$p(\nu=\ell|\ldots)=\frac{1}{z^{\prime\prime}}p(\nu=\ell)\prod_{k\in\mathcal{C}^{(\ell)}}p(Y_{k}^{(\ell)}|{\bf
a}_{y}),$ (81)

where $z^{\prime\prime}$ is a normalising constant to ensure that
probabilities sum to 1.

Note that, when sampling $\nu$, we are effectively comparing different models
for partitioning the $y$’s, and selecting from amongst these according to
probabilities given by Equation (81) above.

## References

  * Aldous et al. (1985) Aldous, D. J., I. A. Ibragimov, J. Jacod, and P. L. Hennequin, eds.
1985\.  Exchangeability and related topics, Berlin, Heidelberg. Springer
Berlin Heidelberg.

  * Alon et al. (1999) Alon, U., N. Barkai, D. A. Notterman, K. Gish, S. Ybarra, D. Mack, and A. J. Levine
1999\.  Broad patterns of gene expression revealed by clustering analysis of
tumor and normal colon tissues probed by oligonucleotide arrays.  Proceedings
Of The National Academy Of Sciences Of The United States Of America,
96(12):6745–6750.

  * Ben-Dor et al. (1999) Ben-Dor, A., R. Shamir, and Z. Yakhini
1999\.  Clustering gene expression patterns.  Journal of computational biology
: a journal of computational molecular cell biology, 6(3-4):281–297.

  * Blackwell and MacQueen (1973) Blackwell, D. and J. B. MacQueen
1973\.  Ferguson Distributions Via Pólya Urn Schemes.  The Annals of
Statistics, 1(2):353–355.

  * Cancer Genome Atlas Network (2012) Cancer Genome Atlas Network
2012\.  Comprehensive molecular portraits of human breast tumours.  Nature,
490(7418):61–70.

  * Crook et al. (2018) Crook, O. M., L. Gatto, and P. D. W. Kirk
2018\.  Fast approximate inference for variable selection in Dirichlet process
mixtures, with an application to pan-cancer proteomics.  arXiv.org.

  * Cui et al. (2007) Cui, Y., X. Z. Fern, and J. G. Dy
2007\.  Non-redundant Multi-view Clustering via Orthogonalization.  In Seventh
IEEE International Conference on Data Mining (ICDM 2007), Pp. 133–142. IEEE.

  * Curtis et al. (2012) Curtis, C., S. P. Shah, S.-F. Chin, G. Turashvili, O. M. Rueda, M. J. Dunning, D. Speed, A. G. Lynch, S. Samarajiwa, Y. Yuan, S. Gräf, G. Ha, G. Haffari, A. Bashashati, R. Russell, S. McKinney, METABRIC Group, A. Langerød, A. Green, E. Provenzano, G. Wishart, S. Pinder, P. Watson, F. Markowetz, L. Murphy, I. Ellis, A. Purushotham, A.-L. Børresen-Dale, J. D. Brenton, S. Tavare, C. Caldas, and S. Aparicio
2012\.  The genomic and transcriptomic architecture of 2,000 breast tumours
reveals novel subgroups.  Nature, 486(7403):346–352.

  * Eisen et al. (1998) Eisen, M. B., P. T. Spellman, P. O. Brown, and D. Botstein
1998\.  Cluster analysis and display of genome-wide expression patterns.
Proceedings Of The National Academy Of Sciences Of The United States Of
America, 95(25):14863–14868.

  * Escobar and West (1995) Escobar, M. and M. West
1995\.  Bayesian density estimation and inference using mixtures.  Journal of
the American Statistical Association, 90(430):577–588.

  * Fop and Murphy (2018) Fop, M. and T. B. Murphy
2018\.  Variable selection methods for model-based clustering.  Statistics
Surveys, 12:18–65.

  * Fritsch and Ickstadt (2009) Fritsch, A. and K. Ickstadt
2009\.  Improved Criteria for Clustering Based on the Posterior Similarity
Matrix.  Bayesian Analysis, 4(2):367–391.

  * Golub et al. (1999) Golub, T. R., D. K. Slonim, P. Tamayo, C. Huard, M. Gaasenbeek, J. P. Mesirov, H. Coller, M. L. Loh, J. R. Downing, M. A. Caligiuri, C. D. Bloomfield, and E. S. Lander
1999\.  Molecular classification of cancer: class discovery and class
prediction by gene expression monitoring.  Science (New York, NY),
286(5439):531–537.

  * Green and Richardson (2001) Green, P. J. and S. Richardson
2001\.  Modelling Heterogeneity with and without the Dirichlet Process.
Scandinavian Journal Of Statistics, 28(2):355–375.

  * Guan et al. (2010) Guan, Y., J. G. Dy, D. Niu, and Z. Ghahramani
2010\.  Variational inference for nonparametric multiple clustering.
MultiClust Workshop.

  * Hastie et al. (2014) Hastie, D. I., S. Liverani, and S. Richardson
2014\.  Sampling from Dirichlet process mixture models with unknown
concentration parameter: mixing issues in large data implementations.
Statistics And Computing, Pp. 1–15.

  * Heyer et al. (1999) Heyer, L. J., S. Kruglyak, and S. Yooseph
1999\.  Exploring expression data: identification and analysis of coexpressed
genes.  Genome research, 9(11):1106–1115.

  * Ishwaran and James (2001) Ishwaran, H. and L. F. James
2001\.  Gibbs Sampling Methods for Stick-Breaking Priors.  Journal Of The
American Statistical Association, 96(453):161–173.

  * Ishwaran and Zarepour (2002) Ishwaran, H. and M. Zarepour
2002\.  Exact and approximate representations for the sum Dirichlet process.
Can J Stat, 30(2):269–283.

  * Jain and Neal (2004) Jain, S. and R. M. Neal
2004\.  A split-merge Markov chain Monte Carlo procedure for the dirichlet
process mixture model.  Journal Of Computational And Graphical Statistics,
13(1):158–182.

  * Jain and Neal (2007) Jain, S. and R. M. Neal
2007\.  Splitting and Merging Components of a Nonconjugate Dirichlet Process
Mixture Model.  Bayesian Analysis, 2(3):445–472.

  * Kalli et al. (2011) Kalli, M., J. E. Griffin, and S. G. Walker
2011\.  Slice sampling mixture models.  Statistics And Computing,
21(1):93–105.

  * Kirk et al. (2012) Kirk, P., J. E. Griffin, R. S. Savage, Z. Ghahramani, and D. L. Wild
2012\.  Bayesian correlated clustering to integrate multiple datasets.
Bioinformatics, 28(24):3290–3297.

  * Kristensen et al. (2014) Kristensen, V. N., O. C. Lingjærde, H. G. Russnes, H. K. M. Vollan, A. Frigessi, and A.-L. Børresen-Dale
2014\.  Principles and methods of integrative genomic analyses in cancer.
Nature reviews. Cancer, 14(5):299–313.

  * Kuijjer et al. (2018) Kuijjer, M. L., J. N. Paulson, P. Salzman, W. Ding, and J. Quackenbush
2018\.  Cancer subtype identification using somatic mutation data.  British
journal of cancer, 118(11):1492–1501.

  * Law et al. (2003) Law, M. H., A. K. Jain, and M. Figueiredo
2003\.  Feature Selection in Mixture-Based Clustering.  In Advances in Neural
Information Processing Systems 15, S. Becker, S. Thrun, and K. Obermayer,
eds., Pp. 641–648.  MIT Press.

  * Law et al. (2004) Law, M. H. C., M. A. T. Figueiredo, and A. K. Jain
2004\.  Simultaneous feature selection and clustering using mixture models.
Pattern Analysis and Machine Intelligence, IEEE Transactions on,
26(9):1154–1166.

  * Li and Shafto (2011) Li, D. and P. Shafto
2011\.  Bayesian Hierarchical Cross-Clustering.  In Proceedings of the
Fourteenth International Conference on Artificial Intelligence and Statistics,
Pp. 443–451, Fort Lauderdale, FL, USA. PMLR.

  * Liverani et al. (2015) Liverani, S., D. I. Hastie, L. Azizi, M. Papathomas, and S. Richardson
2015\.  PReMiuM: An R Package for Profile Regression Mixture Models Using
Dirichlet Processes.  Journal Of Statistical Software, 64(7).

  * Lock and Dunson (2013) Lock, E. F. and D. B. Dunson
2013\.  Bayesian consensus clustering.  Bioinformatics, P. btt425.

  * McLachlan et al. (2002) McLachlan, G. J., R. W. Bean, and D. Peel
2002\.  A mixture model-based approach to the clustering of microarray
expression data.  Bioinformatics (Oxford, England), 18(3):413–422.

  * Miller and Harrison (2017) Miller, J. W. and M. T. Harrison
2017\.  Mixture Models With a Prior on the Number of Components.  Journal Of
The American Statistical Association, 113(521):340–356.

  * Molitor et al. (2010) Molitor, J., M. Papathomas, M. Jerrett, and S. Richardson
2010\.  Bayesian profile regression with an application to the National Survey
of Children’s Health.  Biostatistics (Oxford, England), 11(3):484–498.

  * Neal (2000) Neal, R. M.
2000\.  Markov Chain Sampling Methods for Dirichlet Process Mixture Models.
Journal of Computational and Graphical Statistics, 9(2):249–265.

  * Niu et al. (2010) Niu, D., J. G. Dy, and M. I. Jordan
2010\.  Multiple Non-redundant Spectral Clustering Views.  In Proceedings of
the 27th International Conference on International Conference on Machine
Learning, Pp. 831–838, USA. Omnipress.

  * Niu et al. (2014) Niu, D., J. G. Dy, and M. I. Jordan
2014\.  Iterative Discovery of Multiple Alternative Clustering Views.  Pattern
Analysis and Machine Intelligence, IEEE Transactions on, 36(7):1340–1353.

  * Papathomas et al. (2012) Papathomas, M., J. Molitor, C. Hoggart, D. Hastie, and S. Richardson
2012\.  Exploring data from genetic association studies using Bayesian
variable selection and the Dirichlet process: application to searching for
gene-gene patterns.  Genetic epidemiology, 36(6):663–674.

  * Parker et al. (2009) Parker, J. S., M. Mullins, M. C. U. Cheang, S. Leung, D. Voduc, T. Vickery, S. Davies, C. Fauron, X. He, Z. Hu, J. F. Quackenbush, I. J. Stijleman, J. Palazzo, J. S. Marron, A. B. Nobel, E. Mardis, T. O. Nielsen, M. J. Ellis, C. M. Perou, and P. S. Bernard
2009\.  Supervised risk predictor of breast cancer based on intrinsic
subtypes.  Journal of clinical oncology : official journal of the American
Society of Clinical Oncology, 27(8):1160–1167.

  * Perou et al. (2000) Perou, C. M., T. Sørlie, M. B. Eisen, M. van de Rijn, S. S. Jeffrey, C. A. Rees, J. R. Pollack, D. T. Ross, H. Johnsen, L. A. Akslen, O. Fluge, A. Pergamenschikov, C. Williams, S. X. Zhu, P. E. Lønning, A. L. Børresen-Dale, P. O. Brown, and D. Botstein
2000\.  Molecular portraits of human breast tumours.  Nature,
406(6797):747–752.

  * Rasmussen (2000) Rasmussen, C. E.
2000\.  The Infinite Gaussian Mixture Model.  In Advances in Neural
Information Processing Systems, Volume 12, Pp. 554–560. MIT Press.

  * Richardson and Green (1997) Richardson, S. and P. J. Green
1997\.  On Bayesian Analysis of Mixtures with an Unknown Number of Components
(with discussion).  Journal Of The Royal Statistical Society Series
B-Statistical Methodology, 59(4):731–792.

  * Savage et al. (2013) Savage, R. S., Z. Ghahramani, J. E. Griffin, P. Kirk, and D. L. Wild
2013\.  Identifying cancer subtypes in glioblastoma by combining genomic,
transcriptomic and epigenomic data.  arXiv.org.

  * Shahbaba and Neal (2009) Shahbaba, B. and R. Neal
2009\.  Nonlinear models using Dirichlet process mixtures.  Journal of Machine
Learning Research (JMLR), 10:1829–1850.

  * Shen et al. (2009) Shen, R., A. B. Olshen, and M. Ladanyi
2009\.  Integrative clustering of multiple genomic data types using a joint
latent variable model with application to breast and lung cancer subtype
analysis.  Bioinformatics, 25(22):2906–2912.

  * Son et al. (2005) Son, C. G., S. Bilke, S. Davis, B. T. Greer, J. S. Wei, C. C. Whiteford, Q.-R. Chen, N. Cenacchi, and J. Khan
2005\.  Database of mRNA gene expression profiles of multiple human organs.
Genome research, 15(3):443–450.

  * Sørlie et al. (2001) Sørlie, T., C. M. Perou, R. Tibshirani, T. Aas, S. Geisler, H. Johnsen, T. Hastie, M. B. Eisen, M. van de Rijn, S. S. Jeffrey, T. Thorsen, H. Quist, J. C. Matese, P. O. Brown, D. Botstein, P. E. Lønning, and A. L. Børresen-Dale
2001\.  Gene expression patterns of breast carcinomas distinguish tumor
subclasses with clinical implications.  Proceedings Of The National Academy Of
Sciences Of The United States Of America, 98(19):10869–10874.

  * Suchard et al. (2010) Suchard, M. A., Q. Wang, C. Chan, J. Frelinger, A. Cron, and M. West
2010\.  Understanding GPU Programming for Statistical Computation: Studies in
Massively Parallel Massive Mixtures.  J Comput Graph Stat, 19(2):419–438.

  * Tadesse et al. (2005) Tadesse, M. G., N. Sha, and M. Vannucci
2005\.  Bayesian Variable Selection in Clustering High-Dimensional Data.
Journal Of The American Statistical Association, 100(470):602–617.

  * Taschler et al. (2019) Taschler, B., F. Dondelinger, and S. Mukherjee
2019\.  Model-based clustering in very high dimensions via adaptive
projections.  arXiv.org.

  * The Cancer Genome Atlas (2012) The Cancer Genome Atlas
2012\.  Comprehensive molecular portraits of human breast tumours.  Nature,
487(7407):61–70.

  * Walker (2007) Walker, S. G.
2007\.  Sampling the Dirichlet mixture model with slices.  Communications In
Statistics-Simulation And Computation, 36(1):45–54.

  * Witten and Tibshirani (2010) Witten, D. M. and R. Tibshirani
2010\.  A framework for feature selection in clustering.  Journal of the
American Statistical Association, 105(490):713–726.  PMID: 20811510.

  * Yeung and Ruzzo (2001) Yeung, K. Y. and W. L. Ruzzo
2001\.  Principal component analysis for clustering gene expression data.
Bioinformatics (Oxford, England), 17(9):763–774.