自举分层/多级数据(重采样群集)

Question

我正在生成一个脚本,用于从cats数据集(从-MASS-包中)创建bootstrap样本.

根据Davidson和Hinkley教科书[1],我进行了一个简单的线性回归,并采用了一个基本的非参数程序,用于从iid观察引导,即对重新采样.

原始样本的形式如下:

Bwt   Hwt

2.0   7.0
2.1   7.2

...

1.9    6.8

通过单变量线性模型,我们想通过他们的大脑重量来解释猫的重量.

代码是:

library(MASS)
library(boot)


##################
#   CATS MODEL   #
##################

cats.lm <- glm(Hwt ~ Bwt, data=cats)
cats.diag <- glm.diag.plots(cats.lm, ret=T)


#######################
#   CASE resampling   #
#######################

cats.fit <- function(data) coef(glm(data$Hwt ~ data$Bwt)) 
statistic.coef <- function(data, i) cats.fit(data[i,]) 

bootl <- boot(data=cats, statistic=statistic.coef, R=999)

现在假设存在一个聚类变量cluster = 1, 2,..., 24(例如,每只猫属于一个给定的垃圾).为简单起见,假设数据是平衡的:我们对每个簇有6个观察值.因此,24窝中的每一只由6只猫组成(即n_cluster = 6和n = 144).

可以通过以下方式创建伪cluster变量:

q <- rep(1:24, times=6)
cluster <- sample(q)
c.data <- cbind(cats, cluster)

我有两个相关的问题:

如何根据(聚集的)数据集结构模拟样本？也就是说,如何在集群级别重新采样？我想用替换对集群进行采样,并在原始数据集中设置每个选定集群内的观测值(即,使用替换集群进行采样,而不替换每个集群内的观测值).

这是戴维森提出的策略(第100页).假设我们绘制B = 100样本.它们中的每一个应由24个可能的递归聚类(例如cluster = 3, 3, 1, 4, 12, 11, 12, 5, 6, 8, 17, 19, 10, 9, 7, 7, 16, 18, 24, 23, 11, 15, 20, 1)组成,并且每个聚类应包含与原始数据集相同的6个观察结果.怎么做R？(有或没有-boot-包裹.)你有其他建议继续吗？

第二个问题涉及初始回归模型.假设我采用固定效应模型,具有集群级拦截.它是否改变了采用的重采样程序？

[1] Davidson,AC,Hinkley,DV(1997).Bootstrap方法及其应用.剑桥大学出版社.

Answer 1

如果我理解正确,那么您正在尝试将其c.data作为输入:

• 使用替换重新取样群集
• 保持随机样本中每个群集与原始数据集(即数据集)中的点之间的关联
• 使用采样的集群创建引导程序

这是一个实现此目的的脚本,您可以将其包装到函数中重复R次,其中R是引导复制的数量

q <- rep(1:24, times=6)
cluster <- sample(q)
c.data <- cbind(cats, cluster)

# get a vector with all clusters
c <- sort(unique(c.data$cluster))

# group the data points per cluster
clust.group <- function(c) {
    c.data[c.data$cluster==c,]
}

clust.list <- lapply(c,clust.group)

# resample clusters with replacement
c.sample <- sample(c, replace=T)

clust.sample <- clust.list[c.sample]

clust.size <- 6

# combine the cluster list back to a single data matrix
clust.bind <- function(c) {
    matrix(unlist(c),nrow=clust.size)
}

c.boot <- do.call(rbind,lapply(clust.sample,clust.bind))

# Just to maintain columns name
colnames(c.boot) <- names(c.data)

# the new data set (single bootstrap replicate)
c.boot