在秒 layer/attribute 上提取具有缓冲条件的栅格的平均值

Question

我有这个星星对象（也可以格式化为光栅）：

stars object with 2 dimensions and 2 attributes
attribute(s):
   LST_mean        elevation    
 Min.   :14.98   Min.   :296.0  
 1st Qu.:16.89   1st Qu.:346.9  
 Median :17.64   Median :389.3  
 Mean   :17.52   Mean   :389.2  
 3rd Qu.:18.18   3rd Qu.:428.3  
 Max.   :20.11   Max.   :521.6  
dimension(s):
  from to  offset    delta                       refsys point values    
x   71 83 4387654  860.241 DHDN / 3-degree Gauss-Kru... FALSE   NULL [x]
y   33 41 5598885 -860.241 DHDN / 3-degree Gauss-Kru... FALSE   NULL [y]

它有 2 个属性（在栅格情况下为图层）：温度和海拔。使用温度，我想 select 落在缓冲区内的像素和 return 平均值，仅针对每次与所考虑的像素的海拔差异小于 90 米的像素。

有什么办法吗？计算缓冲区内像素的平均值非常容易，但我找不到对它们设置任何条件的方法。

非常感谢您的帮助和建议。也非常欢迎使用 satrs 以外的其他包的方法:)

Answer 1

请参阅下面使用 terra 的解决方案。代码使用terra::extract创建两个对应的list：

像素值
周围缓冲区值

随后使用 mapply 对值进行成对处理，其功能与您建议的功能相似。

这是我第一次使用 terra，但似乎 terra::extract 比 raster::extract 快得多，因此即使对于大型栅格，此解决方案也是可行的。

正在创建示例数据：

library(sf)
library(terra)

r = rast(ncol = ncol(volcano), nrow = nrow(volcano), xmin = 0, xmax = ncol(volcano), ymin = 0, ymax = nrow(volcano))
values(r) = volcano
s = r
s[] = rnorm(ncell(s))
r = c(r, s)
crs(r) = ""
plot(r)

正在计算缓冲区：

pnt = as.points(r, values = FALSE)
pol = buffer(pnt, 10)

从点中提取栅格值：

x = extract(r, pnt)
head(x)

##      ID lyr.1       lyr.1
## [1,]  1   100 -0.03525223
## [2,]  2   100  0.31525467
## [3,]  3   101  0.94054608
## [4,]  4   101  0.37209238
## [5,]  5   101 -0.38388234
## [6,]  6   101 -0.03120593

正在从缓冲区中提取栅格值：

y = extract(r, pol)
head(y)

##      ID lyr.1       lyr.1
## [1,]  1   100 -0.03525223
## [2,]  1   100  0.31525467
## [3,]  1   101  0.94054608
## [4,]  1   101  0.37209238
## [5,]  1   101 -0.38388234
## [6,]  1   101 -0.03120593

现在，可以使用 mapply 顺序处理提取的值。首先，我们将对象转换为 lists:

x = as.data.frame(x)
x = split(x, x$ID)
y = as.data.frame(y)
y = split(y, y$ID)

接下来，我们使用mapply进行必要的计算，每次考虑当前焦点值 x 和周围缓冲区值 y:

f = function(x, y) {
    d = abs(x[, 2] - y[, 2])  ## differences
    values = y[, 3]  ## values
    mean(values[d < 5], na.rm = TRUE)  ## Mean of subset
}
result = mapply(f, x, y)

最后，将结果放回光栅模板：

u = r[[1]]
values(u) = result
plot(u)

在秒 layer/attribute 上提取具有缓冲条件的栅格的平均值

extract mean value of raster with buffer condition on second layer/attribute

buffer

r

raster

spatial

r-stars