查找 R 中数据点一定半径内的点数

Find the Number of Points Within a Certain Radius of a Data Point in R

我有 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
}

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

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

第三,最后,您子集 procedures 然后测量行数,而不是仅仅使用子集本身的长度。

最后(但不太重要),您将结果一次写入一个值到一个新列中。您可以使用 sapply 一次性 gulp 完成所有这些操作。

因此您的代码可以替换为更简单的代码,例如:

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))
  })

当你提出问题时,你应该总是包括一些示例数据。像这样

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)

你的数据在longitude/latitude吗?如果是这样,您需要使用适当的方法来计算距离。例如

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

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

两者都计算从所有程序到所有医院的距离。对于非常大的数据集,这将失败,但根据您的描述,它应该可以正常工作。 现在你可以使用一个阈值(这里是 400,000 米)来找出每个医院的那个距离内有哪些程序

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

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

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

如果你的数据不是longitude/latitude,你可以使用

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