R语言dggridR包如何实现六边形网格相邻单元查询功能

内置函数直接调用 #

dggridR自带直接查询网格邻接关系的内置函数dgneighbors()，完全不需要手动编写六边形顶点计算逻辑，该函数原生适配DGGS全球网格的拓扑规则，能自动处理极地单元、五边形变形单元、跨180度日期变更线单元的邻接判断，结果准确且运行效率高。

基础调用方法非常简单，传入你构造的网格对象和目标单元的SEQNUM编号，就能直接返回该单元所有相邻单元的编号：

# 示例：查询编号为10的网格单元的所有相邻单元
target_cell <- 10
neighbor_cells <- dgneighbors(dggs, target_cell)

如果需要获取全量网格的邻接关系，可以批量调用后整理成标准邻接表，方便后续做空间分析：

library(dplyr)
# 获取网格内所有单元的编号
all_cell_ids <- 1:dgmaxcell(dggs)
# 批量查询邻接关系
adj_table <- lapply(all_cell_ids, function(cid){
  nb_ids <- dgneighbors(dggs, cid)
  data.frame(
    cell = cid,
    neighbor = nb_ids
  )
}) %>% bind_rows()
# 如果需要无向不重复的邻接对，执行去重即可
adj_table <- adj_table %>%
  rowwise() %>%
  mutate(pair_id = paste(min(cell, neighbor), max(cell, neighbor), sep = "_")) %>%
  ungroup() %>%
  distinct(pair_id, .keep_all = TRUE) %>%
  select(-pair_id)

自定义实现逻辑（无内置函数时用） #

如果因为版本限制无法使用内置函数，也不需要做复杂的六边形顶点相交计算，用中心距离阈值法就能实现，逻辑更简单还不容易出错：

• 第一步：通过dgSEQNUM_to_GEO()获取所有网格单元的中心点经纬度
• 第二步：对每个目标单元，计算它和其余所有单元中心点的球面距离（必须用球面距离算法，不能直接对经纬度算欧氏距离，避免高纬、跨日期变更线区域出错）
• 第三步：筛选出和目标单元中心距离小于1.5倍网格标称间距的单元，就是该单元的相邻单元。六边形网格中直接相邻单元的中心距等于网格间距，间隔一个单元的中心距为网格间距的√3≈1.732倍，1.5倍阈值可以100%准确筛出直接相邻单元，不会误选间隔单元。

现有代码的小修正 #

你贴的代码里地震统计部分有两处小问题，会导致运行报错：

1. 重复定义了dggs对象，两次构造完全相同的1000英里网格属于冗余代码
2. 第二次统计时summarise(total=sum())没有指定求和列，运行会报错

修正后的统计部分代码：

# 统计每个单元的地震数量、震级总和（以dgquakes自带的mag列为例）
quakecounts <- dgquakes %>% 
  group_by(cell) %>% 
  summarise(
    quake_num = n(),
    total_mag = sum(mag, na.rm = TRUE)
  )

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