[Rd] Scope question concerning calls within a user defined function

Duncan Murdoch murdoch at stats.uwo.ca
Wed May 6 18:47:33 CEST 2009


On 5/6/2009 12:18 PM, Terrence Ireland wrote:
> The following is a simple example with a poor solution that shows my
> difficulties with scope.  The function /logit.test /has 3 arguments:  
> /model.start,/
> an initial model; /model.finish/, an all-inclusive model, /my.data/, a 
> dataset, in
> this case trivial. 
> 
> There are 2 function calls in l/ogit.test,/ first to /glm/ to get an 
> initial fit (local variable
> /logit/) using /model.start;/ then a call to /stepAIC/ using the initial 
> fit.
> 
> I had thought that l/ogit/ would carry along the dataset, /test.data/ 
> for use in
> the call to /stepAIC/.  Instead it complains that it cannot find 
> /my.data/.  It seems
> to have found the name, not the value in the /call/ attribute. 
> 
> I fixed the problem by creating a global variable /my.data/ within the 
> function,
> naming it l/ogit.test.global,/ not a very good solution.
> 
> My question is:  How do I handle scope correctly?

The general principle is this:

Variables in formulas are looked up in the data= argument of a modelling 
function; if that is missing (or the variable was not found), then in 
the environment where the formula was created.

I haven't traced through the code, but I believe stepAIC also uses the 
environment of the formula to find the dataset (i.e. it assumes the 
dataset and formula were created in the same place.

Since you didn't do that, you need some trickery.  You want my.data to 
be found in the environment of the formula.  my.data is local to your 
logit.test function, so changing logit.test to look like this:

logit.test <- function(model.start,model.finish,my.data)
   {
     environment(model.start) <- environment()
     environment(model.finish) <- environment()
     logit <- glm(model.start,binomial(link =
"logit"),my.data,control=glm.control(epsilon=0.0001),x=TRUE,y=TRUE)
     Table <- stepAIC(logit,scope=model.finish,direction="both",trace=1,k=2)
   }

should work.  The only potential problem is if your formulas have 
variables named model.start, model.finish, my.data, logit, or Table: 
those will be found before global variables of the same name.

I would avoid <<- in logit.test.global, and use the same technique there.

Duncan Murdoch

> 
> Thanks,
> 
> Terry Ireland
> 
> Script of Experiment:
> 
>  > library(MASS)
> 
>  > test.data
>   Approve Score OldScore
> 1     Yes     1       12
> 2     Yes     3       10
> 3     Yes     5        8
> 4      No     4        9
> 5      No     6        7
> 6      No     8        5
>  > test.start
> Approve ~ 1
>  > test.finish
> Approve ~ Score + OldScore
>  > logit.test
> function(model.start,model.finish,my.data)
>   {
>     logit <- glm(model.start,binomial(link = 
> "logit"),my.data,control=glm.control(epsilon=0.0001),x=TRUE,y=TRUE);
>     Table <- stepAIC(logit,scope=model.finish,direction="both",trace=1,k=2);
>   }
>  > logit.test.global
> function(model.start,model.finish,my.data)
>   {
>     my.data <<- my.data
>     logit <- glm(model.start,binomial(link = 
> "logit"),my.data,control=glm.control(epsilon=0.0001),x=TRUE,y=TRUE);
>     Table <- stepAIC(logit,scope=model.finish,direction="both",trace=1,k=2);
>   }
> 
>  > logit.test(test.start,test.finish,test.data)

In this call, test.start gets assigned to model.start in the function, 
test.finish to model.finish, and test.data to my.data.  Within 
logit.test, the glm() call looks fine, because you pass my.data to it.  But

> Start:  AIC=10.32
> Approve ~ 1
> 
> Error in inherits(x, "data.frame") : object "my.data" not found
>  > logit.test.global(test.start,test.finish,test.data)
> Start:  AIC=10.32
> Approve ~ 1
> 
>            Df Deviance     AIC
> + Score     1   4.8089  8.8089
> + OldScore  1   4.8089  8.8089
> <none>          8.3178 10.3178
> 
> Step:  AIC=8.81
> Approve ~ Score
> 
>         Df Deviance     AIC
> <none>       4.8089  8.8089
> - Score  1   8.3178 10.3178
>  >
> 
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel



More information about the R-devel mailing list