Two different ways to write a mixture model. Question about pm.NormalMixture


#1

Hello!

I am rewrite a following mixture model:

model = pm.Model()
with model:
    # cluster sizes
    p = pm.Dirichlet('p', a=np.array([1., 1., 1.]), shape=k)
    # ensure all clusters have some points
    p_min_potential = pm.Potential('p_min_potential', tt.switch(tt.min(p) < .1, -np.inf, 0))


    # cluster centers
    means = pm.Normal('means', mu=[0, 0, 0], sd=15, shape=k)
    # break symmetry
    order_means_potential = pm.Potential('order_means_potential',
                                         tt.switch(means[1]-means[0] < 0, -np.inf, 0)
                                         + tt.switch(means[2]-means[1] < 0, -np.inf, 0))

    # measurement error
    sd = pm.Uniform('sd', lower=0, upper=20)

    # latent cluster of each observation
    category = pm.Categorical('category',
                              p=p,
                              shape=ndata)

    # likelihood for each observed value
    points = pm.Normal('obs',
                       mu=means[category],
                       sd=sd,
                       observed=data)

I want to rewrite it to a version that doesn’t involve categorical random variable, i.e. using normal mixture.

with pm.Model() as model2:
    # cluster sizes
    p = pm.Dirichlet('p', a=np.array([1., 1., 1.]), shape=k)

    # cluster centers
    means = pm.Normal('means', mu=[0, 0, 0], sd=15, shape=k)

    # measurement error
    sd = pm.Uniform('sd', lower=0, upper=20)

    points = pm.NormalMixture('mix', w=p, mu=means, sd=sd, observed=data, comp_shape=k, shape=ndata)

I keep getting following errors

static_cast<int>( ). 6 errors generated.. ", '[Elemwise{sub,no_inplace}(<TensorType(float64, col)>, <TensorType(float64, row)>)]')

By the way, the data used are:

import numpy as np
import pandas as pd
import matplotlib
matplotlib.use('TkAgg')
import matplotlib.pyplot as plt
import seaborn as sns


sns.set_context('paper')
sns.set_style('darkgrid')

import pymc3 as pm
import theano.tensor as tt


# simulate data
np.random.seed(12345)
k = 3
ndata = 500
spread = 5
centers = np.array([-spread, 0, spread])

# simulate data from mixture distribution
v = np.random.randint(0, k, ndata)
data = centers[v] + np.random.randn(ndata)

#2

The set up to NormalMixture is not quite correct, you should do:

with pm.Model() as model2:
    # cluster sizes
    p = pm.Dirichlet('p', a=np.array([1., 1., 1.]), shape=k)
    # cluster centers
    means = pm.Normal('means', mu=[0, 0, 0], sd=15, shape=k)
    # measurement error
    sd = pm.Uniform('sd', lower=0, upper=20)
    points = pm.NormalMixture('mix', w=p, mu=means, sd=sd, observed=data)
    trace = pm.sample()

(The comp_shape is more for multi-dimensional mixture).

Also, to help model convergence (mixture model are difficult to sample due to label switching!), you should try something like:

with pm.Model() as model2:
    # cluster sizes
    p = pm.Dirichlet('p', a=np.array([1., 1., 1.]), shape=k)
    # cluster centers
    means = pm.Normal('means', mu=[0, 0, 0], sd=15, shape=k, 
                      transform=pm.distributions.transforms.ordered,
                      testval=np.asarray([0., 1., 2.]))
    # measurement error
    sd = pm.Uniform('sd', lower=0, upper=20)
    points = pm.NormalMixture('mix', w=p, mu=means, sd=sd, observed=data)
    trace = pm.sample()

#3

Thanks for your reply @junpenglao

Unfortunately, the two snippets in your answer don’t compile…
Gives a similar error that I originally posted.


#4

Are you on the newest pymc3 release? if so try to upgrade to master branch…


#5

So the issue was that I had to put

import theano
theano.config.gcc.cxxflags = "-Wno-c++11-narrowing"

With above, your code compiles. Would upgrading to master branch solve that problem?
FWIW, I am using pymc3==3.6 using conda.


#6

hmmmm usually conda should handle these quite effortlessly - do you see any warning when you import theano?


#7

Yes. I see following error

In [1]: import theano
WARNING (theano.configdefaults): install mkl with `conda install mkl-service`: No module named 'mkl'

#8

And installing mkl-service doesn’t resolve the problem. I still need to add

import theano
theano.config.gcc.cxxflags = "-Wno-c++11-narrowing"

in order for the code to work :frowning:


#9

Regrettably, once theano has been installed, installing mkl-service later does not cause theano to use it. You have to uninstall theano, install mkl-service and then reinstall theano to make it work. If you are using Windows, I recommend that you install run conda install mkl-service libpython m2w64-toolchain, before installing theano (or just use the entire command they suggest, if you want to start on a fresh environment).