Splitting modes (gmix_kurt.m)

In a method proposed by N. Vlassis and A. Likas [62], the number of modes in a Gaussian mixture is determined by monitoring the weighted kurtosis for each mode. Putting their equation for one-dimensional ${\bf z}$ in our notation, Vlassis et al define

$\displaystyle \kappa_i = \frac{
\sum_{n=1}^N w_{n,i} \left(
\frac{ {\bf z}_n ...
...\sqrt{\mbox{\boldmath$\Sigma$}_i} }
\right)^4
}{
\sum_{n=1}^N w_{n,i}
} - 3
$

where

$\displaystyle w_{n,i} = \frac{
{\cal N}({\bf z}_n,\mbox{\boldmath$\mu$}_{i},\m...
...^N {\cal N}({\bf z}_n,\mbox{\boldmath$\mu$}_{i},\mbox{\boldmath$\Sigma$}_i)
}
$

If $\vert\kappa_i\vert$ is too high for any mode $i$, they split the mode into two. We modify this for higher dimension and use the skew in addition to the kurtosis. Extending to higher dimension is done by projecting each data sample ${\bf z}_n$ onto the $j$-th principal axis of $\Sigma$$_i$ in turn. Let $z_{n,i}^j \stackrel{\mbox{\tiny $\Delta$}}{=}({\bf z}_n-\mbox{\boldmath $\mu$}_i)^\prime {\bf v}_{ij}$ where ${\bf v}_{ij}$ is the $j$-th column of ${\bf V}$, obtained from the SVD of $\Sigma$$_i$ (see discussion in section 13.2.5). Thus, for each $j$,
  1. Let

    $\displaystyle \kappa_{i,j} = \frac{
\sum_{n=1}^N w_{n,i} \left(
\frac{ z_{n,i}^j } {s_i}
\right)^4 }
{ \sum_{n=1}^N w_{n,i} } - 3
$

  2. Let

    $\displaystyle \psi_{i,j} = \frac{
\sum_{n=1}^N w_{n,i} \left(
\frac{ z_{n,i}^j } {s_i}
\right)^3
}
{ \sum_{n=1}^N w_{n,i} }
$

  3. Let

    $\displaystyle m_{i,j} = \vert\kappa_{i,j}\vert + \vert\psi_{i,j}\vert
$

where

$\displaystyle s_i^2 = \frac{
\sum_{n=1}^N w_{n,i} \left(
z_{n,i}^j
\right)^2
}
{ \sum_{n=1}^N w_{n,i} }
$

Now, if $m_{i,j} > \tau$, for any $j$, split mode $i$. Split the mode by creating modes at

   $\displaystyle \mbox{\boldmath$\mu$}$$\displaystyle =$   $\displaystyle \mbox{\boldmath$\mu$}$$\displaystyle _i + {\bf v}_{ij} S_{i,j}
$

and

   $\displaystyle \mbox{\boldmath$\mu$}$$\displaystyle =$   $\displaystyle \mbox{\boldmath$\mu$}$$\displaystyle _i - {\bf v}_{ij} S_{i,j}
$

where $S_{i,j}$ is the $j$-th singular value of $\Sigma$$_i$. The same covariance $\Sigma$$_i$ is used for each new mode. Of course, the decision of whether to split or not depends on the mixing proportion $\alpha_i$ as well. No splitting occurs if $\alpha_i$ is too small.

In the following example, we create data with a gap in it. We begin iterating with a single mode. The kurtosis/skew algorithm above is able to assign modes until it is finally happy after 8 modes (Figure 13.5).

Figure 13.5: Results of bottom up PDF estimation. One mode (left), two modes (center), and after convergence at 8 modes (right).
\includegraphics[width=2.1in,height=2.1in, clip]{ku1.eps} \includegraphics[width=2.1in,height=2.1in, clip]{ku2.eps} \includegraphics[width=2.1in,height=2.1in, clip]{ku8.eps}
The calling syntax for software/gmix_kurt.m is
	gparm = gmix_kurt(gparm,x,[kurt_thresh],[debug]);
The optional threshold parameter (default=1.0) allows control over splitting. A higher threshold is less likely to split. The optional debug parameter, if set to 1, will print out kurtosis and skew information.