2014-03-04 3 views
4

Я хочу рассчитать расстояние между двумя связанными наборами пространственных координат (program и admin в моем поддельном наборе данных). Данные находятся в широком формате, поэтому обе пары координат находятся в одной строке.рассчитать расстояние между каждой парой координат в широком информационном кадре

library(sp) 
set.seed(1) 
n <- 100 
program.id <- seq(1, n) 
c1 <- cbind(runif(n, -90, 90), runif(n, -180, 180)) 
c2 <- cbind(runif(n, -90, 90), runif(n, -180, 180)) 
dat <- data.frame(cbind(program.id, c1, c2)) 
names(dat) <- c("program.id", "program.lat", "program.long", "admin.lat", "admin.long") 
head(dat) 
#  program.id program.lat program.long admin.lat admin.long 
# 1    1 -42.20844  55.70061 -41.848523 62.536404 
# 2    2 -23.01770 -52.84898 -50.643849 -145.851172 
# 3    3 13.11361 -82.70635 3.023431 -2.665397 
# 4    4 73.47740 177.36626 -41.588893 -13.841337 
# 5    5 -53.69725  48.05758 -57.389701 -44.922049 
# 6    6 71.71014 -103.24507 3.343705 176.795719 

Я знаю, как создать матрицу расстояний между program или admin используя sp пакет:

ll <- c("program.lat", "program.long") 
coords <- dat[ll] 
dist <- apply(coords, 1, 
       function(eachPoint) spDistsN1(as.matrix(coords), 
              eachPoint, longlat=TRUE)) 

Но то, что я хочу сделать, это создать NX1 вектор расстояний (dist.km) между каждым пару координат и добавьте его в dat.

#  program.id program.lat program.long admin.lat admin.long dist.km 
# 1    1 -42.20844  55.70061 -41.848523 62.536404 567.35 
# 2    2 -23.01770 -52.84898 -50.643849 -145.851172 8267.86 
# ... 

Любые предложения? Я потратил некоторое время на рассмотрение старых вопросов, но ничего не кажется совершенно правильным. Счастлив быть доказанным.

Update

@ решение Amit работает для моих игрушечного набора данных:

apply(dat,1,function(x) spDistsN1(matrix(x[2:3],nrow=1),x[3:4],longlat=TRUE)) 

Но я думаю, что нужно поменять порядок латы, длинный порядок лат длинных столбцов так долго доходит до лат. От ?spDistsN1:

pts: A matrix of 2D points, first column x/longitude, second column y/latitude, or a SpatialPoints or SpatialPointsDataFrame object 

Кроме того, если я не понял логики, я думаю, что решение Amit должно захватить COLS [2: 3] и [4: 5], а не [2: 3] и [3: 4 ].

Моя задача теперь применяется к моим фактическим данным. Я воспроизвел часть ниже.

library(sp) 
dat <- structure(list(ID = 1:4, 
         subcounty = c("a", "b", "c", "d"), 
         pro.long = c(33.47627919, 31.73605491, 31.54073482, 31.51748984), 
         pro.lat = c(2.73996953, 3.26530095, 3.21327597, 3.17784981), 
         sub.long = c(33.47552, 31.78307, 31.53083, 31.53083), 
         sub.lat = c(2.740362, 3.391209, 3.208736, 3.208736)), 
       .Names = c("ID", "subcounty", "pro.long", "pro.lat", "sub.long", "sub.lat"),  
       row.names = c(NA, 4L), class = "data.frame") 
head(dat) 
#  ID subcounty pro.long pro.lat sub.long sub.lat 
# 1 1   a 33.47628 2.739970 33.47552 2.740362 
# 2 2   b 31.73605 3.265301 31.78307 3.391209 
# 3 3   c 31.54073 3.213276 31.53083 3.208736 
# 4 4   d 31.51749 3.177850 31.53083 3.208736 
apply(dat, 1, function(x) spDistsN1(matrix(x[3:4], nrow=1), 
            x[5:6], 
            longlat=TRUE)) 

Я получаю ошибку: Error in spDistsN1(matrix(x[3:4], nrow = 1), x[5:6], longlat = TRUE) : pts must be numeric

Я запутался, потому что эти столбцы являются числовыми:

> is.numeric(dat$pro.long) 
[1] TRUE 
> is.numeric(dat$pro.lat) 
[1] TRUE 
> is.numeric(dat$sub.long) 
[1] TRUE 
> is.numeric(dat$sub.lat) 
[1] TRUE 
+1

вы пробовали: применить (Дат, 1, функция (х) spDistsN1 (матрица (х [2: 3], nrow = 1), x [3: 4], longlat = TRUE))? – amit

+0

@amit, я этого не делал. Я полагал, что ответ, вероятно, будет включать одну из прикладных функций, но я не знал правильной спецификации матрицы. Это похоже на решение. Я был бы рад принять его, если вы хотите добавить ответ. –

+0

пока это работает и полезно - я счастлив. Меня не волнует репутация, но спасибо за предложение. – amit

ответ

5

Проблема у Вас есть что apply(...) принуждают первый аргумент матрицы , По определению, матрица должна иметь все элементы одного и того же типа данных. Поскольку один из столбцов в dat (dat$subcounty) является символом char, apply(...) забирает все на char. В вашем тестовом наборе данных все было числовым, поэтому у вас не было этой проблемы.

Это должно работать:

dat$dist.km <- sapply(1:nrow(dat),function(i) 
       spDistsN1(as.matrix(dat[i,3:4]),as.matrix(dat[i,5:6]),longlat=T)) 
+0

спасибо за объяснение, @jlhoward. это работает. очень признателен. –

+2

Я столкнулся с этим решением сегодня, так как у меня была аналогичная ситуация. Мне нравится эта идея. Интересно, можем ли мы сделать работу еще лучше. У меня большой набор данных, такой как 2 ГБ, и пробовал этот код с data.table. Обработка действительно на какое-то время. Для каждой строки мы просим R создать две матрицы и обработать вычисление. Я скорее думаю, что создаю SPDF и обрабатываю ту же работу. По крайней мере для каждой строки нам не нужно преобразовывать DF в матрицу. Любая мысль? Я также задаюсь вопросом, есть ли еще одна функция, обрабатывающая ту же работу быстрее. – jazzurro

+0

@jazzurro, я считаю, что существует более быстрое решение, использующее 'data.table' и' geosphere' http://stackoverflow.com/questions/36817423/how-to-efficiently-calculate-distance-between-pair-of-cordinate -using-data-tab –

3

Существует гораздо быстрее, используя решение data.table и geosphere.

library(data.table) 
library(geosphere) 

setDT(dat)[ , dist_km := distGeo(matrix(c(pro.long, pro.lat), ncol = 2), 
            matrix(c(sub.long, sub.lat), ncol = 2))/1000] 

Benchmark:

library(sp) 

jlhoward <- function(dat) { dat$dist.km <- sapply(1:nrow(dat),function(i) 
          spDistsN1(as.matrix(dat[i,3:4]),as.matrix(dat[i,5:6]),longlat=T)) } 

rafa.pereira <- function(dat2) { setDT(dat2)[ , dist_km := distGeo(matrix(c(pro.long, pro.lat), ncol = 2), 
                   matrix(c(sub.long, sub.lat), ncol = 2))/1000] } 


> system.time(jlhoward(dat)) 
    user system elapsed 
    8.94 0.00 8.94 

> system.time(rafa.pereira(dat)) 
    user system elapsed 
    0.07 0.00 0.08 

данных

dat <- structure(list(ID = 1:4, 
         subcounty = c("a", "b", "c", "d"), 
         pro.long = c(33.47627919, 31.73605491, 31.54073482, 31.51748984), 
         pro.lat = c(2.73996953, 3.26530095, 3.21327597, 3.17784981), 
         sub.long = c(33.47552, 31.78307, 31.53083, 31.53083), 
         sub.lat = c(2.740362, 3.391209, 3.208736, 3.208736)), 
       .Names = c("ID", "subcounty", "pro.long", "pro.lat", "sub.long", "sub.lat"),  
       row.names = c(NA, 4L), class = "data.frame") 

# enlarge dataset to 40,000 pairs 
dat <- dat[rep(seq_len(nrow(dat)), 10000), ] 
+1

Рафа, спасибо за ваше сообщение и ваш ответ. Ваше решение, безусловно, быстрее! – jazzurro

Смежные вопросы