A Mixture/Clustering problem

[Giovanni Petris]

A big thank you to all the development and support team of NIMBLE for the great job!!

I am trying to implement an MCMC sampler for a mixture-of-Normals model, but I am getting some errors that I don't quite understand.

The mixture parameters are all known in my case (they come from a mixture approximation to the log-chiSquare_1  which is routinely used in Stochastic Volatility models), so the only unknown 'parameters' are the indicator latent variables.  The following is the error message I am getting. Code is attached, if anybody is kind enough to try and help me!

Thank you in advance!

Giovanni

> ####

> #### 10-component mixture of normals approximation to log chi^2_1

> ####

> w0 <- c(0.00609, 0.04775, 0.13057, 0.20674, 0.22715,

+         0.18842, 0.12047, 0.05591, 0.01575, 0.00115)

> mn0 <- c(1.92677, 1.34744, 0.73504, 0.02266, -0.85173,

+          -1.97278, -3.46788, -5.55246, -8.68384, -14.65000)

> var0 <- c(0.11265, 0.17788, 0.26768, 0.40611, 0.62699,

+           0.98583, 1.57469, 2.54498, 4.16591, 7.33342)

> 

> ####

> #### Simulate log Chi^2_1 data

> ####

> set.seed(62)

> n <- 10L

> y <- log(rnorm(n)^2)

> 

> ####

> #### Set up NIMBLE model

> ####

> nimbleOptions(allowDynamicIndexing = TRUE)

> logChi2Code <- nimbleCode({

+     for (i in 1:n) {

+         Ind[i] ~ dcat(w0[1:10])

+         y[i] ~ dnorm(mn0[Ind[i]], var = var0[Ind[i]])

+     }

+ })

> 

> logChi2Model <- nimbleModel(logChi2Code, constants = list(n = n, w0 = w0, mn0 = mn0, var0 = var0),

+                             data = list(y = y), inits = list(Ind = rep(5, n)))

defining model...

Adding w0 as data for building model.

building model...

setting data and initial values...

running calculate on model (any error reports that follow may simply reflect missing values in model variables) ... Error in envRefInferField(x, what, getClass(class(x)), selfEnv) : 

  ‘var0’ is not a valid field or method name for reference class “code_modelClass_UID_355_UID_356”

checking model sizes and dimensions...Error in envRefInferField(x, what, getClass(class(x)), selfEnv) : 

  ‘mn0’ is not a valid field or method name for reference class “code_modelClass_UID_355_UID_356”

Error in envRefInferField(x, what, getClass(class(x)), selfEnv) : 

  ‘var0’ is not a valid field or method name for reference class “code_modelClass_UID_355_UID_356”

In addition: Warning message:

In model$checkBasics() :

  Unable to calculate parameter 'mean'; this may simply reflect that there are missing values in model variables.

 note that missing values (NAs) or non-finite values were found in model variables: lifted_sqrt_oPvar0_oBInd_oBi_cB_cB_cP_L3. This is not an error, but some or all variables may need to be initialized for certain algorithms to operate properly.

model building finished.

Warning message:

In model$checkBasics() :

  Unable to calculate parameter 'var'; this may simply reflect that there are missing values in model variables.

> 

[Christopher Paciorek]

Hi Giovanni,

The problem is that we take anything that is a 'constant' and plug in the actual numerical values in the model. So 'var0' (and also 'mn0') no longer exist in the model as variables once this substitution happens. This is fine when 'Ind' is also provided as constants but in your case NIMBLE is trying to use a dynamic index to index a constant, which won't work. We overlooked this possibility and hence don't provide a useful error message here.

All you need to do is the following. First, provide 'var0' and 'mn0' via 'inits', so that we treat 'var0' and 'mn0' as variables in the model. Second, we need some help here with the dimensions since there is no easy way for us to know that mn0 and var0 are vectors of length 10 (this is because we don't require initial values for model variables when a model is created and because mu0 and var0 are not created via a deterministic or stochastic assignment). The use of the 'dimensions' argument allows NIMBLE to know what size these two variables are.

logChi2Model <- nimbleModel(logChi2Code, constants = list(n = n, w0 = w0), # , mn0 = mn0, var0 = var0),

                            data = list(y = y), inits = list(var0 = var0, mn0 = mn0, Ind = rep(5, n)),

                            dimensions = list(mn0 = 10, var0 = 10))

-Chris

[Giovanni Petris]

Hi Chris,

Thanks a lot: it works perfectly now, with the changes you suggested.

Best,

Giovanni