每个网格单元区域的 r 栅格统计

r raster statistics per grid cell zone

我想计算格线内区域中的网格单元统计信息(见图表)。

获取显示的网格单元格的 cellstats 并不困难,但是如何为所有网格单元格完成此操作?

library(raster)

filename <- system.file("external/test.grd", package="raster")
r=raster(filename)

plot(r)
e=extent(180000,181000,330000,331000)
plot(e,add=T)
grid()

x=crop(r,e)
cellStats(x,mean)

你的可重现示例对我不起作用,因为我不确定为什么使用 r[r>500]=1r[r<=500]=NA 而你的原始情节不同?!另外,我不确定 rx 是什么?!一种方法如下:

library(raster)

filename <- system.file("external/test.grd", package="raster")
r=raster(filename)

# r[r>500]=1
# r[r<=500]=NA

plot(r)
xmin <- seq(178000, 181000, 1000)
xmax <- seq(179000, 182000, 1000)
ymin <- seq(329000, 333000, 1000)
ymax <- seq(330000, 334000, 1000)

zonal.mtx <- matrix(nrow=length(xmin), ncol=length(ymin))

library(spatialEco)
for (i in 1:length(xmin)){
  for (j in 1:length(ymin)){
    e=extent(xmin[i],xmax[i],ymin[j],ymax[j])
    plot(e,add=T)
    grid()
    p <- as(e, 'SpatialPolygons')  
    crs(p) <- crs(r)

    p$ID=j
    p <- SpatialPolygonsDataFrame(p,p@data) 

    zonal.mtx[i,j] <- zonal.stats(x=p, y=r, stat=sum, trace=TRUE, plot=TRUE)
    #print(zonal.mtx[i,j])
  }

}

注意: 这会在计算区域统计数据时绘制数据的每个部分。我相信您已经知道您还可以定义自己的函数以在 zonal.stats 而不是 summean

中使用