вычислить расстояние между выборками в dataframe с применением apply() и пользовательской функцией с несколькими аргументами

Question

У меня очень большой набор данных образцов и их географических местоположений (> 500 выборок). Я хотел бы вычислить расстояние между всеми образцами (т. Е. Сделать матрицу расстояний), используя два разных метода: большое расстояние между изолятами и большое расстояние по кругу, используя путевые точки человеческих миграций. Первое прямое с пакетом геосферы, и я достиг его с кодом, указанным после набора данных образца (n = 10).

test <- structure(list(sample_lon = c(85.1, 101.65, 101.52, 100.5, 77.67, 
78.01, 41.8376999, 136.3070068, 43.0671997, 33.6925011), sample_lat = c(20.95, 
3.11, 3.07, 13.76, 27.49, 27.18, 56.2234001, 49.0937004, 52.3828011, 
-12.3848), Cairo = c(5482676.53004203, 7979572.75185404, 7969372.29645241, 
7300548.40307534, 4535773.29985337, 4577259.89703268, 3035531.48539288, 
8559356.13799981, 2675638.0088102, 4698651.89376045), Istanbul = c(5775286.47479145, 
8414621.15676231, 8406733.79230612, 7555850.43260897, 4744788.34257317, 
4791494.45027787, 1968286.36800993, 7814352.55494176, 1704767.09047958, 
5939243.14729151), PhnomPenh = c(2300013.05256412, 910323.277871997, 
918740.411068992, 487913.521378381, 3302731.94513236, 3256808.88281153, 
7291948.53468698, 5171445.2592504, 7089458.06789339, 8189260.33147658 
), AddisAbada = c(5137690.28464086, 6987142.79928228, 6973508.88041575, 
6744246.81620126, 4567281.84480851, 4588137.16563698, 5241878.016939, 
9802623.85641921, 4823403.30909928, 2433156.16935115), UN = c("Southern-Asia", 
"South-Eastern-Asia", "South-Eastern-Asia", "South-Eastern-Asia", 
"Southern-Asia", "Southern-Asia", "Eastern-Europe", "Eastern-Europe", 
"Eastern-Europe", "Eastern-Africa"), continent = c("Asia", "Asia", 
"Asia", "Asia", "Asia", "Asia", "Europe", "Europe", "Europe", 
"Africa")), .Names = c("sample_lon", "sample_lat", "Cairo", "Istanbul", 
"PhnomPenh", "AddisAbada", "UN", "continent"), row.names = c(NA, 
10L), class = "data.frame") 

require(geosphere) 
test.dvse <- apply(test[c("sample_lon", "sample_lat")], 1, FUN=function(X) distVincentyEllipsoid(X, test[c("sample_lon", "sample_lat")])) 

Это успешно возвращает матрицу расстояний между всеми изолятами.

Для анализа путевых точек я создал пользовательскую функцию, которая смотрит на континенты, из которых берутся образцы, а затем вычисляет расстояние между изолятами с помощью ряда правил.

dCI <- distVincentyEllipsoid(c(31,30), c(28,41)) #Cairo to Istanbul 
dCP <- distVincentyEllipsoid(c(31,30), c(104,11)) #Cairo to Phnom Penh 
dIP <- distVincentyEllipsoid(c(28,41), c(104,11)) #Istanbul Phnom Penh 

calcWayDist <- function(lon1, lat1, con1, lon2, lat2, con2) { 
    if (con1 == con2) { 
    dd = distVincentyEllipsoid(c(lon1,lat1), c(lon2,lat2)) 
    } 
    else{ 
    if (setequal(c("Africa", "Europe"), c(con1, con2)) == TRUE) { 
     if (which("Africa" == c(con1, con2)) == 1) { 
     afr = distVincentyEllipsoid(c(lon1,lat1), c(31,30)) #dist from isolate to Cairo 
     eur = distVincentyEllipsoid(c(lon2,lat2), c(28,41)) #dist from isolate to Istanbul 
     dd = afr + dCI + eur 
     } 
     else { 
     afr = distVincentyEllipsoid(c(lon2,lat2), c(31,30)) #dist from isolate to Cairo 
     eur = distVincentyEllipsoid(c(lon1,lat1), c(28,41)) #dist from isolate to Istanbul 
     dd = afr + dCI + eur 
     } 
    } 
    else if (setequal(c("Africa", "Asia"), c(con1, con2)) == TRUE) { 
     s1 = distVincentyEllipsoid(c(lon1,lat1), c(31,30)) #dist from isolate to Cairo 
     s2 = distVincentyEllipsoid(c(lon2,lat2), c(31,30)) #dist from isolate to Cairo 
     dd = s1 + s2 
    } 
    else if (setequal(c("Africa", "Melanesia"), c(con1, con2)) == TRUE) { 
     if (which("Africa" == c(con1, con2)) == 1) { 
     afr = distVincentyEllipsoid(c(lon1,lat1), c(31,30)) #dist from isolate to Cairo 
     mel = distVincentyEllipsoid(c(lon2,lat2), c(104,11)) #dist from isolate to Phnom Penh 
     dd = afr + dCP + mel 
     } 
     else { 
     afr = distVincentyEllipsoid(c(lon2,lat2), c(31,30)) #dist from isolate to Cairo 
     mel = distVincentyEllipsoid(c(lon1,lat1), c(104,11)) #dist from isolate to Phnom Penh 
     dd = afr + dCP + mel 
     } 
    } 
    else if (setequal(c("Europe", "Asia"), c(con1, con2)) == TRUE) { 
     s1 = distVincentyEllipsoid(c(lon1,lat1), c(28,41)) #dist from isolate to Cairo 
     s2 = distVincentyEllipsoid(c(lon2,lat2), c(28,41)) #dist from isolate to Cairo 
     dd = s1 + s2 
    } 
    else if (setequal(c("Europe", "Melanesia"), c(con1, con2)) == TRUE) { 
     if (which("Europe" == c(con1, con2)) == 1) { 
      eur = distVincentyEllipsoid(c(lon1,lat1), c(28,41)) #dist from isolate to Istanbul 
      mel = distVincentyEllipsoid(c(lon2,lat2), c(104,11)) #dist from isolate to Phnom Penh 
      dd = eur + dIP + mel 
     } 
     else { 
      eur = distVincentyEllipsoid(c(lon2,lat2), c(28,41)) #dist from isolate to Istanbul 
      mel = distVincentyEllipsoid(c(lon1,lat1), c(104,11)) #dist from isolate to Phnom Penh 
      dd = eur + dIP + mel 
     } 
    } 
    else if (setequal(c("Asia", "Melanesia"), c(con1, con2)) == TRUE) { 
     s1 = distVincentyEllipsoid(c(lon1,lat1), c(104,11)) #dist from isolate to Phnom Penh 
     s2 = distVincentyEllipsoid(c(lon2,lat2), c(104,11)) #dist from isolate to Phnom Penh 
     dd = s1 + s2 
    } 
    } 
return(dd) 
    } 

Это работает при желании при вводе вручную значений. Например,

calcWayDist(3.6,15.6,"Africa",51.9,31.9,"Asia") 

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

test.waypoint <- apply(test, 1, FUN=function(X) calcWayDist(X["sample_lon"], X["sample_lat"], X["continent"], test["sample_lon"], test["sample_lat"], test["continent"])) 

Ошибка в .pointsToMatrix (р1): долгота> 360 Кроме того: Предупреждение Сообщение: В случае (con1 == CON2) {: условие имеет длину> 1 и только первый элемент будет использоваться Вызывается из: .pointsToMatrix (p1)

следующий раз я пытался «векторизации» входные аргументы для моей определенной функции пользователя:

calcWayDistbyVec <- function(s1, s2) { 
    if (s1[3] == s2[3]) { 
    dd = distVincentyEllipsoid(c(as.numeric(s1[1]),as.numeric(s1[2])), c(as.numeric(s2[1]),as.numeric(s2[2]))) 
    } 
    else{ 
    if (setequal(c("Africa", "Europe"), c(s1[3], s2[3])) == TRUE) { 
     if (which("Africa" == c(s1[3], s2[3])) == 1) { 
     afr = distVincentyEllipsoid(c(as.numeric(s1[1]),as.numeric(s1[2])), c(31,30)) #dist from isolate to Cairo 
     eur = distVincentyEllipsoid(c(as.numeric(s2[1]),as.numeric(s2[2])), c(28,41)) #dist from isolate to Istanbul 
     dd = afr + dCI + eur 
     } 
     else { 
     afr = distVincentyEllipsoid(c(as.numeric(s2[1]),as.numeric(s2[2])), c(31,30)) #dist from isolate to Cairo 
     eur = distVincentyEllipsoid(c(as.numeric(s1[1]),as.numeric(s1[2])), c(28,41)) #dist from isolate to Istanbul 
     dd = afr + dCI + eur 
     } 
    } 
    else if (setequal(c("Africa", "Asia"), c(s1[3], s2[3])) == TRUE) { 
     d1 = distVincentyEllipsoid(c(as.numeric(s1[1]),as.numeric(s1[2])), c(31,30)) #dist from isolate to Cairo 
     d2 = distVincentyEllipsoid(c(as.numeric(s2[1]),as.numeric(s2[2])), c(31,30)) #dist from isolate to Cairo 
     dd = d1 + d2 
    } 
    else if (setequal(c("Africa", "Melanesia"), c(s1[3], s2[3])) == TRUE) { 
     if (which("Africa" == c(s1[3], s2[3])) == 1) { 
     afr = distVincentyEllipsoid(c(as.numeric(s1[1]),as.numeric(s1[2])), c(31,30)) #dist from isolate to Cairo 
     mel = distVincentyEllipsoid(c(as.numeric(s2[1]),as.numeric(s2[2])), c(104,11)) #dist from isolate to Phnom Penh 
     dd = afr + dCP + mel 
     } 
     else { 
     afr = distVincentyEllipsoid(c(as.numeric(s2[1]),as.numeric(s2[2])), c(31,30)) #dist from isolate to Cairo 
     mel = distVincentyEllipsoid(c(as.numeric(s1[1]),as.numeric(s1[2])), c(104,11)) #dist from isolate to Phnom Penh 
     dd = afr + dCP + mel 
     } 
    } 
    else if (setequal(c("Europe", "Asia"), c(s1[3], s2[3])) == TRUE) { 
     d1 = distVincentyEllipsoid(c(as.numeric(s1[1]),as.numeric(s1[2])), c(28,41)) #dist from isolate to Cairo 
     d2 = distVincentyEllipsoid(c(as.numeric(s2[1]),as.numeric(s2[2])), c(28,41)) #dist from isolate to Cairo 
     dd = d1 + d2 
    } 
    else if (setequal(c("Europe", "Melanesia"), c(s1[3], s2[3])) == TRUE) { 
     if (which("Europe" == c(s1[3], s2[3])) == 1) { 
      eur = distVincentyEllipsoid(c(as.numeric(s1[1]),as.numeric(s1[2])), c(28,41)) #dist from isolate to Istanbul 
      mel = distVincentyEllipsoid(c(as.numeric(s2[1]),as.numeric(s2[2])), c(104,11)) #dist from isolate to Phnom Penh 
      dd = eur + dIP + mel 
     } 
     else { 
      eur = distVincentyEllipsoid(c(as.numeric(s2[1]),as.numeric(s2[2])), c(28,41)) #dist from isolate to Istanbul 
      mel = distVincentyEllipsoid(c(as.numeric(s1[1]),as.numeric(s1[2])), c(104,11)) #dist from isolate to Phnom Penh 
      dd = eur + dIP + mel 
     } 
    } 
    else if (setequal(c("Asia", "Melanesia"), c(s1[3], s2[3])) == TRUE) { 
     d1 = distVincentyEllipsoid(c(as.numeric(s1[1]),as.numeric(s1[2])), c(104,11)) #dist from isolate to Phnom Penh 
     d2 = distVincentyEllipsoid(c(as.numeric(s2[1]),as.numeric(s2[2])), c(104,11)) #dist from isolate to Phnom Penh 
     dd = d1 + d2 
    } 
    } 
return(dd) 
    } 

Это работает для индивидуального расчета:

calcWayDistbyVec(c(3.6,15.6,"Africa"), c(41.8,56.2,"Europe")) 

[1] 6435307

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

test.waypoint <- apply(test[c("sample_lon", "sample_lat", "continent")], 1, FUN=function(X) calcWayDistbyVec(X, test[c("sample_lon", "sample_lat", "continent")])) 

Ошибка наследует (р, "SpatialPoints"): (список) объект не может быть принужден к типу 'двойной' В дополнение: предупреждающее сообщение: В случае (s1 [3] == с2 [3]) {: условие имеет длину> 1 и только первый элемент будет использоваться Вызывается из: наследует (р, «SpatialPoints»)

Если я это сделать для одной точки, что делает работа хотя.

test.waypoint.ind <- apply(test[c("sample_lon", "sample_lat", "continent")], 1, FUN=function(X) calcWayDistbyVec(X, c(3.6,15.6,"Africa"))) 

Если кто-то может определить, что не работает, я бы очень признателен! Я посмотрел на сообщения Apply User Defined Functions, Apply User-Defined Function to Specific Dataframe Columns, которые помогли мне прийти к моему расчету первой расчётной матрицы, но, похоже, здесь не совсем применимы.

Answer 1

Я передал этот вопрос в исследовательскую комиссию R в своем университете, а кто-то, у кого нет учетной записи переполнения стека, предложил решение. Он не использует apply().

calcWayDistMODIFIED <- function(lon1, lat1, con1, lon2, lat2, con2) { 
    if (con1 == con2) { 
    dd = distVincentyEllipsoid(c(lon1,lat1), c(lon2,lat2)) 
    } else if (con1 == "Africa" & con2 == "Europe") { 
    afr = distVincentyEllipsoid(c(lon1,lat1), c(31,30)) #dist from isolate to Cairo 
    eur = distVincentyEllipsoid(c(lon2,lat2), c(28,41)) #dist from isolate to Istanbul 
    dd = afr + dCI + eur 
    } else if (con1 == "Europe" & con2 == "Africa") { 
    afr = distVincentyEllipsoid(c(lon2,lat2), c(31,30)) #dist from isolate to Cairo 
    eur = distVincentyEllipsoid(c(lon1,lat1), c(28,41)) #dist from isolate to Istanbul 
    dd = afr + dCI + eur 
    } else if (setequal(c("Africa", "Asia"), c(con1, con2))) { 
    s1 = distVincentyEllipsoid(c(lon1,lat1), c(31,30)) #dist from isolate to Cairo 
    s2 = distVincentyEllipsoid(c(lon2,lat2), c(31,30)) #dist from isolate to Cairo 
    dd = s1 + s2 
    } else if (con1 == "Africa" & con2 == "Melanesia") { 
    afr = distVincentyEllipsoid(c(lon1,lat1), c(31,30)) #dist from isolate to Cairo 
    mel = distVincentyEllipsoid(c(lon2,lat2), c(104,11)) #dist from isolate to Phnom Penh 
    dd = afr + dCP + mel 
    } else if (con1 == "Melanesia" & con2 == "Africa") { 
    afr = distVincentyEllipsoid(c(lon2,lat2), c(31,30)) #dist from isolate to Cairo 
    mel = distVincentyEllipsoid(c(lon1,lat1), c(104,11)) #dist from isolate to Phnom Penh 
    dd = afr + dCP + mel 
    } else if (setequal(c("Europe", "Asia"), c(con1, con2))) { 
    s1 = distVincentyEllipsoid(c(lon1,lat1), c(28,41)) #dist from isolate to Cairo 
    s2 = distVincentyEllipsoid(c(lon2,lat2), c(28,41)) #dist from isolate to Cairo 
    dd = s1 + s2 
    } else if (con1 == "Europe" & con2 == "Melanesia") { 
    eur = distVincentyEllipsoid(c(lon1,lat1), c(28,41)) #dist from isolate to Istanbul 
    mel = distVincentyEllipsoid(c(lon2,lat2), c(104,11)) #dist from isolate to Phnom Penh 
    dd = eur + dIP + mel 
    } else if (con1 == "Melanesia" & con2 == "Europe") { 
    eur = distVincentyEllipsoid(c(lon2,lat2), c(28,41)) #dist from isolate to Istanbul 
    mel = distVincentyEllipsoid(c(lon1,lat1), c(104,11)) #dist from isolate to Phnom Penh 
    dd = eur + dIP + mel 
    } else if (setequal(c("Asia", "Melanesia"), c(con1, con2))) { 
    s1 = distVincentyEllipsoid(c(lon1,lat1), c(104,11)) #dist from isolate to Phnom Penh 
    s2 = distVincentyEllipsoid(c(lon2,lat2), c(104,11)) #dist from isolate to Phnom Penh 
    dd = s1 + s2 
    } 
    return(dd) 
} 

distance_function = function(d) { 
    mat = matrix(0, ncol = nrow(d), nrow = nrow(d)) 
    for(i in 1:nrow(mat)) { 
    for(j in 1:nrow(mat)) { 
     mat[i,j] = calcWayDistMODIFIED(lon1 = d[i,'sample_lon'], lat1 = d[i,'sample_lat'], con1 = d[i,'continent'], lon2 = d[j,'sample_lon'], lat2 = d[j,'sample_lat'], con2 = d[j,'continent']) 
    } 
    } 
    return(mat) 
}