根据面板数据中的索引进行块采样

Question

我有一个面板数据，即n 的每个观察值 (nxt) 的 t 行，例如

data("Grunfeld", package="plm")
head(Grunfeld)
firm year   inv  value capital
   1 1935 317.6 3078.5     2.8
   1 1936 391.8 4661.7    52.6
   1 1937 410.6 5387.1   156.9
   2 1935 257.7 2792.2   209.2
   2 1936 330.8 4313.2   203.4
   2 1937 461.2 4643.9   207.2

我想进行块引导，即我想用替换重新采样，在观察到它的所有年份取一个公司 [i]。例如，如果year=1935:1937 和公司 1 是随机抽取的，我希望公司 [1] 将在新样本中出现 3 次，对应于year=1935:1937。如果重新绘制，则必须再次绘制 3 次。此外，我需要将我自己的函数应用到新的引导样本中，并且我需要这样做 500 次。 我当前的代码是这样的：

library(boot)
boot.fun <- function(data) {
   est.boot = myfunction(y=Grunfeld$v1, x=Grunfeld$v2, other parameters)
   return(est.boot)
}
boot.sim <- function(data, mle) {
data =  sample(data, ?? ) #
return(data)
}

start.time = Sys.time()
result.boot <- boot(Grunfeld, myfunction( ... ), R=500, sim = "parametric",  
               ran.gen = boot.sim)
Sys.time() - start.time

我正在考虑通过以正确的方式指定 data = sample(data, ?? ) 来重新采样，因为它工作平稳且干净，使用列 firm 作为索引。我怎么能那样做？还有其他更有效的选择吗？

编辑。 我不一定需要新的 boot.function。我只需要一个（可能是快速的）代码，它允许通过替换重新采样，然后我将它作为ran.gen=code.which.works 放在boot 参数中。 输出应该是与原始样本相同维度的样本，即使公司可以被随机挑选两次或更多次（或不被挑选）。例如结果可能是

head(GrunfeldResampled)
firm year   inv  value capital
   2 1935 257.7 2792.2   209.2
   2 1936 330.8 4313.2   203.4
   2 1937 461.2 4643.9   207.2
   1 1935 317.6 3078.5    2.8
   1 1936 391.8 4661.7    52.6
   1 1937 410.6 5387.1   156.9
   2 1935 257.7 2792.2   209.2
   2 1936 330.8 4313.2   203.4
   2 1937 461.2 4643.9   207.2
   9 1935 317.6 3078.5   122.8
   9 1936 391.8 4661.7   342.6
   9 1937 410.6 5387.1   156.9

基本上，我需要将每个公司视为block，因此重采样应适用于整个区块。希望这可以澄清

Answer 1

显然，在这个答案中，每家公司都被观察了 20 年，所以我不会有任何问题来证明：

data("Grunfeld", package="plm") #load data

解决方案

#n is the the firms column, df is the dataframe
myfunc <- function(n,df) {      #define function
 unique_firms <- unique(n)      #unique firms
 sample_firms <- sample(unique_firms, size=length(unique_firms), replace=T ) #choose from unique firms randomly with replacement
 new_df <- do.call(rbind, lapply(sample_firms, function(x)  df[df$firm==x,] ))  #fetch all years for each randomly picked firm and rbind
}

a <- myfunc(Grunfeld$firm, Grunfeld) #run function 

输出

> str(a)
'data.frame':   200 obs. of  5 variables:
 $ firm   : int  4 4 4 4 4 4 4 4 4 4 ...
 $ year   : int  1935 1936 1937 1938 1939 1940 1941 1942 1943 1944 ...
 $ inv    : num  40.3 72.8 66.3 51.6 52.4 ...
 $ value  : num  418 838 884 438 680 ...
 $ capital: num  10.5 10.2 34.7 51.8 64.3 67.1 75.2 71.4 67.1 60.5 ...

如您所见，dim 与输入的 data.frame 完全相同

对于您的数据，解决方案将是：

myfunc <- function(n,df) {      #define function
  unique_firms <- unique(n)      #unique firms
  print(unique_firms)
  sample_firms <- sample(unique_firms, size=length(unique_firms), replace=T ) #choose from unique firms randomly with replacement
  new_df <- do.call(rbind, lapply(sample_firms, function(x)  df[df$country==x,] ))  #fetch all years for each randomly picked firm and rbind
}

和输出：

> str(a)
'data.frame':   848 obs. of  18 variables:
 $ isocode  : Factor w/ 106 levels "AGO","ALB","ARG",..: 82 82 82 82 82 82 82 82 61 61 ...
 $ time     : int  2 3 4 5 6 7 8 9 2 3 ...
 $ country  : num  80 80 80 80 80 80 80 80 59 59 ...
 $ year     : int  1975 1980 1985 1990 1995 2000 2005 2010 1975 1980 ...
 $ gdp      : num  184619 210169 199343 268870 305255 ...
 $ pop      : num  33.4 34.9 36.6 37.8 38.3 ...
 $ gdp_k    : num  5526 6022 5443 7117 7969 ...
 $ co2      : num  340353 431436 426881 431052 350874 ...
 $ co2_k    : num  10191 12333 11674 11407 9128 ...
 $ oecd     : int  1 1 1 1 1 1 1 1 1 1 ...
 $ LI       : int  0 0 0 0 0 0 0 0 0 0 ...
 $ LMI      : int  0 0 0 0 0 0 0 0 0 0 ...
 $ UMI      : int  0 0 0 0 0 0 0 0 0 0 ...
 $ HI       : int  1 1 1 1 1 1 1 1 1 1 ...
 $ gdpk     : num  5531 6018 5449 7118 7971 ...
 $ co2k     : num  10196 12355 11668 11412 9162 ...
 $ co2_k.lag: num  8595 10191 12333 11674 11407 ...
 $ gdp_k.lag: num  4730 5526 6022 5443 7117 ...

Answer 2

您可以使用引导函数的“strata”参数来执行此操作。这称为分层引导。 更改代码的最后一行：

result.boot <- boot(Grunfeld, boot.fun, R=500, sim = "ordinary",  
                strata = Grunfeld$firm)

我抑制了参数ran.gen & sim

我建议对引导功能进行这些更改，使其正常工作：

boot.fun <- function(d, i) { # d being your data, i the set of indices)
   est.boot = myfunction(y=d[i ,]$v1, x=d[i, ]$v2, other parameters)
   return(est.boot)
}