Ken,
This is expected because of Jensen's inequality; for any convex
function f() and a random variable Y, you have E(f(Y)) >= f(E(Y))
where in your case f() = exp() and Y = XB.
Kosuke
--
Department of Politics
Princeton University
Hi Zelig listers - I have discovered something
that is vexing me
about how zelig's creates point predictions. This happens for the
neg bin and Poisson models. E(Y) = exp(XB) but when I take the
mean of the s.out$qi$ev, it is consistently larger than exp(XB).
Here is an example, with the code at the bottom in case anyone
wants to recreate the output. Any thoughts appreciated!
Ken
### Test Zelig's predictions for neg
binomial and Poisson
### mean of ev seems to consistently OVERESTIMATE exp(XB)
library(foreign)
library(MASS)
library(Zelig)
## note: dataset is
http://www.kenbenoit.net/courses/quant2/
weede.dta
weede <- read.dta("weede.dta")
## replicate part of Table 3 from Benoit (1996), negative binomial
z.out <- zelig(butterw ~ poldem65, model="negbin", data=weede)
x.out <- setx(z.out, poldem65=c(0,20,55,85,100))
s.out <- sim(z.out, x=x.out, num=50000)
coef(z.out)
(Intercept) poldem65
0.506841 -0.007094
data.frame(expXB=exp(cbind(1,c(0,20,55,85,100))
%*% coef(z.out)),
+ zeligMean=apply(s.out$qi$ev, 2, mean),
+ zeligSE=apply(s.out$qi$ev, 2, sd))
expXB zeligMean zeligSE
1 1.6600 1.7281 0.5012
2 1.4405 1.4718 0.3090
3 1.1237 1.1361 0.1681
4 0.9083 0.9305 0.2061
5 0.8166 0.8486 0.2392
predict(z.out,
new=data.frame(poldem65=c(0,20,55,85,100)))
1 2 3 4
5
0.50684 0.36496 0.11666 -0.09617 -0.20258
## now with Poisson
z.out <- zelig(butterw ~ poldem65, model="poisson", data=weede)
x.out <- setx(z.out, poldem65=c(0,20,55,85,100))
s.out <- sim(z.out, x=x.out, num=50000)
coef(z.out)
(Intercept) poldem65
0.51672 -0.00729
data.frame(expXB=exp(cbind(1,c(0,20,55,85,100))
%*% coef(z.out)),
+ zeligMean=apply(s.out$qi$ev, 2, mean),
+ zeligSE=apply(s.out$qi$ev, 2, sd))
expXB zeligMean zeligSE
1 1.6765 1.7002 0.2942
2 1.4491 1.4597 0.1824
3 1.1227 1.1281 0.1088
4 0.9022 0.9129 0.1375
5 0.8087 0.8238 0.1565
..........................................................................
The code file is here:
### Test Zelig's predictions for neg binomial and Poisson
### mean of ev seems to consistently OVERESTIMATE exp(XB)
library(foreign)
library(MASS)
library(Zelig)
## note: dataset is
http://www.kenbenoit.net/courses/quant2/weede.dta
weede <- read.dta("weede.dta")
## replicate part of Table 3 from Benoit (1996), negative binomial
z.out <- zelig(butterw ~ poldem65, model="negbin", data=weede)
x.out <- setx(z.out, poldem65=c(0,20,55,85,100))
s.out <- sim(z.out, x=x.out, num=50000)
coef(z.out)
data.frame(expXB=exp(cbind(1,c(0,20,55,85,100)) %*% coef(z.out)),
zeligMean=apply(s.out$qi$ev, 2, mean),
zeligSE=apply(s.out$qi$ev, 2, sd))
predict(z.out, new=data.frame(poldem65=c(0,20,55,85,100)))
## now with Poisson
z.out <- zelig(butterw ~ poldem65, model="poisson", data=weede)
x.out <- setx(z.out, poldem65=c(0,20,55,85,100))
s.out <- sim(z.out, x=x.out, num=50000)
coef(z.out)
data.frame(expXB=exp(cbind(1,c(0,20,55,85,100)) %*% coef(z.out)),
zeligMean=apply(s.out$qi$ev, 2, mean),
zeligSE=apply(s.out$qi$ev, 2, sd))
-
Zelig Mailing List, served by Harvard-MIT Data Center
Send messages: zelig(a)lists.gking.harvard.edu
[un]subscribe Options:
http://lists.gking.harvard.edu/?info=zelig
Zelig program information:
http://gking.harvard.edu/zelig/
-
Zelig Mailing List, served by Harvard-MIT Data Center
Send messages: zelig(a)lists.gking.harvard.edu
[un]subscribe Options: