Hi,
I have come across the following problem which seems to be a scoping issue but I'm at a loss to see why this is so or to find a good workaround.
Suppose I have a function to get a prediction after model selection using the step function.
step.pred <- function(dat, x0) { fit.model <- step(lm(y~., data=dat), trace=F) predict(fit.model, x0, se.fit=T) }
This function works sometimes for example
set.seed(1) X.1 <- data.frame(x1=rnorm(20), x2=rnorm(20), x3=rnorm(20)) y.1 <- 5+as.matrix(X.1[,1:2])%*%matrix(c(1,1))+rnorm(20) Xy.1 <- data.frame(X.1,y=y.1) x0.1 <- data.frame(x1=-1,x2=-1, x3=-1) step.pred(Xy.1, x0.1)
$fit [1] 3.359540
$se.fit [1] 0.523629
$df [1] 16
$residual.scale [1] 1.093526
but most often it crashes as in
It does not crash. Your outdated version reports an error. R-1.8.0 does not. So please upgrade!
Uwe Ligges
set.seed(2)>
X.2 <- data.frame(x1=rnorm(20), x2=rnorm(20), x3=rnorm(20))
y.2 <- 5+as.matrix(X.2[,1:2])%*%matrix(c(1,1))+rnorm(20)
Xy.2 <- data.frame(X.2,y=y.2)
x0.2 <- data.frame(x1=-1,x2=-1, x3=-1)
step.pred(Xy.2, x0.2)
Error in model.frame.default(formula = y ~ x1 + x2, data = dat,
drop.unused.levels = TRUE) : Object "dat" not found
The difference seems to be that for the first dataset, step retains all three variables whereas for the second it drops one of them.
step(lm(y~.,data=Xy.1), trace=F)
Call: lm(formula = y ~ x1 + x2 + x3, data = Xy.1)
Coefficients:
(Intercept) x1 x2 x3 4.8347 0.8937 1.0331 -0.4516
step(lm(y~.,data=Xy.2), trace=F)
Call: lm(formula = y ~ x1 + x2, data = Xy.2)
Coefficients:
(Intercept) x1 x2 5.0802 0.9763 1.1369
One possible workaround is to explicitely assign the local variable dat in the .GlobalEnv as in
step.pred1 <- function(dat, x0) { assign("dat",dat, envir=.GlobalEnv) fit.model <- step(lm(y~., data=dat), trace=F) predict(fit.model, x0, se.fit=T) }
I don't like this method since it would overwrite anything else called dat in .GlobalEnv. I realize that I could give it an obscure name but the potential for damage still remains. Am I missing something obvious here? If not, is it possible to work around this problem in such a way that .GlobalEnv does not need to be touched?
In S-Plus I would use assign("dat",dat, frame=1)
which works but that is not available (AFAIK) in R. Is there
something similar that I can use?
I am using R 1.6.1 for Unix on a Sun Workstation. I know that I need
to upgrade but our sysadmin doesn't regard it as priority!
Thanks for any help you can give for this. Angelo
------------------------------------------------------------------ | Angelo J. Canty Email: [EMAIL PROTECTED] | | Mathematics and Statistics Phone: (905) 525-9140 x 27079 | | McMaster University Fax : (905) 522-0935 | | 1280 Main St. W. | | Hamilton ON L8S 4K1 |
______________________________________________ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
______________________________________________ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
