是什么原因导致R中栅格计算中calc和cellStats之间的差异？

Question

我正在处理一个由20个图层组成的数据集，这些图层堆叠在一个RasterBrick中（源于一个数组）。我研究了层的总和，用 "calc "和 "cellStats "计算。我使用calc计算总数值的总和，使用cellStats查看每层数值的平均值（对时间序列有用）。然而，当我把每层的平均值相加时，它是另一个计算值相加的一半。是什么原因造成这种差异？我忽略了什么？

代码看起来像这样。

testarray <- runif(54214776,0,1) 
# Although testarray should contain a raster of 127x147 with 2904 time layers. 
# Not really sure how to create that yet. 

for (i in 1830:1849){
  slice<-array2[,,i]
  r <- raster(nrow=(127*5), ncol=(147*5), resolution =5, ext=ext1, vals=slice)
  x <- stack(x , r)
}

brickhp2 <- brick(x)

r_sumhp2 <- calc(brickhp2, sum, na.rm=TRUE)
r_sumhp2[r_sumhp2<= 0] <- NA


SWEavgpertimestepM <- cellStats(brickhp2, stat='mean', na.rm=TRUE)

目标是比较用 "calc(x, sum) "计算的层数之和 和用 "cellStats(x, mean) "计算的平均值之和。

Rasterbrick是这样的(600kb, GTiff) 。http:/www.filedropper.combrickhp2*如果有更好的分享方式，请告诉我。

Answer 1

当你使用 "calc "和 "cellStats "时，就会产生困惑。calc 在砖块上按像素操作（即在每个像素上执行20个值的计算，并返回一个单一的光栅层）和 cellStats 它对每个光栅层单独执行计算，并返回每个层的单一值。你可以看到，如果你使用这段代码，结果是相当的。

library(raster)
##set seed so you get the same runif vals
set.seed(999)

##create example rasters
ls=list()
for (i in 1:20){
  r <- raster(nrow=(127*5), ncol=(147*5), vals=runif(127*5*147*5))
  ls[[i]] <- r
}
##create raster brick
brickhp2 <- brick(ls)

##calc sum (pixel-wise)
r_sumhp2 <- calc(brickhp2, sum, na.rm=TRUE)
r_sumhp2 ##returns raster layer

##calc mean (layer-wise)
r_meanhp2 <- cellStats(brickhp2, stat='mean', na.rm=TRUE)
r_meanhp2 ##returns vector of length nlayers(brickhp2)

##to get equivalent values you need to divide r_sumhp2 by the number of layers 
##and then calculate the mean
cellStats(r_sumhp2/nlayers(brickhp2),stat="mean")

[1] 0.4999381

##and for r_meanhp2 you need to calculate the mean of the means
mean(r_meanhp2)

[1] 0.4999381

你需要自己决定你的应用是要使用像素还是图层的结果。