如何将SpatialGridDataFrame对象导入ggplot2进行地图绘制？

我用ggplot2绘制地图时，用于填充的预测梯度值存储在SpatialGridDataFrame对象中，运行绘图代码时出现以下错误：

Error in fortify():! data must be a <data.frame>, or an object coercible by fortify(), not an S4 object with class

我的R脚本如下：

ama <- read.csv('Pr.csv')

#rainfall
ama %>% as.data.frame %>% ggplot(aes(LON, LAT)) +
geom_point(aes(size= rainfall), color="blue", alpha=3/4) + 
  ggtitle("RF") + coord_equal() + theme_bw()
coordinates(ama) <- ~ LON + LAT
proj4string(ama) <- ('+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs')

#fitting variogram
Annual.vgm <- variogram(log(rainfall) ~ LON + LAT, ama, width=0.1)
Annual.fit = fit.variogram(Annual.vgm, vgm(c("Gau", "Sph", "Mat", "Exp")), fit.kappa = FALSE)
plot(Annual.vgm, Annual.fit)

#newGride creation
names(grd)       <- c("LON", "LAT")
coordinates(grd) <- c("LON", "LAT")
gridded(grd)     <- TRUE  # Create SpatialPixel object
fullgrid(grd)    <- TRUE
crs(grd)= crs('+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs')

#kriging interpolation
pred <- krige(rainfall~1, ama, grd, Annual.fit)

#plot the prediction
plot(pred["var1.pred"], crs=TRUE)

#study area map
mystudyarea <- read_sf('mystudyarea.shp')

# ploting using ggplot
ggplot() +
 geom_sf(data=mystudyarea) + geom_sf(pred, 
mapping=aes(fill=pred$var1.pred))+
coord_sf(datum = st_crs(mystudyarea))

错误原因是geom_sf()仅支持sf包的空间对象，而pred是sp包的SpatialGridDataFrame（S4类对象），无法被ggplot直接解析。提供两种可行的解决方法：

方法1：转换为sf对象（推荐） #

将SpatialGridDataFrame转换为sf格式后，即可直接用geom_sf()绘制：

1. 确保已安装并加载sf包（未安装先执行install.packages("sf")）
2. 使用st_as_sf()完成格式转换
3. 修改ggplot代码，直接引用转换后的对象及列名

修改后的绘图代码片段：

# 转换SpatialGridDataFrame为sf对象
pred_sf <- st_as_sf(pred)

# ggplot绘图
ggplot() +
  geom_sf(data = mystudyarea, fill = NA, color = "black") + # 绘制研究区边界
  geom_sf(data = pred_sf, aes(fill = var1.pred)) + # 绘制插值结果
  coord_sf(datum = st_crs(mystudyarea)) +
  scale_fill_viridis_c(name = "Rainfall") + # 可选：添加渐变配色
  theme_bw()

方法2：转换为数据框用栅格绘制 #

如果你的插值网格是规则的，也可以将对象转换为普通数据框，用geom_raster()绘制：

# 转换为带坐标的数据框
pred_df <- as.data.frame(pred, xy = TRUE)

# ggplot绘图
ggplot() +
  geom_sf(data = mystudyarea, fill = NA, color = "black") +
  geom_raster(data = pred_df, aes(x = x, y = y, fill = var1.pred), alpha = 0.8) +
  coord_sf(datum = st_crs(mystudyarea)) +
  scale_fill_viridis_c(name = "Rainfall") +
  theme_bw()

library(sf)
library(ggplot2)
library(gstat)
library(dplyr)

ama <- read.csv('Pr.csv')

#rainfall
ama %>% as.data.frame %>% ggplot(aes(LON, LAT)) +
  geom_point(aes(size= rainfall), color="blue", alpha=3/4) + 
  ggtitle("RF") + coord_equal() + theme_bw()
coordinates(ama) <- ~ LON + LAT
proj4string(ama) <- ('+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs')

#fitting variogram
Annual.vgm <- variogram(log(rainfall) ~ LON + LAT, ama, width=0.1)
Annual.fit = fit.variogram(Annual.vgm, vgm(c("Gau", "Sph", "Mat", "Exp")), fit.kappa = FALSE)
plot(Annual.vgm, Annual.fit)

#newGride creation
names(grd)       <- c("LON", "LAT")
coordinates(grd) <- c("LON", "LAT")
gridded(grd)     <- TRUE  # Create SpatialPixel object
fullgrid(grd)    <- TRUE
crs(grd)= crs('+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs')

#kriging interpolation
pred <- krige(rainfall~1, ama, grd, Annual.fit)

#plot the prediction
plot(pred["var1.pred"], crs=TRUE)

#study area map
mystudyarea <- read_sf('mystudyarea.shp')

# 转换为sf对象并绘图
pred_sf <- st_as_sf(pred)

ggplot() +
  geom_sf(data = mystudyarea, fill = NA, color = "black") +
  geom_sf(data = pred_sf, aes(fill = var1.pred)) +
  coord_sf(datum = st_crs(mystudyarea)) +
  ggtitle("Kriging Prediction of Rainfall") +
  scale_fill_viridis_c(name = "Rainfall") +
  theme_bw()

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