How does one extract the final equation from glm poisson model?

Question

I have a Poisson model that is performing well. Now we need to put it into Java code and release it to the world. What is the equation that I plug the Poisson coefficients into?

Similar to this question: Find the equation from generalized linear model output

Here is some sample code from R that shows what I am doing:

d = data.frame(y=(10:100)^-1*100, x=10:100)
m = glm(y~x, data=d, family=poisson(link="log"))
plot(d$x, predict.glm(m, type="response"))
points(d$x, d$y, col="blue")
summary(m)

output:

Call:
glm(formula = y ~ x, family = poisson(link = "log"), data = d)

Deviance Residuals: 
     Min        1Q    Median        3Q       Max  
-0.34781  -0.25933  -0.05075   0.21286   1.20265  

Coefficients:
             Estimate Std. Error z value Pr(>|z|)    
(Intercept)  2.155676   0.126778  17.004   <2e-16 ***
x           -0.025912   0.002819  -9.191   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for poisson family taken to be 1)

    Null deviance: 104.5820  on 90  degrees of freedom
Residual deviance:   8.5623  on 89  degrees of freedom
AIC: Inf

Number of Fisher Scoring iterations: 4

> cor(d$y, predict.glm(m, type="response"))
[1] 0.9461024

That is pretty good! :)

Answer 1

The equation is

$$\log(\mu_i) = \beta_0 + \beta_1 x_i$$

where $\mu_i$ is the conditional expectation of $y_i$, $E(y | x)$, $\beta_0$ is the coefficient marked Intercept and $\beta_1$ the coefficient marked x. The $\log$ bit is the link function you specified. Hence to get actual predictions on the scale of your response data $y$, you need to apply the inverse of the link function (anti-log) to the both sides of the equation:

$$\mu_i = \exp(\beta_0 + \beta_1 x_i)$$

$\mu_i$ is then the predicted mean count given the value of x.

Print out the coefficients if you want them:

coef(m)

They may be a little more precise than the ones in the summary() output and so your Java code will be closer to the predictions given by predict().

The model doesn't look fantastic, despite the high linear correlation. There is bias throughout the range of x. Without knowing anything about the data, did you consider a model with $x$ and $x^2$?

Answer 2

This is great, and I have used it to verify part of a plot in which I am modeling the number of trees counted at certain elevations. However, I am using a zero-inflated Poisson model through zeroinfl() in pscl and the predicted values of zeroinfl() begin to deviate from the equation provided by Gavin after about 1000m in elevation. I assume this is because the probability of zero counts increases after about an elevation of 1000m, and so the zeroinfl() model begins to account for this increased zero count probability. Thats great and it does a good job of reflecting the actual observations, but I also need to know the equation for this line. What would be the zero-inflated Poisson version of this model? Something involving probabilities? Below is a plot in which the dots are the actual observation, the red is the predicted values based on zeroinfl() and the blue is the predicted valued based on the equation provided by Gavin and the coefficients of the zeroinfl() model. I am using a logarithmic representation of the y-axis to clarify the deviation. Included is a general representation of the code.

LogMega[,1] <- yjPower(LogMega[,1],0) #This is a log transform of the average density
m1 <- zeroinfl(Avg_dens~Z, data = LogMega) #the zero-inflated model
m2 <- 2.718^(coef(m1)[1]+(coef(m1)[2]*Z), data=LogMega) #the Poisson equation model
newdata1 <- expand.grid(LogMega[,n+3]) #create a table of the zeroinfl() predicted values
colnames(newdata1)<-"pred"
newdata1$resp <- predict(m1,newdata1)
p <- ggplot(LogMega,aes_string(x=names(LogMega)[n+3],y=LogMega[1]))
theme(
  geom_line(data=newdata1, aes(x=pred,y=resp), size=1, color="red")+
  geom_line(data=LogMega, aes(x=LogMega[,n+3], y=m2), size=1, color="blue"))