在R中为多国家不规则边界的插值数据绘制掩膜地图

针对你需要在指定经纬度范围（x:20.375-31.375，y:52.375-61.375）内为多国家边界添加掩膜的需求，提供两种可行方案：

方法1：基于raster包结合rnaturalearth数据实现掩膜 #

该方法延续你原有的raster工作流，适配多国家边界场景：

# 加载所需包
library(ggplot2)
library(gstat)
library(sp)
library(raster)
library(rnaturalearth)
library(rnaturalearthdata)

# 1. 获取目标区域的国家边界数据
world_sf <- ne_countries(scale = "medium", returnclass = "sf")
# 裁剪到指定范围
target_extent <- st_bbox(c(xmin=20.375, xmax=31.375, ymin=52.375, ymax=61.375), crs=st_crs(world_sf))
target_countries <- st_crop(world_sf, target_extent)
# 转成sp格式适配raster包（新版raster可直接支持sf，此步可省略）
target_countries_sp <- as(target_countries, "Spatial")

# 2. 对插值结果做裁剪和掩膜
idw.r <- rasterFromXYZ(idw.output2[, c("lon", "lat", "var1.pred")])
idw.crp <- crop(idw.r, target_countries_sp)
idw.msk <- mask(idw.crp, target_countries_sp)
idw.msk.dfr <- as.data.frame(rasterToPoints(idw.msk))
names(idw.msk.dfr)[1:2] <- c("lon", "lat")

# 3. 绘制结果
ggplot() +
  geom_tile(data = idw.msk.dfr, alpha = 0.8, aes(x = lon, y = lat, fill = round(var1.pred, 0))) +
  scale_fill_gradient(low = "cyan", high = "orange") +
  geom_sf(data = target_countries, color = "grey", fill = NA) +
  geom_point(data = df2, aes(x = lon, y = lat), shape = 18, colour = "red") +
  coord_sf(xlim = c(20.375, 31.375), ylim = c(52.375, 61.375), expand = FALSE)

方法2：基于sf+ggplot2直接绘制（无需转raster） #

该方法利用sf空间交集功能直接实现掩膜，避免栅格转换，代码更简洁：

# 加载所需包
library(ggplot2)
library(sf)
library(rnaturalearth)
library(rnaturalearthdata)

# 1. 将插值结果转成栅格化sf面要素
# 创建匹配插值分辨率的栅格
grd_sf <- st_make_grid(
  st_bbox(c(xmin=20.375, xmax=31.375, ymin=52.375, ymax=61.375), crs=4326),
  cellsize = 0.01, # 和你的插值分辨率一致
  what = "polygons"
)
# 将插值结果匹配到栅格面
idw_grid_sf <- st_sf(geometry = grd_sf)
idw_sf <- st_as_sf(idw.output2, coords = c("lon", "lat"), crs = 4326)
idw_grid_sf <- st_join(idw_grid_sf, idw_sf, join = st_contains)

# 2. 获取目标区域国家边界并做空间交集
world_sf <- ne_countries(scale = "medium", returnclass = "sf")
target_countries <- st_crop(world_sf, st_bbox(c(xmin=20.375, xmax=31.375, ymin=52.375, ymax=61.375), crs=4326))
idw_masked_sf <- st_intersection(idw_grid_sf, target_countries)

# 3. 绘制结果
ggplot() +
  geom_sf(data = idw_masked_sf, aes(fill = round(var1.pred, 0)), alpha = 0.8) +
  scale_fill_gradient(low = "cyan", high = "orange") +
  geom_sf(data = target_countries, color = "grey", fill = NA) +
  geom_point(data = df2, aes(x = lon, y = lat), shape = 18, colour = "red") +
  coord_sf(xlim = c(20.375, 31.375), ylim = c(52.375, 61.375), expand = FALSE)

关键说明 #

• 两种方案均无需额外下载全球ESRI shp，直接通过rnaturalearth获取标准化国家边界数据
• 方法1适配你原有工作流，学习成本低；方法2基于现代sf框架，避免栅格转换的性能损耗

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