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

[Marc Schwartz]

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