Использование простого цикла for для пространственных данных

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

locations <-read.csv("distances.csv")

Location возвращает следующую таблицу:

       City Type       long      lat
1 Sheffield  EUR  -1.470085 53.38113
2        HK WRLD 114.109497 22.39643
3    Venice  EUR  12.315515 45.44085
4  New York WRLD -74.005941 40.71278

Моя цель в этой конкретной части задачи - создать таблицу расстояний (в километрах) между каждым из городов в виде корреляционной матрицы с диагональю 0 (т. Е. Все города находятся на нулевом расстоянии от самих себя).

Для этого я использую пакет sp, для которого требуется матрица значений long-lat, поэтому я могу удалить текст следующим образом:

datmax <- data.matrix(locations)
datmax2 <- datmax[,-1:-2]

Инструмент spDistsN1 позволяет мне получить эту информацию, сравнивая расстояние между всеми городами в матрице от одного отдельного города. Ясно, что я могу использовать следующее выражение для получения расстояний до всех городов от Шеффилда (город или строка №1):

km <- spDistsN1(datmax2, datmax2[1,], longlat=TRUE)

Это правильно дает:

[1]    0.000 9591.009 1329.882 5436.133

Однако, чтобы добиться желаемого результата в стиле корреляционной матрицы, я хочу добиться этого для каждого из городов, поэтому я попытался написать цикл for:

for (i in 1:nrow(datmax2)){
  kmnew <- spDistsN1(datmax2, datmax2[i,], longlat=TRUE)
}

Это дает мне правильные значения для NY:

[1]  5436.133 12967.023  6697.541     0.000

Итак, я предполагаю, что в цикле я перезаписал один город другим. Я ценю помощь, показывающую мне, в чем я ошибаюсь. Большое спасибо.


person RichS    schedule 03.04.2015    source источник


Ответы (3)


Сначала объявите матрицу и используйте свой итератор i, чтобы указать строку, которую нужно заполнить:

kmnew <- matrix(NA, nrow=4, ncol=4)
for (i in 1:nrow(datmax2)){
  kmnew[i,] <- spDistsN1(datmax2, datmax2[i,], longlat=TRUE)
}

colnames(kmnew) <- locations$City
rownames(kmnew) <- locations$City

Результаты

> kmnew

          Sheffield        HK   Venice  New York
Sheffield     0.000  9591.009 1329.882  5436.134
HK         9591.009     0.000 9134.698 12967.024
Venice     1329.882  9134.698    0.000  6697.541
New York   5436.134 12967.024 6697.541     0.000
person Dominic Comtois    schedule 03.04.2015

Я не уверен, что это то, что вы ищете

library(sp)

# Provide data for reproducibility
locations <- data.frame(City=c("Sheffield", "HK", "Venice", "New York"),
                    Type=c("EUR", "WRLD", "EUR", "WRLD"),
                    long=c(-1.470085, 114.109497, 12.315515, -74.005941),
                    lat=c(53.38113, 22.39643, 45.44085, 40.71278))

km <- apply(as.matrix(locations[, c(-1, -2)]), 1, function(x){
  spDistsN1(as.matrix(locations[, c(-1, -2)]), x, longlat=TRUE)
})

km <- data.frame(locations[, 1],  km)
names(km) <- c("City", as.character(locations[, 1]))
km

Результаты

       City Sheffield        HK   Venice  New York
1 Sheffield     0.000  9591.009 1329.882  5436.134
2        HK  9591.009     0.000 9134.698 12967.024
3    Venice  1329.882  9134.698    0.000  6697.541
4  New York  5436.134 12967.024 6697.541     0.000
person dimitris_ps    schedule 03.04.2015
comment
Все ответы здесь действительно полезны. С Домиником мне легче всего справиться, основываясь на моем существующем мыслительном процессе. «apply» - элегантное решение, спасибо Димитрису и Николаю. Я проверю этот пакет, хотя, как вы, вероятно, можете сказать по моим попыткам, выход за рамки 2D-мышления может занять у меня некоторое время! - person RichS; 03.04.2015

Вы можете попробовать функцию distm из пакета geosphere:

 distm(datmax2)
 #        [,1]     [,2]    [,3]     [,4]
 #[1,]       0  9586671 1329405  5427956
 #[2,] 9586671        0 9130036 12962132
 #[3,] 1329405  9130036       0  6687416
 #[4,] 5427956 12962132 6687416        0

Он возвращает расстояние в метрах и учитывает геометрию земли.

person nicola    schedule 03.04.2015
comment
да, но предполагается, что Земля - ​​это сфера; sp::spDistsN1 предполагает, что это эллипсоид (WGS84). - person Edzer Pebesma; 03.04.2015