代码之家 › 专栏 › 技术社区 › RInatM

R大数据分组引导统计

statistics-bootstrap data-manipulation data.table r

RInatM · 技术社区 · 11 年前

我想引导一个包含组的数据集。一个简单的场景是引导简单的方法:

data <- as.data.table(list(x1 = runif(200), x2 = runif(200), group = runif(200)>0.5))
stat <- function(x, i) {x[i, c(m1 = mean(x1), m2 = mean(x2)), by = "group"]}
boot(data, stat, R = 10)

这给了我错误 incorrect number of subscripts on matrix 因为 by = "group" 部分我设法用子集法解决了这个问题,但我不喜欢这个解决方案。有没有更简单的方法来完成这类任务?

特别是,我想在统计函数中引入一个额外的参数,比如 stat(x, i, groupvar) 并将其传递给引导函数,如 boot(data, stat(groupvar = group), R = 100) ?

3 回复 | 直到 9 年前

swihart 8 年前

使用

 boot       * 1.3-18  2016-02-23 CRAN (R 3.2.3)                        
 data.table * 1.9.7   2015-10-05 Github (Rdatatable/data.table@d607425)

我在使用OP代码时收到一个错误,答案由@eddi提供:

data <- as.data.table(list(x1 = runif(200), x2 = runif(200), group = runif(200)>0.5))
stat <- function(x, i) {x[i, c(m1 = mean(x1), m2 = mean(x2)), by = "group"]}
data[, list(list(boot(.SD, stat, R = 10))), by = group]$V1

生成错误消息:

Error in eval(expr, envir, enclos) : object 'group' not found

通过删除来修复错误 by=group 从函数 stat :

set.seed(1000)
data <- as.data.table(list(x1 = runif(200), x2 = runif(200), group = runif(200)>0.5))
stat <- function(x, i) {x[i, c(m1 = mean(x1), m2 = mean(x2))]}
data[, list(list(boot(.SD, stat, R = 10))), by = group]$V1

它会产生以下Bootstrap统计结果:

[[1]]

ORDINARY NONPARAMETRIC BOOTSTRAP


Call:
boot(data = .SD, statistic = stat, R = 10)


Bootstrap Statistics :
     original       bias    std. error
t1* 0.5158232  0.004930451  0.01576641
t2* 0.5240713 -0.001851889  0.02851483

[[2]]

ORDINARY NONPARAMETRIC BOOTSTRAP


Call:
boot(data = .SD, statistic = stat, R = 10)


Bootstrap Statistics :
     original        bias    std. error
t1* 0.5142383 -0.0072475030  0.02568692
t2* 0.5291694 -0.0001509404  0.02378447

下面,我修改了示例数据集,以突出显示哪个Bootstrap Statistics与哪个组列组合:

考虑组1和组2,组1的x1的平均值为10,组2的x2的平均值分别为10000和2000:

data2 <- as.data.table(list(x1 = c(runif(100, 9,11),runif(100, 1999,2001)), x2 = c(runif(100, 9999,10001),runif(100, 7999,8001)), group = rep(c(1,2), each=100)))
stat <- function(x, i) {x[i, c(m1 = mean(x1), m2 = mean(x2))]}
data2[, list(list(boot(.SD, stat, R = 10))), by = group]$V1

它给出了:

[[1]]

ORDINARY NONPARAMETRIC BOOTSTRAP


Call:
boot(data = .SD, statistic = stat, R = 10)


Bootstrap Statistics :
      original       bias    std. error
t1*   10.00907  0.007115938  0.04349184
t2* 9999.90176 -0.019569568  0.06160653

[[2]]

ORDINARY NONPARAMETRIC BOOTSTRAP


Call:
boot(data = .SD, statistic = stat, R = 10)


Bootstrap Statistics :
    original       bias    std. error
t1* 1999.965  0.031694179  0.06561209
t2* 8000.110 -0.006569872  0.03992401

eddi 11 年前

这应该做到:

data[, list(list(boot(.SD, stat, R = 10))), by = group]$V1

Community Reversed Engineer 7 年前

在进入分组部分之前,代码中有很多问题。

你的意思是这样的吗?

data <- as.data.frame(list(x1 = runif(200), x2 = runif(200), group = factor(sample(letters[1:2]))))
stat <- function(x, i)  c(m1 = mean(x$x1[i]), m2 = mean(x$x2[i]))

> stat(x,1:10)
       m1        m2 
0.4465738 0.5522221

然后你就可以担心了 doing it by group however you choose to .

例如:

library(plyr)
dlply( data, .(group), function( dat ) boot(dat, stat, R=10) )

对于更大的数据集,请尝试 data.table :

by( seq(nrow(data)), data$group, function(idx) myboot(data[idx,]))

我去了 by() 而不是 数据表 的 ,by= 参数,因为您希望输出为列表。可能有一些功能我不知道,但我找不到(请参阅编辑历史记录以了解它引起的问题)。

子设置仍然通过 数据表 的 [] 方法,所以它应该足够快。