[R-sig-ME] nlme model not working but lme models are fine [SEC=UNCLASSIFIED]
Douglas Bates
bates at stat.wisc.edu
Thu Jan 12 20:38:23 CET 2012
On Wed, Jan 11, 2012 at 4:35 PM, Gosse, Michelle
<Michelle.Gosse at foodstandards.gov.au> wrote:
> Hi all,
>
> I've replied to my message so hopefully the archive will stay tracking this as a single question.
>
> Having examined the log likelihood formulae I was given for the SAS code, read Chapter 7 of Pinheiro & Bates, and searching for help, I have got as far as (note, I am not using a grouped data frame as I get an error in nlme.fomula message saying that the "starting values for the fixed component are the wrong length" when I do this):
>
> male.lme3 <- lme(BoxCoxXY ~ ordered(AgeFactor) + IntakeDay,
> data=Male.Data,
> random= ~1|RespondentID)
>
> and then I am trying to run the following, and yes I wish the covariates to be additive the model I am basing this on (from SAS) has additive covariates:
> Male.nlme <- nlme(BoxCoxXY ~ A + B*factor(AgeFactor) + C*factor(IntakeDay),
> data=Male.Data,
> fixed= A + B + C ~ 1,
> random=A ~1,
> group=RespondentID,
> start=fixef(male.lme3)
> )
>
> I get the error " Error in eval(expr, envir, enclos) : object 'A' not found"
>
> I don't have a more complicated model as the REML is taking care of it. Using the Wafer example in Chapter 8 of Pinheiro & Bates, and also this thread:
> https://stat.ethz.ch/pipermail/r-sig-mixed-models/2011q1/005473.html
>
> I thought I just had to show how the coefficients relate to the covariates and intercept. Clearly I have made an error of judgement, and I'm still not sure where I have gone wrong with the model. The Wafer nlme model uses A, B, and C in its model statement and the linked thread uses b0,...,b5.
>
> Would someone mind pointing out where I have gone wrong this time?
The names of your starting values must correspond to the names used in
the formula (A, B and C). I don't think that formula is what you want
because ordered(AgeFactor) will be converted to integer values 1, 2,
..., # of levels of AgeFactor and the same with factor(IntakeDay).
It is highly unlikely that this model does what you expect it to do.
In the end you would just end up with an awkward nonlinear model
formula representing a linear model.
You said that you were asked to fit a nonlinear mixed-effects model.
I recommend that you go back to the person who suggested this and get
clarification on what model they intended. Converting a linear
mixed-effects model to a nonlinear model formula, in which all the
parameters occur linearly, and fitting that is not a nonlinear
mixed-effects model. It is the same linear mixed-effects model fit
inefficiently.
> -----Original Message-----
> From: r-sig-mixed-models-bounces at r-project.org [mailto:r-sig-mixed-models-bounces at r-project.org] On Behalf Of Gosse, Michelle
> Sent: Wednesday, January 11, 2012 4:05 PM
> To: 'Douglas Bates'
> Cc: r-sig-mixed-models at r-project.org
> Subject: Re: [R-sig-ME] nlme model not working but lme models are fine [SEC=UNCLASSIFIED]
>
> Hi,
>
> Thanks for the help. I'm now trying to figure out (1) the function and (2) how to specify it in nlme.
>
> cheers
> Michelle
>
> -----Original Message-----
> From: dmbates at gmail.com [mailto:dmbates at gmail.com] On Behalf Of Douglas Bates
> Sent: Wednesday, January 11, 2012 8:56 AM
> To: Gosse, Michelle
> Cc: r-sig-mixed-models at r-project.org
> Subject: Re: [R-sig-ME] nlme model not working but lme models are fine [SEC=UNCLASSIFIED]
>
> On Tue, Jan 10, 2012 at 1:05 PM, Gosse, Michelle <Michelle.Gosse at foodstandards.gov.au> wrote:
>> Hi all,
>>
>> I have a nicely working lme model from the lme4 package, however I have been requested to produce a nonlinear model and I am having problems with the model specification.
>
> So what is the nonlinear model? You haven't specified it in your call to nlme.
>
>> My working lme4 code is:
>> Male.lme1 <- lmer(BoxCoxXY ~ AgeFactor + IntakeDay + (1|RespondentID),
>> data=Male.Data,
>> weights = SampleWeight)
>>
>> Where:
>> BoxCoxXY are transformed nutrient intake values AgeFactor is the age band for the subject, this has been converted to a factor. There are 4 agebands in the current dataset, and these are in ascending order. The values for AgeFactor are 1,2,3,4.
>> IntakeDay is either "IntakeDay1" or IntakeDay2", depending on whether it is the first 24-hour recall period or the second, for food intake (used to calculate the raw nutrient intake) RespondentID is the subject number, this has been converted to a factor I'm doing the analysis separately by gender, so all my data relate to males, hence the data frame names.
>>
>> I am now trying to fit this data using nlme in the nlme package. I've been taking hints from the Machines and Wafers data, as these seem closest to my data frame. I just can't get the nlme model to work.
>>
>> In my R session, at this point I've detached lme4 and attached nlme. These two bits of code work fine, and the lme coefficient estimates from the nlme package are exactly the same as those from the lme in the lme4 package, so I believe I have the linear model specified correctly in nlme, the issue seems to be getting the nlme model to work. The working code in nlme is (I am using lme to generate the fixed effect coefficients to put into the nlme model as the starting parameters):
>>
>> Male.Group <- groupedData(BoxCoxXY ~ RespondentID|AgeFactor, data=Male.Data) # This seems to be the most sensible grouping as the subjects are grouped within age, with repeated intake measures grouped within subject.
>>
>> male.lme2 <- lme(BoxCoxXY ~ ordered(AgeFactor) + IntakeDay,
>> data = Male.Group,
>> random = ~ 1 | RespondentID)
>>
>>
>> When I run the next bit,I get an error. Error and traceback() provided after the syntax:
>>
>> Male.nlme <- nlme(BoxCoxXY ~ ordered(AgeFactor) + IntakeDay,
>> fixed = ordered(AgeFactor) + IntakeDay ~ 1,
>> random = RespondentID ~ 1,
>> data = Male.Group,
>> start= fixef(male.lme2)
>> )
>
> Your formula is not an nlme specification. The right hand side of the formula should be a function call using nonlinear model parameters and covariates. You are using a linear model formula on the right hand side and this will not give the result you are expecting.
>
>> I get the error message:
>> Error in eval(expr, envir, enclos) : object 'AgeFactor' not found
>>
>> The results of traceback() are:
>> 8: eval(expr, envir, enclos)
>> 7: eval(x[[length(x)]], dat)
>> 6: FUN(X[[1L]], ...)
>> 5: lapply(form, function(x, dat, N) {
>> val <- eval(x[[length(x)]], dat)
>> if (length(val) == 1) {
>> return(as.factor(rep(val, N)))
>> }
>> else {
>> return(as.factor(val)[drop = TRUE])
>> }
>> }, dat = object, N = nrow(object))
>> 4: getGroups.data.frame(dataMix, eval(parse(text = paste("~1",
>> deparse(groups[[2]]),
>> sep = "|"))))
>> 3: getGroups(dataMix, eval(parse(text = paste("~1",
>> deparse(groups[[2]]),
>> sep = "|"))))
>> 2: nlme.formula(BoxCoxXY ~ ordered(AgeFactor) + IntakeDay, fixed =
>> ordered(AgeFactor) +
>> IntakeDay ~ 1, random = RespondentID ~ 1, data = Male.Group,
>> start = fixef(male.lme2))
>> 1: nlme(BoxCoxXY ~ ordered(AgeFactor) + IntakeDay, fixed =
>> ordered(AgeFactor) +
>> IntakeDay ~ 1, random = RespondentID ~ 1, data = Male.Group,
>> start = fixef(male.lme2))
>>
>> Could someone advise me where I have gone wrong? I don't understand how nlme cannot find AgeFactor where the other syntax didn't have an issue.
>
> UNCLASSIFIED
>
> **********************************************************************
> This email and any files transmitted with it are confide...{{dropped:9}}
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
> UNCLASSIFIED
>
> **********************************************************************
> This email and any files transmitted with it are confide...{{dropped:9}}
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
More information about the R-sig-mixed-models
mailing list