Fitting the same models in nlme and lme4

Question

The data are from here

library(nlme)
dat0 <- read.table("aids.dat2",head=T)
dat1 <- dat0[dat0$day<=90, ]   # use only first 90-day data
dat2 <- dat1[!apply(is.na(dat1),1,any),]  # remove missing data 

# Next, let's treat the data as longitudinal (or grouped) data 
aids.dat <- groupedData(lgcopy ~ day | patid, data=dat2)

# A NLME model fit, with random effects on all 4 parameters 
start <- c(10,0.5,6,0.005)  # starting value 

aids.dat$log10copy = log10(aids.dat$lgcopy)

nlme.fit <- nlme(log10copy ~ exp(p1-b1*day) + exp(p2-b2*day + 1),
                 fixed = list(p1 ~ 1, b1 ~ 1, p2 ~ 1, b2 ~ 1),
                 random = list(patid = pdDiag(list(p1 ~ 1, b1 ~ 1, p2 ~ 1, b2 ~ 1))),
                 data =aids.dat, start=c(start)) 
summary(nlme.fit)

Here I fit a nonlinear mixed effects model using nlme in the nlme package. The model has 4 fixed effects and 4 random effects. I specified a diagonal structure on the variance-covariance matrix, and each patid forms a group.

library(lme4)
deriv_mod <- deriv( ~ exp(p1 - b1*t) + exp(p2 - b2*t + 1), 
                    c("p1", "b1", "p2", "b2"), function(t, p1, b1, p2, b2){})
nlmer.fit <- nlmer(deriv_mod ~ list(p1 ~ 1, b1 ~ 1, p2 ~ 1, b2 ~ 1) + 
                     list(p1 ~ 1, b1 ~ 1, p2 ~ 1, b2 ~ 1), data = aids.dat, start = c(start))

Here, I would like to fit the same model using the lme4 package. From the documentation it seems that the formula for nlmer must also have a gradient component, thus I used the deriv function first. However, I am not sure how to specify the rest of the parameters? The

deriv_mod ~ list(p1 ~ 1, b1 ~ 1, p2 ~ 1, b2 ~ 1) + 
                     list(p1 ~ 1, b1 ~ 1, p2 ~ 1, b2 ~ 1)

is to specify 4 fixed effects (in the first list object) and their corresponding 4 random effects (in the second list object). However, I am not quite sure how to specify a diagonal variance-covariance structure and make sure that the observations are grouped by patid, like I had specified in random = list(patid = pdDiag(list(p1 ~ 1, b1 ~ 1, p2 ~ 1, b2 ~ 1))) with nlme.

`lme4` doesn't allow you as much flexibility with the variance-covariance structure as `nlme`... are your RE's cross-classified, nested, etc? Can you give more details, please? — alexwhitworth, Aug 10 '17 at 19:00
The REs are neither cross-classified nor nested. They are assumed to be independent from one another since I specified a diagonal structure on the variance covariance matrix. — Adrian, Aug 12 '17 at 04:07

alexwhitworth · Accepted Answer · 2018-02-13T18:26:04.750

3

The standard way to specify fixed effects FE1 ... FE4 and independent random effects RE1 ... RE4 is shown below

mod_fit <- lme4::nlmer(Y ~ FE1 + FE2 + FE3 + FE4 + 
  (1|RE1) + (1|RE2) + (1|RE3) + (1|RE4), data= dat)

The nlme package has slightly different syntax than lme4 package.

mod_fit <- nlme::nlme(Y ~ FE1 + FE2 + FE3 + FE4 + 
      (1|RE1) + (1|RE2) + (1|RE3) + (1|RE4), 
  fixed= FE1 + FE2 + FE3 + FE4 ~ Y,
  groups= 1 ~ RE1 + RE2 + RE3 + RE4,
  data= dat)

That said, I'm not sure I completely understand the nuances of your question, so it's possible your situation means this needs to be slighly modified. If you provide comments, I'm happy to revise my answer as needed

edited Feb 13 '18 at 18:26

answered Aug 12 '17 at 14:08

alexwhitworth

4,839
5
32
59

I would like to fit the same exact model using `nlme` and `nlmer` (i.e. I'd like to see similar if not the same parameter estimates, SEs, etc) from both outputs. I was able to fit the model in `nlme` but not in `nlmer`. – Adrian Feb 13 '18 at 16:38
@Adrian updating answer now – alexwhitworth Feb 13 '18 at 18:17
`lme4.fit <- lme4::nlmer(log10copy ~ exp(p1-b1*day) + exp(p2-b2*day + 1) + (1|p1) + (1|b1) + (1|p2) + (1|b2), data= aids.dat)` I tried this but it says `log10copy` is not found? But the variable does exist in my data.frame? – Adrian Feb 13 '18 at 20:34
1

Please run my first chunk of code in my original post... it defines `log10copy` as a variable. This is a reproducible example unless I missed something...... – Adrian Feb 14 '18 at 17:24
@Adrian Thanks! You're getting an inaccurate error message. The problem appears to be with passing `deriv_mod` into the call stack. specifically, `nlformula` is not able to parse your model specification. It's not quite clear to me what your intent is; but you can try specifying the model more direclty in `nlmer` – alexwhitworth Feb 14 '18 at 17:42

Fitting the same models in nlme and lme4

1 Answers1