[R] scoping problem

Christos Hatzis christos at nuverabio.com
Wed Oct 24 22:20:04 CEST 2007


Another way to do this without messing with environments is to update the
data and formula locally within the function and re-run the regression on
the updated model/data:

tukey.test <- function(m) {
    ud <- data.frame( m$model, pred2=m$fitted.values^2 )
    uf <- update.formula(formula(m$terms), ~ . + pred2)
    summary(lm(uf, ud))$coef
}

> data(BOD)
> m1 <- lm(demand~Time,BOD)
> tukey.test(m1)
              Estimate Std. Error    t value  Pr(>|t|)
(Intercept) 12.5854086  5.0344850  2.4998403 0.0877190
Time         7.5390634  6.1271204  1.2304415 0.3062126
pred2       -0.1096744  0.1148654 -0.9548081 0.4101142
 

> -----Original Message-----
> From: r-help-bounces at r-project.org 
> [mailto:r-help-bounces at r-project.org] On Behalf Of Sandy Weisberg
> Sent: Wednesday, October 24, 2007 3:44 PM
> To: r-help at stat.math.ethz.ch
> Subject: [R] scoping problem
> 
> I would like to write a function that computes Tukey's 1 df 
> for nonadditivity.  Here is a simplified version of the 
> function I'd like to
> write:  (m is an object created by lm):
> 
> tukey.test <- function(m) {
>   m1 <- update(m, ~.+I(predict(m)^2))
>   summary(m1)$coef
>   }
> 
> The t-test for the added variable is Tukey's test.  This won't work:
> 
> data(BOD)
> m1 <- lm(demand~Time,BOD)
> tukey.test(m1)
> 
> Error in predict(m) : object "m" not found
> 
> This function doesn't work for two reasons:
> 	1. The statement m1 <- update(m, ~.+I(predict(m)^2)) 
> can't see 'm' in the call to predict.
> 	2. If in creating m missing values had been present, 
> then predict(m), even if it could be computed, could be of 
> the wrong length.
> 
> Can anyone help?
> 
> 
> 
> --
> Sanford Weisberg, sandy at stat.umn.edu
> Office and mailing address:
> University of Minnesota, School of Statistics
> 312 Ford Hall, 224 Church St. SE, Minneapolis, MN  55455 
> 612-625-8355, FAX 612-624-8868
> 
> St. Paul office:
> 146 Classroom-Office Building, 612-625-8777
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
>



More information about the R-help mailing list