4

在我的统计课中，我们使用 Stata，因为我是 R 用户，所以我想在 R 中做同样的事情。我得到了正确的结果，但获得像置信区间这样简单的东西似乎有点尴尬。

这是我的粗略解决方案：

library(quantreg)
na = round(runif(100, min=127, max=144))
f <- rq(na~1, tau=.5, data=ds)
s <- summary.rq(f, se="boot", R=1000)
coef(s)[1]
coef(s)[1]+ c(-1,1)*1.96*coef(s)[2]

我还在引导包上做了一些实验，但我还没有让它工作：

library(boot)
b <- boot(na, function(w, i){ 
        rand_bootstrap_sample = w[i]
        f <- rq(rand_bootstrap_sample~1, tau=.5)
        return(coef(f)) 
        }, R=100)
boot.ci(b)

给出一个错误：

bca.ci 中的错误（boot.out，conf，index[1L]，L = L，t = to，t0 = t0.o，：估计调整“a”为 NA

我的问题：

  • 我不想知道是否有另一种更好的方法来获得置信区间
  • 为什么引导代码抱怨？
4

1 回答 1

4

您的示例没有为我提供错误消息（Windows 7/64，R 2.14.2），因此可能是随机种子的问题。因此，如果您使用某种随机方法发布示例，最好添加一行 set.seed; 见例子。

注意错误信息是指boot.ci的bca类型；因为这个经常抱怨，所以通过明确给出类型来取消选择它。

我不知道你为什么在引导程序中使用相当复杂的 rq 。如果您真的想分析 rq，请忘记下面的简单示例，但请提供更多详细信息。

library(boot)
set.seed(4711)
na = round(runif(100, min=127, max=144))

b <- boot(na, function(w, i) median(w[i]), R=1000)
boot.ci(b,type=c("norm","basic","perc"))
于 2012-04-17T07:19:01.987 回答