Re: [stan-users] Rstan gaussian mixture model 'Scale parameter is 0, but must be > 0'

[Mike Lawrence]

"If the number in the 'count' column is small, do not ask about this message on stan-users."

Now, regarding slow sampling, you should tell us how long it takes, what iteration value you specified, and what the final effective sample size estimates are.

On Friday, 21 October 2016, Shawn Lee <ls199...@gmail.com> wrote:

At first I generate synthetic data as follows:

K=10

N=5000

theta=runif(K,0,1)

theta=theta/sum(theta)

mu=runif(K,0,1)

sigma=runif(K,0,2)

prob=function(y,theta,mu,sigma)

{

  p=c(1:K)

  for(k1 in c(1:K))

  {

    p[k1]=theta[k1]*dnorm(y,mu[k1],sqrt(sigma[k1]))

  }

  return(sum(p))

}

y=rep(NA,N)  # generate y based on rejection-method

n1=1

while(n1<(N+1))

{

  temp1=runif(1,-5,5)

  p=prob(temp1,theta,mu,sigma)

  temp2=runif(1,0,0.2)

  if(temp2<p)

  {

    y[n1]=temp1

    n1=n1+1

  }

}

Then build the model: 

data{

int<lower=1> K;

int<lower=1> N;

real y[N];

}

parameters{

simplex[K] theta;

real mu[K];

real<lower=0> sigma[K];

}

model {

real ps[K];

sigma ~ cauchy(0,0.5);

mu~normal(0,1);

for(n in 1:N){

for(k in 1:K){

ps[k]=log(theta[k])+normal_lpdf(y[n]|mu[k],sigma[k]);

}

target+=log_sum_exp(ps);

}

}

The sampling is very slow and some errors appears:

SAMPLING FOR MODEL 'gmm' NOW (CHAIN 1).

WARNING: The initial buffer, adaptation window, and terminal buffer

         overflow the total number of warmup iterations.

         Defaulting to a 15%/75%/10% partition,

           init_buffer = 3

           adapt_window = 20

           term_buffer = 2

Chain 1, Iteration:  1 / 50 [  2%]  (Warmup)

Chain 1, Iteration:  5 / 50 [ 10%]  (Warmup)

Chain 1, Iteration: 10 / 50 [ 20%]  (Warmup)

Chain 1, Iteration: 15 / 50 [ 30%]  (Warmup)

Chain 1, Iteration: 20 / 50 [ 40%]  (Warmup)

Chain 1, Iteration: 25 / 50 [ 50%]  (Warmup)

Chain 1, Iteration: 26 / 50 [ 52%]  (Sampling)

Chain 1, Iteration: 30 / 50 [ 60%]  (Sampling)

Chain 1, Iteration: 35 / 50 [ 70%]  (Sampling)

Chain 1, Iteration: 40 / 50 [ 80%]  (Sampling)

Chain 1, Iteration: 45 / 50 [ 90%]  (Sampling)

Chain 1, Iteration: 50 / 50 [100%]  (Sampling)

 Elapsed Time: 25.5386 seconds (Warm-up)

               28.9616 seconds (Sampling)

               54.5001 seconds (Total)

The following numerical problems occured the indicated number of times after warmup on chain 1

                                                                                count

Exception thrown at line 19: normal_log: Scale parameter is 0, but must be > 0!     3

When a numerical problem occurs, the Hamiltonian proposal gets rejected.

See http://mc-stan.org/misc/warnings.html#exception-hamiltonian-proposal-rejected

If the number in the 'count' column is small, do not ask about this message on stan-users.

Could anyone help me with this problem? I think there might be something wrong with the prior. 

Thank you.

--
You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+unsubscribe@googlegroups.com.
To post to this group, send email to stan-...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

--

--
Mike Lawrence
Graduate Student
Department of Psychology & Neuroscience
Dalhousie University

~ Certainty is (possibly) folly ~

[Bob Carpenter]

sigma is occassionally hitting 0 due to underflow in the transforms.

There's a problem with mixture models in general that they don't have
maximum likelihoods. If you put one mixture component be responsible
for one data item, as the scale (std dev param) goes to zero and the location
(mean param) goes to the data item's value, the likelihood grows.

One way to avoid this is to force the components to have more responsibility.
You can do this with a prior on the simplex. Something like Dirichlet(alpha)
where alpha is a vector values >> 1 can help. Also putting a prior on the
scales that explicitly avoids zero (like a lognormal) can help.

- Bob

> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.

> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.
>
>
> --
>
> --
> Mike Lawrence
> Graduate Student
> Department of Psychology & Neuroscience
> Dalhousie University
>
> ~ Certainty is (possibly) folly ~
>

> --
> You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.

> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.

[Shawn Lee]

Thank you for your response, I have solved the 0 problem.

After sampling, I find all the elements of theta equals to 0.1(I defined there are 10 components). 

Is it because the label switching problem? 

在 2016年10月22日星期六 UTC+8上午2:47:31，Bob Carpenter写道：

[Bob Carpenter]

Maybe. Try running a single chain and see if you have
a problem. Or look at the traceplots of multiple parameters
and see if you can detect the label switching that way.

- Bob

[Stephen Martin]

Just to add to this: Yes, it's probably due to label switching. You're marginalizing over multiple modes that are equivalent in inference, but labels are switched, so everything winds up being the same when averaged.
This is *not* an issue if you have between-chain label switching, meaning that chain 1 found and explored 1 mode, and chain 2 found and explored another mode, but neither chain found the other modes.
If this is the case, your trace plots will tell you pretty quickly --- The traces of the individual chains will just be flipped from one another. In this case, just extract out individual chains or chains that converged to similar modes and look at the posterior summaries.

However, if your chains are switching modes during simulation (chain 1 explores modes 1, 2, 3, ...), you have an issue that I have no idea how to solve from MCMC output itself. I've only stopped this by applying some identifiability constraints that either creates one mode, or pushes the modes further apart so any one given chain is unlikely to 'find' the other one throughout its simulation. It's not ideal, but it works well enough.

[Bob Carpenter]

There is no way to stop this if the model allows mixing between
the modes/labels. What you need to do is define inferences that
aren't sensitive to the labels. These will go through.

We still worry about convergence running a single model, but we
really don't have great advice (nor does anyone else) about this
label switching issue. It gets in the way of our trying to evaluate
converegence more than anything else.

- Bob

[Stephen Martin]

Just curious, what sort of inferences aren't sensitive to the labels?
Say you're wanting to fit a latent group model from data. You want to
know what probability each person has of belonging to each of K groups,
then the parameters of said K groups. Is there a way of obtaining such
parameters in a way that doesn't have a label switching problem?

[Bob Carpenter]

On Nov 1, 2016, at 8:47 PM, Stephen Martin <hwki...@gmail.com> wrote:
>
> Just curious, what sort of inferences aren't sensitive to the labels?

You can do prediction, i.e., likelihood for a new observation x', where
you estimate p(y' | x', theta) where theta is all the parameters, x is
predictors and y is observed data. It's the same likelihood calculation
as in the model.

You can do similarity, i.e., Pr[items n and n' belong to same group | x, y].

> Say you're wanting to fit a latent group model from data. You want to
> know what probability each person has of belonging to each of K groups,
> then the parameters of said K groups.

Exactly because of label switching this isn't a valid inference.

> Is there a way of obtaining such
> parameters in a way that doesn't have a label switching problem?

No. You can't evaluate what the probability is that item n
belongs to mixture component k. You can get the marginal for a parameter,
such as mu, given item i. It'll be multimodal.

Another thing I should put in the manual!

- Bob

[Stephen Martin]

In my experience, you *can* do this so long as you specify enough constraints to avoid within-chain label switching.
If each chain settles on one mode without swapping modes, extracting similar chains or unimodal chains seems to work just fine.
I am working on a publication that uses mixtures extensively, using stan, and it appears to work fine on simulated data. A few constraints (valid constraints) keeps it either unimodal or the modes far enough apart that the chains never swap modes during simulation. Some parameters are invariant to label switching, but the latent group membership component of the model, of course, is not invariant and requires the constraint.

If you're going to update the manual, I'd at least be fair to the issue.
The full marginal posterior may not be interpretable if labels switch, but the type of label switching makes a big difference.
If each chain explores separate modes but never switch during simulation, this is ok for inference on EACH of the chains. You can examine the convergence of similar chains, or manually compute Rhat on the MCMC output using the correct columns --- If chain 1 and chain 2 label-switched, you can compute Rhat on, e.g., theta[1] for chain 1 and theta[2] for chain 2; if those converged, you're golden, just make sure you do your posterior inference on either one of the chains or adjust the MCMC output manually and assess convergence across both.

If chains label-switched during simulation, you're a bit screwed, but constraints can avoid this problem.

[Stephen Martin]

As a followup, I thought I'd mention the constraints that helped me. One is not preferred, but more generally applicable; the other is preferred, but only applicable if you're adding a certain prediction to the model.

1) Ordering constraints; this has been mentioned extensively. Place an ordering constraint on, e.g., gaussian means if that's the more important aspect of the different gaussian distributions to you. Yes, this can lead to bias in some cases, and can create a 'trough' between modes that a chain can get stuck in. If a chain gets stuck in the trough, you can usually tell it in the trace plots (chains 1 and 2 are label-swapped, but chain 3 is exactly between them bouncing endlessly between them). In this case, I've found that discarding chain 3 is 'fine', as it has nothing to do with my inference, per se, and more to do with a sampling problem --- The expectactions indicated by chain 3 are only applicable if I think label switching occurs in the population, which makes no intuitive sense to me, so it's out.

2) If you're predicting latent group membership via, e.g., logistic or softmax regression, impose a positivity constraint on one of the predictors for one of the groups. In the case of softmax, I recommend keeping group 1's log-odds at 0, to set the scale, then group 2's first predictor to be positive. That also seems to work; it works better if the first predictor is actually useful though, of course (otherwise it'll bunch up around 0); if need be, you could run it without the constraint, see which predictors appear to be multimodal (+2, or -2), then set the constraint so that it will find the mode with the beta = +2.

[Michael Betancourt]

I want to point out that this conversation has been focused on degenerate modes

only, where label switching is a consequence of non-identified latent labels.  But

once you add a measurement model and data in complex mixture models you often

get true multi modality in the posterior, where the modes are not equivalent (LDA

being a classic example of this).  In this case none of the strategies discussed are 

valid and one has to adopt a more sophisticated approach.  In particular, in practice 

it is critical to verify that any modes being considered are actually degenerate.

[Bob Carpenter]

> On Nov 2, 2016, at 12:17 PM, Stephen Martin <hwki...@gmail.com> wrote:
>
> In my experience, you *can* do this so long as you specify enough constraints to avoid within-chain label switching.

Sure, you can pick a mode and go with it. I'm just
saying that it's not a Bayesian predictive inference
because the identity of the mixture component isn't
identified.

> If each chain settles on one mode without swapping modes, extracting similar chains or unimodal chains seems to work just fine.
> I am working on a publication that uses mixtures extensively, using stan, and it appears to work fine on simulated data. A few constraints (valid constraints) keeps it either unimodal or the modes far enough apart that the chains never swap modes during simulation. Some parameters are invariant to label switching, but the latent group membership component of the model, of course, is not invariant and requires the constraint.

With constaints, such as orderings, you can sometimes
identify the model and make these inferences because
now there's something to break the symmetry on the
label switching. You may still miss a lot of probability
mass, though if you get stuck in a mode that isn't representative
of the full posterior.

> If you're going to update the manual, I'd at least be fair to the issue.

Will do.

- Bob