[R-sig-ME] Some Basic lmer Questions

Andrew Robinson A.Robinson at ms.unimelb.edu.au
Mon Jan 14 21:50:24 CET 2008


Hi Kevin,

welcome to the wonderful world of mixed-models!  I heartily recommend
that you read the following books: Gelman and Hill, and Pinehiro and
Bates. 

On Mon, Jan 14, 2008 at 03:21:20PM -0500, Kevin E. Thorpe wrote:
> I am taking my first cautious steps into the mixed-models pool and I
> have a few, probably basic questions.
> 
> The data I am faced with are lab values taken at regular time intervals
> (0, 4, 8, 12 and 24 hours) following a surgery.
> 
> > str(trop)
> 'data.frame':	790 obs. of  6 variables:
>  $ pid  : int  0 0 0 0 0 1 1 1 1 1 ...
>  $ ittrx: int  1 1 1 1 1 2 2 2 2 2 ...
>  $ pprx : int  1 1 1 1 1 2 2 2 2 2 ...
>  $ rx3  : Factor w/ 3 levels "On","Off","Converted": 1 1 1 1 1 2 2 2 2 2 ...
>  $ hours: num  0 4 8 12 24 0 4 8 12 24 ...
>  $ trop : num  2.12 9.51 5.79 4.37 1.8 NA NA NA NA NA ...
>  - attr(*, "reshapeLong")=List of 4
>   ..$ varying:List of 1
>   .. ..$ : chr  "Trop0" "Trop4" "Trop8" "Trop12" ...
>   ..$ v.names: chr "trop"
>   ..$ idvar  : chr "pid"
>   ..$ timevar: chr "hours"
> 
> 
> Of interest is whether or not there are differences among the groups
> represented by rx3 above.  If we pretend for the moment that the time
> effect is linear and there is no treatment by time interaction, I
> would be inclined to test for differences as follows (also
> ignoring any correlation structure).
> 
> trop.lme0 <- lmer(trop~hours+(1|pid),data=trop,method="ML")
> trop.lme1 <- lmer(trop~rx3+hours+(1|pid),data=trop,method="ML")
> anova(trop.lme0,trop.lme1)
> 
> I seem to recall hearing/reading that the LRT from anova() is
> appropriate for maximum-likelihood but not REML which is why
> I used method="ML".  So, is this the right approach or have
> I seriously misunderstood something?

It's under discussion.  The current recommendation is to use
mcmcsamp to sample the posterior distribution.  

> Next, assuming I have not done anything egregious, I want to
> turn to the non-linearity of the time effect.  I found that
> I can use ns() in the splines package to include a spline term
> for hours, but is this the right approach?  I could also see
> making hours a factor to allow for non-linearity as well.

It's plausible.  You could also try GAMM in mgcv.
 
> Finally, (and feel free to point me at suitable references)
> how does one determine the appropriate correlation structure
> to use in these models?

See above!

Andrew
-- 
Andrew Robinson  
Department of Mathematics and Statistics            Tel: +61-3-8344-9763
University of Melbourne, VIC 3010 Australia         Fax: +61-3-8344-4599
http://www.ms.unimelb.edu.au/~andrewpr
http://blogs.mbs.edu/fishing-in-the-bay/




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