[R] glm or transformation of the response?

[Greg Snow]

The default link function for the glm poisson family is a log link, which means that it is fitting the model:

log(mu) ~ b0 + b1 * x

But the data that you generate is based on a linear link.  Therefore your glm analysis does not match with how the data was generated (and therefore should not necessarily be the best fit).  Either analyze using glm and a linear link, or generate the data based on a log link (e.g. rpois(40, exp(seq(1,3, length.out=40))) ).

Hope this helps,

--
Gregory (Greg) L. Snow Ph.D.
Statistical Data Center
Intermountain Healthcare
greg.snow at imail.org
801.408.8111


> -----Original Message-----
> From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-
> project.org] On Behalf Of Christoph Scherber
> Sent: Tuesday, November 25, 2008 7:53 AM
> To: r-help at stat.math.ethz.ch
> Subject: [R] glm or transformation of the response?
>
> Dear all,
>
> For an introductory course on glm?s I would like to create an example
> to show the difference between
> glm and transformation of the response. For this, I tried to create a
> dataset where the variance
> increases with the mean (as is the case in many ecological datasets):
>
>         poissondata=data.frame(
>         response=rpois(40,1:40),
>         explanatory=1:40)
>
>         attach(poissondata)
>
> However, I have run into a problem because it looks like the lm model
> (with sqrt-transformation)
> fits the data best:
>
> ##
>
> model1=lm(response~explanatory,poissondata)
> model2=lm(sqrt(response+0.5)~explanatory,poissondata)
> model3=lm(log(response+1)~explanatory,poissondata)
> model4=glm(response~explanatory,poissondata,family=poisson)
> model5=glm(response~explanatory,poissondata,family=quasipoisson)
> model6=glm.nb(response~explanatory,poissondata)
> model7=glm(response~explanatory,quasi(variance="mu",link="identity"))
>
>
> plot(explanatory,response,pch=16)
> lines(explanatory,predict(model1,explanatory=explanatory))
> lines(explanatory,(predict(model2,explanatory=explanatory))^2-
> 0.5,lty=2)
> lines(explanatory,exp(predict(model3,explanatory=explanatory))-1,lty=3)
> lines(explanatory,exp(predict(model5,explanatory=explanatory)),lty=1,co
> l="red")
> lines(explanatory,predict(model6,explanatory=explanatory,type="response
> "),lty=1,col="blue")
> lines(explanatory,predict(model7,explanatory=explanatory,type="response
> "),lty=1,col="green")
>
> ##
>
> The only model that performs equally well is model7.
>
> How would you deal with this kind of analysis? What would be your
> recommendation to the students,
> given the fact that most of the standard glm models obviously don?t
> seem to produce good fits here?
>
> Many thanks and best wishes
> Christoph
>
> (using R 2.8.0 on Windows XP)
>
>
>
>
>
> --
> Dr. rer.nat. Christoph Scherber
> University of Goettingen
> DNPW, Agroecology
> Waldweg 26
> D-37073 Goettingen
> Germany
>
> phone +49 (0)551 39 8807
> fax   +49 (0)551 39 8806
>
> Homepage http://www.gwdg.de/~cscherb1
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-
> guide.html
> and provide commented, minimal, self-contained, reproducible code.