3

Assume the following easy example of a glm regression with an offset:

numberofdrugs <- rpois(84, 10)
healthvalue   <- rpois(84,75)
age           <- rnorm(84,50,5)
test          <- glm(healthvalue~age, family=poisson, offset=log(numberofdrugs))
summary(test)
fitted(test) # How to get one of these values manually?
  • How can I compute the fitted values manually?
  • Also, why is there no estimation of log(numberofdrugs)? (In the book Generalized Linear Models on page 205-207 there is an example where the offset is estimated. It was done to see if the coefficient is close to one. It's 0.903 (see page 207 if you've this classic book) and from this follows, that there is nearly a constant rate in the number of damage incident!)

Previous related questions asked:

gung - Reinstate Monica
  • 132,789
  • 81
  • 357
  • 650
MarkDollar
  • 5,575
  • 14
  • 44
  • 60
  • 1
    Please use code (four-spaces indent, `{}` button) instead of cite (`>` indent, `"` button). –  Jun 08 '11 at 10:48
  • 1
    Possible duplicate of [When to use an offset in a Poisson regression?](https://stats.stackexchange.com/questions/11182/when-to-use-an-offset-in-a-poisson-regression) – kjetil b halvorsen Sep 17 '17 at 13:06
  • "How can I compute the fitted values manually?", eg, is sufficiently distinct that I think this can be left open as a separate question. – gung - Reinstate Monica Sep 18 '17 at 17:15

3 Answers3

3

There should not be an estimate of the offset: this offset is (could be) different for every observation (the whole idea is that you monitor the number of events within a (linear) 'timemeasure' (here apparently numberofdrugs).

There is no one 'population' offset you could estimate: person 1 is going to have 5 drugs administered, person 2 maybe 10, and you assume that the number of events (healthvalue?) is linear with this person's numberofdrugs (as per the answer to your previous question).

I don't have the book at hand, and Amazon won't let me look at the pages you mention, but I suppose something else is happening there (maybe simply the average numberofdrugs in the population)?

Nick Sabbe
  • 12,119
  • 2
  • 35
  • 43
  • Ok that makes sense to me, but how are the fitted values calculated? – MarkDollar Jun 08 '11 at 12:54
  • 2
    Once you have the beta_0 and beta_1 from the model, then for any observation, given its `age` and `numberofdrugs`, you calculate `beta_0+beta_1*age+log(numberofdrugs)`. You can see the offset as simply another variable for which you fix the coefficient (well, of its logtransform) to 1 in the Poisson regression. – Nick Sabbe Jun 08 '11 at 13:13
2

(Much of this has already been addressed in the other answers and in other threads on the site, including your own previous question linked here. I will just address the remaining parts.)

There is no estimate of log(numberofdrugs) because that coefficient is forced to be $1$. That's what an offset is. On the other hand, you can estimate a coefficient for log(numberofdrugs) instead of forcing it to be $1$. (It may help to read my answer here: In a Poisson model, what is the difference between using time as a covariate or an offset?)

Moreover, you can test a fitted coefficient against any null value, not just $0$. Thus, you could estimate a coefficient for log(numberofdrugs) and test it against $1$, if that were the hypothesis of interest for a study. This should be fairly trivial with any statistical software. The Wald test for a coefficient from a GLiM is:
$$ z = \frac{(\hat\beta - {\rm null})}{SE} $$ The resulting $z$-statistic can be checked against a standard normal distribution. Here is a simple example, coded in R using your data: 

set.seed(3917)  # this makes the example exactly reproducible
numberofdrugs <- rpois(84, 10)
healthvalue   <- rpois(84,75)
age           <- rnorm(84,50,5)
test.coef     <- glm(healthvalue~age+log(numberofdrugs), family=poisson)
summary(test.coef)
# ...
# Coefficients:
#                      Estimate Std. Error z value Pr(>|z|)    
# (Intercept)         4.2275043  0.1482237  28.521   <2e-16 ***
# ...
# log(numberofdrugs)  0.0346964  0.0403885   0.859    0.390    
# ...
abs(0.0346964-1)/0.0403885
# [1] 23.90046       # <- the z-value
2*pnorm(23.90046, lower.tail=FALSE)
# [1] 3.029191e-126  # <- the p-value

This example shows that, on the one hand log(numberofdrugs) is not significantly related to healthvalue (which makes sense, given there is no relationship in how the data are generated). On the other hand, the fitted coefficient does significantly differ from $1$. That means that using log(numberofdrugs) as an offset is not sensible; there is not a consistent rate of healthvalue per unit of numberofdrugs.

gung - Reinstate Monica
  • 132,789
  • 81
  • 357
  • 650
1

About the practical part -- outputs of glm or summary are just lists which are pretty-printed for user convenience. You can see their full structure calling unclass on them and extract single values as usual, with a help of $, [[]] and [] operators.