R语言按HEX_Tag_ID分组匹配坐标自动计算移动距离方法

我借助points_in_circle函数按日期计算不同HEX_Tag_ID对应的移动距离，结果存为distance_m字段。现有代码的中心点坐标硬编码为ID是3D6.15341BBB4F在Datum为"2022-04-21"的经纬度，其他ID（比如3D6.15341BC60F）在2022-04-21的坐标和这个值不一致，导致这类ID的距离计算出错。
需要实现自动化计算逻辑：以每个唯一HEX_Tag_ID在Datum="2022-04-21"对应的经纬度作为自身的中心点，计算各ID的移动距离，后续新增ID时自动适配，不需要手动修改坐标参数。核心规则：中心点坐标必须取每个HEX_Tag_ID在Datum为"2022-04-21"的对应记录值。

现有错误代码：

library(spatialrisk)
Distance_ind <- points_in_circle(Distance_ind, lat_center = 51.93349, lon_center = 4.70912, lon = Longitude, lat = Latitude, radius = 1e6)

示例数据框Distance_ind结构如下：

structure(list(Datum_Tijd = structure(c(1650535073, 1653637570, 
1653285342, 1654242563, 1654578739, 1654837567, 1653899310, 1655108033, 
1653900853, 1653286920, 1653639081, 1654838913, 1653032146, 1655110494, 
1654580434, 1654244601, 1650533634), tzone = "", class = c("POSIXct", 
"POSIXt")), Datum = structure(c(19103, 19139, 19135, 19146, 19150, 
19153, 19142, 19156, 19142, 19135, 19139, 19153, 19132, 19156, 
19150, 19146, 19103), class = "Date"), Tijd = c("11:57:53", "09:46:10", 
"07:55:42", "09:49:23", "07:12:19", "07:06:07", "10:28:30", "10:13:53", 
"10:54:13", "08:22:00", "10:11:21", "07:28:33", "09:35:46", "10:54:54", 
"07:40:34", "10:23:21", "11:33:54"), Reader_ID = c("A0", "A0", 
"A0", "A0", "A0", "A0", "A0", "A0", "A0", "A0", "A0", "A0", "A0", 
"A0", "A0", "A0", "A0"), HEX_Tag_ID = c("3D6.15341BBB4F", "3D6.15341BBB4F", 
"3D6.15341BBB4F", "3D6.15341BBB4F", "3D6.15341BBB4F", "3D6.15341BBB4F", 
"3D6.15341BBB4F", "3D6.15341BBB4F", "3D6.15341BC60F", "3D6.15341BC60F", 
"3D6.15341BC60F", "3D6.15341BC60F", "3D6.15341BC60F", "3D6.15341BC60F", 
"3D6.15341BC60F", "3D6.15341BC60F", "3D6.15341BC60F"), Longitude = c(4.70912, 
4.70917, 4.70918, 4.70918, 4.70914, 4.70914, 4.70927, 4.70921, 
4.70904, 4.70903, 4.70906, 4.709, 4.70901, 4.70903, 4.70902, 
4.70902, 4.70925), Latitude = c(51.933491, 51.934189, 51.9342, 
51.934269, 51.93428, 51.934292, 51.934341, 51.934441, 51.932499, 
51.932491, 51.932468, 51.932369, 51.93235, 51.932339, 51.932331, 
51.932308, 51.934891), x = c(108366.11, 108370.273, 108370.972, 
108371.043, 108368.304, 108368.317, 108377.308, 108373.285, 108359.578, 
108358.882, 108360.921, 108356.692, 108357.36, 108358.724, 108358.028, 
108358.004, 108376.503), y = c(438553.585, 438631.208, 438632.425, 
438640.102, 438641.351, 438642.686, 438648.054, 438659.218, 438443.273, 
438442.39, 438439.812, 438428.836, 438426.716, 438425.479, 438424.595, 
438422.037, 438709.256), `Lengte_(cm)` = c(9.7, 9.7, 9.7, 9.7, 
9.7, 9.7, 9.7, 9.7, 11.2, 11.2, 11.2, 11.2, 11.2, 11.2, 11.2, 
11.2, 11.2), Geslacht = c("vrouw", "vrouw", "vrouw", "vrouw", 
"vrouw", "vrouw", "vrouw", "vrouw", "vrouw", "vrouw", "vrouw", 
"vrouw", "vrouw", "vrouw", "vrouw", "vrouw", "vrouw"), Sloot = c("22", 
"22", "22", "22", "22", "22", "22", "22", "22", "22", "22", "22", 
"22", "22", "22", "22", "22"), Lengte_8e_lichting = c(NA, NA, 
NA, NA, NA, NA, NA, 10.3, NA, NA, NA, NA, NA, 12.5, NA, NA, NA
), Lengteklasse = structure(c(4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L), .Label = c("6", "7", "8", 
"9", "10", "11", "12", "13"), class = "factor"), distance_m = c(0.111319490512219, 
77.8879654012696, 79.1440538153195, 86.8156131287763, 87.953110772786, 
89.2887843834125, 95.2906913528804, 106.044905287657, 110.454187269549, 
111.379609951795, 113.843032832766, 125.060674080157, 127.12861901825, 
128.277561285803, 129.201736097579, 131.75853920649, 156.213638336975
)), row.names = c(NA, -17L), class = c("tbl_df", "tbl", "data.frame"
))

points_in_circle本身针对全局固定中心点设计，要实现按ID匹配对应中心点，两种方案都支持新增ID自动适配，无需手动修改参数：

方法1：分组计算（代码最简洁） #

按HEX_Tag_ID分组，每组内先提取该ID在2022-04-21的经纬度作为中心点，再调用points_in_circle计算组内所有记录到中心点的距离，最后合并结果。

library(spatialrisk)
library(dplyr)

Distance_ind <- Distance_ind %>%
  group_by(HEX_Tag_ID) %>%
  group_modify(~{
    # 提取当前ID在2022-04-21的中心点坐标
    center_coords <- .x %>% filter(Datum == as.Date("2022-04-21")) %>% slice(1)
    # 组内计算距离
    points_in_circle(
      .x,
      lat_center = center_coords$Latitude,
      lon_center = center_coords$Longitude,
      lon = Longitude,
      lat = Latitude,
      radius = 1e6
    )
  }) %>%
  ungroup()

如果某个ID在2022-04-21有多条记录，上述代码默认取第一条作为中心点，需要取均值/其他聚合规则的话，修改center_coords的提取逻辑即可。

方法2：关联中心点逐行计算（大数据量下效率更高） #

如果数据量较大，分组操作效率低，可以先提取所有ID的中心点做成查找表，关联回原表后直接用球面距离公式逐行计算，不需要循环调用points_in_circle。spatialrisk内置的haversine函数计算逻辑和points_in_circle完全一致，输出距离单位为米。

library(spatialrisk)
library(dplyr)

# 生成所有ID的中心点查找表
id_center <- Distance_ind %>%
  filter(Datum == as.Date("2022-04-21")) %>%
  select(HEX_Tag_ID, center_lon = Longitude, center_lat = Latitude) %>%
  distinct(HEX_Tag_ID, .keep_all = T)

# 关联中心点后逐行计算距离
Distance_ind <- Distance_ind %>%
  left_join(id_center, by = "HEX_Tag_ID") %>%
  mutate(
    distance_m = haversine(
      lon1 = Longitude, lat1 = Latitude,
      lon2 = center_lon, lat2 = center_lat
    )
  ) %>%
  select(-center_lon, -center_lat)

两种方法计算结果完全一致，后续新增ID时只要数据中包含该ID在2022-04-21的坐标记录，就会自动匹配对应中心点计算，不需要修改代码参数。

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