ARTICLE DETAIL

资讯详情

深耕网站建设与运营推广的一线实战洞察。

物种分布模型SDM全流程:R语言从数据清洗到集成建模实操指南

物种分布模型SDM全流程:R语言从数据清洗到集成建模实操指南 第一次接到物种分布模型的需求时我以为只要把一堆经纬度丢进R语言跑一下dismo包出一张预测图任务就算完成了。结果真正跑起来才发现数据清洗、环境变量选择、模型评估每个环节都藏着能让你返工一星期的坑。后来做得多了才慢慢把流程理顺。这篇博文就把我理解的物种分布模型Species Distribution Model简称SDM和一套可以复现的R语言实操流程完整写出来。内容不只讲模型代码更多是讲“为什么这样做”以及各种文档里不会告诉你的实战判断。适合刚接触生态位建模的研究生、做生物多样性评估的从业者以及想用免费数据预测某个物种潜在分布区的爱好者。1. SDM到底做了什么从分布点到一张有意义的概率图1.1 一次建模请求背后的真实需求有人拿一份红树林分布点数据来找我说想预测2050年哪些海岸可能适合红树林扩展。这种需求在生态保护、入侵物种防控、病虫害传播评估里非常常见。本质上他想做的事是已知这个物种在哪里出现想推算它还会在哪里出现。物种分布模型解决的就是这个问题——利用已知的物种出现记录结合环境变量温度、降水、海拔、土壤、土地利用等建立环境条件与物种出现之间的关系然后把这个关系投影到整个研究区域输出一张每个格点“环境适宜程度”的栅格图。刚开始做SDM的人最容易把它当成一次“点击运行”的机械操作。但实际项目里一个模型能否使用取决于你如何定义研究区域、选择环境变量、处理采样偏差以及最终如何解读输出。我曾见过有人把MaxEnt默认输出的相对适生指数直接解释成“出现概率”后续保护决策完全跑偏这是很危险的事。1.2 核心假设物种分布是环境筛过的结果SDM的理论基础来自生态位理论。一个物种能在某个地方生存必须满足它的生理耐受范围比如很多两栖动物需要潮湿环境不耐高温干旱高山植物依赖低温与特定降水模式。环境条件像一层筛子把不适宜的区域过滤掉剩余的地方才是潜在分布区。这个假设听起来简单实际使用时要清楚三个前提物种处于平衡分布即它已经扩散到了所有能到达的适宜区域。入侵物种刚引入时往往未达到平衡模型会严重低估潜在分布。采样记录能真实反映物种的环境偏好而不是偏向于道路沿线、保护区或城市附近。环境变量与物种分布的关系在时空外推时保持稳定。把今天的模型投影到2050年依赖的正是这个假设。R语言在SDM领域成为主流原因很实际生态学数据预处理、栅格计算、模型拟合、结果绘图全都能在一个工作流里完成不需要来回切换软件。dismo、raster/terra、biomod2、sdm、ENMeval等一系列包覆盖了从数据清洗到模型集成的完整链路。下面我就按自己跑项目的顺序来展开。2. 建模前先问三件事数据、坐标系和采样偏差2.1 分布点数据的来源与清洗分布点数据最常见的来源是GBIF、地方生物多样性监测、标本馆数字化记录和自己团队的野外调查。很多人辛辛苦苦整理好几十页Excel导入R后第一件事就是直接建模结果模型的图出来连自己都不信。我现在的习惯是拿到分布点后先执行一套固定清洗流程删除经纬度缺失、经纬度都是0、经纬度明显落在海洋或研究区域之外的记录。删除重复记录包括同一网格内的重复点和完全重复的坐标。检查是否有离谱坐标比如某个昆虫“出现”在南极大陆这种基本是录入错误。做空间稀疏化如果环境变量分辨率是1公里那么同一个1公里格子里保留一条记录就够了。否则模型会被密集采样区带偏相当于给某个地方加了重量。用R做简单稀疏化可以自己写个网格筛选也可以用spThin包。我自己常用方法是把坐标转成栅格单元编号每个栅格单元随机保留一个点library(raster) # 假设occ是数据框包含lon和lat两列res是你的栅格分辨率 occ$cell - cellFromXY(baseline_raster, occ[, c(lon, lat)]) occ_sparse - occ[!duplicated(occ$cell), ]2.2 环境变量怎么选、怎么踩点环境变量是SDM的“营养原料”。常用来源包括WorldClim的19个生物气候变量、CHELSA的高分辨率气候数据、土壤数据库以及一些专门的地形因子。这里最关键的不是变量越多越好而是变量之间不要高度共线。我在初学阶段犯过的错误是把全部19个bio变量都丢进模型结果AUC高得吓人但预测图只在已知点位周围画了几个圆斑。为什么因为bio-variables之间强相关模型把噪声也背下来了。现在我的做法是先根据物种生物学知识初筛变量。例如热带物种不适合把“最冷季度降水”当作核心变量。用Pearson相关系数或方差膨胀因子VIF检查共线性通常保留相关系数绝对值低于0.7的变量组合。优先选择有明确生理意义、测量误差小的变量比如温度季节性和年降水比你更“易解释”。背景点backgroud points的取样同样重要。MaxEnt这类presence-background模型不需要真正的不存在点它会在研究区域内随机抽取指定数量的背景点用来表示环境空间中的“可获得性”。如果研究区域范围太大背景点里面混入了大量明显不适合栖息的环境模型会变得很保守范围太小背景点与分布点环境差异不够模型又会过度自信。常规做法是让背景点数达到分布点数的10倍以上并限制在研究区域的实际可到达范围内。2.3 坐标系混乱是SDM第一杀手这是我踩过最深、最痛、最不想让别人再踩的坑。SDM建模过程涉及分布点坐标和环境栅格坐标两个必须保持一致。最常见的问题是分布点是WGS84经纬度环境栅格却是某个投影坐标系或者反过来。extract时如果坐标系对不上轻则取出的环境值是NA重则每个点都错位到几百公里以外而模型还以为世界本来就该这样。建议每次建模前都执行一个检查crs_distribution - CRS(projlonglat datumWGS84) crs_env - crs(predictors) # 如果两者不一致用sp或sf转成统一坐标系另外需要注意国内部分数据源会使用GCJ-02等加密坐标。如果你手里是加密坐标不转回WGS84就进入SDM位置会偏移数百米到数公里。对于遥感分辨率几公里的模型这可能只是平移一两个像元但也会影响提取变量的一致性。这个问题我建议在拿到数据的第一时间就确认而不是等模型跑完一张分布图挂在展示墙上后再被质疑。3. 用dismo跑通你的第一个MaxEnt模型3.1 准备数据对象dismo包里的maxent()函数是R语言里面调用MaxEnt的经典方式。它依赖MaxEnt的Java版本所以先把Java Runtime装上确保rJava::.jinit()能正常执行。这一步经常有人卡住报错信息通常是“Java not found”或者“rJava package error”。解决方法不是改模型参数而是先解决系统Java环境。我们直接用dismo自带的懒猴Bradypus示例数据跑通全流程。library(dismo) library(raster) # 读取dismo自带的懒猴分布点和环境栅格 occurrence_file - system.file(ex/bradypus.csv, package dismo) occ - read.csv(occurrence_file) # 第一列是species后面才是经度、纬度 occ - occ[, 2:3] names(occ) - c(lon, lat) raster_file - system.file(ex/bradypus.grd, package dismo) predictors - raster(raster_file)这个bradypus.grd文件是一个多波段的RasterBrick大概包含温度、降水等几个环境变量。如果你想用一组自己的环境栅格可以用stack()把所有tif文件叠到一起predictors - stack(bio1.tif, bio4.tif, bio12.tif, bio15.tif)3.2 运行模型并生成预测图接下来就是一句代码的事me - maxent(predictors, occ) p - predict(me, predictors) plot(p)maxent()把分布点坐标和整个背景栅格环境值都收集起来重采样出约一万个背景点然后训练MaxEnt模型。predict()将模型映射回栅格空间得到0到1之间的相对适生指数图。注意这里的值是一个连续性得分不是严格的概率。如果你需要把结果导出成tif可以这样writeRaster(p, suitability.tif, overwrite TRUE)运行之后第一件事不是看预测图而是看me这个对象里保存的两个信息变量贡献率和响应曲线。用plot(me)能输出MaxEnt自带的那几张诊断图。3.3 MaxEnt核心结果怎么看贡献率、响应曲线和预测图MaxEnt输出里有一张“Percent contribution”的表显示每个环境变量对模型拟合的贡献百分比。很多文章只放这张表就算解释变量重要性实际上不够严谨。变量贡献率来自算法寻优过程中的路径累积当两个变量高度相关时贡献率会在二者之间随机分配不能完全代表生态重要性。更可靠的做法是看刀切法Jackknife结果单独用某个变量建模以及排除某个变量建模综合判断变量是否携带了独立信息。响应曲线是另一项关键输出。它告诉你物种适生概率如何随环境变量变化是单峰曲线还是单调上升、单调下降比如某物种在年均温20°C附近出现最大适生值那么这条曲线就能用于讨论气候变暖对物种迁移方向的影响。如果响应曲线出现特别剧烈的锯齿或边界处断崖通常说明变量范围没有覆盖完整生态位或模型过拟合。预测图需要结合生态常识检查。例如一个陆地哺乳动物的预测分布图大面积落在水体里肯定是背景区域或环境数据出了问题一个热带物种的预测分布图跑到寒带大概率是环境变量选择或外推范围出了问题。模型输出只是数学结果最终是否可信需要人来判断。4. 别满足于单模型biomod2集成建模让预测更稳4.1 为什么单模型不够MaxEnt在presence-background模型里非常好用但我不建议把一个项目押在单个算法上。不同算法隐含的数据假设差别很大GLM偏线性拟合GAM能拟合非线性随机森林擅长捕捉复杂交互MaxEnt用的是最大熵原理。小样本场景下不同算法的排序和预测图经常出现明显差异。集成建模的思路很简单用同一套分布点和环境变量同时训练多个模型再根据各自表现分配权重最终把所有模型预测图加权平均。这样做可以降低单一算法带来的结构不确定性结果也更稳健。biomod2是R语言里最常用的集成建模框架它把数据格式化、伪缺样本生成、模型训练、评估、投影和集成封装在一起。4.2 biomod2建模流程与代码骨架biomod2的函数随着版本迭代改过几次名字这里给你一个当前常用思路的代码骨架。如果你的包版本不同以help(package biomod2)里的函数说明为准。library(biomod2) library(raster) myRespName - Bradypus myResp - rep(1, nrow(occ)) # 已知分布点全部标为1 myRespXY - occ # 坐标列 myExpl - stack(predictors) # 环境栅格stack # 第一步格式化数据同时生成伪缺样本(Pseudo-absences) myBiomodData - BIOMOD_FormatingData( resp.var myResp, expl.var myExpl, resp.xy myRespXY, resp.name myRespName, PA.nb.rep 2, # 重复生成2套背景点 PA.nb.absences 1000, # 每套背景点数量 PA.strategy random # 随机背景点 ) # 第二步训练多个模型 myBiomodModelOut - BIOMOD_Modeling( bm.format myBiomodData, modeling.id SDM_example, models c(GLM, RF, MAXENT.Phillips), NbRunEval 3, # 3次重复抽样 DataSplit 80, # 80%训练20%验证 Prevalence 0.5, VarImport 2, # 计算变量重要性重复2次 models.eval.meth c(TSS, ROC) )运行后get_evaluations(myBiomodModelOut)会给出每个模型在每次重复中的ROC和TSS分数。我只看这些分数能不能稳定在可接受范围而不会只挑最好的一次来写报告。4.3 集成预测与模型权重biomod2的集成模型会从已经训练好的单模型中挑出表现达标的模型做加权平均生成最终的集成预测图。myBiomodEM - BIOMOD_EnsembleModeling( bm.mod myBiomodModelOut, models.chosen myBiomodModelOutmodels.computed, em.by all, metric.select.thresh c(0.8, 0.85) ) myBiomodProj - BIOMOD_Projection( bm.mod myBiomodModelOut, new.env myExpl, proj.name current, selected.models all ) myBiomodEF - BIOMOD_EnsembleForecasting( bm.em myBiomodEM, bm.proj myBiomodProj )集成预测图会比单模型更平滑个别算法在环境外推时出现的极端值也会被稀释。但集成不是万能的如果所有单模型都共享同一个错误假设集成模型也会继承这个错误。比如所有模型都用了有偏的背景点那么集成之后只是把偏差平均了一下而已。5. 模型评估不是走流程AUC、TSS和阈值选择5.1 别拿训练数据评估模型MaxEnt在训练数据上的AUC经常能到0.95以上甚至接近1。为什么因为模型就是在这些数据上拟合出来的它当然认识这些点。评估必须用独立验证数据。最简单的做法是把分布点随机分成两组比如80%训练20%验证。dismo里这样操作set.seed(123) group - kfold(occ, k 5) occ_train - occ[group ! 1, ] occ_test - occ[group 1, ]然后在训练数据上重新跑模型用验证点加背景点评估。这里有一个更严格的方案是空间交叉验证比如使用blockCV包按空间格子和环境梯度划分数据。随机划分验证仍然会受到空间自相关的影响因为训练点附近同一个环境区的验证点可能过于相似导致AUC被高估。空间交叉验证更贴近真实应用场景尤其当你要外推到别的区域时。5.2 阈值选择把连续适生图变成有/无分布图MaxEnt输出的连续适生指数在管理上往往需要转换成二值图哪些区域属于“适生区”哪些属于“不适生区”。阈值的选择直接影响面积估算结果。常用方法有三种固定阈值比如0.5简单但缺乏依据。最小训练存在阈值Minimum Training PresenceMTP保证所有训练分布点都被划为适生区适合分布点很少且漏检风险高的场景。最大化敏感度与特异度之和MaxSSS它试图在“漏检”和“误报”之间找到平衡是目前论文里最常用的阈值之一。使用dismo的evaluate()与threshold()可以快速得到这些阈值# 用背景点作为假缺数据做模型评估 bg - randomPoints(predictors, n 1000) model_predictions - predict(me, predictors) test_pres - extract(model_predictions, occ_test) test_bg - extract(model_predictions, bg) ev - evaluate(p test_pres, a test_bg) evauc th - threshold(ev, spec_sens) th之后把这个阈值套用到预测图上就能生成有/无分布二值图。5.3 评估指标的实际解读局限AUC反映的是模型把出现点排序在背景点之前的概率。如果背景点不是真实缺数据只是在环境空间里随机抽的“可用点”那么AUC的意义就是排序能力而不是真实的判别能力。一个模型AUC0.8可以解释为“在随机出现点与随机背景点的比较中有80%的概率模型给出现点的得分更高”。这并不等于它能准确预测物种是否存在。TSS灵敏度特异度-1取值范围在-1到1之间比AUC更容易被管理者接受。一般来说TSS0.4时模型效果较弱0.4-0.6中等0.6较强。但所有阈值类指标都依赖验证数据的质量。如果你手里的出现点在环境空间里高度聚集或者背景点范围不合适指标数值都会失真。所以我始终建议把评估指标和预测图一并汇报不要只放一个孤零零的数字。6. 实战中踩过的坑空间自相关、未来气候和过度拟合6.1 采样偏差是SDM的隐形元凶很多物种分布数据库的采样坐标沿着公路、河流、城市和保护区密集分布偏远地区几乎没有记录。模型会把这种“采样努力程度”当成环境偏好结果在公路附近大幅预测出高适生区在真正偏远的适宜区反而预测得很低。这个问题在志愿者上报类数据里尤其突出。我个人的处理步骤是空间稀疏化让分布点在空间上尽量均匀。如果知道采样努力图把它纳入背景点生成策略比如在采样密集区多抽背景点。使用“目标组背景”方法即选择与目标物种采样方式相似的其他物种出现点作为背景代表“这个区域确实被采样过”。这比随机背景更贴近“可探测性”。空间自相关的另一个杀伤力是让验证AUC虚高。如果训练点和验证点来自同一片高强度采样区域模型在验证集上表现得很好迁移到新区域就崩盘。这就是为什么越来越多人推荐用空间分块交叉验证而不是随机交叉验证来报告模型性能。6.2 未来气候情景和数据外推的边界预测未来分布是SDM高频应用场景也是误用最多的场景。CMIP6提供的不同气候模型GCM在某些区域给出的未来温度和降水趋势甚至可能是相反的。你只挑一个GCM来预测相当于把所有赌注押在一个“未来故事”上。可靠做法是同时下载多个GCM、多个共享社会经济路径SSPs的数据分别投影然后输出预测结果的中位数和不确定性区间。还有一点常被忽略未来气候环境可能超出当前训练数据中的环境范围这种“新型环境组合”会让MaxEnt等算法外推到一个从未见过的空间。输出结果看起来高大上实际只是数学外推甚至完全违背物种生理限制。应对方法有两种一是使用predict()时设置clampTRUE把超出训练范围的变量截断在训练边界二是明确把外推区域标注为“模型外推区”不做精细解读。生态位保守性问题同样影响未来预测。模型假设物种的环境偏好不随时间改变但物种完全可能通过表型可塑性、进化或迁移到新微环境来适应变化。短期预测相对可信向2050年甚至2100年外推时不确定性会急剧放大。最稳妥的呈现方式是给出不同情景下一系列预测图而不是一张看起来非常确定的单幅图。6.3 过度拟合与参数简化MaxEnt虽然好用但默认参数配置并不是万能的尤其在小样本情况下特别容易过拟合。它的feature classes包括了Linear、Quadratic、Product、Threshold、Hinge等组合数据量小的时候模型会学出一堆不必要的高频细节预测图呈现碎片化变量响应曲线出现极端尖峰。解决思路是使用ENMeval包做模型调参比较不同feature class组合和正则化倍率下的AUC与AICc选择既稳又简的模型。不要只看AUC最高因为AUC最高往往伴随过度拟合。一个我习惯的检查方法是看预测图是否光滑、响应曲线是否符合生态常识如果一幅预测图只把已知分布点画成几个孤立的红点那基本就是过拟合了。变量数量也需要克制。我见过有人用20多个环境变量做一套几十个分布点的模型即使模型评估指标很高也不值得信任。有效的经验法则是模型复杂度不要超过样本量所能支撑的范围。一个样本量只有50左右的物种3到5个环境变量就已经很多了再多就是在背数据。再补充一个容易被忽略的细节将来建立模型后一定要保留一套“不参与任何训练和调参”的独立分布点数据只用于最终验证。如果能组织野外团队去模型预测的高适生但没有记录的空白区域做实地调查那比任何统计指标都有说服力。项目汇报时一套实地抽样照片和记录往往比几十页模型诊断图更能让决策者信任。我在实际项目中最后悔的往往不是模型选得不对而是在数据前期处理上太急躁。有一次把坐标系问题养到画图阶段才发现整整两天的工作全部作废。现在每接到一个新物种、新区域我都会先在R里画一遍分布点和环境边界确认它们落在同一片空间再开始碰模型参数。这个习惯帮我省下的时间比我学任何算法都快得多。如果你也想从零开始跑通SDM我的建议是别一上来就上集成模型和未来气候。先用dismo自带数据把MaxEnt全流程跑通理解预测图长什么样、响应曲线怎么读然后再加伪缺样本、模型评估、模型集成一步一步来。物种分布模型不是生成一张图的按钮它是由问题定义、数据质量、模型假设和现实约束共同决定的生态学工具。用得好它会成为保护和决策里的有力依据用得糙它只是装饰在PPT上的彩色地图。
返回列表