在R中计算光栅文件中像素数量（count）

Question

我有大量的栅格文件和一个在栅格文件范围内的多边形。我想要获得多边形内每个栅格文件的像素数（计数）。 此外，我想创建一个表格，其中包含栅格文件的名称和每个栅格文件的像素数（列出）。 我已经尝试过堆叠，但是我无法跟踪名称。在R中是否有执行此任务的其他方法？

- whj

---

请添加有关您的数据的更多信息（栅格文件？多边形？），并可能提供一个可重现的示例。 - Val

Answer 1

请务必包含示例数据

library(raster)
fn <-system.file("external/rlogo.grd", package="raster")
s <- stack(fn)
s[[1]][1:5000] <- NA
s[[2]][5001:ncell(s)] <- NA
names(s)
#[1] "red"   "green" "blue" 
p <- rbind(c(5,20), c(25,55), c(50, 20), c(20,6), c(5,20))
pol <- spPolygons(p)

plot(s, addfun=function() lines(pol, lwd=2))

我不太确定您想要什么。如果您可以堆叠栅格，则所有栅格的单元格（像素）数量将相同（您说您可以）。我假设您希望获得不是NA的单元格的总和。如果您实际上具有不同起点/分辨率的栅格，则可以重复这些步骤，但无需将它们堆叠到RasterStack中，而是需要调整方法以计数NA单元格。

较小对象的简单方法：

m <- mask(s, pol)
cellStats(m, function(i, ...) sum(!is.na(i)))
# red green  blue 
# 600   506  1106 

如果内存不足，你可以执行以下操作：

m <- mask(s, pol)
x <- reclassify(m, cbind(-Inf, Inf, 1))
names(x) <- names(m)
cellStats(x, 'sum')
#red green  blue 
#600   506  1106 

您也可以尝试以下方法：

extract(s, pol, fun=function(x,...)length(na.omit(x)))
#     red green blue
#[1,] 600   506 1106

如果您想计算所有单元格（包括NA），可以进行如下操作：

# example RasterLayer
r <- s[[1]]
# this step may help in speed if your polygon is small relative to the raster
r <- crop(r, pol)

x <- rasterize(pol, r, 1)
cellStats(x, 'sum')
#[1] 1106

- Robert Hijmans