Non-Identified Random Effects

[Michael Betancourt]

I’m building a relatively simple linear model with multiple group-specific
intercepts (I’ll fall back to the statistics notation and call these “random
effects”). For simplicity let’s say the model is

y_{n} ~ N(alpha0 + alpha1[group1[n]] + alpha2[group2[n]], sigma),

with no hierarchical structure.

In this case the alpha2 are all consistent with zero, but the model
cannot identify them all uniquely so instead there are extremely
strong linear correlations between the various alpha2 and alpha0,
leading to slow sampling.

I know this has come up before but I don’t remember there being
any cute workarounds. Have there been any developments since?
Anyone have a good intuition for how the distribution of the data
across the various groups can improve or worsen the non-identifiability
(for example, in this case the alpha1 are identified just fine)?

[Mick Cooney]

Do you have any idea for why the alpha2 are not identifying?

Has it something to do with the system you are modelling?

--
Mick Cooney
mickc...@gmail.com

[Krzysztof Sakrejda]

I usually drop one of the alpha2's and be done with it, the only cute 

trick I can remember comes from BUGS and uses transformed

alpha2's with their mean subtracted out.  That's still not identified

and doesn't work with HMC, do you get to add priors? Mean-centering

to de-couple it from the intercept, plus a zero-centered weak prior would

make the whole thing weakly identified.  

....if there is a trick to this I'm curious too. 

Krzysztof

[Michael Betancourt]

Yeah, I have reasonable priors that are the only thing keeping

the model reasonable.  It looks like all of the alpha2 are positive

correlated with each other but negatively correlated with alpha0

which would make some sense (and suggest dropping alpha0

altogether), but then why aren’t the alpha1 also non-identified?

Too little data to constrain the random effects onto the 

non-identified submanifold?

--
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.
To post to this group, send email to stan-...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Bob Carpenter]

This comes up in IRT (Rasch) models all the time, which
are built up from exactly that structure (with an added
multiplicative invariance).

Given the vectorized multilevel density

y ~ N(alpha0 + alpha1[group1] + alpha2[group2], sigma)

I'm assuming you have something like these hierarchical priors:

alpha0 ~ normal(0, const);

alpha1 ~ normal(0, tau1);

alpha2 ~ normal(0, tau2);

The priors being centered at 0 will already weakly identify
everything, but not as strongly as if you literally force
the alpha1 and alpha2 to be centered:

alpha1[1:(K - 1)] <- alpha1_raw;
alpha1[K] <- -sum(alpha1_raw);

I wasn't sure whether Krzysztof was suggesting that, or
another common approach, which sets

alpha1[K] <- 0
alpha2[K] <- 0

The latter will force alpha0 (albeit weakly) to be the right thing
so that alpha0 + alpha1[K] + alpha2[K] is correct. This approach is also
recommended in Gelman and Hill around the discussion of ideal point models
(a kind of IRT-like voting model). But I've found it doesn't work so
well in practice if K is large.

Gelman and Hill also recommend subtracting the mean, but as
Krzysztof points out, that produces an improper posterior, so won't
work in HMC (I doubt it actually works in BUGS for the same reason, though
it'll get close because Gibbs doesn't move very fast with correlations).

The problem I often run into if the alpha2 really all should be
roughly zero is numerical, because then you're way down in the neck
of the funnel and need small movements for alpha2 and large movements
for sigma (on the unconstrained scale), thus building in a kind of
ill-conditioning for adaptation.

- Bob

[Michael Betancourt]

I don’t have any hierarchical structure (yet), all of the intercepts
are independent and coupled only by the data which wants
almost everything to be zero. K is quite large here (~100
for both groups) so having only K - 1 intercepts and forcing
the Kth to the sum of the first (K - 1) to compensate is the
best option?

[Bob Carpenter]

Yes, I think that's the best option, but Ben and Jonah may
have better suggestions.

I'd also strongly recommend a hierarchical prior. It usually
fits better than guessing the prior scale ahead of time for these
models (I guess that'a corollary of the folk theorem).

- Bob

[Michael Betancourt]

Hierarchy comes next — I’m restricted by an existing model.

[Krzysztof Sakrejda]

Bob, is this dicsussion in the manual? Seems like it's a super common problem

and I like how clear your presentation is! Krzysztof

On Wednesday, June 15, 2016 at 2:21:53 PM UTC-4, Bob Carpenter wrote:

This comes up in IRT (Rasch) models all the time, which
are built up from exactly that structure (with an added
multiplicative invariance).

Given the vectorized multilevel density

  y ~ N(alpha0 + alpha1[group1] + alpha2[group2], sigma)

I'm assuming you have something like these hierarchical priors:

  alpha0 ~ normal(0, const);

  alpha1 ~ normal(0, tau1);

  alpha2 ~ normal(0, tau2);

The priors being centered at 0 will already weakly identify
everything, but not as strongly as if you literally force
the alpha1 and alpha2 to be centered:

  alpha1[1:(K - 1)] <- alpha1_raw;
  alpha1[K] <- -sum(alpha1_raw);

I was suggesting (potentially from Gelman and Hill) the following in transformed parameters:

alpha1 <- alpha1 - mean(alpha1)

But that's just the backwards less-functional version of what you are suggesting.  Does your

version end up with an odd definition for alpha0?  Is that fixable?

I wasn't sure whether Krzysztof was suggesting that, or
another common approach, which sets

  alpha1[K] <- 0
  alpha2[K] <- 0

Right, this is the most common approach as far as I can tell. 

 

>> 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.
>
>
> --
> 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.

[Bob Carpenter]

The problem with the improper posterior case is
is discussed in the problematic posteriors chapter.

The parameter vector centering transform is defined
(it used to be a terrible hack of mine with a simplex,
but Ben suggested the approach below and I updated the
manual).

- 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.
> >
> >
> > --
> > 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.

> > To post to this group, send email to stan-...@googlegroups.com.
> > For more options, visit https://groups.google.com/d/optout.
>
>
> --
> 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.

[Andrew Gelman]

Hi, Mike. Can you give the whole model? If it's written as is, with flat priors on everything, it's all nonidentifiable, so I assume there's more going on here that I'm mssing.
A

[Michael Betancourt]

There are weakly-informative priors and some slope terms, but that’s about it.
Can’t share the model publicly, but I’ll send it to you privately.