Question

I'm trying to write a function that collects some calls I use often in scripts
I use the sleepstudy data of the lme4 package in my examples
Here's (a simplified version of) the function I started with:

trimModel1 <- function(frm, df) {
  require(LMERConvenienceFunctions)
  require(lme4)

  lm<-lmer(frm,data=df)
  lm.trimmed = romr.fnc(lm, df)
  df = lm.trimmed$data
  # update initial model on trimmed data
  lm<-lmer(frm,data=df)
#   lm@call$formula<-frm
  mcp.fnc(lm)
  lm
}

When I call this function like below:

(fm1<-trimModel1(Reaction ~ Days + (Days|Subject),sleepstudy))

The first three lines of the output look like this:

Linear mixed model fit by REML 
Formula: frm    
Data: df

If I had called the commands of the trimModel1 function in the console the first three lines of the summary of the model look like this:

Linear mixed model fit by REML 
Formula: Reaction ~ Days + (Days | Subject) 
   Data: sleepstudy 

The difference is a problem because several packages that use the lme4 package make use of the formula and data fields. For instance the effects package uses these fields and a command like below will not work when I use the trimModel1 function above:

library(effects)
plot(allEffects(fm1))

I looked around on stackoverflow and R discussion groups for a solution and saw that you could change the formula field of the model. If you uncomment the lm@call$formula<-frm line in the trimModel1 function the formula field in the summary is displayed correctly. Unfortunately when I run a function from the effects package now I still get the error:

Error in terms.formula(formula, data = data) : 
  'data' argument is of the wrong type

This is because the data field is still incorrect.
Another possible solution I found is this function:

trimModel2 <- function(frm, df) {
  require(LMERConvenienceFunctions)
  require(lme4)

  lm<-do.call("lmer",list(frm,data=df))
  lm.trimmed = romr.fnc(lm, df)
  df = lm.trimmed$data
  # update initial model on trimmed data
  lm<-do.call("lmer",list(frm,data=df))
  mcp.fnc(lm)
  lm
}

When I now type the following commands in the console I get no errors:

(fm2<-trimModel2(Reaction ~ Days + (Days|Subject),sleepstudy))
plot(allEffects(fm2))

The allEffects function works but now the problem is that the the summary of the fm2 model displays the raw sleepstudy data. That is not a big problem with the sleepstudy data but with very large datasets sometimes Rstudio crashed when displaying a model.
Does anyone know how to make one (or both) of these functions work correctly?
I think I have to change the fm1@call$data field but I don't know how.

Was it helpful?

Solution

Do it like this:

trimModel1 <- function(frm, df) {
  require(LMERConvenienceFunctions)
  require(lme4)
  dfname <- as.name(deparse(substitute(df)))

  lm<-lmer(frm,data=df)
  lm.trimmed = romr.fnc(lm, df)
  df = lm.trimmed$data
  # update initial model on trimmed data
  lm<-lmer(frm,data=df)
  lm@call$formula <- frm
  lm@call$data <- dfname
  mcp.fnc(lm)
  lm
}

It's the "deparse-substitute trick" to get an object name from the object itself.

Licensed under: CC-BY-SA with attribution
Not affiliated with StackOverflow
scroll top