基于时间匹配计算叉角羚位置与在建风机的最近距离

首先确保两个sf对象的时间列均为POSIXct格式，坐标系保持一致（若不一致用st_transform转换为同一投影坐标系，比如UTM，确保距离计算单位为米）。以下是两种高效的实现方法，替代仅处理奇数行的错误循环：

方法1：用purrr逐行处理（中小数据集友好） #

通过定义自定义函数，遍历每个叉角羚位置，筛选对应时间的在建风机，计算最近距离：

library(sf)
library(purrr)

# 标准化时间格式（若尚未转换）
pronghorn_sf$DateTime <- as.POSIXct(pronghorn_sf$DateTime)
turbine_sf$construction_start <- as.POSIXct(turbine_sf$construction_start)
turbine_sf$construction_end <- as.POSIXct(turbine_sf$construction_end)

# 转换为投影坐标系（示例用UTM 12N，根据你的研究区域调整EPSG代码）
pronghorn_sf <- st_transform(pronghorn_sf, 32612)
turbine_sf <- st_transform(turbine_sf, 32612)

# 定义计算单条叉角羚记录最近在建风机距离的函数
calc_nearest_turbine <- function(row) {
  # 筛选当前时间处于施工时段内的风机
  active_turbines <- turbine_sf[
    row$DateTime >= turbine_sf$construction_start & 
    row$DateTime <= turbine_sf$construction_end, 
  ]
  
  # 无匹配风机时返回NA，否则计算最小距离
  if (nrow(active_turbines) == 0) {
    return(NA_real_)
  } else {
    st_distance(row, active_turbines) |> as.numeric() |> min()
  }
}

# 应用到所有叉角羚记录
pronghorn_sf$nearest_active_turbine_m <- map_dbl(pronghorn_sf, calc_nearest_turbine)

方法2：用data.table非等连接（大数据集高效） #

利用非等时间匹配+空间距离计算，避免循环，大幅提升处理速度：

library(sf)
library(data.table)

# 标准化时间与坐标系（同方法1）
pronghorn_sf$DateTime <- as.POSIXct(pronghorn_sf$DateTime)
turbine_sf[, `:=`(construction_start = as.POSIXct(construction_start), 
                   construction_end = as.POSIXct(construction_end))]
pronghorn_sf <- st_transform(pronghorn_sf, 32612)
turbine_sf <- st_transform(turbine_sf, 32612)

# 转为data.table并保留几何列
setDT(pronghorn_sf, keep.rownames = "pronghorn_id")
setDT(turbine_sf)

# 非等连接：匹配叉角羚时间落在风机施工时段内的所有组合
matched_pairs <- pronghorn_sf[turbine_sf, 
                              on = .(DateTime >= construction_start, DateTime <= construction_end),
                              allow.cartesian = TRUE,
                              .(pronghorn_id, pronghorn_geom = i.geometry, turbine_geom = geometry)]

# 计算每对的距离（米）
matched_pairs[, distance_m := st_distance(pronghorn_geom, turbine_geom) |> as.numeric()]

# 按叉角羚ID分组，取最小距离
min_distances <- matched_pairs[, .(nearest_active_turbine_m = min(distance_m)), by = pronghorn_id]

# 合并回原数据集，无匹配的填充NA
pronghorn_sf <- pronghorn_sf[min_distances, on = "pronghorn_id"]
pronghorn_sf[is.na(nearest_active_turbine_m), nearest_active_turbine_m := NA_real_]

原代码问题排查 #

如果你的原代码使用了for循环且仅处理奇数行，大概率是循环步长设置错误（比如seq(1, nrow(pronghorn_sf), 2)），将步长改为1即可遍历所有行，但循环效率远低于上述两种方法，不推荐。

关键注意事项 #

1. 坐标系一致性：必须确保两个sf对象使用同一投影坐标系（非地理坐标系，如WGS84），否则距离计算结果为度而非米。
2. 施工时段定义：若风机有多个施工阶段，需明确“在建”对应的时段列，调整筛选条件即可。
3. 时间格式：所有时间列必须转换为POSIXct，否则无法进行时间区间匹配。

内容的提问来源于stack exchange，提问作者c.rob4