Skip to content

HPD intervals for fixed parameters in model

4 messages · John.Morrongiello at csiro.au, Ben Bolker

Original
John.Morrongiello at csiro.au
# 
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 think if you want to compute HPD intervals you are more or less
stuck with MCMCglmm or glmmADMB ....
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
Ben Bolker
#
On 12-06-04 03:03 AM, John.Morrongiello at csiro.au wrote:
As a quick possible workaround, try defining data$logAge <-
log(data$Age) and using (logAge|FishID) in your model formula; I'm not
sure that glmmADMB can handle defining variable transformations on the
fly in that context.

  Also be aware that the default for lme4 is to use the full correlation
matrix while glmmADMB uses a diagonal matrix, so (e.g.) (logAge|FishID)
would estimate three parameters (variation in intercept across fish,
variation in slope across fish, covariance between intercept and slope)
while glmmADMB would only estimate the first two.  (I think you have to
specify this choice explicitly in MCMCglmm.)
John.Morrongiello at csiro.au
# 
Log transforming the Age variable prior to putting it in the model solved the problem of a random slope in glmmADMB. 

Thanks also to Ben for the note about the different covariance matrices in use. I did a bit more reading on this, especially for MCMCglmm, and came up with the following model structures. Are the following lmer, MCMCglmm and glmmADMB roughly comparable?

*Variation in fish intercepts, variation in Age slopes, correlation between intercept and slope:

M1lmer<-lmer(log(Increment)~logAge+temp+(logAge|FishID)+(1|fYear),data,REML=T)
M1mcmc<-MCMCglmm(log(Increment)~logAge+temp,random=~us(1+logAge):FishID + fYear,rcov= ~units, data=data, prior=prior)##us allows for intercept and slope covariance
##not possible with glmmADMB??

* Variation in fish intercepts and age slopes, but no correlation between intercept and slope:

M2lmer<-lmer(log(Increment)~logAge+temp+(1+FishID)+(0+logAge|FishID)+(1|Year),data,REML=F)##REML= FALSE is equivalent to glmmADMB?
M2mcmc<-MCMCglmm(log(Increment)~logAge+temp,random=~idh(1+logAge):FishID + fYear,rcov= ~units, data=data, prior=prior)##idh sets covariance between intercept and slope to zero (also could use idh(logAge):FishID+FishID)
M2glmmADMB<-glmmadmb(log(Increment)~logAge+temp+(logAge|FishID)+(1|fYear),NSWtiger,family='gaussian')

Cheers

John

-----Original Message-----
From: Ben Bolker [mailto:bbolker at gmail.com] 
Sent: Monday, 4 June 2012 5:12 PM
To: Morrongiello, John (CMAR, Hobart)
Cc: r-sig-mixed-models at r-project.org
Subject: Re: [R-sig-ME] HPD intervals for fixed parameters in model
On 12-06-04 03:03 AM, John.Morrongiello at csiro.au wrote:
As a quick possible workaround, try defining data$logAge <-
log(data$Age) and using (logAge|FishID) in your model formula; I'm not
sure that glmmADMB can handle defining variable transformations on the
fly in that context.

  Also be aware that the default for lme4 is to use the full correlation
matrix while glmmADMB uses a diagonal matrix, so (e.g.) (logAge|FishID)
would estimate three parameters (variation in intercept across fish,
variation in slope across fish, covariance between intercept and slope)
while glmmADMB would only estimate the first two.  (I think you have to
specify this choice explicitly in MCMCglmm.)
Ben Bolker
# 
<John.Morrongiello at ...> writes:
[snip]
intercept and slope:
random=~us(1+logAge):FishID + fYear,rcov= ~units,
          data=data, prior=prior)
I believe the glmmadmb syntax here would be

M1glmmadmb <- glmmadmb(log(Increment)~logAge+temp+(logAge|FishID)+(1|fYear),
        data,corStruct="full",family="gaussian")

A couple of notes:

* glmmadmb() uses _Poisson_ family by default (it was originally
designed for GLMMs, not LMMs)
* I believe that neither MCMCglmm nor glmmadmb is doing REML here
 (although technically MCMCglmm isn't doing likelihood estimation
  in any case, and at the moment I can't think about the correspondence
 between marginal posterior means and REML estimates of variances ...)
I think so.
random=~idh(1+logAge):FishID + fYear,rcov=~units, 
     data=data, prior=prior)
 ##idh sets covariance between intercept and slope to zero (also could
##  use idh(logAge):FishID+FishID)

 M2glmmADMB<-glmmadmb(log(Increment)~logAge+temp+(logAge|FishID)+(1|fYear),
  NSWtiger,family='gaussian')

  I think that's all right.

  You know, the easiest way to check these correspondences is just
to try them (especially on simple simulated data where you know
what the true answer should be).