4

我已经尝试了这两种方法,但我发现两种方法都有困难。我试图更好地解释这是我的问题,然后再告诉你我对这两种方法的问题。

我有数据集“acceptances”,其中我有医院每天的接受次数,其中包含前面描述的独立变量。医院有三个地方可以进行访问。所以在我的数据集中,我每天有 3 行,每个地方一个。数据集看起来像:

Date        Place    NumerAccept    weekday month   NoConvention    Rain

2008-01-02  Place1        203       wed     Gen         0             1
2008-01-02  Place2         70       wed     Gen         0             1
2008-01-02  Place3          9       wed     Gen         0             1
2008-01-03  Place1        345       thu     Gen         0             1
2008-01-03  Place2         24       thu     Gen         0             1
2008-01-03  Place3         99       thu     Gen         0             1
2008-01-04  Place1        339       fri     Gen         0             0
2008-01-04  Place2         36       fri     Gen         0             0
2008-01-04  Place3        101       fri     Gen         0             0

.... 依此类推...直到昨天我才拥有数据集,所以最后三行是 2013 年 7 月 29 日昨天的接受情况。现在我进行泊松回归:

poisson_reg=glm(NumeberAccept ~ 1 + weekday + month + place + NoConvention + Rain, 
                    family = poisson(link = log), data = acceptances)

现在对于我的预测,我创建了一个新的数据集 Acceptances_2,我想从中计算未来 2 个月的接受数的预测间隔!!所以第一行是今天的录取数,最后一行是 9 月 29 日的录取数。


我不知道这个问题是否已经有了答案,但我无法找到它。我正在尝试在 R 中进行泊松回归,并且我想获得预测区间。我看到 predict 函数为lm它提供了写作'interval="prediction"',但它不适用于predict.glm!

有人知道是否有办法获得这些预测区间?如果你有一些例子,你能输入代码吗?

所以我必须计算医院每天接受的人数,我有以下代码:

poisson_reg=glm(NumeberAccept ~ 1 + weekday + month + place + NoConvention + Rain, 
                family = poisson(link = log), data = dataset)
summary(poisson_reg)

现在,如果我输入 R,predict(poisson_reg, newdata, type="responce")我可以预测每天的接受次数,但我也需要预测间隔!我看到,对于"lm"predict 调用中的类对象,您可以编写:predict(poisson_reg, newdata, interval="prediction")它给出了 95% 的预测区间。有没有办法用类的对象获得相同的结果"glm"

4

2 回答 2

5

这可能更像是一个统计问题而不是编程问题,但是:

从上一个问题中窃取示例数据:

ex <- read.table(
  header=TRUE, text='
Number.Accepted  Weekday    Month   Place
  20    6   8   1
  16    7   8   1
  12    4   8   2
  11    7   1   1
  12    1   4   1
  12    7   10  2
  13    5   6   2
')
ex.glm <- glm(Number.Accepted ~ Weekday + Month + Place,
              family = poisson, data = ex)

我们想要预测间隔的数据框:

newdata <- data.frame(Weekday=c(5,6),Month=c(9,9),Place=c(1,1))

像这样的东西:

bootSimFun <- function(preddata,fit,data) {
    bdat <- data[sample(seq(nrow(data)),size=nrow(data),replace=TRUE),]
    bfit <- update(fit,data=bdat)
    bpred <- predict(bfit,type="response",newdata=preddata)
    rpois(length(bpred),lambda=bpred)
}

您也可以replicate()从基础 R 使用,但plyr::raply()很方便:

library(plyr)
set.seed(101)
simvals <- raply(500,bootSimFun(preddata=newdata,fit=ex.glm,data=ex))
t(apply(simvals,2,quantile,c(0.025,0.975)))
##    2.5% 97.5%
## 1 7.000    40
## 2 7.475    36
于 2013-07-29T14:08:14.243 回答
3

考虑 Zelig 包。请参阅此处的泊松小插图 - http://rss.acs.unt.edu/Rdoc/library/Zelig/doc/poisson.pdf

Zelig 有一个统一的方法,不仅可以建模(为此,glm() 及其各种链接函数就足够了),还可以提取和绘制感兴趣的数量。特别是,要对预测范围(而不仅仅是预期范围)进行建模,您必须同时模拟您的系数(系统分量)和误差项(随机分量)。简单地平均误差项,我认为这是 predict.glm() 所做的,将为您提供更窄的预期范围。

Zelig 有一个函数 sim() 可以模拟系统和随机分量,并输出可用于绘制预测范围和预期范围的内存对象。它还有一个函数 setx(),如果你想模拟你的解释变量给定值的预测不确定性,你可以在 sim() 之前使用它。见这里——http://rss.acs.unt.edu/Rdoc/library/Zelig/html/setx.html

这一切都始于这篇论文: http: //gking.harvard.edu/files/abs/making-abs.shtml。Zelig 基本上是 Clarify 长大的。

于 2013-07-29T15:05:11.660 回答