[Rd] Curious finding in MASS:::confint.glm() tied to eval()

[Marc Schwartz]

Hi all,

Has anyone had a chance to look at this and either validate my finding
or tell me that my brain has turned to mush?

Either would be welcome...  :-)

Thanks,

Marc


On Wed, 2006-12-13 at 13:53 -0600, Marc Schwartz wrote:
> Greetings all,
> 
> I was in the process of creating a function to generate profile
> likelihood confidence intervals for a proportion using a binomial glm.
> This is a component of a larger function to generate and plot confidence
> intervals for proportions using the above, along with bootstrap (BCa),
> Wilson and Exact to visually demonstrate the variation across the
> methods to some folks.
> 
> I had initially tried this using:
> 
>   R version 2.4.0 Patched (2006-11-19 r39944)
> 
> However, I just verified that it still occurs in:
> 
>   R version 2.4.1 RC (2006-12-11 r40157)
> 
> In both cases, I started R using '--vanilla' and this is under FC6,
> fully updated.
> 
> 
> Using the following code, it runs fine:
> 
> x <- 14
> n <- 35
> conf <- 0.95
> DF <- data.frame(y = x / n, weights = n)
> mod <- glm(y ~ 1, weights = weights, family = binomial(), data = DF)
> 
> pl.ci <- plogis(confint(mod, level = conf))
> Waiting for profiling to be done...
> 
> > pl.ci
>     2.5 %    97.5 % 
> 0.2490412 0.5651094 
> 
> 
> However, when I encapsulate the above in a function:
> 
> PropCI <- function(x, n, conf = 0.95)
> {
>   DF <- data.frame(y = x / n, weights = n)
>   mod <- glm(y ~ 1, weights = weights, family = binomial(), data = DF)
>   plogis(confint(mod, level = conf))
> }
>   
> 
> I get the following:
> 
> # Nothing else in the current global environment
> > ls()
> [1] "PropCI"
> 
> > PropCI(14, 35)
> Waiting for profiling to be done...
> Error in inherits(x, "data.frame") : object "DF" not found
> 
> 
> If however, I create a 'DF' in the global environment (which may NOT be
> the same 'DF' values as within the function!!):
> 
> DF <- data.frame(y = 14 / 35, weights = 35)
> 
> > DF
>     y weights
> 1 0.4      35
> 
> > ls()
> [1] "DF"     "PropCI"
> 
> > PropCI(14, 35)
> Waiting for profiling to be done...
>     2.5 %    97.5 % 
> 0.2490412 0.5651094 
> 
> 
> To my point above about the 'DF' in the global environment not being the
> same as the 'DF' within the function body:
> 
> > DF <- data.frame(y = 5 / 40, weights = 40)
> > DF
>       y weights
> 1 0.125      40
> 
> > PropCI(14, 35)
> Waiting for profiling to be done...
>     2.5 %    97.5 % 
> 0.3752233 0.4217556 
> 
> 
> Rather scary I think....
> 
> 
> 
> So, this suggested a problem in where confint.glm() was looking for
> 'DF'.
> 
> I subsequently traced through the code using debug(), having removed
> 'DF' from the global environment:
> 
> > debug(MASS:::confint.glm)
> > PropCI(14, 35)
> debugging in: confint.glm(mod, level = conf)
> 
> ...
> 
> debug: object <- profile(object, which = parm, alpha = (1 - level)/4, 
>     trace = trace)
> Browse[1]> 
> Error in inherits(x, "data.frame") : object "DF" not found
> 
> 
> which brought me to profile.glm():
> 
> > debug(MASS:::profile.glm)
> > PropCI(14, 35)
> Waiting for profiling to be done...
> debugging in: profile.glm(object, which = parm, alpha = (1 - level)/4,
> trace = trace)
> 
> ...
> 
> debug: mf <- update(fitted, method = "model.frame")
> Browse[1]> 
> Error in inherits(x, "data.frame") : object "DF" not found
> 
> 
> which brought me to update.default():
> 
> > debug(update.default)
> > PropCI(14, 35)
> Waiting for profiling to be done...
> debugging in: update.default(fitted, method = "model.frame")
> 
> ...
> 
> debug: if (evaluate) eval(call, parent.frame()) else call
> Browse[1]> 
> Error in inherits(x, "data.frame") : object "DF" not found
> 
> 
> which brought me to eval():
> 
> > debug(eval)
> > PropCI(14, 35)
> debugging in: eval(mf, parent.frame())
> 
> ...
> 
> debugging in: eval(mf, parent.frame())
> debug: .Internal(eval(expr, envir, enclos))
> Browse[1]> 
> Error in inherits(x, "data.frame") : object "DF" not found
> 
> 
> 
> Unfortunately, at this point, largely due to lack of sleep of late, my
> eyes are sufficiently glazed over and my brain sufficiently vapor locked
> that the resolution is not immediately clear and I have not yet dug into
> the C code for eval().
> 
> Presumably, either I am missing something fundamental here, or there is
> a bug where, environment-wise, these respective functions are looking in
> the wrong place for 'DF', probably based upon the default environment
> arguments noted above.
> 
> Any comments from a fresh pair of eyes would be most welcome.
> 
> Regards and thanks,
> 
> Marc Schwartz