Работа с координатами и огромными наборами данных в R

Я немного борюсь с двумя наборами данных, содержащими координаты людей и вышек сотовой связи:

  • Первый набор данных о 9459 человек с 1214 переменными, включая их широту и долготу в градусах.
  • второй набор данных о 31 176 вышках сотовой связи с 4 переменными, включая их широту и долготу в градусах и диапазон в метрах.

Я хотел бы определить, находится ли человек в зоне действия хотя бы одной из вышек сотовой связи, и создать манекен, равный 1, если это так.

Однако из-за размера наборов данных я не могу объединить их с помощью команды cross-join. Я попытался использовать пакет geosphere со следующей командой:

distm(c(df1$longitude, df2$latitude), c(df2$longitude, df2$latitude), fun= distHaversine)

К сожалению, это не работает, так как два набора данных имеют разный размер. Любая идея о том, как решить эту проблему?


person William L.    schedule 02.03.2019    source источник
comment
Возможно, вы можете сделать оба набора данных равными, добавляя значения 0 к меньшему, пока оба не станут равными.   -  person LocoGris    schedule 02.03.2019
comment
Я попробовал то, что вы предложили, но это не работает даже с функцией gc() и циклом. Я всегда получаю одно и то же сообщение об ошибке о размере памяти.   -  person William L.    schedule 04.03.2019


Ответы (2)


Как правило, это можно сделать гораздо эффективнее, чтобы максимизировать использование ОЗУ и процессора и сократить накладные расходы. Однако, если то, что вы пытаетесь сделать, является одноразовой операцией, приведенного ниже подхода должно быть достаточно (занимает около 5 минут на текущем ноутбуке).

Вспомогательная функция

# More info: https://github.com/RomanAbashin/distGeo_v
distGeo_v <- function(x, y, xx, yy) { 
    if(!"geosphere" %in% installed.packages())  {
        stop("The 'geosphere' package needs to be installed for this function to work.")
    }
    matrix(.Call("_inversegeodesic", 
                 as.double(x), as.double(y), as.double(xx), as.double(yy), 
                 as.double(6378137), 1/298.257223563, PACKAGE='geosphere'), 
           ncol = 3, byrow = TRUE)[,1]
}

Данные

library(geosphere)
library(tidyverse)
set.seed(1702)

users <- tibble(userid = 1:10000,
                x = rnorm(10000, 16.3738, 5),
                y = rnorm(10000, 48.2082, 5))

towers <- tibble(lon = rnorm(35000, 16.3738, 10),
                 lat = rnorm(35000, 48.2082, 10),
                 range = runif(35000, 50, 10000))

Код

result <- NULL
for(i in 1:nrow(users)) {

    is_match <- users[i, 1:3] %>%
        tidyr::crossing(towers[, 1:3]) %>%
        filter(distGeo_v(x, y, lon, lat) <= range) %>%
        nrow() > 0

    result <- bind_rows(result, tibble(userid = users$userid[i],
                                       match = is_match))

}

Результат

> head(result)
# A tibble: 6 x 2
  userid match
   <int> <lgl>
1      1 TRUE 
2      2 FALSE
3      3 FALSE
4      4 TRUE 
5      5 FALSE
6      6 FALSE

Теперь вы можете left_join преобразовать результат в исходные данные.

person Roman    schedule 20.08.2019

Ниже я добавляю решение с использованием пакета пространственного риска. Основные функции этого пакета написаны на C++ (Rcpp) и поэтому работают очень быстро.

Функция пространственного риска::points_in_circle() вычисляет наблюдения в пределах радиуса от центральной точки. Обратите внимание, что расстояния рассчитываются по формуле Хаверсина. Поскольку каждый элемент вывода представляет собой фрейм данных, purrr::map_dfr используется для связывания их вместе:

library(tibble)
library(spatialrisk)
library(dplyr)

set.seed(1702)
users <- tibble(userid = as.character(1:10000),
                lon = rnorm(10000, 16.3738, 1),
                lat = rnorm(10000, 48.2082, 1))

towers <- tibble(lon = rnorm(35000, 16.3738, 1),
                 lat = rnorm(35000, 48.2082, 1))

# Users with tower within 200 meters
purrr::map2_dfr(users$lon, users$lat, 
                   ~points_in_circle(towers, .x, .y, radius = 200)[1,], 
                   .id = "userid") %>%
     mutate(inrange = ifelse(is.na(distance_m), FALSE, TRUE))
person mharinga    schedule 01.11.2019