r - 如何确保 spatstat::owin(poly=<polygon>) 中的多边形没有 "negative area"

我是 R spatstat 包的新用户，在使用 owin() 创建多边形观察窗口时遇到问题。代码如下:

library("maps")
library ("sp")`
library("spatstat")
mass.map <- map("state", "massachusetts:main", fill=T) # This returns a data frame includding x and y components that form a polygon of massachusetts mainland`

mass.win <- owin(poly=data.frame(x=mass.map$x, y=mass.map$y)

我尝试了诸如颠倒多边形顺序之类的方法，但得到了同样的错误。

 mass.win <- owin(poly=data.frame(x=rev(mass.map$x), y=rev(mass.map$y)))

多边形包含重复的顶点



然后我想也许 map() 返回的多边形并不意味着要馈送给 owin()。所以我尝试加载马萨诸塞州的形状文件(此时我完全猜测)。:

x <- readShapePoly("../Geog/OUTLINE25K_POLY") ## The shape file for MASS, loaded from MassGIS website
mass.poly <- x <- readShapePoly("../Geog/OUTLINE25K_POLY", force_ring=T, delete_null_obj=T) ## I got following error whether or not I used force_ring

我收到了提示相交顶点等的消息，但它无法构建多边形。

关于问题的一些背景:我正在尝试使用 spatstat 中的函数进行空间相对风险计算，即案例与控件密度的空间比率。为此，我需要一个观察窗口和该窗口内的点图。我可以作弊，使观察窗口成为围绕马萨诸塞州的矩形，但这可能会扭曲海岸附近的值。无论如何，我想学习如何为我使用此软件包所做的任何 future 工作正确地做到这一点。感谢您的任何帮助，您可以提供。

最佳答案

编辑 2021-02-04:我再次遇到了这个顺时针/逆时针问题，并找到了我自己的答案，这是唯一一个解决 spatstat 问题的答案。答案不是很有帮助。对其进行了改进，使其可用于更广泛的应用。spatstat 包需要逆时针指定的 owin 对象坐标。引用自 1.36-0 版文档:

library("maps")
library ("sp")
library("spatstat")
mass.map <- map("state", "massachusetts:main", fill=T)

首先，您需要确定多边形是顺时针运行还是逆时针运行。方法 here 可用于找出答案。将其公式化为函数:

#' @title Check whether points for an owin are clockwise
#' @param x a dataframe with x coordinates in the first column and y coordinates in the second.
#' @details Similarly to owin, the polygon should not be closed
#' @return A logical telling whether the polygon is arranged clockwise.
#' @author The idea has been scavenged from https://stackoverflow.com/a/1165943/1082004

clockwise <- function(x) {

  x.coords <- c(x[[1]], x[[1]][1])
  y.coords <- c(x[[2]], x[[2]][1])

  double.area <- sum(sapply(2:length(x.coords), function(i) {
    (x.coords[i] - x.coords[i-1])*(y.coords[i] + y.coords[i-1])
  }))

  double.area > 0
}

现在，查看 mass.map 中的坐标是否为顺时针:

clockwise(data.frame(x=mass.map$x, y=mass.map$y))
#> TRUE

他们是。使用 rev 函数逆时针旋转坐标，该函数反转 x 和 y 的向量:

mass.win <- owin(poly=data.frame(x=rev(mass.map$x), y=rev(mass.map$y)))
plot(mass.win)

关于r - 如何确保 spatstat::owin(poly=<polygon>) 中的多边形没有 "negative area"，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/18885777/