如何计算SpatRasterDataset中子数据集的平均栅格数据

问题背景 #

我创建了一个SpatRasterDataset对象，代码如下：

library(terra)

r = rast(system.file("ex/logo.tif", package="terra"))   

s = terra::sds(r, r/2, r/3)

names(s) = letters[1:3]

该对象的信息：

> s
class       : SpatRasterDataset 
subdatasets : 3 
dimensions  : 77, 101 (nrow, ncol)
nlyr        : 3, 3, 3 
resolution  : 1, 1  (x, y)
extent      : 0, 101, 0, 77  (xmin, xmax, ymin, ymax)
coord. ref. : Cartesian (Meter) 
source(s)   : logo.tif, memory 
names       : a, b, c 

我需要计算这个对象中3个SpatRaster的平均值，结果要求是一个维度为77×101的单个SpatRaster（单图层），每个单元格的值是三个子数据集对应位置所有单元格值的平均值。

我尝试了terra::app(s, mean)，但得到的是一个3图层的SpatRaster，不符合预期：

> terra::app(s, mean)
class       : SpatRaster 
dimensions  : 77, 101, 3  (nrow, ncol, nlyr)
resolution  : 1, 1  (x, y)
extent      : 0, 101, 0, 77  (xmin, xmax, ymin, ymax)
coord. ref. : Cartesian (Meter) 
source(s)   : memory
varname     : logo 
names       :      red,    green,     blue 
min values  :   0.0000,   0.0000,   0.0000 
max values  : 155.8333, 155.8333, 155.8333 

解决方案 #

terra::app()对SpatRasterDataset操作时，默认是按对应图层分别计算平均值（比如三个子数据集的红色通道取平均，绿色通道取平均，以此类推），所以得到3图层结果。要实现所有值的全局单元格平均，需要先将SpatRasterDataset转换为一个包含所有图层的SpatRaster，再计算每个单元格的平均值：

# 将SpatRasterDataset转换为包含所有子图层的SpatRaster
combined_raster <- rast(s)

# 计算每个单元格的平均值（所有9个图层的均值）
mean_raster <- app(combined_raster, mean, na.rm = TRUE)

查看结果：

> mean_raster
class       : SpatRaster 
dimensions  : 77, 101, 1  (nrow, ncol, nlyr)
resolution  : 1, 1  (x, y)
extent      : 0, 101, 0, 77  (xmin, xmax, ymin, ymax)
coord. ref. : Cartesian (Meter) 
source(s)   : memory
names       :     mean 
min values  :   0.0000 
max values  : 116.8750 

如果你的需求是保留三个通道（每个通道取三个子数据集对应通道的平均），那terra::app(s, mean)的结果其实符合预期，只需明确需求即可。

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