为 r 中的核密度估计图 (density.ppp) 选择像素大小

Question

我有一个多边形 shapefile 和一个分布在多边形内的点的 shapefile。我已经使用 spatstat 包中的 density.ppp 函数基于每个多边形包含的点为每个多边形创建了核密度估计 (KDE) 地图。我现在希望创建具有不同像素大小的不同 kde，以便选择最适合我的模型的那个。我尝试在 as.mask 函数中使用 eps 参数，但这只改变了 window 的像素大小而不是内核映射本身，所以结果没有改变。在抛出 density.ppp 函数的整个手册之后，我能找到的所有看起来相关的是 spatstat.geom 包中的 pixellate.ppp 函数，但我不确定如何使用它与 density.ppp.

关于如何更改 kde 的像素大小有什么建议吗？

library(sf)
library(spatstat)

buffer <- st_read("gis/layers/buffers.shp")
pbb<- st_read("gis/layers/points_by_buffer.shp")

for (p in 1:10) {
 if(p %in% pbb$field_id) {
   
   poly123 <- pbb[pbb$field_id == p,]
   
   C <- as.owin(buffer$geometry[p])

   W<- as.mask(C, eps = 100)

   point<- ppp(poly123$X,poly123$Y, window = W)

   sigma <- bw.diggle(point)
   
   d <- density.ppp(point, kernel = "gaussian", sigma=sigma, positive = TRUE, at="pixels", )
   plot(d) 
}

Answer 1

更改像素大小只会更改结果的分辨率，但基本数学原理是相同的。这只是关于生成的网格有多精细的问题。您可以通过在对 density.ppp() 的调用中使用参数 eps 或 dimyx 来更改此设置。这当然是相关的，但在很大程度上答案只是：更精细的分辨率更好（或者至少通常不会更差，除非你走到极端并遇到数值不稳定性）。更精细的分辨率使估计强度的图看起来更好，但当你说你想“选择最适合我的模型的那个”时，这可能不是你想要的。查看 sigma 指定的所选平滑带宽的效果通常更有趣。您刚刚使用 Diggle 规则解决了这个问题，但这绝不是既定的“最佳实践”。实际上选择带宽是一个没有简单通用解决方案的问题。您可以尝试一系列不同的带宽，看看哪一个最适合您的目的。像素分辨率应该足够高以使绘图看起来漂亮，并且足够低以使计算不会花费太长时间。

Answer 2

我认为您混淆了“分辨率”和“带宽”。

像素分辨率是指输出中像素的大小 - 实际上是计算强度估计的网格点的间距。对于大多数用途，此间距并不重要。

“带宽”是应用于数据的平滑程度。大带宽将产生相对平坦的表面，而小带宽将产生高尖峰表面。带宽确实有很大的影响，需要根据数据选择合适的带宽。

像素分辨率由参数 eps 和 dimyx 控制，可以传递给 density.ppp。（诚然，很难从帮助文件中读到这一点，但是 help(density.ppp) 提到不匹配的参数 ... 被传递给 pixellate.ppp，而 help(pixellate.ppp) 提到不匹配的参数 ... 被传递给 as.mask，最后 as.mask 的帮助描述了参数 eps 和 dimyx。这就是参数传递在 R 中的工作方式。）

平滑带宽由参数 sigma 控制。有关如何 select 和控制平滑带宽的信息，请参阅 help(density.ppp)。另请阅读“另请参阅”部分。

Answer 3

三个不同的东西：

分辨率：计算密度的网格点间距。
带宽：将应用于输入数据的平滑量。
长度单位：物理长度单位（如米、公里）；给出空间坐标的单位；密度值使用相同的单位表示（密度是每平方单位的点数）。

如果您更改 分辨率（您可以在调用 density.ppp 时使用参数 eps 或 dimyx），密度值不会有太大变化；您只是在更改计算密度的位置的间距。

如果您更改带宽（您可以在调用 density.ppp 时使用参数 sigma 和 adjust） , 密度值将发生显着变化。当您增加带宽时，每个数据点有效地分布在更远的距离上，并且密度估计变得更加平滑和平坦。

如果你改变长度单位，一切都会改变一个比例因子。密度值是“每单位面积的点数”。如果密度值为 3，这意味着我们预计每平方单位平均看到 3 个点。这里的'unit'是表示空间坐标的长度单位。如果原点坐标以米为单位，则密度值为'number of points per square metre'。如果你想改变这个，有两种方法：

重新缩放输入数据（推荐）：例如，如果原始数据坐标以米为单位，则可以通过将坐标除以 1000 将其转换为公里。我们建议使用函数 rescale.ppp 来执行此操作：Xnew <- rescale(Xold, 1000)。您可以使用函数 unitname（如 unitname(Xnew) <- "km"）或使用 rescale.ppp（如 Xnew <- rescale(Xold, 1000, "km")).重新缩放点模式后，您调用 density.ppp（如 Dnew <- density(Xnew, ...) ，带宽以千米为单位 ），结果密度值以每平方公里点数表示.
重新缩放结果（不推荐）：如果您的坐标以米为单位并且您后来决定想知道以每平方公里的点数表示的密度，那么您可以只需将比例因子应用于图像（通过键入 Dnew <- D * 1000^2，其中 D 是原始密度图像）并重新缩放图像的空间坐标（通过键入 Dnew <- rescale(Dnew, 1000, "km")）。我们不推荐这样做，因为它很容易搞砸。

有关详细信息，请参阅 the spatstat book。

为 r 中的核密度估计图 (density.ppp) 选择像素大小

choosing pixel size for kernel density estimate map (density.ppp) in r

r

spatial

kernel-density

spatstat