R语言绘制物种分布栅格时手动设置scale_fill_viridis最大值

问题描述 #

我有若干物种分布范围的栅格数据，需要在R语言中绘制地图。这些栅格存在区域重叠，希望图例能显示重叠物种的数量。由于该图例要用于多张地图，需要适配所有地图中的最高物种数（当前示例中最大物种数为3，但我想手动将图例基准设为更高值，比如4）。

我的处理代码如下：

set.seed(123)
library(rnaturalearth)
library(terra)
library(sf)

world <- ne_countries(returnclass = "sv")
# 转换为sf格式用于绘图
world_st <- ne_countries(returnclass = "sf")

# 读取IUCN物种分布矢量数据
fish <- vect("../../Databases/IUCN/HAGFISH/HAGFISH.shp")
sp_list <- unique(fish$sci_name)

# 按物种拆分矢量数据
mylist <- list()
i = 1

for(sp in sp_list){
  
  mylist[[i]] <- fish[fish$sci_name == sp,]
  i = i + 1
  
}

# 创建空栅格模板
r <- rast(extent = ext(world),
          crs = "EPSG:4326",
          resolution = 0.1)

# 栅格化世界范围作为基准
world_r <- terra::rasterize(world, r)

# 栅格化物种分布（取前40个物种）
mylist_r <- list()

for(i in 1:40){
  
  mylist_r[[i]] <- terra::rasterize(mylist[[i]], r)
  # 重采样匹配基准栅格
  mylist_r[[i]] <- terra::resample(mylist_r[[i]], world_r, method = "near") 
  i = i + 1
  
}

# 创建栅格集合并求和得到物种重叠数
mylistsprc <- sprc(mylist_r)
mylistsum <- mosaic(mylistsprc, fun = "sum")
print(paste0("Max number of species: ", minmax(mylistsum)[2]))
# 输出结果：[1] "Max number of species: 3"

后续绘图时，我尝试用以下代码设置图例，但结果不符合预期，不清楚是操作有误还是需要其他方法手动设置最大值：

scale_fill_viridis(direction = -1, discrete = T, breaks = seq(1, max(maxsp), 1), limits = c(1, max(maxsp))) 

解决方法 #

核心问题分析 #

你当前的代码问题在于：

• 当设置discrete = TRUE时，limits参数需要包含所有要显示的级别（从1到你设定的最大值，比如4），而非仅设置首尾；
• 代码中max(maxsp)未明确定义，需替换为你手动指定的目标最大值（比如4）。

修改后的ggplot绘图代码示例 #

结合ggplot2与terra的兼容格式，调整后的图例设置如下：

library(ggplot2)

# 将求和后的栅格转为ggplot兼容的数据框
sum_df <- as.data.frame(mylistsum, xy = TRUE)
colnames(sum_df)[3] <- "species_count"

# 绘制地图
ggplot() +
  geom_sf(data = world_st, fill = "white", color = "gray50") +
  geom_raster(data = sum_df, aes(x = x, y = y, fill = factor(species_count))) +
  # 统一图例级别为1-4，适配所有地图
  scale_fill_viridis(
    direction = -1,
    discrete = TRUE,
    name = "重叠物种数",
    breaks = seq(1, 4, 1),
    limits = as.character(seq(1, 4, 1)),  # 离散型需指定所有显示级别
    na.value = NA  # 隐藏无物种分布的区域
  ) +
  theme_bw() +
  coord_sf(crs = "EPSG:4326")

关键调整点 #

1. 离散图例的limits设置：
   当discrete = TRUE时，limits必须包含所有想要显示的图例级别（如c("1","2","3","4")），即使当前地图中没有达到4的区域，图例也会保留该级别，实现多张地图的图例统一。
2. 栅格数据类型转换：
   将求和后的栅格数值转为因子类型，确保ggplot将其识别为离散变量，避免图例自动省略不存在的级别。
3. 手动指定最大值：
   直接将seq(1,4,1)中的4替换为全局最高物种数，所有地图复用该设置即可保持图例一致。

替代方案（terra原生绘图） #

若不想使用ggplot，可直接用terra的plot函数手动配置图例：

library(viridis)

# 定义对应1-4级别的颜色映射
cols <- viridis(4, direction = -1)
names(cols) <- seq(1,4)

# 绘图并设置统一图例
plot(mylistsum, 
     col = cols,
     legend = TRUE,
     breaks = seq(0.5, 4.5, 1),  # 连续转离散的断点设置
     lab.breaks = seq(1,4,1),
     main = "物种分布重叠数")
plot(world, add = TRUE, border = "gray50")

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