Infinite Dirichlet mixture models learning via expectation propagation

Abstract

In this article, we propose a novel Bayesian nonparametric clustering algorithm based on a Dirichlet process mixture of Dirichlet distributions which have been shown to be very flexible for modeling proportional data. The idea is to let the number of mixture components increases as new data to cluster arrive in such a manner that the model selection problem (i.e. determination of the number of clusters) can be answered without recourse to classic selection criteria. Thus, the proposed model can be considered as an infinite Dirichlet mixture model. An expectation propagation inference framework is developed to learn this model by obtaining a full posterior distribution on its parameters. Within this learning framework, the model complexity and all the involved parameters are evaluated simultaneously. To show the practical relevance and efficiency of our model, we perform a detailed analysis using extensive simulations based on both synthetic and real data. In particular, real data are generated from three challenging applications namely images categorization, anomaly intrusion detection and videos summarization.

The calculation of \(Z_i\) in Eq. (17)

The normalized constant $Z_i$ in Eq. (17) can be calculated as

$$\begin{aligned} Z_i \!=\! \int f_i(\varTheta )q^{\setminus i}(\varTheta )d\varTheta \!=\!\sum _{j=1}^M \bar{\lambda }_j\prod _{s=1}^{j-1}(1-\bar{\lambda }_s) \!\int \! \mathrm{Dir }\left( \mathbf X _i|{\varvec{\alpha }}_j\right) N\left( {\varvec{\alpha }}_j|{\varvec{\mu }}_j^{\setminus i},A_j^{\setminus i}\right) \mathrm{d }{\varvec{\alpha }}_j\nonumber \\ \end{aligned}$$

(27)

where $\bar{\lambda }_j$ is the expected value of $\lambda _j$. Since the integration involved in Eq. (27) is analytically intractable, we tackle this problem by adopting the Laplace approximation to approximate the integrand with a Gaussian distribution as suggested in Ma and Leijon (2010).

First, we define $h({\varvec{\alpha }}_j)$ as the integrand in Eq. (27):

$$\begin{aligned} h({\varvec{\alpha }}_j) =\mathrm{Dir }(\mathbf X _i|{\varvec{\alpha }}_j)\mathcal{N }\left( {\varvec{\alpha }}_j|\varvec{\mu }_j^{\setminus i},A^{\setminus i}_{j}\right) \end{aligned}$$

(28)

Then, the normalized distribution for this integrand which is indeed a product of a Dirichlet distribution and a Gaussian distribution is given by

$$\begin{aligned} \mathcal{H }({\varvec{\alpha }}_j) =\frac{h({\varvec{\alpha }}_j)}{\int h({\varvec{\alpha }}_j)d{\varvec{\alpha }}_j} \end{aligned}$$

(29)

Our goal for the Laplace method the goal is to find a Gaussian approximation which is centered on the mode of the distribution $\mathcal{H }({\varvec{\alpha }}_j)$. We may obtain the mode ${\varvec{\alpha }}_j^*$ numerically by setting the first derivative of $\ln h({\varvec{\alpha }}_j)$ to 0, where $\ln h({\varvec{\alpha }}_j)$ can be calculated by

$$\begin{aligned} \ln h({\varvec{\alpha }}_j)&= \ln \frac{\sum _{l=1}^D\alpha _{jl}}{\prod _{l=1}^D\varGamma (\alpha _{jl})}\nonumber \\&+ \sum _{l=1}^D(\alpha _{jl}-1)\ln X_{il}\!-\! \frac{1}{2}\left( {\varvec{\alpha }}_j \!-\! {\varvec{\mu }}^{\setminus i}_{j}\right) ^T A^{\setminus i}_j \left( {\varvec{\alpha }}_j \!-\! {\varvec{\mu }}^{\setminus i}_{j}\right) \!+\!\text{ const. }\nonumber \\ \end{aligned}$$

(30)

Subsequently, we can calculate the first and second derivatives with respect to ${\varvec{\alpha }}_j$ as

$$\begin{aligned} \frac{\partial \ln h({\varvec{\alpha }}_j)}{\partial {\varvec{\alpha }}_j} =\left[ \begin{array}{c} \varPsi \left( \displaystyle \sum \limits _{l=1}^D\alpha _{jl}\right) - \varPsi (\alpha _{j1}) + \ln X_{i1}\\ \vdots \\ \varPsi \left( \displaystyle \sum \limits _{l=1}^D\alpha _{jl}\right) - \varPsi (\alpha _{jD}) + \ln X_{iD} \end{array}\right] -A_j^{\setminus i}\left( {\varvec{\alpha }}_j- {\varvec{\mu }}^{\setminus i}_j\right) \end{aligned}$$

(31)

and

$$\begin{aligned} \frac{\partial ^2\ln h({\varvec{\alpha }}_j)}{\partial {\varvec{\alpha }}_j^2} = \left[ \begin{array}{c@{\quad }c@{\quad }c} \varPsi '\left( \displaystyle \sum \limits _{l=1}^D\alpha _{jl}\right) - \varPsi '(\alpha _{j1}) &{} \cdots &{} \varPsi '\left( \displaystyle \sum \limits _{l=1}^D\alpha _{jl}\right) \\ \vdots &{} \ddots &{}\vdots \\ \varPsi '\left( \displaystyle \sum \limits _{l=1}^D\alpha _{jl}\right) &{} \cdots &{}\varPsi '\left( \displaystyle \sum \limits _{l=1}^D\alpha _{jl}\right) - \varPsi '(\alpha _{jD}) \end{array}\right] -A^{\setminus i}_{j}\nonumber \\ \end{aligned}$$

(32)

where $\varPsi (\cdot )$ is the digamma function. Then, we can approximate $h({\varvec{\alpha }}_j)$ using the obtained mode as

$$\begin{aligned} h({\varvec{\alpha }}_j)\simeq h({\varvec{\alpha }}_j^*)\exp \bigg (-\frac{1}{2}({\varvec{\alpha }}_j- {\varvec{\alpha }}_j^*)\widehat{A}_{j}({\varvec{\alpha }}_j- {\varvec{\alpha }}_j^*)\bigg )\quad \end{aligned}$$

(33)

where the precision matrix $\widehat{A}_{j}$ is given by

$$\begin{aligned} \widehat{A}_{j} = - \left. \frac{\partial ^2\ln h({\varvec{\alpha }}_j)}{\partial {\varvec{\alpha }}_j^2} \right| _{{\varvec{\alpha }}_j ={\varvec{\alpha }}_j^*} \end{aligned}$$

(34)

Therefore, the integration of $h({\varvec{\alpha }}_j)$ can be approximated by using Eq. (33) as

$$\begin{aligned} \int h({\varvec{\alpha }}_j)d{\varvec{\alpha }}_j \simeq h({\varvec{\alpha }}_j^*)\int \exp \left( \!-\!\frac{1}{2}({\varvec{\alpha }}_j\!-\!{\varvec{\alpha }}_j^*)\widehat{A}_{j}({\varvec{\alpha }}_j\!-v{\varvec{\alpha }}_j^*)\right) d{\varvec{\alpha }}_j\!=\! h({\varvec{\alpha }}_j^*) \frac{(2\pi )^{D/2}}{|\widehat{A}_j|^{1/2}}\nonumber \\ \end{aligned}$$

(35)

Finally, we can rewrite Eq. (27) as following:

$$\begin{aligned} Z_i=\sum _{j=1}^M \bar{\lambda }_j\prod _{s=1}^{j-1}(1-\bar{\lambda }_s)h\left( {\varvec{\alpha }}_j^*\right) \frac{(2\pi )^{D/2}}{|\widehat{A}_j|^{1/2}} \end{aligned}$$

(36)