如何将rasterBrick拆分为多个rasterStack并计算分组均值汇总

解决方案 #

方法1：使用raster包自带的stackApply函数（推荐，代码简洁不易出错） #

该函数专门用于栅格堆栈/砖的图层分组统计，完全匹配需求，代码如下：

library(raster)

# 定义分组规则
block <- 7
n_layers <- nlayers(sstd)
# 给每个图层分配分组编号：1-7为第1组，8-14为第2组，最后5个为第26组
group_index <- rep(1:ceiling(n_layers/block), each = block)[1:n_layers]

# 按分组直接计算均值，一步得到目标MASTER堆栈
MASTER <- stackApply(sstd, indices = group_index, fun = mean, na.rm = TRUE)

运行后直接得到所需结果，无需手动拆分循环。

方法2：修正原有循环写法 #

你原来的循环存在两个问题：

1. 索引i:i+7受运算优先级影响，实际解析为(i:i)+7，只会返回单个值
2. 没有处理最后一组不足7个图层的边界情况

修正后的代码如下：

library(raster)
outstk <- stack()
block <- 7
n_layers <- nlayers(sstd)

for(i in seq(1, n_layers, by = block)){
  # 计算当前组的结束索引，最大不超过总图层数
  end_idx <- min(i + block -1, n_layers)
  # 提取当前组图层
  stk <- sstd[[i:end_idx]]
  # 计算组内均值
  mean_layer <- calc(stk, mean, na.rm = TRUE)
  # 加入结果堆栈
  outstk <- stack(outstk, mean_layer)
}
MASTER <- outstk

结果验证 #

执行nlayers(MASTER)可看到图层数为26，和ceiling(181/7)=26的分组数完全匹配。

内容的提问来源于stack exchange，提问作者jjulip