Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

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

Tags:

r

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.

like image 611
Erik Avatar asked Aug 02 '12 13:08

Erik


1 Answers

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.

like image 154
Josh O'Brien Avatar answered Oct 31 '22 03:10

Josh O'Brien