Найдите количество точек в пределах определенного радиуса точки данных в 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
}

r geospatial raster geosphere
person aport550    schedule 04.06.2020    source источник
comment
Здесь могут помочь ответы: stackoverflow.com/q/21977720/12265198. Я бы предложил использовать fields пакетную функцию rdist.earth. Вы получаете расстояние в километрах или милях между двумя матрицами координат долгота / широта.   -  person bzki    schedule 04.06.2020

Ответы (2)


arrow_upward
2
arrow_downward

Когда вы задаете вопрос, вы всегда должны включать некоторые примеры данных. Нравится

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)

Or 

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

Оба вычисляют расстояние от всех процедур до всех больниц. Это не удастся с очень большими наборами данных, но, судя по тому, что вы описываете, должно работать нормально. Теперь вы можете использовать порог (здесь 400000 м), чтобы узнать, какие процедуры находятся на таком расстоянии от каждой больницы.

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

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

Таким образом, процедуры 1, 2 и 3 находятся на таком же расстоянии от больницы 2.

Если ваши данные не являются долготой / широтой, вы можете использовать

 d <- pointDistance(procedures, hosp, lonlat=FALSE)
person Robert Hijmans    schedule 05.06.2020

arrow_upward
1
arrow_downward

Здесь есть несколько вещей, которые можно улучшить. Во-первых, вы фактически не рассчитываете процедуры, выполняемые в радиусе 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))
  })
person Allan Cameron    schedule 04.06.2020
comment
Спасибо, очень признателен! есть ли способ улучшить это с точки зрения времени выполнения? - person aport550; 05.06.2020