如何统一栅格堆栈中非 NA 像元的数量？

Question

我的光栅文件具有相同的分辨率和范围，但 NA 的数量不同。我想统一它们之间的 NA 数量。如果在所有栅格文件中都不是 NA，是否可以通过将单元格视为非 NA 来做到这一点？

这里是一个例子：

   library(raster)
library(terra)

f <- system.file("external/test.grd", package="raster")
r1 <- raster(f)
r2 <- calc(r1, fun=function(x){ x[x < 500] <- NA; return(x)} )
r1 <- calc(r1, fun=function(x){ x[x > 1200] <- NA; return(x)} )
raste <-  rast(r1)
rNA <- terra:: global(!(is.na(raste)), sum, na.rm=TRUE)
print(paste0("Non-NA of r1", rNA))
raste <-  rast(r2)
rNA <- terra:: global(!(is.na(raste)), sum, na.rm=TRUE)
print(paste0("Non-NA of r2", rNA))

我希望 r1 和 r2 具有相同数量的非 NA 单元格。我有两个以上的栅格，所以我想知道是否可以对大量文件进行此操作。

Answer 1

同时使用 raster 和 terra 可能会有点混乱，所以我将只使用 terra（但您可以对 raster 执行相同操作，使用 stack 代替 c 和 cellStats而不是global。

您的示例数据

library(terra)
f <- system.file("external/test.grd", package="raster")
r <- rast(f)
r1 <- clamp(r, upper=1200, values=FALSE)
r2 <- clamp(r, lower=500, values=FALSE)

global(!(is.na(r1)), sum)
#       sum
#lyr.1 3145
global(!(is.na(r2)), sum)
#      sum
#lyr.1 802

解决方案：

r <- c(r1, r2)
names(r) <- c("r1", "r2")

m <- any(is.na(r))
x <- mask(r, m, maskvalue=1)

global(!(is.na(x)), sum, na.rm=TRUE)
#   sum
#r1 769
#r2 769

我喜欢使用any(is.na())，因为它清楚地表明了意图是什么。

但是您可以通过许多其他方式之一来组合这些层。只要您不使用na.rm=TRUE，其中一层中带有NA 的单元格将在输出中为NA。例如sum、diff、prod、mean 或 app。

m <- sum(r)
x <- mask(r, m)
global(!(is.na(x)), sum, na.rm=TRUE)