如何在R中用公式计算数据集间的坐标差值与欧氏距离

现有两个R数据集i和k，结构如下：

i=structure(list(State = c("x1", "x2", "x3", "x4", "x5", "x6", 
"x7", "x8", "x9", "x10"), centroid_lat = c(25.0567526488506, 
24.5506057345197, 25.108046481826, 25.0238523142178, 24.994141005348, 
25.0811976922289, 24.9630616206443, 24.9924361621629, 25.2103128046815, 
25.2173817704152), centroid_lon = c(76.5385899060845, 76.7828927388028, 
76.4028924220128, 76.4542417907297, 76.4114079037918, 76.4363750139148, 
76.4794068845473, 76.4328442796457, 76.5591513892293, 76.600503986722
), fin = c(3219L, 1519L, 2184L, 1919L, 3003L, 2040L, 1843L, 3502L, 
1085L, 2215L), fin2 = c(3789L, 1452L, 2420L, 1971L, 1719L, 2006L, 
1463L, 1278L, 1351L, 3163L), fin3 = c(1513L, 1740L, 2115L, 2431L, 
1148L, 2251L, 1370L, 2122L, 2679L, 1999L), fin4 = c(1903L, 1891L, 
1687L, 1301L, 3080L, 2713L, 2006L, 3891L, 2537L, 1145L)), class = "data.frame", row.names = c(NA, 
-10L))

k=structure(list(State = c("x11", "x12", "x13", "x14", "x15", "x16", 
"x17", "x18", "x19", "x20"), centroid_lat = c(24.05675265, 23.55060573, 
24.10804648, 24.02385231, 23.99414101, 35.08119769, 34.96306162, 
34.99243616, 35.2103128, 35.21738177), centroid_lon = c(75.53858991, 
75.78289274, 75.40289242, 75.45424179, 75.4114079, 86.43637501, 
86.47940688, 86.43284428, 86.55915139, 86.60050399), fin1 = c(NA, 
NA, NA, NA, NA, NA, NA, NA, NA, NA), fin2 = c(NA, NA, NA, NA, 
NA, NA, NA, NA, NA, NA), fin3 = c(NA, NA, NA, NA, NA, NA, NA, 
NA, NA, NA), fin4 = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA)), class = "data.frame", row.names = c(NA, 
-10L))

需求：将i中每组坐标（centroid_lat、centroid_lon）与k中所有坐标进行差值计算（diff.lat = i$centroid_lat - k$centroid_lat，diff.lon = i$centroid_lon - k$centroid_lon），并按公式L = SQRT((X1-X2)^2 + (Y1-Y2)^2)计算欧氏距离（X1/Y1为i的坐标，X2/Y2为k的坐标），最终输出包含两个数据集的State、原始坐标、差值及欧氏距离的结果。

解决方案 #

方法1：基础R实现 #

无需额外包，通过生成所有坐标组合、合并数据、计算差值和距离得到结果：

# 生成i和k的所有State组合
all_pairs <- expand.grid(i_State = i$State, k_State = k$State, stringsAsFactors = FALSE)

# 合并i和k的坐标信息
all_pairs <- merge(all_pairs, i[, c("State", "centroid_lat", "centroid_lon")], by.x = "i_State", by.y = "State")
all_pairs <- merge(all_pairs, k[, c("State", "centroid_lat", "centroid_lon")], by.x = "k_State", by.y = "State", suffixes = c("_i", "_k"))

# 计算纬度和经度差值
all_pairs$diff.lat <- all_pairs$centroid_lat_i - all_pairs$centroid_lat_k
all_pairs$diff.lon <- all_pairs$centroid_lon_i - all_pairs$centroid_lon_k

# 计算欧氏距离
all_pairs$euclidean_distance <- sqrt((all_pairs$centroid_lat_i - all_pairs$centroid_lat_k)^2 + (all_pairs$centroid_lon_i - all_pairs$centroid_lon_k)^2)

# 调整列顺序，让结果更清晰
result <- all_pairs[, c("i_State", "k_State", "centroid_lat_i", "centroid_lon_i", "centroid_lat_k", "centroid_lon_k", "diff.lat", "diff.lon", "euclidean_distance")]

# 查看前6行结果
head(result)

方法2：tidyverse框架实现 #

使用dplyr和tidyr简化代码，逻辑更直观：

library(dplyr)
library(tidyr)

result_tidy <- i %>%
  # 重命名i的列，避免合并冲突
  select(State_i = State, centroid_lat_i = centroid_lat, centroid_lon_i = centroid_lon) %>%
  # 生成i和k的所有坐标组合
  cross_join(k %>% select(State_k = State, centroid_lat_k = centroid_lat, centroid_lon_k = centroid_lon)) %>%
  # 计算差值和欧氏距离
  mutate(
    diff.lat = centroid_lat_i - centroid_lat_k,
    diff.lon = centroid_lon_i - centroid_lon_k,
    euclidean_distance = sqrt((centroid_lat_i - centroid_lat_k)^2 + (centroid_lon_i - centroid_lon_k)^2)
  ) %>%
  # 调整列顺序
  select(State_i, State_k, centroid_lat_i, centroid_lon_i, centroid_lat_k, centroid_lon_k, diff.lat, diff.lon, euclidean_distance)

# 查看前6行结果
head(result_tidy)

输出结果示例 #

执行上述代码后，结果的前6行如下：

> head(result)
  i_State k_State centroid_lat_i centroid_lon_i centroid_lat_k centroid_lon_k  diff.lat  diff.lon euclidean_distance
1      x1     x11       25.05675       76.53859        24.05675        75.53859 1.0000000 1.0000000          1.4142136
2      x2     x11       24.55061       76.78289        24.05675        75.53859 0.4938500 1.0000000          1.1152792
3      x3     x11       25.10805       76.40289        24.05675        75.53859 1.0513000 1.0000000          1.4530499
4      x4     x11       25.02385       76.45424        24.05675        75.53859 0.9671000 1.0000000          1.3903237
5      x5     x11       24.99414       76.41141        24.05675        75.53859 0.9374000 1.0000000          1.3700365
6      x6     x11       25.08120       76.43638        24.05675        75.53859 1.0244500 1.0000000          1.4357750

内容的提问来源于stack exchange，提问作者psysky