在R中查找数据点的某个半径内的点数答案

【问题标题】：Find the Number of Points Within a Certain Radius of a Data Point in R在R中查找数据点的某个半径内的点数
【发布时间】：2020-09-23 21:28:05
【问题描述】：

我有 2 个数据集，一个用于医院，另一个用于程序。每个数据集都有纬度和经度坐标。程序要么在医院内进行，要么在医院外进行，尽管如果在医院提供坐标，则不一定精确。我试图在每家医院周围形成一定大小的半径，并确定平均有多少手术点落在该半径内。因此，例如，如果我有 100 家医院和 3000 个程序，我想在所有医院周围形成一个半径，然后查看平均有多少家医院落入该指定半径内。我的初始代码如下，但我知道这可以更快地完成。用 R 编码。谢谢！

for(i in 1:NROW(hospitals)){
  hospital <- hospitals[i,]
  radius <- .016

  # find all the procedures that lie in the .016 sized radius from this hospital

  hospital$latitude_low <- hospital$lat - radius
  hospital$longitude_low <- hospital$long - radius
  hospital$latitude_high <- hospital$lat + radius
  hospital$longitude_high <- hospital$long + radius

  in_rad <- procedures[(procedures$long >= hospital$longitude_low & procedures$long <= 
  hospital$longitude_high & procedures$lat <= hospital$latitude_high & procedures$lat >= 
  hospital$latitude_low),]

  num <- NROW(in_rad)
  hospitals[i,]$number_of_procedures <- num
}

【问题讨论】：

这里的答案可能会有所帮助：stackoverflow.com/q/21977720/12265198。我建议使用fields 包函数rdist.earth。您可以得到两个 lon/lat 坐标矩阵之间的公里或英里距离。

标签： r geospatial raster geosphere

【解决方案1】：

这里有几件事可以改进。首先，您实际上不是在计算距离医院 0.16 个单位半径内完成的程序，而是在以医院为中心的 0.32 * 0.32 个单位正方形内完成的程序。对于特定问题可能不是什么大问题，但实际上它可以更快地计算出特定距离内的点，正如您实际想要的那样。

其次，即使您只打算使用一次，您也倾向于存储您计算过的任何变量。这有助于理解代码，但有时效率较低，而且肯定会使您的代码更长，特别是如果您喜欢使用long_descriptive_variable_names。

最后，您对procedures 进行子集化，然后测量行数，而不是仅使用子集本身的长度。

最后（但不太重要），您将结果一次写入一个值到新列中。您可以使用sapply 一口气完成所有这些操作。

所以你的代码可以用更简单的东西代替，比如：

hospitals$number_of_procedures <- sapply(1:NROW(hospitals), function(i)
  {
    d <- (procedures$long - hospitals[i,]$long)^2 + (procedures$lat - hospitals[i,]$lat)^2
    length(which(d < 0.16^2))
  })

【讨论】：

谢谢，非常感谢！有没有办法在运行时改善这一点？

【解决方案2】：

当您提出问题时，您应该始终包含一些示例数据。像这样

lat <- c(-23.8, -25.8)
lon <- c(-49.6, -44.6)
hosp <- cbind(lon, lat)


lat <- c(-22.8, -24.8, -29.1, -28, -20)
lon <- c(-46.4, -46.3, -45.3, -40, -30)
procedures <- cbind(lon, lat)

您的数据是经度/纬度吗？如果是这样，您需要使用适当的方法来计算距离。例如

 library(geosphere)
 dm <- distm(procedures, hosp)

或者

 library(raster)
 d <- pointDistance(procedures, hosp, lonlat=TRUE)

两者都计算从所有程序到所有医院的距离。对于非常大的数据集，这将失败，但根据您的描述，它应该可以正常工作。现在您可以使用阈值（此处为 400,000 m）来找出每个医院的该距离内有哪些程序

apply(d < 400000, 2, which)
#[[1]]
#[1] 1 2

#[[2]]
#[1] 1 2 3

所以程序 1、2 和 3 都在到医院 2 的距离内

如果你的数据不是经度/纬度，可以使用

 d <- pointDistance(procedures, hosp, lonlat=FALSE)

【讨论】：