glm.nb confidence bands around regression line

[Rob Sparks]

Hi I would like to create 95% confidence bands around a prediction response curve.  When I plot the confidence bands it doesn’t encompass the whole lcl and ucl of the confidence band.  Also the prediction line is about 0.5 units from the actual prediction response value.

 

I am not sure where I went wrong?

 

Here is the code I used and my data is attached (text,comma delimited).  The MASS package needs to be loaded for glm.nb.

 

 

vespnb5<-glm.nb(SumOftamano~grassNew+grassNew2+OtherNew,data=VESP)

 

grid<-with(VESP,expand.grid(grassNew=c(0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8),grassNew2=c(0.01,0.04,0.09,0.16,0.25,0.36,0.49,0.64),OtherNew=mean(OtherNew)))

 

grid$SumOftamano<-stats::predict(vespnb5,type="response",newdata=grid)

err<-stats::predict(vespnb5,newdata=grid,type="response",se=TRUE)

grid$ucl<-err$fit + 1.96*err$se.fit

grid$lcl<-err$fit - 1.96*err$se.fit

 

p<-ggplot(VESP,aes(grassNew,SumOftamano))

p+geom_smooth(aes(ymin=lcl,y=ucl),data=grid)

 

I would appreciate any ideas on how to correctly graph the predicted response curve with confidence bands.

 

 

Rob Sparks

Research Biologist/GIS Manager

Rocky Mountain Bird Observatory

230 Cherry St. Suite 150

Fort Collins, CO 80521

970-482-1707 ext: 19

www.rmbo.org

 

 

 

 

[Rob Sparks]

Ok I should have been more careful. Here is the correct code for the glm.nb
part.

vespnb5<-glm.nb(SumOftamano~grassNew+I(grassNew^2)+OtherNew,data=VESP)
##The second variable is the quadratic of grassNew.

p+geom_smooth(aes(min=lcl,max=ucl),data=grid,stat="identity") ##I replaced
the ymin with min and y with max and now I get the

#correct response curve at the right value.

Now I just need to adjust the predicted response 95% confidence band to
include the lcl and ucl values.

I am using R 2.8.1

Any ideas,
Rob


-----Original Message-----
From: s...@z107.de [mailto:s...@z107.de]
Sent: Tuesday, June 02, 2009 2:00 PM
To: Rob Sparks
Subject: Re: glm.nb confidence bands around regression line

hi Rob,

I don't think I can help with your problem, but I saw, that grassNew2
doesn't exist in your attached data file, so will have a problem to
reproduce your example.

regards,
stefan

> <http://www.rmbo.org> www.rmbo.org
>
>
>
>
>
>
>
>
>
>
> >

> "SumOftamano","bareNew","grassNew","OtherNew","shrub"
> 1,0.30,0.53,0.15,0.01
> 0,0.21,0.58,0.19,0.04
> 0,0.74,0.06,0.19,0.12
> 0,0.74,0.14,0.10,0.07
> 0,0.58,0.00,0.41,0.17
> 0,0.71,0.03,0.24,0.11
> 0,0.51,0.28,0.20,0.06
> 0,0.44,0.30,0.24,0.07
> 0,0.64,0.22,0.13,0.04
> 0,0.50,0.11,0.37,0.00
> 0,0.63,0.19,0.16,0.04
> 0,0.63,0.19,0.17,0.02
> 2,0.30,0.29,0.39,0.02
> 1,0.56,0.22,0.20,0.04
> 0,0.25,0.30,0.44,0.03
> 0,0.50,0.11,0.38,0.06
> 0,0.69,0.10,0.19,0.03
> 0,0.15,0.59,0.25,0.03
> 1,0.41,0.20,0.38,0.03
> 0,0.65,0.04,0.30,0.10
> 0,0.07,0.80,0.11,0.02
> 0,0.55,0.07,0.36,0.12
> 0,0.72,0.10,0.17,0.14
> 1,0.16,0.61,0.21,0.10
> 0,0.42,0.07,0.50,0.03
> 0,0.59,0.20,0.19,0.01
> 0,0.78,0.17,0.03,0.03
> 2,0.50,0.19,0.29,0.08
> 0,0.55,0.18,0.25,0.19
> 3,0.38,0.21,0.40,0.08
> 6,0.82,0.09,0.07,0.08
> 1,0.45,0.16,0.38,0.03
> 0,0.40,0.08,0.51,0.05
> 1,0.56,0.24,0.19,0.01
> 5,0.10,0.48,0.40,0.00
> 5,0.31,0.59,0.08,0.00
> 0,0.23,0.29,0.46,0.01
> 16,0.47,0.29,0.22,0.09
> 0,0.41,0.30,0.27,0.18
> 0,0.27,0.40,0.31,0.10
> 4,0.33,0.38,0.27,0.01
> 0,0.25,0.18,0.56,0.04
> 0,0.36,0.12,0.50,0.03
> 0,0.20,0.17,0.61,0.08
> 0,0.16,0.75,0.07,0.03
> 0,0.10,0.78,0.11,0.07
> 1,0.09,0.78,0.12,0.04
> 0,0.13,0.64,0.22,0.05
> 0,0.12,0.65,0.22,0.05
> 0,0.14,0.64,0.21,0.07
> 0,0.30,0.39,0.30,6.25e-03
> 0,0.52,0.11,0.36,0.08
> 0,0.49,0.10,0.40,0.04
> 0,0.43,0.09,0.46,0.07
> 0,0.43,0.10,0.45,0.03
> 0,0.59,0.10,0.30,0.05
> 0,0.50,0.35,0.13,0.07
> 0,0.64,0.24,0.10,0.02
> 0,0.46,0.33,0.20,0.05
> 0,0.93,0.02,0.04,0.00
> 0,0.34,0.49,0.16,0.03
> 0,0.85,0.09,0.04,0.01
> 0,0.81,0.05,0.12,0.00
> 0,0.80,0.11,0.08,0.00
> 0,0.40,0.42,0.16,0.00
> 0,0.97,0.01,0.01,0.06
> 0,0.76,0.16,0.07,0.01
> 0,0.95,0.00,0.04,0.03
> 0,0.84,0.03,0.11,0.00
> 0,0.76,0.17,0.05,0.00
> 0,0.77,0.11,0.10,1.38888888888889e-03
> 0,0.94,0.01,0.03,0.01
> 0,0.66,0.09,0.23,9.02777777777778e-03
> 0,0.88,0.03,0.08,9.02777777777778e-03
> 0,0.88,0.01,0.09,0.00
> 0,0.75,0.23,0.01,0.00
> 0,0.67,0.26,0.06,6.25e-03
> 0,0.65,0.21,0.12,0.00
> 0,0.72,0.27,4.16088765603329e-03,0.00
> 0,0.68,0.23,0.08,0.00
> 0,0.52,0.07,0.39,0.02
> 2,0.39,0.20,0.40,0.04
> 0,0.77,0.10,0.11,4.16666666666667e-03
> 0,0.84,0.09,0.06,2.77777777777778e-03
> 0,0.82,0.04,0.13,1.38888888888889e-03
> 0,0.83,0.13,0.02,0.04
> 0,0.78,0.10,0.11,5.55555555555556e-03
> 0,0.59,0.09,0.30,6.94444444444444e-04
> 0,0.55,0.19,0.24,0.00
> 0,0.77,0.14,0.08,0.00
> 0,0.75,0.06,0.18,0.00
> 0,0.29,0.39,0.31,0.00
> 0,0.70,0.15,0.14,0.00
> 0,0.62,0.24,0.13,0.00
> 0,0.78,0.03,0.18,0.00
> 0,0.92,0.03,0.03,2.08333333333333e-03
> 0,0.72,0.09,0.17,0.00
> 0,0.54,0.24,0.21,0.00
> 0,0.45,0.25,0.28,0.00
> 0,0.66,0.11,0.21,2.08333333333333e-03
> 0,0.73,0.07,0.19,0.00
> 0,0.81,0.08,0.10,6.94444444444444e-03
> 0,0.81,0.07,0.10,0.00
> 0,0.68,0.07,0.24,0.00
> 0,0.76,0.03,0.20,0.00
> 0,0.62,0.16,0.20,0.00
> 5,0.13,0.65,0.20,0.02
> 9,0.11,0.56,0.32,0.03
> 3,0.17,0.36,0.45,0.02
> 3,0.17,0.36,0.45,0.02
> 8,0.17,0.44,0.38,0.00
> 9,0.16,0.40,0.42,2.08333333333333e-03
> 7,0.10,0.48,0.40,0.00
> 17,0.11,0.52,0.35,0.03
> 10,0.22,0.42,0.35,0.02
> 6,0.09,0.49,0.40,0.04
> 3,0.24,0.22,0.53,7.63888888888889e-03
> 5,0.25,0.24,0.50,0.02
> 8,0.02,0.65,0.32,0.02
> 7,0.17,0.61,0.20,1.38888888888889e-03
> 10,0.07,0.62,0.30,0.00
> 5,0.17,0.24,0.57,0.01
> 2,0.28,0.21,0.49,0.02
> 0,0.19,0.25,0.55,1.38888888888889e-03
> 3,0.21,0.22,0.55,0.01
> 1,0.31,0.25,0.42,0.00
> 6,0.05,0.53,0.41,0.00
> 2,0.18,0.42,0.39,0.02
> 2,0.12,0.47,0.40,0.00
> 3,0.05,0.40,0.53,0.00
> 6,0.11,0.41,0.46,0.01
> 5,0.08,0.35,0.55,0.00
> 3,0.22,0.32,0.45,0.08
> 0,0.11,0.44,0.44,2.08333333333333e-03
> 3,0.14,0.38,0.46,6.25e-03
> 4,0.40,0.34,0.25,0.17
> 5,0.13,0.51,0.35,0.03
> 2,0.25,0.55,0.19,0.00
> 0,0.47,0.22,0.29,0.13
> 0,0.43,0.10,0.45,0.23
> 2,0.21,0.54,0.23,0.05
> 0,0.34,0.19,0.46,4.16666666666667e-03
> 28,0.13,0.48,0.38,0.06
> 6,0.23,0.32,0.43,0.07
> 0,0.26,0.32,0.40,0.07
> 0,0.27,0.49,0.23,0.11
> 15,0.08,0.27,0.64,0.05
> 1,0.26,0.33,0.40,0.06
> 3,0.09,0.66,0.23,0.03
> 5,0.14,0.36,0.48,0.15
> 2,0.01,0.72,0.25,0.11
> 3,0.13,0.48,0.37,0.05
> 0,0.32,0.17,0.49,2.77777777777778e-03
> 0,0.36,0.30,0.32,0.00
> 1,0.11,0.48,0.40,0.06
> 0,0.17,0.51,0.30,0.05
> 0,0.08,0.45,0.45,0.01
> 1,0.18,0.53,0.27,9.72222222222222e-03
> 0,0.25,0.53,0.21,0.06
> 0,0.36,0.30,0.32,0.07
> 0,0.31,0.30,0.38,0.00
> 0,0.23,0.27,0.49,2.77777777777778e-03
> 4,0.21,0.42,0.35,0.06
> 1,0.53,0.41,0.04,0.01
> 1,0.28,0.40,0.30,0.05
> 8,0.40,0.42,0.17,0.01
> 0,0.60,0.39,-2.82885431400283e-03,0.01
> 0,0.39,0.48,0.11,0.04
> 3,0.29,0.66,0.04,0.15
> 1,0.27,0.27,0.45,0.16
> 0,0.49,0.37,0.13,0.03
> 0,0.17,0.24,0.58,0.03
> 1,0.45,0.40,0.13,0.03
> 2,0.46,0.57,-0.03,0.15
> 0,0.44,0.51,0.04,0.04
> 4,0.64,0.32,0.03,0.02
> 2,0.50,0.38,0.11,3.47222222222222e-03
> 6,0.37,0.37,0.25,0.00
> 0,0.54,0.37,0.07,0.07
> 3,0.50,0.29,0.19,0.05
> 0,0.70,0.33,-0.03,0.05
> 1,0.78,0.16,0.04,0.02
> 0,0.37,0.46,0.15,6.25e-03
> 9,0.57,0.36,0.05,4.16666666666667e-03
> 0,0.37,0.57,0.05,0.06
> 2,0.33,0.57,0.08,8.33333333333333e-03
> 2,0.73,0.23,0.03,9.72222222222222e-03
> 1,0.45,0.29,0.24,0.01
> 0,0.51,0.45,0.03,4.16666666666667e-03
> 0,0.43,0.44,0.12,3.47222222222222e-03
> 7,0.40,0.45,0.13,0.01
> 10,0.55,0.21,0.23,0.04
> 0,0.69,0.03,0.27,0.00
> 0,0.79,0.03,0.16,0.00
> 0,0.74,0.11,0.14,0.00
> 0,0.77,0.08,0.13,0.00
> 0,0.71,0.07,0.20,0.00
> 0,0.59,0.20,0.19,0.00
> 0,0.76,0.05,0.17,0.00
> 0,0.73,0.04,0.21,6.94444444444444e-04
> 0,0.81,0.05,0.13,6.94444444444444e-04
> 0,0.74,0.09,0.15,0.00
> 2,0.27,0.40,0.31,0.00
> 0,0.20,0.43,0.36,0.00
> 0,0.11,0.62,0.26,0.00
> 2,0.28,0.53,0.17,1.38888888888889e-03
> 0,0.13,0.61,0.24,8.33333333333333e-03
> 2,0.30,0.22,0.46,2.77777777777778e-03
> 0,0.14,0.63,0.22,0.00
> 0,0.18,0.81,2.08768267223382e-03,2.08333333333333e-03
> 0,0.29,0.63,0.07,0.04
> 2,0.46,0.34,0.19,0.00
> 3,0.22,0.53,0.24,0.02
> 0,0.38,0.36,0.24,0.01
> 0,0.52,0.14,0.33,2.08333333333333e-03
> 0,0.76,0.06,0.17,0.00
> 0,0.79,0.06,0.14,0.00
> 0,0.69,0.14,0.15,0.00
> 0,0.64,0.09,0.25,0.00
> 0,0.70,0.10,0.18,6.94444444444444e-03
> 4,0.21,0.45,0.32,8.33333333333333e-03
> 2,0.10,0.36,0.53,0.01
> 0,0.27,0.05,0.66,0.07
> 0,0.36,0.11,0.51,0.04
> 0,0.32,0.23,0.43,0.03
> 2,0.16,0.37,0.45,4.86111111111111e-03
> 1,0.08,0.64,0.27,7.63888888888889e-03
> 1,0.32,0.36,0.31,0.05
> 0,0.06,0.51,0.41,6.25e-03
> 1,0.10,0.33,0.55,0.01
> 2,0.13,0.49,0.37,0.00
> 2,0.16,0.34,0.49,0.01
> 0,0.47,0.32,0.20,0.02
> 6,0.27,0.10,0.61,0.00
> 1,0.38,0.28,0.32,0.04
> 2,0.40,0.16,0.43,0.05
> 9,0.16,0.36,0.47,6.94444444444444e-04
> 0,0.43,0.16,0.40,0.04
> 0,0.06,0.19,0.73,0.02
> 0,0.17,0.17,0.64,0.06
> 0,0.11,0.28,0.59,0.04
> 0,0.33,0.29,0.37,0.14
> 2,0.17,0.45,0.37,0.02
> 1,0.20,0.41,0.37,0.09
> 0,0.17,0.55,0.27,0.05
> 0,0.20,0.35,0.44,0.02
> 2,0.18,0.31,0.49,0.01
> 8,0.31,0.42,0.25,0.01
> 0,0.34,0.37,0.27,0.06
> 0,0.06,0.38,0.55,0.02
> 0,0.20,0.18,0.60,0.03
> 0,0.33,0.19,0.46,2.77777777777778e-03
> 3,0.08,0.55,0.35,6.25e-03
> 0,0.04,0.72,0.23,0.02
> 0,0.17,0.39,0.42,0.02
> 6,0.26,0.35,0.37,0.01
> 0,0.02,0.48,0.49,0.02
> 1,0.10,0.62,0.27,7.63888888888889e-03
> 0,0.20,0.46,0.32,0.07
> 2,0.05,0.65,0.29,4.86111111111111e-03
> 0,0.04,0.56,0.38,0.04
> 0,0.12,0.43,0.44,0.02
> 1,0.50,0.31,0.17,0.07
> 0,0.36,0.23,0.40,0.10
> 5,0.51,0.30,0.18,0.05
> 6,0.38,0.40,0.21,0.02
> 4,0.56,0.23,0.19,9.02777777777778e-03
> 0,0.37,0.37,0.25,0.02
> 0,0.23,0.48,0.28,0.02
> 3,0.30,0.21,0.48,0.07
> 5,0.20,0.22,0.57,0.07
> 7,0.21,0.36,0.41,9.02777777777778e-03
> 2,0.14,0.39,0.46,0.06
> 4,0.16,0.34,0.49,0.04
> 2,0.01,0.78,0.20,0.07
> 5,0.06,0.63,0.30,0.04
> 2,0.06,0.68,0.25,0.06
> 3,0.02,0.59,0.38,0.01
> 1,0.08,0.57,0.33,0.09
> 3,0.02,0.58,0.38,0.04
> 6,0.07,0.65,0.26,0.00
> 0,0.11,0.21,0.66,0.00
> 0,0.14,0.23,0.62,0.00
> 4,0.28,0.29,0.41,0.01
> 4,0.24,0.41,0.34,0.04
> 7,0.09,0.34,0.55,0.00
> 3,0.11,0.74,0.14,1.38888888888889e-03
> 3,0.04,0.66,0.28,7.63888888888889e-03
> 0,0.17,0.45,0.37,1.38888888888889e-03
> 5,0.13,0.39,0.46,9.72222222222222e-03
> 6,0.06,0.58,0.35,6.25e-03
> 0,0.12,0.63,0.24,0.01
> 5,0.04,0.49,0.46,0.00
> 3,0.11,0.39,0.48,0.03
> 5,0.22,0.35,0.42,0.04
> 3,0.18,0.38,0.42,0.02
> 2,0.12,0.50,0.36,0.08
> 1,0.11,0.43,0.45,0.02
> 3,0.10,0.65,0.23,0.02
> 2,0.16,0.50,0.33,0.01
> 3,0.07,0.44,0.47,0.01
> 4,0.10,0.47,0.41,0.00
> 3,0.11,0.51,0.36,0.04
> 7,0.02,0.37,0.60,0.00
> 1,0.04,0.47,0.48,0.02
> 6,0.05,0.64,0.29,0.02
> 4,0.05,0.43,0.50,0.03
> 6,0.16,0.37,0.46,0.09
> 2,0.08,0.27,0.63,0.05
> 7,0.09,0.49,0.41,0.03
> 9,0.03,0.53,0.43,0.03
> 0,0.08,0.65,0.26,4.16666666666667e-03
> 0,0.17,0.48,0.34,0.09
> 0,0.08,0.31,0.60,0.13
> 0,0.19,0.44,0.36,0.06
> 15,0.20,0.27,0.52,6.25e-03
> 6,0.06,0.54,0.38,0.02
> 0,0.09,0.56,0.34,5.55555555555556e-03
> 5,0.23,0.43,0.32,5.55555555555556e-03
> 3,0.18,0.25,0.55,0.05
> 3,0.09,0.46,0.44,0.05
> 7,0.03,0.79,0.17,0.01
> 1,0.09,0.31,0.58,0.10
> 0,0.12,0.23,0.63,0.03
> 1,0.09,0.39,0.51,0.05
> 0,0.09,0.46,0.44,0.01
> 1,0.11,0.43,0.44,0.02
> 1,0.16,0.33,0.50,0.01
> 4,0.12,0.33,0.53,0.01
> 6,0.03,0.38,0.57,0.01
> 0,0.05,0.43,0.51,9.72222222222222e-03
> 6,0.10,0.31,0.58,0.06
> 3,0.17,0.29,0.53,0.05
> 1,0.18,0.36,0.44,9.72222222222222e-03
> 3,0.08,0.32,0.58,0.11
> 0,0.38,0.27,0.34,0.06
> 5,0.10,0.27,0.62,6.25e-03
> 1,0.18,0.17,0.63,0.09
> 3,0.12,0.36,0.51,5.55555555555556e-03
> 0,0.43,0.54,0.01,0.17
> 1,0.08,0.39,0.52,0.10
> 0,0.20,0.19,0.59,0.09
> 0,0.21,0.14,0.63,0.04
> 1,0.42,0.06,0.51,0.04
> 2,0.16,0.15,0.67,0.06
> 3,0.20,0.10,0.68,0.02
> 6,0.08,0.44,0.46,0.14
> 0,0.49,0.30,0.19,0.07
> 3,0.31,0.08,0.60,0.13
> 0,0.26,0.15,0.57,0.02
> 0,0.39,0.22,0.37,0.03
> 2,0.14,0.15,0.69,0.06
> 1,0.34,0.40,0.24,0.00
> 0,0.64,0.02,0.32,0.01
> 0,0.67,0.14,0.17,0.14
> 0,0.53,0.15,0.31,0.13
> 1,0.75,0.16,0.08,0.10
> 1,0.64,0.05,0.29,2.77777777777778e-03
> 0,0.34,0.10,0.55,1.38888888888889e-03
> 1,0.25,0.37,0.37,1.38888888888889e-03
> 0,0.92,0.08,-4.18410041841004e-03,4.16666666666667e-03
> 2,0.17,0.56,0.25,0.01
> 0,0.55,0.43,0.01,5.55555555555556e-03
> 0,0.62,0.23,0.13,0.05
> 6,0.20,0.31,0.47,2.77777777777778e-03
> 0,0.22,0.47,0.29,0.02
> 1,0.18,0.31,0.49,0.01
> 1,0.02,0.62,0.35,0.00
> 0,0.37,0.28,0.33,0.08
> 0,0.40,0.28,0.30,0.07

[hadley wickham]

Hi Rob,

Here's how I'd do it. I think the main thing you were missing is that
you need to calculate the confidence bands on the raw predictions and
then back transform to the original scale.

Hadley

VESP <- read.csv("VESP.txt")

library(MASS)
vespnb5 <- glm.nb(SumOftamano ~ poly(grassNew, 2) + OtherNew, data = VESP)

grid <- with(VESP, expand.grid(
grassNew = seq(0.1, 0.8, by = 0.1),
OtherNew = mean(OtherNew)
))

err <- predict(vespnb5, newdata=grid, se=TRUE)

grid$SumOftamano <- exp(err$fit)
grid$ucl <- with(err, exp(fit + 1.96 * se.fit))
grid$lcl <- with(err, exp(fit - 1.96 * se.fit))


library(ggplot2)
ggplot(VESP, aes(grassNew, SumOftamano)) +
geom_point() +
geom_smooth(aes(min = lcl, max = ucl), data = grid,stat = "identity")

--
http://had.co.nz/