[R-sig-ME] HPD intervals for fixed parameters in model

John.Morrongiello at csiro.au John.Morrongiello at csiro.au
Mon Jun 4 09:03:41 CEST 2012


Hi Ben

Thanks very much for the response. I've tried both MCMCglmm and glmmADMB and am happy with the parameter estimates and the HPD intervals they produce. The only issue was specifying the random age slope in glmmADMB (as per M1admb below). It seems to only accept random factors- is it possible to specify a random slope? 

M1admb<-glmmadmb(log(Increment)~log(Age)+temp+(log(Age)|FishID)+(1|fYear),data,family='gaussian')

Warning message:
In Ops.factor(log(Age), FishID) : | not meaningful for factors
Error in model.matrix(as.formula(paste("~", lbit)), mdata) : 
  error in evaluating the argument 'object' in selecting a method for function 'model.matrix': Error in parse(text = x) : <text>:2:0: unexpected end of input
1: ~ 
  ^

Cheers

John


Message: 4
Date: Fri, 1 Jun 2012 16:44:51 +0000 (UTC)
From: Ben Bolker <bbolker at gmail.com>
To: r-sig-mixed-models at r-project.org
Subject: Re: [R-sig-ME] HPD intervals for fixed parameters in model
	with	crossed and correlated random effects
Message-ID: <loom.20120601T183904-68 at post.gmane.org>
Content-Type: text/plain; charset=us-ascii

 <John.Morrongiello at ...> writes:


> I'd very much like to calculate some HPD intervals for the fixed
> parameters in a lmer model containing crossed and correlated random
> effects. However, I have not been able to find a solution for this
> type of model after a pretty thorough search of discussion threats
> or worked examples. Is this possible at the moment with lme4
> (function or work around), or do I have to change packages and try
> something like mcmcglmm?

  I think if you want to compute HPD intervals you are more or less
stuck with MCMCglmm or glmmADMB ....
 
> My data contains 1917 growth measurements (Increment) from 704 fish
> (FishID) and the data is sampled from 29 years (Year). I'm
> interested in the effects of age (2-18 years) and water temperature
> (temp) on growth.  My model is:
 
> M1<-lmer(log(Increment)~log(Age)+temp+(log(Age)|FishID)+(1|Year),
>   data,REML=T)
> 
> HPDinterval requires an mcmc object, but when I try
> mcmcsamp(M1, n=1000)
> 
> I receive the following error:
> Error in .local(object, n, verbose, ...) :
>   Code for non-trivial theta_T not yet written


> I understand this is because MCMCsamp doesn't yet cope with
> correlated random effects? On the glmm wikidot page there is code to
> generate confidence and prediction intervals on predictions- could I
> use some of N this to get a what I'm after (if so, which bits)? I
> haven't used mcmcglmm before, so ideally would like to stay with the
> more familiar lme4 if possible to expedite the analyses! 

   I don't think so, I think you'll need MCMCglmm or glmmADMB.

In principle glmmADMB should work with your model more or less as-is
(although I don't remember whether it can do REML-y stuff or not),
and you can use mcmc=TRUE to get it to run an MCMC chain on the
parameters for you.  However, based on my limited experience, if
you do want to get MCMC working it will be easier using MCMCglmm
(which is intrinsically based on MCMC, and does a bunch of clever
stuff to make the chains mix well).

  good luck,
    Ben Bolker



More information about the R-sig-mixed-models mailing list