2
$\begingroup$

I have obtained some strange results for the following data when conducting Poisson regression in R.

> RHT1b
   f x1 y1 y2
1 f1  0 35  1
2 f2  2 70  4
3 f3  0  5  1
4 f4  9 37  4
5 f5  0  3  0

> summary(amod2b)

Call:
glm(formula = x1 ~ y2, family = poisson, data = RHT1b)

Deviance Residuals: 
       1         2         3         4         5  
-0.00008  -1.71860  -0.00008   1.36550   0.00000  

Coefficients:
            Estimate Std. Error z value Pr(>|z|)
(Intercept)   -26.61    9977.85  -0.003    0.998
y2              7.08    2494.46   0.003    0.998

(Dispersion parameter for poisson family taken to be 1)

    Null deviance: 24.9766  on 4  degrees of freedom
Residual deviance:  4.8182  on 3  degrees of freedom
AIC: 15.485

Number of Fisher Scoring iterations: 18

I had an issue with the data earlier where x1 row 4 read "8", and row 5 read "1". When I ran the regression the results seemed reliable and were nearly significant. After correcting the data to how it is shown above and updating the model, I have results which literally do not make sense to me.

My question is, why has this happened? Is is possible the value of 1 changing to 0 exceeded some sort of threshold of zero counts that the Poisson model can't handle? Would a negative binomial model be more appropriate?

Improve this question
$\endgroup$
5
  • 1
    $\begingroup$ well, its either the effect of row 4 being different or row 5 being different... $\endgroup$
    – Spacedman
    Commented Nov 27, 2012 at 10:40
  • $\begingroup$ I understand that is what caused the change. What I am more concerned about though is how changing x1, row 4 and 5 by only 1, has changed the significance of the slope from nearly significant to not significant at all. It seems to me that the significance should have improved. $\endgroup$
    – J_L
    Commented Nov 27, 2012 at 10:44
  • $\begingroup$ Yes,and if you narrow it down to exactly which of those two changes is affecting it you might get some insight... $\endgroup$
    – Spacedman
    Commented Nov 27, 2012 at 10:49
  • 3
    $\begingroup$ You only have 5 data points. Changing 2 of them is changing 40% of your data! It looks like you have zero inflation, but with only 5 data points, it's impossible to tell. And pretty much any model will have overfitting issues. $\endgroup$
    – Peter Flom
    Commented Nov 27, 2012 at 12:25
  • 2
    $\begingroup$ The link function is the logarithm. In the first case, the three zeros for $x1$ reflect extraordinary uncertainty: they represent essentially all the possible values less than $\log(1)$, which is all negative numbers. The enormous standard errors in the output reflect that. Changing a single $0$ to $1$ is informative: it suggests the true mean at $y=0$ is around $(0+0+1)/3$, whose logarithm is around $-1$. In effect, an infinite band of uncertainty from $-\infty$ to $0$ shrinks to a rather small, finite range centered around $-1$. (Changing the $9$ to $8$ hardly matters.) $\endgroup$
    – whuber ♦
    Commented Nov 27, 2012 at 20:57

1 Answer 1

Reset to default
9
$\begingroup$

Plot the profile likelihood functions for the intercept & coefficient to see what's going on. Basically you can rule out high values for the intercept, & low values for the coefficient; but for any low value of the intercept you can find a high value of the coefficient with which your data will be consistent. If you want better estimates, you don't need to fit a different model to these data, but to collect more data.

In reply to your comment: Say for the y2 coefficient, take fixed values (something like 1, 2, ... , 8) & fit a GLM for each (you'll just be fitting the intercept). Take the log-likelihood for each GLM and plot against the value of the coefficient.

What you'll see is that there's no maximum to the log-likelihood - it plateaus. The tell-tale signs in your results are ...

(1) 18 iterations to get to the "maximum" likelihood fit - R is actually just calculating that the likelihood hardly changes between the last iterations.

(2) Ridiculously large standard errors for the coefficient estimates. This is because they're estimated from the rate of change of the slope at the "maximum", & are only any use when there is a maximum & the log-likelihood function is more or less a parabola. (Look up Wald tests.)

You still have a difference of 20 between the residual & null deviances, which is a lot on one degree of freedom - though it would probably be better to simulate than to assume a chi-square distribution. At least small or negative values of the y2 coefficient don't look plausible.

Note that zero-inflated or negative binomial models will add parameters, when your data aren't enough to pin down a Poisson !

Improve this answer
$\endgroup$
2
  • 1
    $\begingroup$ Thanks. Can you provide instructions on how to plot the "profile likelihood functions for the intercept and coefficients? $\endgroup$
    – J_L
    Commented Nov 27, 2012 at 18:32
  • 1
    $\begingroup$ +1. The end of the first paragraph--"collect more data"--is the best advice. The little "sensitivity analysis" of changing a single data value from $0$ to $1$ already reveals much about the limitations of this dataset. $\endgroup$
    – whuber ♦
    Commented Nov 27, 2012 at 21:03

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.