NLME fixed and random effects syntax

1,236 views
Skip to first unread message

Catherine Hulshof

unread,
Dec 6, 2013, 3:52:21 AM12/6/13
to davi...@googlegroups.com
Querido D-RUG:

I'm trying to fit a non-linear mixed effects model with 'division' as a fixed effect, 'plot' as a random effect and curve parameters a and b as mixed effects (varying by division). Hopefully it's just the nlme syntax that's throwing me off but I can't figure out how to include the division and plot effects (even after reading documentation, googling, r-sig-mixed-model archive). Does it matter that plots are nested within divisions? 

Simplified example follows. What am I not seeing?

Mucho thanks, 
Catherine


#######################
DIVISION=c("Prairie","Prairie","Savanna","Savanna","Savanna",
           "Savanna","Savanna","Subtropical","Subtropical",
           "Subtropical","Subtropical","Subtropical","Subtropical",
           "Subtropical","Subtropical","Subtropical","Subtropical","Desert",
           "Desert","Desert","Desert")
PLOT=c(1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21)
HT=c(10,20,13,14,15,18,20,25,28,29,18,14,18,23,28,29,12,13,16,18,19)
DIA=c(20,40,30,30,30,45,50,53,54,56,29,49,59,39,39,49,29,49,29,29,20)
LOCO=data.frame(DIVISION,PLOT,HT,DIA)
grouped=groupedData(HT~DIA|DIVISION,data=LOCO)

power=nlme(HT~a*DIA^b,
           data=grouped,
           fixed=a+b ~ 1,
           random=a +b ~ 1|DIVISION,
           start=c(a=3,b=0.5))

--
Catherine M. Hulshof, PhD
Environmental Science and Policy Department
University of California - Davis, Davis, CA 95616
http://catherinehulshof.wordpress.com/

Jens Stevens

unread,
Dec 7, 2013, 4:30:23 PM12/7/13
to davi...@googlegroups.com
Although I'm not familiar with nlme syntax, from your example, it doesn't look like you need plot as a random effect, since you appear to have only one data point from each plot. In that case, it looks like your syntax gives you what you want: coef(power) returns a value of a and b for each division (in essence, it sums the fixed effects estimates of a and b with the random effects estimates of division). If you actually have multiple data points from a given plot, then it would make sense to include plot as a random effect (I think), but I'm not sure what the syntax would be. I think the Pinhiero and Bates book (2000?) would be worth looking at for this.

Possible alternatives to this kind of analysis are discussed here: https://groups.nceas.ucsb.edu/non-linear-modeling/projects/OrangeTree/WRITEUP/OrangeTree.pdf

Daniel Fulop

unread,
Dec 7, 2013, 8:26:09 PM12/7/13
to davi...@googlegroups.com
Hi Catherine,

You definitely have a syntax problem with respect to the random effects formula, and additionally I think you're going beyond nlme's random effects capabilities.  In nlme's random formula (as far as I know) there are never any parameters to the left of the tilde.  Moreover, I think nlme cannot handle crossed random effects in a complex mixed effect model such as yours.  In any case, you would need to use pdIdent and pdBlocked objects to specify crossed random effects in nlme.

You're better off sticking with lme4, which has better capabilities for complex random effects specifications.  However, you won't be able to fit a multivariate (more than 1 response variable) model in lme4.  To include both the multivariate aspect and the crossed and nested random/mixed effects you could try MCMCglmm.

HTH,
Dan.

December 6, 2013 at 12:52 AM
--
Check out our R resources at http://www.noamross.net/davis-r-users-group.html
---
You received this message because you are subscribed to the Google Groups "Davis R Users' Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email to davis-rug+...@googlegroups.com.
Visit this group at http://groups.google.com/group/davis-rug.
For more options, visit https://groups.google.com/groups/opt_out.

--
Daniel Fulop, Ph.D.
Postdoctoral Scholar
Dept. Plant Biology, UC Davis
Maloof Lab, Rm. 2220
Life Sciences Addition, One Shields Ave.
Davis, CA 95616

510-253-7462
dfu...@ucdavis.edu

Daniel Fulop

unread,
Dec 7, 2013, 8:28:58 PM12/7/13
to davi...@googlegroups.com
Sorry, I just realized you're dealing with a non-linear model.  In any case, it still seems like you might have to look beyond nlme, and the resource Jens provided seems very appropriate.

Catherine Hulshof

unread,
Dec 8, 2013, 8:06:13 PM12/8/13
to davi...@googlegroups.com
Jens, Dan,

Thanks for the suggestions, definitely a few good leads!

Saludos,
CMH


On Sat, Dec 7, 2013 at 5:28 PM, Daniel Fulop <dfulo...@gmail.com> wrote:
Sorry, I just realized you're dealing with a non-linear model.  In any case, it still seems like you might have to look beyond nlme, and the resource Jens provided seems very appropriate.

--
Check out our R resources at http://www.noamross.net/davis-r-users-group.html
---
You received this message because you are subscribed to the Google Groups "Davis R Users' Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email to davis-rug+...@googlegroups.com.
Visit this group at http://groups.google.com/group/davis-rug.
For more options, visit https://groups.google.com/groups/opt_out.
Reply all
Reply to author
Forward
0 new messages