《Methods in Ecology and Evolution》：OccuGAMs: Non-linear occupancy and abundance modelling with imperfect detection
编辑推荐：
本综述系统介绍了将广义加性模型（GAMs）整合到分层占有率和丰度模型（HOAMs）框架中的方法学创新（OccuGAMs）。通过模拟数据和东南亚相机陷阱实证研究，证明OccuGAMs能更精准地捕捉物种对环境因子的非线性响应（如阈值效应和超丰度现象），相比传统多项式模型具有更好的曲线拟合能力和稳定性。文章为生态学家提供了在JAGS、Stan及R包mvgam/flocker中的实现代码，推动复杂生态关系的探索。
引言：生态监测技术的革新与统计模型挑战
近年来，相机陷阱和声学传感器等野生动物监测技术的飞速发展（Burton等，2015；Lahoz-Monfort与Magrath，2021），革命性地改变了在自然生境中研究动物种群的能力。这些工具能够经济高效地收集大量观测数据，尤其适用于偏远或难以进入的区域（Bruce等，2025）。生态学家日益依赖先进的统计模型来分析这些设备产生的大型检测历史数据，例如物种分布模型、分层占有率和丰度模型（HOAMs）、动态种群监测以及结构方程模型（Goldstein等，2024；Rozylowicz等，2024）。其中，分层模型因其能够明确地将不完美的检测过程与产生观测值的生态过程（即占有率和丰度等状态变量）分离开来而显得尤为强大和常见（MacKenzie等，2002；Royle，2004）。这一点至关重要，因为物种的可探测性会随着环境协变量（例如不同植被类型的可见度）和局部行为而变化，这可能使状态变量和栖息地关联的估计产生偏差（Gu与Swihart，2004；Kéry与Schmidt，2008）。
传统的分层占有率和丰度模型（HOAMs）通常在检测和状态变量方程中使用广义线性模型（GLMs）（Kéry与Schaub，2011；MacKenzie等，2002）。然而，GLMs相对刚性，在捕捉非线性栖息地关联或时间趋势的形状方面并非最优。解决此限制的一个常见策略是引入关注协变量的多项式项（Briscoe等，2021；Dehaudt等，2022；Ruiz-Gutiérrez等，2010）。但是，多项式近似存在若干众所周知的局限性。多项式项引入了全局基函数，这限制了响应曲线只能呈现平滑但本质上是非局部的行为（Harrell，2015；Magee，1998）。这通常导致不合理的对称性、边界不稳定性以及对异常值的高度敏感性（Aho等，2014；Gelman与Imbens，2019）。此外，多项式模型通常没有正则化机制来强制执行局部平滑。因此，协变量空间某一区域的数据可能对遥远区域的拟合产生不成比例的影响（Harrell，2015）。尽管存在这些严重限制，多项式仍然是模拟非线性关系时广泛使用的方法。
野生动物常常表现出非线性的时间趋势或对环境条件的类阈值响应。例如，在多个生态系统中观察到了阈值响应，包括大型动物丰度与人类基础设施距离的关系（Potvin等，2005）、美国森林鸟类与栖息地破碎化的关系（Morante-Filho等，2015），以及热带广布性哺乳动物与油棕种植园食物补贴距离的关系（Luskin, Albert,与Tobler，2017；Luskin, Brashares等，2017）。更复杂的非线性模式也可能通过生态和人为过程的相互作用而产生。然而，对162项动物生态学研究的回顾发现，只有14.2%的研究报告了对非线性的检验，这表明在拟合GLMs时，线性假设很少得到充分关注（Heit等，2024）。因此，需要更灵活的方法来将非线性关系纳入HOAMs，以更好地表征生态复杂性。
广义加性模型（GAMs）为非线性HOAMs提供了一种未被充分利用但强大的方法（Bled等，2013；Kéry与Royle，2020）。这类“OccuGAMs”通过惩罚样条允许对非线性效应进行平滑的、数据驱动的估计，从而扩展了HOAM框架。这些平滑函数表示为基函数的线性组合。在贝叶斯框架中，系数向量通常被赋予一个以零为中心的多元正态分布，其精度矩阵与惩罚矩阵的加权和成正比（Wood，2016）。惩罚矩阵定义了函数如何因粗糙度而受到惩罚，而平滑参数则控制着这种惩罚的强度。惩罚矩阵和平滑参数共同通过惩罚过度的“摆动性”来控制拟合函数的灵活性（Wood，2012017）。与多项式不同，多项式的每个系数是针对协变量的固定变换独立估计的，而惩罚样条系数是联合估计的（Wood，2012017）。因此，模型可以平滑地适应数据中的非线性模式，同时通过正则化限制过拟合（Wood，2012017）。此外，由于GAMs中使用的惩罚样条基函数具有局部支持性（它们仅在协变量范围的有限部分内活跃），单个数据点的影响仍然是局部的（Wood，2012017）。这在目标是理解栖息地或检测关联的功能形式时非常有价值（Strebel等，2014），例如当生态阈值引发急剧增加、减少或平台状形状时（Morse等，2003；Rhodes等，2008）。
尽管OccuGAMs具有优势，但迄今为止在野生动物生态学文献中仍很少出现，其在分层建模中的使用主要局限于模拟空间自相关（Bled等，2013；Rushing等，2019，2020）。这部分原因可能是拟合这些模型需要高级技能，即需要在贝叶斯建模框架（如JAGS、Stan或Nimble）中编写自定义代码。编码能力一直是生态学家更普遍面临的障碍（Bolker等，2013）。重要的是，近期的R包（mvgam和flocker）现在使得能够在熟悉的界面中使用“开箱即用”的功能来拟合带有平滑函数的HOAMs类型（Clark与Wells，2023；Socolar与Mills，2023）。具体来说，mvgam支持丰度（N-mixture）模型，而flocker支持占有率模型。然而，目前缺乏对OccuGAMs与其GLM对应模型之间的直接比较，这代表了理解它们相对性能和生态适用性的一个关键空白。
研究方法：模拟与实证分析框架
本研究通过模拟实验和实证案例研究，系统比较了传统HOAMs（包含线性和多项式项）与OccuGAMs的模型准确性和可解释性。
广义加性分层模型的一般表述
对于占有率模型，潜在占有率状态Zi被建模为以概率ψi为参数的伯努利随机变量。在传统模型中，占有率概率ψi通过logit连接函数与环境协变量关联，通常使用线性预测变量或多项式扩展。为了放松这些约束，可以将平滑项纳入线性预测变量中，用惩罚样条平滑函数s(xi)替代GLMs。检测概率pij同样可以使用GAM而非GLM来模拟检测协变量的非线性效应。观测到的检测数据yij则以真实占有率状态Zi为条件，建模为伯努利随机变量。
对于N-mixture模型（丰度模型），潜在状态是站点i的丰度Ni，通常建模为以平均预期丰度λi为参数的泊松随机变量。协变量效应可以像占有率模型一样，使用线性预测变量或平滑项来表达。观测计数yij则以Ni为条件进行建模，检测概率pij同样可以使用GAM或GLM建模。
模拟研究设计
模拟研究使用mvgam和flockerR包进行，两者均使用Stan作为后端。mvgam用于模拟和拟合N-mixture模型，而flocker用于占有率模型。模拟的关系包括简单（线性）和复杂（例如分段或非单调）趋势，以反映生态上真实的模式。
关键模拟要素包括：为四个假设物种模拟物种计数历史或存在-缺失数据，这些物种在平均丰度（或占有率）、基础检测概率以及对检测协变量的响应（线性或非线性）方面有所不同。潜在丰度（或占有率）从泊松分布（或伯努利分布）中抽取，其均值由物种特异性平均丰度（或占有率）、协变量效应和附加高斯噪声定义。协变量-响应关系使用六种功能形式之一建模：（1）单调阈值，（2）线性，（3）三次多项式，（4）分段阈值，（5）超丰度，（6）非单调阈值。研究设计在站点数量和重复次数上有所不同，代表小、中、大和超大样本量。每种情景重复50次，共产生4800个独特的计数历史（或占有率数据集）。
对每个模拟数据集，使用四种不同的N-mixture和占有率模型公式进行估计，这些公式在状态公式中是单变量的：（1）带有线性项的传统模型；（2）带有二次项的传统模型；（3）带有三次项的传统模型；（4）带有薄板回归样条的OccuGAM。每种模型类型还使用两种替代的检测子模型进行拟合：一种指定检测协变量的线性效应，另一种使用GAM。这导致每个数据集总共有八个模型变体，共估计76,800个模型。
模型评估使用两种适当的评分规则：能量评分和变差函数评分，两者均基于完整的后验预测分布计算。能量评分量化预测分布与真实值之间的整体差异，捕捉位置和分布上的差异，并对有偏或低估/高估不确定性的预测进行惩罚。变差函数评分评估预测之间成对距离的差异，对于评估模型是否正确捕捉底层响应表面的平滑度或相关结构特别有用。还使用归一化均方根误差（NRMSE）评估了替代检测模型对潜在丰度和占有率概率估计的影响。
东南亚热带哺乳动物对人为干扰的响应的案例研究
案例研究分析了来自东南亚10个低地热带森林景观的21次相机陷阱调查的数据。将相机位置聚合到3平方公里2的六边形采样单元中，以生成空间标准化的数据集。研究了四种哺乳动物（豚尾猴、野猪、水鹿和麂鹿）在异质环境中的占有率和丰度。包含了四个人为干扰度量作为协变量：森林完整性、人类足迹指数、油棕覆盖率和森林覆盖率。
对于占有率分析，估计了物种在网格单元（3平方公里2六边形）内存在的概率。潜在状态被建模为伯努利随机变量。占有率概率被建模为单个环境协变量的平滑函数，包含每个景观-年份的随机截距和具有K个基函数的薄板回归样条平滑。为了评估协变量效应的替代表示，还估计了递增阶数的多项式变体（线性、二次、三次）。检测概率包含了相机努力度和一个过度离散随机效应（ODRE）。
对于丰度分析，估计了站点i的潜在丰度Ni。潜在状态被建模为泊松随机变量，其对数均值包含一个景观-年份随机截距和一个具有K个基函数的薄板回归样条。与占有率模型类似，也估计了每个干扰协变量的三种多项式公式（线性、二次、三次）。观测模型与占有率模型相同，包含了相机努力度和ODRE的效应。假设在调查期间种群是封闭的。
模型评估包括计算GAM与非GAM模型估计关系之间的NRMSE。还通过计算贝叶斯p值和C-hat值（通过后验预测检查）来检查模型拟合优度和过度离散性。使用留一法（LOO）交叉验证比较模型预测性能，计算每个模型的预期对数点wise预测密度（elpd），并使用∣Δelpd∣ > 2σ作为评估预测准确性差异强度的标准。
研究结果：模型性能与生态启示
模拟数据中的功能形状恢复
OccuGAMs在恢复真实的丰度和占有率关系方面优于线性、二次和三次公式，且在不同样本量下表现稳定。多项式模型的性能随着数据量增加而改善，反映了其准确估计曲率对更大样本的依赖性。基于变差函数评分的排名大致相似。
在所有情景中，多项式模型有时能捕捉非线性响应的广泛特征，但通常缺乏灵活性来恢复更复杂或不规则的形状。在某些特定情况下存在匹配，例如，二次模型能很好地近似单调阈值关系。但当真实关系是线性时，多项式模型表现更差。相比之下，OccuGAMs具有相对鲁棒性。在线性情景中，它们排名前两位，仅次于线性模型，并在所有响应形状下保持强劲性能。然而，OccuGAMs的性能受物种平均丰度的影响。在N-mixture分析中，OccuGAMs对于平均丰度为2的物种（物种二）表现较差。
所有模型的绝对性能随着样本量减少而下降，这与在有限数据下重建功能响应能力降低一致。在10个站点的样本量下，非线性函数的估计不可靠，具有高方差且明显偏离真实响应。性能在样本量为25及以上时显著改善。
检测子模型对占有率和丰度估计的影响
对于N-mixture模型，当真实检测概率随检测协变量非线性变化时，将GAMs纳入检测子模型改善了潜在丰度的估计。非线性检测函数显著影响N-mixture模型估计，使用线性检测子模型通常产生更高的NRMSE和较差的真实丰度恢复。通过整合GAM使检测子模型匹配非线性协变量，降低了NRMSE，提高了跨状态模型公式的丰度估计保真度。使用线性检测子模型估计非线性协变量的N-mixture模型也更容易产生不可接受的R-hat值。
对于占有率模型分析，各种子模型公式之间没有明显差异。NRMSE基本上不受检测协变量形式或检测子模型的影响。与N-mixture模型相比，占有率模型出现参数R-hat偏高的可能性较小。
真实数据中传统HOAMs与OccuGAMs的拟合
由于真实关系未知，比较性质不同。OccuGAM和多项式模型在模型拟合优度或过度离散方面没有明显差异。所有N-mixture模型的贝叶斯p值介于0.25和0.75之间，C-hat值小于1.1，表明模型拟合可接受且没有过度离散的证据。
大多数模型在物种和变量间显示出高度一致性，表现为传统模型和OccuGAMs趋势方向的一致性。总体而言，占有率与四个干扰协变量之间的栖息地关联不如相对丰度明显——大多数是平坦或略微增加的。然而，显著的例外是水鹿~油棕和麂鹿~油棕的占有率存在强烈的负相关关系。
NRMSE表明，最接近GAM产生的曲线的栖息地关联曲线因物种和协变量而异。没有观察到系统性的偏差偏向更简单或更复杂的模型。关键的是，有三个例子中GAM明显偏离所有模型类型：豚尾猴~森林覆盖率、麂鹿~森林覆盖率和豚尾猴~人类足迹。在前两个例子中，GAM产生了相对丰度与森林覆盖率之间初始的平坦关系，随后是一个 punctuated increase，这被二次和三次模型都近似得不好。在后者中，相对丰度随人类足迹出现初始的急剧下降，随后是平台期，这被三次和二次模型捕捉得很差。
对相对丰度N-mixture结果的定性解释表明，有八种组合（50%）中，GAM模拟的栖息地关联明显是非线性的。例如，水鹿对油棕覆盖率的响应，如GAMs所模拟的，是非线性的，反映了相对丰度的初始增加，随后在油棕覆盖率大于30%–40%的区域急剧下降。豚尾猴对油棕的存在表现出积极的响应，但其丰度是逐渐线性增加的，缺乏明确的拐点。野猪OccuGAM显示对附近油棕覆盖率的响应初始平坦，随后当油棕覆盖率超过约60%时，相对丰度呈近指数增长，这反映了先前的研究。野猪二次模型与GAM最相似。
N-mixture模型预测性能
N-mixture模型的预测性能在所有物种-协变量组合中相似，差异很少超过2σ。OccuGAMs在16个案例中的3个（18.75%）是表现最好的模型。然而，OccuGAMs在16个案例中的14个（87.50%）Δ elpd < 2σ，表明预测性能通常与最佳表现模型没有显著差异。简单线性模型在6个案例中（37.50%）具有最高的预测准确性，二次模型在2个案例中（12.50%），三次模型在5个案例中（31.25%）。
当三次模型表现最佳时，GAM的功能形状通常暗示了高度非线性的关系，尽管豚尾猴~油棕是一个明显的例外。不出所料，GAMs比它们的三次等效模型灵活得多，当最佳表现模型简单或复杂时，通常能实现Δ elpd < 2σ。
GAMs与多项式N-mixture模型预测响应之间的差异
GAMs与多项式对应模型之间拟合关系的NRMSE概述表明，当NRMSE随着多项式次数增加而增加时，GAM暗示了一个较不复杂（线性）的关系。或者，NRMSE可以随着多项式次数增加而减少，这意味着GAM暗示了一个更复杂的关系，甚至三次多项式也无法充分捕捉。
有4个（25%）案例中线性模型的NRMSE低于二次和三次模型，这意味着GAM暗示了一个“简单”的线性关系。二次模型在6个（37.5%）案例中最接近地近似了GAM，而三次模型在其余6个（37.5%）案例中具有最低的NRMSE。值得注意的是，这四种情况中有三种包含单个协变量（森林覆盖率）。物种间的最低NRMSE没有明显趋势。然而，一些协变量特异性模式是明显的：对于四个物种中的三个，二次模型对森林完整性具有最低的NRMSE。相比之下，对于四个物种中的三个，三次模型最接近地反映了GAM对森林覆盖率的预测。
讨论：OccuGAMs的生态学应用与展望
本研究证明了将GAMs整合到分层占有率和丰度模型（HOAM）框架中的实用性，可以在考虑不完美检测的同时，灵活地模拟物种-栖息地关系。与假设一致，OccuGAMs在模拟情景中大多数情况下优于多项式模型公式，能够更准确地恢复复杂的、非线性的物种-栖息地关系。此外，在真实关系是线性的情景下，OccuGAMs始终排名第二，仅次于线性模型。这展示了它们适应更简单功能形式的能力——这是相对于多项式公式的一个优势，多项式公式常常由于虚假曲率而表现更差。然而，对于低丰度物种，准确性下降，反映了
