Show Sidebar Hide Sidebar

Concentration Prior Type Analysis of Variation Bayesian Gaussian Mixture in Scikit-learn

This example plots the ellipsoids obtained from a toy dataset (mixture of three Gaussians) fitted by the BayesianGaussianMixture class models with a Dirichlet distribution prior (weight_concentration_prior_type='dirichlet_distribution') and a Dirichlet process prior (weight_concentration_prior_type='dirichlet_process'). On each figure, we plot the results for three different values of the weight concentration prior.

The BayesianGaussianMixture class can adapt its number of mixture componentsautomatically. The parameter weight_concentration_prior has a direct link with the resulting number of components with non-zero weights.

Specifying a low value for the concentration prior will make the model put most of the weight on few components set the remaining components weights very close to zero. High values of the concentration prior will allow a larger number of components to be active in the mixture.

The Dirichlet process prior allows to define an infinite number of components and automatically selects the correct number of components: it activates a component only if it is necessary. On the contrary the classical finite mixture model with a Dirichlet distribution prior will favor more uniformly weighted components and therefore tends to divide natural clusters into unnecessary sub-components.

New to Plotly?

Plotly's Python library is free and open source! Get started by downloading the client and reading the primer.
You can set up Plotly to work in online or offline mode, or in jupyter notebooks.
We also have a quick-reference cheatsheet (new!) to help you get started!

Version

In [1]:
import sklearn
sklearn.__version__
Out[1]:
'0.18.1'

Imports

This tutorial imports BayesianGaussianMixture.

In [2]:
import plotly.plotly as py
import plotly.graph_objs as go
from plotly import tools

import math
import numpy as np
from sklearn.mixture import BayesianGaussianMixture

Calculations

In [3]:
def plot_ellipses(weights, means, covars):
    data = []
    for n in range(means.shape[0]):
        eig_vals, eig_vecs = np.linalg.eigh(covars[n])
        unit_eig_vec = eig_vecs[0] / np.linalg.norm(eig_vecs[0])
        angle = np.arctan2(unit_eig_vec[1], unit_eig_vec[0])
        
        # eigenvector normalization
        eig_vals = 2 * np.sqrt(2) * np.sqrt(eig_vals)
        a =  eig_vals[0]
        b =  eig_vals[1]
        x_origin = means[n][0]
        y_origin = means[n][1]
        x_ = [ ]
        y_ = [ ]
        
        for t in range(0, 360,10):
            x = a*(math.cos(math.radians(t))) + x_origin
            x_.append(x)
            y = b*(math.sin(math.radians(t))) + y_origin
            y_.append(y)
        
        
        elle = go.Scatter(x=x_ , y=y_, mode='lines',
                          
                          line=dict(color='#56B4E9', width=1))
        data.append(elle)
    
    return data


def plot_results(estimator, X, y):
    data = [[], []]
    trace = go.Scatter(x=X[:, 0], y=X[:, 1],
                       mode='markers',
                       marker=dict(color=colors[y]))
    data[0].append(trace)
    k = plot_ellipses(estimator.weights_, estimator.means_,
                      estimator.covariances_)
    
    for i in range(0, len(k)):
        data[0].append(k[i])
    
    for k, w in enumerate(estimator.weights_):
       
        trace = go.Bar(x=[k - .45], y=[w], 
                       marker=dict(color='#56B4E9', 
                                   line=dict(color='black', width=1)))
        data[1].append(trace)
    
    return data
In [4]:
random_state, n_components, n_features = 2, 3, 2
colors = np.array(['#0072B2', '#F0E442', '#D55E00'])

covars = np.array([[[.7, .0], [.0, .1]],
                   [[.5, .0], [.0, .1]],
                   [[.5, .0], [.0, .1]]])
samples = np.array([200, 500, 200])
means = np.array([[.0, -.70],
                  [.0, .0],
                  [.0, .70]])
In [5]:
# mean_precision_prior= 0.8 to minimize the influence of the prior
estimators = [
    ("Finite mixture <br>with a Dirichlet distribution<br>prior and "
      "gamma=", BayesianGaussianMixture(
        weight_concentration_prior_type="dirichlet_distribution",
        n_components=2 * n_components, reg_covar=0, init_params='random',
        max_iter=1500, mean_precision_prior=.8,
        random_state=random_state), [0.001, 1, 1000]),
    ("Infinite mixture <br>with a Dirichlet process<br> prior and gamma=",
     BayesianGaussianMixture(
        weight_concentration_prior_type="dirichlet_process",
        n_components=2 * n_components, reg_covar=0, init_params='random',
        max_iter=1500, mean_precision_prior=.8,
        random_state=random_state), [1, 1000, 100000])]

# Generate data
rng = np.random.RandomState(random_state)
X = np.vstack([
    rng.multivariate_normal(means[j], covars[j], samples[j])
    for j in range(n_components)])
y = np.concatenate([j * np.ones(samples[j], dtype=int)
                    for j in range(n_components)])
In [6]:
# Plot results in two different figures
titles = []
data = []
i = 0
for (title, estimator, concentrations_prior) in estimators:
    
    data.append([[], []])
    for k, concentration in enumerate(concentrations_prior):
        estimator.weight_concentration_prior = concentration
        estimator.fit(X)
        
        titles.append("%s %.1e " % (title, concentration))
        k = plot_results(estimator, X, y)
        data[i][0].append(k[0])
        data[i][1].append(k[1])
        
    i+=1

Finite Mixture

In [7]:
fig = tools.make_subplots(rows=2, cols=3, 
                          subplot_titles=tuple(titles[: 3]))

for i in range(0, 2):
    for j in range(0, len(data[0][i])):
        for k in range(0, len(data[0][i][j])):
            fig.append_trace(data[0][i][j][k], i+1, j+1)

for i in map(str, range(1,6)):
    y = 'yaxis'+i
    x = 'xaxis'+i
    fig['layout'][y].update(showgrid=False)
    fig['layout'][x].update(showgrid=False)
    
fig['layout'].update(height=700, hovermode='closest',
                     showlegend=False)

fig['layout']['yaxis1'].update(title='Estimated Mixtures')
fig['layout']['yaxis4'].update(title='Weight of each component')
This is the format of your plot grid:
[ (1,1) x1,y1 ]  [ (1,2) x2,y2 ]  [ (1,3) x3,y3 ]
[ (2,1) x4,y4 ]  [ (2,2) x5,y5 ]  [ (2,3) x6,y6 ]

In [8]:
py.iplot(fig)
Out[8]:

Infinite Mixture

In [9]:
fig = tools.make_subplots(rows=2, cols=3, 
                          subplot_titles=tuple(titles[3: 6]))

for i in range(0, 2):
    for j in range(0, len(data[1][i])):
        for k in range(0, len(data[1][i][j])):
            fig.append_trace(data[1][i][j][k], i+1, j+1)

for i in map(str, range(1,6)):
    y = 'yaxis'+i
    x = 'xaxis'+i
    fig['layout'][y].update(showgrid=False)
    fig['layout'][x].update(showgrid=False)
    
fig['layout'].update(height=700, hovermode='closest',
                     showlegend=False)

fig['layout']['yaxis1'].update(title='Estimated Mixtures')
fig['layout']['yaxis4'].update(title='Weight of each component')
This is the format of your plot grid:
[ (1,1) x1,y1 ]  [ (1,2) x2,y2 ]  [ (1,3) x3,y3 ]
[ (2,1) x4,y4 ]  [ (2,2) x5,y5 ]  [ (2,3) x6,y6 ]

In [10]:
py.iplot(fig)
Out[10]:

License

Author:

    Thierry Guillemot <thierry.guillemot.work@gmail.com>

License:

    BSD 3 clause
Still need help?
Contact Us

For guaranteed 24 hour response turnarounds, upgrade to a Developer Support Plan.