R语言网格重计算：中心格点与8邻域格点值比较技术求助

我有一个包含lat（纬度）、lon（经度）和X（格点数值）三列的空间数据框，X的取值范围是1-25。需要对每个格点执行以下操作：

• 将当前格点的X值与其8个相邻格点的X值逐一比较
• 统计中心格点X值大于邻域格点的数量（最多8个，小于则不计数）
• 该逻辑需要应用到约30万个格点的数据集

坐标格式说明：lat和lon为十进制坐标，示例如下：

24 24 ,
24 24.1 ,
24 24.2 ,

我尝试用嵌套for循环和条件判断写了R代码，但目前代码只运行第一行，后续无输出，代码如下：

lat <- seq(38, 39.5, by = 0.1)

lon <- seq(30, 40, by = 0.1)

#create empty df
Per_yr_beta <- data.frame(matrix(ncol = 4, nrow = 22779))

#provide column names
colnames(Per_yr_beta) <- c('Lon', 'Lat', 'Total', 'Beta')


for (i in lat) {

if(is.na(i+0.1)==F) {ymax1 = (i+0.1)} else {next}

if(is.na(i-0.1)==F) {ymin1 = (i-0.1)} else {next}

sq_grid <- Per_yr_dummy %>% 
  filter (Lat <= ymax1 & Lat >= ymin1)

for (j in lon) {

  if(is.na(j+0.1)==F) {xmax1 = (i+0.1)} else {next}
  
  if(is.na(j-0.1)==F) {xmin1 = (i-0.1)} else {next}  
  

sq_grid <- sq_grid %>% filter (Lon <= xmax1 & Lon >= xmin1);

  x=0;

  if(sq_grid$Total[5] > sq_grid$Total[1]){x=x+1}
  if(sq_grid$Total[5] > sq_grid$Total[2]){x=x+1} 
  if(sq_grid$Total[5] > sq_grid$Total[3]){x=x+1} 
  if(sq_grid$Total[5] > sq_grid$Total[4]){x=x+1}
  if(sq_grid$Total[5] > sq_grid$Total[6]){x=x+1}
  if(sq_grid$Total[5] > sq_grid$Total[7]){x=x+1}
  if(sq_grid$Total[5] > sq_grid$Total[8]){x=x+1}
  if(sq_grid$Total[5] > sq_grid$Total[9]){x=x+1}

  sq_grid$Beta[5] <- x           
  sq_grid_v <- as.data.frame(sq_grid[5, 1:4])

  z=1;
  Per_yr_beta[z,] <- sq_grid_v
  z=z+1
  
}}

原代码存在多个关键错误，导致无法正常运行：

• 坐标过滤错误：经度过滤时错误使用纬度值i计算xmax1和xmin1，应使用经度值j
• 循环赋值覆盖：变量z每次内层循环都重置为1，所有结果都覆盖到结果框第一行，无后续输出
• 邻域索引依赖问题：硬编码假设过滤后的sq_grid一定是9行且中心在第5位，边界格点（无完整8邻域）会报错或结果错误
• 效率极低：嵌套for循环处理30万格点速度极慢，且每次循环重复过滤数据框，性能冗余

针对空间格点邻域计算，推荐用专业空间包提高效率，以下提供两种可行方案：

方案1：使用terra包（推荐，适合规则格点） #

如果数据是间隔固定为0.1的规则格点，转成栅格处理效率远超循环：

library(terra)

# 将数据框转为栅格对象（假设原数据框名为df，列名lat, lon, X）
r <- rast(df, type="xyz")

# 定义3x3邻域窗口，排除中心格点
w <- matrix(1, nrow=3, ncol=3)
w[2,2] <- 0

# 计算每个格点X大于邻域的数量
count_greater <- focal(r, w, fun=function(x) sum(x[!is.na(x)] < x[5], na.rm=TRUE))

# 转回数据框并合并原数据
result <- as.data.frame(count_greater, xy=TRUE)
colnames(result) <- c("lon", "lat", "Beta")
final_result <- merge(df, result, by=c("lat", "lon"))

方案2：使用dplyr（适合非规则格点） #

通过经纬度匹配邻域，避免循环：

library(dplyr)
library(tidyr)

# 假设原数据框名为df，列名lat, lon, X
df <- df %>%
  mutate(
    # 生成8邻域的经纬度组合
    neighbor_lat = list(c(lat-0.1, lat-0.1, lat-0.1, lat, lat, lat+0.1, lat+0.1, lat+0.1)),
    neighbor_lon = list(c(lon-0.1, lon, lon+0.1, lon-0.1, lon+0.1, lon-0.1, lon, lon+0.1))
  ) %>%
  unnest(c(neighbor_lat, neighbor_lon)) %>%
  # 匹配邻域的X值
  left_join(df, by=c("neighbor_lat"="lat", "neighbor_lon"="lon"), suffix=c("", "_neighbor")) %>%
  # 统计每个中心格点的符合条件数量
  group_by(lat, lon, X) %>%
  summarize(Beta = sum(X > X_neighbor, na.rm=TRUE), .groups="drop")

如果一定要修复原代码，需修正核心错误：

lat <- seq(38, 39.5, by = 0.1)
lon <- seq(30, 40, by = 0.1)

# 初始化空结果框，避免固定行数限制
Per_yr_beta <- data.frame(Lon = numeric(), Lat = numeric(), Total = numeric(), Beta = numeric())

z <- 1  # 将z放在循环外，避免每次重置

for (i in lat) {
  ymax1 <- i + 0.1
  ymin1 <- i - 0.1
  
  # 先过滤纬度范围的格点
  lat_filtered <- Per_yr_dummy %>% filter(Lat >= ymin1 & Lat <= ymax1)
  
  for (j in lon) {
    xmax1 <- j + 0.1  # 修正为使用经度j计算
    xmin1 <- j - 0.1
    
    # 过滤3x3范围的格点
    sq_grid <- lat_filtered %>% filter(Lon >= xmin1 & Lon <= xmax1)
    
    # 确认中心格点存在
    center_idx <- which(sq_grid$Lat == i & sq_grid$Lon == j)
    if (length(center_idx) != 1) next
    
    center_val <- sq_grid$Total[center_idx]
    # 统计中心值大于邻域的数量
    x <- sum(sq_grid$Total[-center_idx] < center_val, na.rm = TRUE)
    
    # 追加结果到数据框
    Per_yr_beta[z, ] <- data.frame(Lon = j, Lat = i, Total = center_val, Beta = x)
    z <- z + 1
  }
}

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