observation-level (or unit-level) random effects with INLA

[Ramon Diaz-Uriarte]

Dear All,

When working with GLMMs, e.g., Poisson or binomial, a possible way of
dealing with overdispersion (specially to account for omitted relevant
variables) is to add what is often called a "unit-level" or
"observation-level" random effect: a random effect for each
observation. This can be done explicitly with, say, glmer, and it is done
implicitly in MCMCglmm.


But I am slightly confused about the possible use of "observation-level"
random effects with INLA.


The RSS 2009 paper, in p. 336, uses an individual-level random effect (to
deal with overdispersion in a Poisson GLMM) but it is individual-level
(each individual measured several times), not observation-level.  The
paper by Holand et al., "Animal models and integrated nested Laplace
approximations.", discusses (p. 8) the overdispersion problem, and in its
supplementary file S4, shows how it adds "an individual random
effect". However, their Poisson and binomial examples only have one
observation per individual. Likewise, the paper by Fong et al. ("Bayesian
inference for generalized linear mixed models") has an example where they
use plate-level random effects, but there is a single observation per
plate. The paper by Martino & Rue, "Case studies in Bayesian Computation
using INLA" in p. 3 says " a small random noise, εi , with high precision
is always automatically added to the model". This sounds similar to that
additive overdispersion added to the linear predictor, but since it says
"with high precision", and I can see no information about it returned in
the model fits, I do not think this really plays the role of that
"observation-level" random effect.



So it is unclear to me if, when one has more than one observation per
individual, one can/should still add that observation-level random effect
(in addition to the individual random effect). Looking in the examples and
googling around I have not found anything else.



What am I missing?


Best,

R.

[INLA help]

On Wed, 2014-04-02 at 04:33 -0700, Ramon Diaz-Uriarte wrote:
> Dear All,
>
> When working with GLMMs, e.g., Poisson or binomial, a possible way of
> dealing with overdispersion (specially to account for omitted relevant
> variables) is to add what is often called a "unit-level" or
> "observation-level" random effect: a random effect for each
> observation. This can be done explicitly with, say, glmer, and it is
> done
> implicitly in MCMCglmm.
>

You can put the random effects how you like...

y ~ ... + f(idx, model='iid')

says that the i'th linear predictor, used to predict, y[i], depends on
the idx[i]'th random effect.

with idx = c(1,2,3,4), then each linear predictor has different random
effect, wheras idx = c(1,1,3,4), the first and second linear predictor
share the random effect wheras the 3rd and 4th linear predicctor use the
3rd and 4th linear predictor (here, there is no 2nd...).

if you need an example, let me know.

H

--
Håvard Rue
he...@r-inla.org

[Ramon Diaz-Uriarte]

Dear Havard,

Thansk for your response. A question follows below:

So, for the unit-level or observation-level random effect this would be just

y ~  ... + f(factor(1:nrow(data)), model='iid'), data = data

 
But, does that make sense? And I am also wondering if the issue addressed here

http://www.r-inla.org/faq#TOC-INLA-seems-to-work-great-for-near-all-cases-but-are-there-cases-where-INLA-is-known-to-have-problems-

would make that problematic.


Best,

R.

[INLA help]

On Thu, 2014-04-03 at 04:26 -0700, Ramon Diaz-Uriarte wrote:

>
>
> So, for the unit-level or observation-level random effect this would
> be just
>
> y ~ ... + f(factor(1:nrow(data)), model='iid'), data = data
>
>
> But, does that make sense?

It does what you ask it to do, he he. I am not used to the unit-level
phrase, so if you have n=nrow(data) observations, and want one random
effect for each, then f(idx, model="iid") where idx = 1:nrow(data)
(factor also works), is the right thing.

> And I am also wondering if the issue addressed here
>
> http://www.r-inla.org/faq#TOC-INLA-seems-to-work-great-for-near-all-cases-but-are-there-cases-where-INLA-is-known-to-have-problems-
>
> would make that problematic.

if you have binary data with very little information, then this applies.
except for the precision of the random effect, you still get great
reseults.

Best


--
Håvard Rue
he...@r-inla.org

[Ramon Diaz-Uriarte]

On Thursday, April 3, 2014 1:38:50 PM UTC+2, help help wrote:

On Thu, 2014-04-03 at 04:26 -0700, Ramon Diaz-Uriarte wrote:

>
>
> So, for the unit-level or observation-level random effect this would
> be just
>
> y ~  ... + f(factor(1:nrow(data)), model='iid'), data = data
>
>  
> But, does that make sense?

It does what you ask it to do, he he.

:-)
 

I am not used to the unit-level
phrase, so if you have n=nrow(data) observations, and want one random
effect for each, then f(idx, model="iid")  where idx = 1:nrow(data)
(factor also works), is the right thing.

OK, understood.

 

> And I am also wondering if the issue addressed here
>
> http://www.r-inla.org/faq#TOC-INLA-seems-to-work-great-for-near-all-cases-but-are-there-cases-where-INLA-is-known-to-have-problems-
>
> would make that problematic.

if you have binary data with very little information, then this applies.
except for the precision of the random effect, you still get great
reseults.

Aha, thanks for the clarification.


Best,

R.

 

[KQ]

Hi Ramon and Håvard,

I know I am joining this thread late, but I just wanted to jump off the same question that Ramon initially had. I am facing the same issue using glmer in the lme4 package using a Poisson distribution and facing problems with overdispersion. Lme4 no longer supports "quasi" families, but other options like glm do not incorporate random effects. Therefore, I have looked into two other options: 1) the package mcmcglmm (which I am not familiar with in theory or practice) and 2) adding in observation level random effects into the lme4 model (my preference). 

I have looked around for a workable example but have thus far have not found any. Would either of you be willing to share an example you might have? 

Håvard: Although you have given Ramon a very good explanation of how to solve the issue, I have to admit that I do not have the modelling experience to replicate what you have taught. The explanation and link (http://www.r-inla.org/faq#TOC-INLA-seems-to-work-great-for-near-all-cases-but-are-there-cases-where-INLA-is-known-to-have-problems-) were quite out of my area of expertise. 

Any simplified explanation of how to choose the level of the unit and then incorporate that into the glmer model would be much appreciated if possible!

Thank you,

K

[Håvard Rue]

Hi,

I'm not sure what it eh problem really. INLA uses a very basic interface
where you connect the linear predictor to a sum of model components, and
the spesification like y ~ (a|b) is not supported, but must be coded
manually. If you have a spesific model in mind, please send me the math
of the model (not just words, he he), and I can write a small simulation
code showing how to do this

Best
H

On Mon, 2014-05-26 at 21:33 -0700, KQ wrote:
> Hi Ramon and Håvard,
> I know I am joining this thread late, but I just wanted to jump off
> the same question that Ramon initially had. I am facing the same issue
> using glmer in the lme4 package using a Poisson distribution and
> facing problems with overdispersion. Lme4 no longer supports "quasi"
> families, but other options like glm do not incorporate random
> effects. Therefore, I have looked into two other options: 1) the
> package mcmcglmm (which I am not familiar with in theory or practice)
> and 2) adding in observation level random effects into the lme4 model
> (my preference).
> I have looked around for a workable example but have thus far have not
> found any. Would either of you be willing to share an example you
> might have?
>
>
> Håvard: Although you have given Ramon a very good explanation of how
> to solve the issue, I have to admit that I do not have the modelling
> experience to replicate what you have taught. The explanation and link
> (http://www.r-inla.org/faq#TOC-INLA-seems-to-work-great-for-near-all-cases-but-are-there-cases-where-INLA-is-known-to-have-problems-) were quite out of my area of expertise.
>
>
>
> Any simplified explanation of how to choose the level of the unit and
> then incorporate that into the glmer model would be much appreciated
> if possible!

--
Håvard Rue
Department of Mathematical Sciences
Norwegian University of Science and Technology
N-7491 Trondheim, Norway
Voice: +47-7359-3533 URL : http://www.math.ntnu.no/~hrue
Mobile: +47-9260-0021 Email: havar...@math.ntnu.no

R-INLA: www.r-inla.org

[Ramon Diaz-Uriarte]

Dear K,

On Tue, 27-05-2014, at 06:33, KQ <piensagl...@gmail.com> wrote:
> Hi Ramon and Håvard,
> I know I am joining this thread late, but I just wanted to jump off the
> same question that Ramon initially had. I am facing the same issue using
> glmer in the lme4 package using a Poisson distribution and facing problems
> with overdispersion. Lme4 no longer supports "quasi" families, but other
> options like glm do not incorporate random effects. Therefore, I have
> looked into two other options: 1) the package mcmcglmm (which I am not
> familiar with in theory or practice) and 2) adding in observation level
> random effects into the lme4 model (my preference).
> I have looked around for a workable example but have thus far have not
> found any. Would either of you be willing to share an example you might
> have?

Hummm... I am not sure I see the problem. Here are a few thoughts:


1. With INLA, as Havard explained, you can add an observation-level random
effect:

inla(Y ~ x1 + x2 + x3 + f(z, model = "iid") + f(observation.id, model =
"iid"), )


where z is the random effect, and observation.id is something you can
create as

observation.id <- factor(1:nrow(mydata))

(converting it to factor is not needed, but I do it anyway, to prevent
problems in other places, not INLA related).


By the way, that is the same as in the page you link below:

idx = 1:n
result = inla( y ~ 1 + f(idx, model="iid")


What is not entirely clear to me is how to asses when you really need to
add that observation-level random effect. But that is a different story.


2. With nlme you can do something similar. The syntax will be different, of
course ( (1|observation.id)), but it is immediate.


3. If you do not have huge that, you can try running both INLA and
MCMCglmm, to double check you are fitting similar models. Or try it in a
small sample of your data.

>
> Håvard: Although you have given Ramon a very good explanation of how to
> solve the issue, I have to admit that I do not have the modelling
> experience to replicate what you have taught. The explanation and link
> (http://www.r-inla.org/faq#TOC-INLA-seems-to-work-great-for-near-all-cases-but-are-there-cases-where-INLA-is-known-to-have-problems-)
> were quite out of my area of expertise.
>

But that is for a binomial with N = 1. You said you had Poisson counts.

> Any simplified explanation of how to choose the level of the unit and then
> incorporate that into the glmer model would be much appreciated if
> possible!

That is what I do not understand. You just create the 1:n identifier of
observations, and add it to the model.


Best,


R.

--
Ramon Diaz-Uriarte
Department of Biochemistry, Lab B-25
Facultad de Medicina
Universidad Autónoma de Madrid
Arzobispo Morcillo, 4
28029 Madrid
Spain

Phone: +34-91-497-2412

Email: rdi...@gmail.com
ramon...@iib.uam.es

http://ligarto.org/rdiaz