R中的光栅：创建特定单元格值的分区计数，无需重新分类

Question

我想知道是否有办法为光栅图层对象创建分区统计数据，特别是R中给定单元值（例如土地使用类别）的计数，而无需重新分类整个光栅。解决方案应具有内存效率，以便处理大型光栅文件，即不需要将值提取到R中的矩阵中。

下面是一个到目前为止我是如何处理它的示例。在这种情况下，我将原始光栅重新分类为仅保留1的利息值和所有其他值的缺失值。

我提出的解决方案创建了冗余数据和额外的处理步骤，以实现我的初始目标。我认为类似于分区（r1[r1==6]，r2，“count”）的方法可以奏效，但显然不行（见下文）。

# generate reproducible Raster 
library("raster")

## RASTER 1 (e.g. land-use classes)
r1 <- raster( crs="+proj=utm +zone=31")
extent(r1) <- extent(0, 100, 0, 100)
res(r1) <- c(5, 5)
values(r1) <- sample(10, ncell(r1), replace=TRUE)
plot(r1)

## RASTER 2 (containing zones of interest)
r2 <- raster( crs="+proj=utm +zone=31")
extent(r2) <- extent(0, 100, 0, 100)
res(r2) <- c(5, 5)
values(r2) <- c(rep(1,100),rep(2,100),rep(3,100),rep(4,100))
plot(r2)

# (1) ZONAL STATISTICS
# a. how many cells per zone (independent of specific cell value)
zonal(r1,r2,"count")

# b. how many cells per zone of specific value 6
zonal(r1[r1==6],r2,"count")
# -> fails

# with reclassification
r1.reclass<-
  reclassify(r1,
             matrix(c(1,5,NA,
                    5.5,6.5,1, #class of interest 
                    6.5,10,NA),
                    ncol=3,
                    byrow = T),
             include.lowest=T # include the lowest value from the table.
           )
zonal(r1.reclass,r2,"count")

Answer 1

您可以使用光栅：：匹配。

zonal(match(r1, 6),r2, "count")

正如您在绘图（匹配（r1，6））中所看到的，它只返回包含所需值的光栅单元格。所有其他单元格均为NA。

很遗憾，try中使用的r1==6返回了一个向量，因此不能再用于焦点。