2
$\begingroup$

I am trying to model disease counts (d) by using population (p) as offset to control for exposure. In R, I found two possible ways to go:

m1 <- glm(d ~ 1 + offset(log(n)), family=poisson, data=dat) m2 <- glm(d ~ 1, family=poisson, data=dat, offset=log(n))

The summary of my and m2 shows that summary(m1) = summary(m2) but if I try to calculate the McFadden through the pR2 (pscl package): McFadden(m1) ≠ McFadden(m2).

Does someone have an explanation for that?

Improve this question
$\endgroup$
3
  • 3
    $\begingroup$ Only you have dat - it'd be useful to illustrate the issue using one of the built-in datasets, as well as providing the rest of the code to calculate McFadden's $R^2$. See stackoverflow.com/help/mcve. $\endgroup$ Commented Jun 24, 2015 at 9:25
  • $\begingroup$ See this post for more information: [stackoverflow.com/questions/31022213/… [1]: stackoverflow.com/questions/31022213/… $\endgroup$ Commented Jun 24, 2015 at 10:39
  • 5
    $\begingroup$ As explained in help center, please don't crosspost but pick one best place for your question. It's also no use pointing to another post that makes the same formatting errors in the hope of avoiding having to fix it here. Please fix the same errors in that other post. $\endgroup$ Commented Jun 24, 2015 at 10:50

2 Answers 2

Reset to default
5
$\begingroup$

Here is the source code of pscl:::pR2.glm:

function (object, ...) { llh <- logLik(object) objectNull <- update(object, ~1) llhNull <- logLik(objectNull) n <- dim(object$model)[1] pR2Work(llh, llhNull, n) } <environment: namespace:pscl>

If the offset is specified in the formula, it gets lost in the second line (update to compute the intercept-only model).

See this example:

library("foreign") ceb <- read.dta("http://data.princeton.edu/wws509/datasets/ceb.dta") ceb$y <- round(ceb$mean*ceb$n, 0) ceb$os <- log(ceb$n) m0 <- glm(y ~ res + offset(os), data=ceb, family=poisson) m1 <- glm(y ~ res, offset=os, data=ceb, family=poisson) all.equal(coef(m0), coef(m1)) # [1] TRUE ### compute null models coef(update(m0, ~1)) # wrong, offset not considered # (Intercept) # 5.02 coef(update(m1, ~1)) # (Intercept) # 1.376 coef(update(m0, ~1, offset=os)) # (Intercept) # 1.376
Improve this answer
$\endgroup$
-3
$\begingroup$

I did some tests and the offset seems not to be included in m2. The right way to go to include an offset in a glm is:

m1 <- glm(d ~ 1 + offset(log(n)), family=poisson, data=dat)

The problem seems not to be due to the pscl package but to the formulation of the glm.

Improve this answer
$\endgroup$
2
  • 2
    $\begingroup$ I would say that "the problem" is the mistake you made in writing your R code rather than "the formulation of the glm." $\endgroup$ Commented Jun 24, 2015 at 13:37
  • 1
    $\begingroup$ (-1) The answer by @rcs shows this is the wrong way to include an offset if you want pR2 to take note of it. $\endgroup$ Commented Jun 24, 2015 at 14:45

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.