使用R语言terra包计算点到国家边界距离结果异常的问题

你的核心问题是使用经纬度坐标系（WGS84）直接计算距离，terra::distance在地理坐标系下无法正确将度转换为千米——经纬度的1度对应的实际地面距离随纬度变化（赤道约111km，高纬度会缩小），此时指定unit="km"只是对度数值做了错误的线性转换，导致结果完全不符合实际。

修正步骤 #

• 选择适合加蓬的投影坐标系（比如UTM 33S，EPSG:32733，或等距投影EPSG:4087，确保距离计算准确）
• 将所有空间数据（多边形、栅格、点）转换到该投影坐标系
• 重新执行距离计算

修正后的完整代码 #

library(terra)
library(geodata)
library(ggplot2)
library(RColorBrewer) # 补充加载依赖包

# 获取加蓬边界数据
gabon <- geodata::gadm('GAB', level = 0, path = getwd())    

# 选择适合加蓬的投影坐标系：UTM 33S（EPSG:32733）
target_crs <- "EPSG:32733"

# 转换边界到投影坐标系
gabon_proj <- project(gabon, target_crs)

# 创建投影后的栅格画布（分辨率设为8000米，和原0.08度尺度接近）
canvas_proj <- terra::rast(ext(gabon_proj), 
                           resolution = 8000,
                           crs = target_crs,
                           vals = 0)
pts_proj <- as.points(canvas_proj)    
pts_proj <- terra::crop(pts_proj, gabon_proj) # 裁剪出加蓬境内的点

# 转换多边形为线要素（投影后）
gabon_lines_proj <- terra::as.lines(gabon_proj)

# 计算点到边界的距离，单位设为km
dis_pts <- terra::distance(pts_proj, gabon_lines_proj, pairwise = FALSE, unit="km")
range(dis_pts) # 现在结果会符合加蓬实际国土范围

# 转换回经纬度用于ggplot可视化
pts_wgs84 <- project(pts_proj, "EPSG:4326")
dat <- data.frame(dist = as.vector(dis_pts), 
                  crds(pts_wgs84))

col_dist <- brewer.pal(11, "RdGy")

ggplot(dat, aes(x, y, fill = dist)) + 
  geom_tile() + 
  scale_fill_gradientn(colours = rev(col_dist))+ 
  labs(fill = "到边界的距离 (km)")+ 
  theme_void()+ 
  theme(legend.position = "bottom")

关键改动说明 #

1. 投影转换：所有空间数据都转换到UTM投影，确保距离计算的准确性——UTM投影是等角投影，在小范围区域内距离计算误差极小。
2. 栅格分辨率：原代码用resolution=0.08（度），转换为投影后用8000（米），两者在加蓬纬度区域的实际尺度接近。
3. 可视化适配：把计算后的点再转回经纬度，避免ggplot展示投影坐标系时出现变形问题。

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