简介

| 📖 上期回顾:在上一期推文中,我们学会了如何用 raster::extract(), 一键为物种分布点匹配气候数据。今天,我们进入物种分布模型(SDM) 的核心环节——模型调参与选择。

我们以经典数据“拟南芥(Arabidopsis thaliana)”为例,使用 ENMeval 包,自动寻找最优的 特征组合(fc) 和 正则化倍数(rm),并采用空间分块交叉验证(block partition),避免模型过拟合。

🎯 为什么要用 ENMeval?

经典的 MaxEnt 模型如果不调参(使用默认的 fc = LQHrm = 1),容易产生两种问题:

  1. 过拟合:模型过度拟合噪声,预测到未来或新区域时失真;

  2. 参数组合主观:手动尝试不同 fc 和 rm 非常耗时,且难以保证找到最优组合。

ENMeval 包(Muscarella et al. 2014) 专门解决上述问题,它能够:

  • 自动进行 网格搜索(grid search)遍历多种特征组合与正则化倍数;

  • 内置 空间分块交叉验证block partition),降低空间自相关对模型评价的干扰;

  • 基于 AICc(小样本校正后的赤池信息准则)客观选出最优模型,平衡拟合优度与复杂度。

因此,在发表物种分布模型相关研究时,使用 ENMeval 调参已是主流推荐做法(参见 Morales et al. 2017; Cobos et al. 2019)。

✉️ 关注本公众号,后台发送关键词”R语言做物种分布模型三“,

🎯 即可领取本期完整示例源代码和数据!

📦 一、加载所需R包

library(dplyr)      # 数据清洗
library(tictoc)     # 计时
library(sp)         # 空间数据
library(terra)      # 栅格处理(替代raster)
library(ENMeval)    # MaxEnt调参核心包
library(ggplot2)    # 可视化

🗺️ 二、准备物种分布点

本期的测试数据是从GBIF获取的拟南芥分布数据,大家可以参照第一期的内容,R语言做物种分布模型一:从GBIF获取物种分布数据,获取自己的物种分布数据进行分析。

如果只需要本期测试数据,可以后台回复**”R语言做物种分布模型三“获取源代码和数据。**

以下代码,我们首先:

1)剔除经度或纬度超出地球合理范围(经度 ±180°,纬度 ±90°)的记录,并删除含有缺失值的行,确保坐标数据有效。

2)将数据框转为R空间点对象,并声明坐标系为 WGS84(经纬度,与 WorldClim 等环境数据一致),便于后续提取栅格值或生成背景点。

occ <- read.table("Arabidopsis thaliana/Arabidopsis thaliana-site.xls",
                header = T,
                sep = '\t')

#清洗经纬度
occ <- occ %>%
  filter(
    Long >= -180 & Long <= 180,
    Lat  >= -90  & Lat  <= 90
  ) %>%
  na.omit()

coordinates(occ) <- ~ Long + Lat
proj4string(occ) <- CRS("+proj=longlat +datum=WGS84")

🌎 三、加载环境变量(栅格)

假设你已经下载了WorldClim的19个生物气候因子(2.5分辨率),存放在本地文件夹。没有下载气候数据的话,可以参考我们的第二期内容,R语言做物种分布模型二:根据经纬度进行气候数据提取

#读取环境数据
env_files <- list.files(
  "D:/桌面/jianzhi/code/气候因子抓取/worldclim/climate/wc2.1_2.5m",
  pattern = "\\.tif$",
  full.names = TRUE
)

env <- rast(env_files)

🎲 四、生成随机背景点(10,000个)

背景点用于MaxEnt的“伪不存在”样本,必须在环境栅格的有效范围内随机抽取。

#生成随机背景点
set.seed(123)

tic("生成背景点:")
bg <- spatSample(
  env,
  size = 10000,
  method = "random",
  na.rm = TRUE,
  as.df = TRUE,
  xy = TRUE
)
toc()

为什么要生成背景点?

MaxEnt 是一种存在数据(presence-only)模型,它不需要真实的“不存在”点,而是通过比较物种存在位置与整个研究区的随机背景点来估计适生分布。这些背景点代表了该区域可用的环境条件。

参数解释

  • size = 10000:背景点数量。通常建议不少于 10,000 个,样本量过大会增加计算负担,过小则影响模型稳定性。

  • method = "random":在整个栅格范围内完全随机抽样(非分层或系统抽样)。

  • na.rm = TRUE:跳过栅格值为 NA 的像元(如海洋或区域外),保证背景点都有有效环境数据。

  • as.df = TRUE:返回数据框格式,方便后续处理。

  • xy = TRUE:同时返回背景点的经纬度坐标(x, y)。

⚠️ 随机种子 set.seed(123) 保证结果可重复。tictoc 包会输出耗时,方便评估性能。

⚙️ 五、ENMevaluate 核心调参

1. 准备输入数据格式

#准备 ENMevaluate 输入
occs <- as.data.frame(occ)[, c("Long", "Lat")]
colnames(occs) <- c("x", "y")

bg_mat <- as.matrix(bg[, c("x", "y")])
occs_mat <- as.matrix(occs)

# 强制转为double,避免底层C++报错
storage.mode(occs_mat) <- "double"
storage.mode(bg_mat) <- "double"

2. 执行网格搜索调参

我们测试:

  • 特征组合 (fc)L(线性)、LQ(线性+二次)、LQH(线性+二次+铰链)

  • 💡 MaxNet 支持的 fc 还包括 H(铰链)、P(乘积)和 T(阈值),但组合数过多会显著增加计算量。本文仅演示常用子集,你可以按需扩展。

  • 正则化倍数 (rm):从1到4,步长0.5(即 1, 1.5, 2, 2.5, 3, 3.5, 4)

  • 🔍 partitions = "block" 会将研究区划分为4个空间块,轮流用3块训练、1块测试,比随机k折更符合空间数据的特性。

  • 这里我们设置的核心数numCores=4,大家可以根据自己的cpu核心数自行调节。

tic("ENMevaluate模型调参:")
eval <- ENMevaluate(
  occs = occs_mat, 
  envs = env, 
  bg   = bg_mat, 
  algorithm = "maxnet",
  tune.args = list(
    fc = c("L","LQ","LQH"),
    rm = seq(1, 4, 0.5)
  ),
  partitions = "block",  # ⭐ 空间交叉验证
  parallel = T,
  numCores = 4
)
toc()

我这里用的4核,运行内存128G的情况下,大概花了半小时。峰值内存花了60G。建议大家至少在32G情况下运行。

Package ecospat is not installed, so Continuous Boyce Index (CBI) cannot be calculated.
*** Running initial checks... ***
* Removed 96 occurrence localities that shared the same grid cell.
* Removed 1 occurrence points with NA predictor variable values.
* Clamping predictor variable rasters...
* Model evaluations with spatial block (4-fold) cross validation and lat_lon orientation...
*** Running ENMeval v2.0.5 with maxnet from maxnet package v0.1.4 ***
Of 28 total cores using 4...
Running in parallel ...
Making model prediction rasters...
ENMevaluate completed in 27 minutes 33.7 seconds.
警告信息:
1: In for (i in seq_len(n)) { :
  关闭不再使用的链结6(<-local.raysync.cn:11517)
2: In for (i in seq_len(n)) { :
  关闭不再使用的链结5(<-local.raysync.cn:11517)
3: In for (i in seq_len(n)) { :
  关闭不再使用的链结4(<-local.raysync.cn:11517)
4: In for (i in seq_len(n)) { :
  关闭不再使用的链结3(<-local.raysync.cn:11517)
> toc()
ENMevaluate模型调参:: 1654.34 sec elapsed

💾 六、保存结果并筛选最优模型

eval 作为ENMevaluate() 返回的结果对象,它是一个 R 列表(S4 类),包含了模型调参的全部信息:

  • @results:所有参数组合的评价指标表(AUC、AICc、delta.AICc 等);

  • @models:每个参数组合对应的 maxnet 模型对象(可用于后续预测);

  • @predictions:每个分块交叉验证的预测栅格(如果设置 do.predictions = TRUE);

  • @partition:空间分块的具体分配信息。

保存为 .rds 文件后,下次可以直接用 readRDS() 加载,无需重新跑模型。

delta.AICc = 0 表示该模型在候选参数中AICc最低,是复杂度与拟合优度的最佳平衡。

saveRDS(eval, file = "ENMevaluate_full_result.rds")
res <- eval@results
write.table("ENMevaluate_fc_rm_res.xls",row.names = F,sep = '\t')

best <- res[res$delta.AICc == 0, ]
best

📊 七、可视化:ΔAICc 热图

用热图直观展示不同 fc 和 rm 组合下的ΔAICc值,颜色越深(ΔAICc越大)表示模型越差。

delta_heatmap <- ggplot(res, aes(x=factor(rm), y=fc, fill=delta.AICc)) +
  geom_tile() +
  scale_fill_viridis_c() +
  labs(x="Regularization Multiplier", y="Feature Classes", fill="ΔAICc") +
  theme_bw()+
  theme(
    panel.grid = element_blank(),
    panel.border = element_rect(size = 1.2),
    legend.key.spacing.y = unit(0.2,'cm'),
    legend.title = element_text(size = 20,face = 'bold'),
    axis.title = element_text(size = 20),
    legend.text = element_text(size = 15),
    axis.ticks = element_line(size = 1),
    axis.text = element_text(size = 15,color = 'black'))+
  coord_equal()

delta_heatmap
ggsave(plot = delta_heatmap,filename = 'ENMevaluate_fc_rm_delta_heatmap.png',height = 4,width = 6,dpi = 300)

图片

🔜 下期预告

模型调完了,参数选好了,那怎么知道模型好不好?未来气候下物种会怎么变化?

下期我们用最优模型:

✅ 保存最佳 maxnet 模型,可重复调用✅ 预测当前气候下的适生分布图✅ 绘制 Boyce 曲线,定量评估模型预测能力✅ 加载未来气候数据,绘制未来情景的适生分布图

👉 关注更新,不迷路

📚 参考文献

1. Muscarella, R., et al. (2014). ENMeval: An R package for conducting spatially independent evaluations and estimating optimal model complexity for Maxent ecological niche models. Methods in Ecology and Evolution, 5(11), 1198–1205.

2. Morales, N. S., et al. (2017). MaxEnt’s parameter configuration and small samples: Are we paying attention to recommendations? A systematic review. PeerJ, 5, e3093.

3. Cobos, M. E., et al. (2019). kuenm: An R package for detailed development of ecological niche models using Maxent. PeerJ, 7, e6281.

Logo

AtomGit 是由开放原子开源基金会联合 CSDN 等生态伙伴共同推出的新一代开源与人工智能协作平台。平台坚持“开放、中立、公益”的理念,把代码托管、模型共享、数据集托管、智能体开发体验和算力服务整合在一起,为开发者提供从开发、训练到部署的一站式体验。

更多推荐