【问题标题】:Resampling with replacement替换重采样
【发布时间】:2015-10-14 21:44:37
【问题描述】:

我正在尝试为我的模型获取(引导的)输入数据。

源文件:https://www.dropbox.com/s/dudzxhozr50uhr7/EddyData_2010.csv?dl=0

library("dplyr")
library("readr")
library("reshape2")
library("ggplot2")

sub <- read_csv("EddyData_2010.csv", 
                col_types = list(col_integer(), col_integer(), col_double(),
                                 col_double(), col_double(), col_double(),
                                 col_double(), col_double(), col_double(),
                                 col_double(), col_double(), col_double()),
                col_names = c("Year", "DoY", "Hour", "NEE", "LE", "H", "Rg",
                              "Tair", "Tsoil", "rH", "Ustar", "VPD")) %>%
  filter(DoY == 170) %>%
  mutate(hour = 1:48) %>%
  select(NEE:hour)

# Number of resampling 
n_resempling <- 1000 
result_resampling <- NULL

# Do the resampling
for (i in 1:n_resempling) {
  result_resampling <- sample(48, 48, replace = T) %>%
    slice(sub, .) %>%
    mutate(resempling_number = i) %>%
    bind_rows(. , result_resampling)
}

这会生成一个带有替换的重采样,例如

输出显示在一天的 48 个半小时内重新采样了 1000 个引导程序并进行了替换。

这是我的问题:

带替换的重采样在一天中随机混杂了半个小时,没有任何控制。例如,我不想把晚上的半小时和白天的半小时混为一谈。结果导致事后进行奇怪的计算。是否有可能以我定义允许和不允许的方式来编码?

例子:

  • 只允许从晚上 10 点到下午 5 点重新采样
  • 夜间时间不能与白天时间重新采样(反之亦然)

【问题讨论】:

  • 在重新采样之前过滤您的原始数据。只想在上午 10 点到下午 5 点之间重新采样?过滤您的数据,使其仅包含该窗口中的时间,然后重新采样。只想要白天的时间?仅将您的数据过滤到一天,然后重新采样...

标签: r sampling statistics-bootstrap


【解决方案1】:

我为 CRD 设计做过简单的引导,但对时间数据没有做过。那是时间序列数据吗?据我了解,您希望仅在下午 2 点而不是下午 3 点对下午 2 点进行采样。了解抽样对于进行正确的分析很重要,因为在 R 中很容易出错。

我注意到您使用循环而不是包进行引导。 我使用 'boot' 包进行简单的引导,所以我求助于谷歌来查看时间数据。这是我发现的,我很抱歉这就是我所拥有的(由于缺乏代表我无法发表评论) 使用引导包我敢打赌任何东西都比使用循环更快 您可以使用 system.time( ) 进行检查,尤其是在您有大量数据的情况下。

https://stat.ethz.ch/R-manual/R-devel/library/boot/html/tsboot.html

这是我在幼稚的引导程序上的工作方式:

my.boot.fnx<-function(var, ind,alpha=0.95){
  newdf      <-na.omit(var[ind])
  mymean     <-mean(newdf,na.rm=TRUE)
  mytrimmean <-mean(newdf, trim=0.1, na.rm=TRUE)
  mymedian   <-median(newdf,na.rm=TRUE)
  mysd       <-sd(newdf,na.rm=TRUE)
  nonmiss    <- length(na.omit(newdf))
  semean     <- mysd/sqrt(nonmiss)
  lcl        <- mymean - qt(1-alpha/2,nonmiss-1)*semean
  ucl        <- mymean + qt(1-alpha/2,nonmiss-1)*semean
  mygini     <-
   sum(abs(outer(newdf,newdf,FUN="")))/
   length(newdf)/(length(newdf)-1)*sqrt(pi)/2
   c(mean=mymean,median=mymedian,se=semean, 
   lcl=lcl,ucl=ucl,sd=mysd,gsd=mygini)
#gini coef is a robust measure of SE
}
strap.df <- boot(df$var,statistic=my.boot.fnx, R=1000)

我还找到了这个时间数据来源 http://eranraviv.com/bootstrapping-time-series-r-code/

对于不同设计的适当分析方法也是很好的资源:

http://people.stat.sfu.ca/~cschwarz/CourseNotes/

无论如何,很抱歉我没有提供太多帮助,但想分享一些想法。

【讨论】:

  • 谢谢@zazu。我做“自制”引导程序的原因仅仅是因为从这些输出中,我将获得我试图评估其不确定性的两个模型的输入。然后,我将对每个模型运行 1000 次。一旦我得到输出,我就会对自举输出的 95% 置信区间的分析感兴趣(从每个模型中获得!)。
  • 不用担心,在您开始自举之前,请注意样本收集的方法。如果它不是 CRD,则不能只运行 Naive 引导程序(我希望我是统计专业来解释原因)。这是一个链接:stats.stackexchange.com/questions/9664/…
  • 你能做到吗?您是否尝试过让时间成为具有独特水平的有序因素?
猜你喜欢
  • 2014-08-19
  • 2020-11-21
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2013-09-18
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
相关资源
最近更新 更多