Erik
Erik

Reputation: 308

Using predict in a function call with NLME objects and a formula

I have a problem with the package NLME using the following code:

library(nlme)
x <- rnorm(100)
z <- rep(c("a","b"),each=50)
y <- rnorm(100)
test.data <- data.frame(x,y,z)
test.fun <- function(test.dat)
{
    form <- as.formula("y~x")   
    ran.form <- as.formula("~1|z")
    modell <- lme(fixed = form, random=ran.form, data=test.dat)
    pseudo.newdata <- test.dat[1,]
    predict(modell, newdata= pseudo.newdata) ###THIS CAUSES THE ERROR!
}

test.fun(test.data)

The predict causes an error and I already found what basically causes it.

The modell object saves how it was called and predict seems to use that to make prediction but is unable to find the formula objects form and ran.form becauses it does not look for them in the right namespace. In fact, I can avoid the problem by doing this:

 attach(environment(form), warn.conflicts = FALSE)
 predict(modell, newdata= pseudo.newdata) 
 detach()

My long term goal however is to save the modell to disk and use them later. I suppose I could try saving the formula objects as well, but this strikes me as a very annoying and cumbersome way to deal with the problem.

I work with automatically generated formula objects instead of writing them down explicitly because I create many models with different definitions in a sort of batch process so I can not avoid them. So my ideal solution would be a way to create the lme object so that I can forget about the formula object afterwards and predict "just works". Thanks for any help.

Upvotes: 9

Views: 1670

Answers (1)

Josh O&#39;Brien
Josh O&#39;Brien

Reputation: 162401

Try replacing lme(arg1, arg2, arg3) with do.call(lme, list(arg1, arg2, arg3)).

library(nlme)
x <- rnorm(100)
z <- rep(c("a","b"),each=50)
y <- rnorm(100)
test.data <- data.frame(x,y,z)
test.fun <- function(test.dat)
{
    form <- as.formula("y~x")   
    ran.form <- as.formula("~1|z")
    ## JUST NEED TO CHANGE THE FOLLOWING LINE
    ## modell <- lme(fixed = form, random=ran.form, data=test.dat)
    modell <- do.call(lme, list(fixed=form, random=ran.form, data=test.data))
    pseudo.newdata <- test.dat[1,]
    predict(modell, newdata= pseudo.newdata) ###THIS CAUSES THE ERROR!
}

test.fun(test.data)
#          a 
# 0.07547742 
# attr(,"label")
# [1] "Predicted values"

This works because do.call() evaluates its argument list in the calling frame, before evaluating the call to lme() that it constructs. To see why that helps, type debug(predict), and then run your code and mine, comparing the debugging messages printed when you are popped into the browser.

Upvotes: 6

Related Questions