如何使用 R 中的 3x3 窗口对光栅执行焦点操作（均值）？我有纬度/经度值答案

【问题标题】：How to perform focal operation (mean) on raster using 3x3 window in R? I have lat/long values如何使用 R 中的 3x3 窗口对光栅执行焦点操作（均值）？我有纬度/经度值
【发布时间】：2018-05-17 18:42:03
【问题描述】：

我有 raster 和 lat/long 值，我想使用 3x3 窗口/内核对这些点执行焦点操作。我是 R 新手。

【问题讨论】：

你看过栅格的focal函数吗：focal(r, w=matrix(1/9,nrow=3,ncol=3)) ？
其实我不想移动窗口。我必须使用经纬度表中指定的点。该代码应计算每个点上 9 个单元格（3x3 窗口）的平均值。

标签： r r-raster operation

【解决方案1】：

这是计算以定义的纬度/经度坐标为中心的 3x3 区域中栅格平均值的工作流程。

对点进行栅格化并扩大生成的栅格以创建 3x3 区域

library(mmand)
library(raster)

# rasterize points based on lat/lon coordinates
z <- rasterize(pts[,2:3], r, field = pts$id)

# dilate z using a 3x3 box
kern <- shapeKernel(c(3,3), type="box")
z[,] <- dilate(as.matrix(z), kern)

plot(z)

计算每个区域中 r 值的平均值

# raster::zonal function (zonal statistics) is used 
# -Inf correspond to NA values and should not be taken into account
# you can change "mean" by the stats you would like to compute
zonal(r, z, fun = "mean")

#     zone     mean
#[1,] -Inf 5.563607
#[2,]    1 5.000000
#[3,]    2 3.444444
#[4,]    3 5.222222

样本数据

library(raster)
set.seed(1)
# Generate raster of random values
r <- raster(crs = CRS("+proj=robin +datum=WGS84"), resolution = c(10, 10))
r[] <- round(runif(ncell(r), 1, 10))
# Generate data frame with lat/lon coordinates
pts <- data.frame(id = 1:3, lon = c(-100, 40, 120), lat = c(-45, 5, 35))

plot(r)
points(pts$lon, pts$lat, pch = 20, cex = 2)

【讨论】：

这种方法在栅格化点时增加了相当大的开销。为什么不将raster::extract 与缓冲区参数的单元格大小一起使用，并且small = FALSE？这相当于一个 3x3 窗口（即使给定与单元质心的大点偏移）。另一种方法是使用raster::focal 对窗口进行显式控制（例如，圆形、不规则），然后使用raster::extract 检索光栅点值。 OP 显然不理解有关使用焦点功能的评论。