拥抱不确定性：海虱幼虫扩散的集成模型

《Ecological Modelling》：Embracing uncertainty: Ensemble models of sea lice larval dispersal

　　蒂姆·M·塞夫奇克（Tim M. Szewczyk）|德米特里·阿列尼克（Dmitry Aleynik）|海伦娜·C·雷纳迪（Helena C. Reinardy）|金·S·拉斯特（Kim S. Last）|安德鲁·C·戴尔（Andrew C. Dale）
英国高地与群岛大学（University of the Highlands and Islands），UHI House，Old Perth Road，因弗内斯（Inverness）IV2 3JH

**摘要**
海虱，尤其是Lepeophtheirus salmonis，对鲑鱼养殖业带来了重大的生态和经济挑战。它们的扩散通过浮游幼虫阶段进行，这些幼虫可能因物理和生物过程的相互作用而广泛传播。直接观察海虱非常困难，因此通常使用生物物理粒子追踪模型来估计感染压力。然而，模型参数化存在相当大的不确定性。我们评估了多种集成方法，旨在利用这种不确定性来改进对养殖场海虱分布的预测。利用2021至2024年间苏格兰鲑鱼养殖场的公开数据以及20种参数化方法，我们比较了三种集成方法：简单平均法、具有空间变化混合权重的方法以及机器学习方法。通过使用成年等效感染压力（根据鱼体上的桡足类压力进行滚动计算）来评估每种方法的性能，以预测每条鱼上的平均成年雌性海虱数量。交叉验证显示，机器学习集成方法在预测养殖场海虱动态方面表现更优，尽管需要更直接的数据来预测幼虫感染压力。总体而言，具有空间变化权重的集成方法优于其他方法，特别是对于空间模式的预测。重新采样实验验证了这些结果在不同方法和集成规模（3-15个）下的稳健性。因此，集成建模是提高海虱感染压力预测的一种有前景的方法，同时具备根据特定目标定制集成的灵活性。更好的海虱扩散预测有助于更好地管理感染情况，最终改善养殖鱼类的福祉和产量，并减少野生鲑鱼的寄生虫负担。

1. **引言**
海虱，尤其是Lepeophtheirus salmonis，是鲑鱼的体外寄生虫，对鲑鱼养殖业和野生鱼类种群构成了严重威胁（Dempster等，2021；Ives等，2023；Overton等，2019）。其扩散过程包括卵黄营养的浮游幼虫阶段，这些幼虫包括两个无感染性的溞状幼虫阶段和一个有感染性的桡足类阶段。桡足类能够寻找并结合到宿主鱼身上，在那里它们会经历四个进一步的发展阶段（chalimus 1-2、pre-adult 1-2），最终成为成虫（Piasecki等，2023；Wootten等，1982）。在挪威和苏格兰，为管理海虱感染采取了主动措施，包括定期监测养殖场中的海虱数量、规定必需或最佳实践的抗虱处理措施，以及根据预测的海虱压力对野生鱼类种群的影响来调整养殖场的规模和位置（SEPA，2025；Vormedal，2024）。抗虱处理可能会增加养殖鲑鱼的压力和死亡率，因此仅在虱子造成的危害超过这些副作用时才会使用（Gillson等，2022；Meyer等，2019；Overton等，2019；Szewczyk等，2024b）。然而，限制成年海虱的数量对于减轻疫情严重程度和减少对区域浮游幼虫种群的影响至关重要（Kragesteen等，2019）。因此，做出有效决策依赖于对早期阶段的准确了解。不幸的是，准确检测早期海虱阶段十分困难。由于检测率低，计数早期阶段的海虱容易出错（Moriarty等，2024）。浮游幼虫体积小（<500 μm），密度低（<1个/立方米），且与其他浮游桡足类难以区分，这限制了目前原位监测的范围和可靠性（Moriarty等，2024；Murphy等，2024；Rabe等，2024；Thevar等，2023；Tryggestad，2023）。因此，监管和政策决策通常依赖于容易识别的养殖场成虫和pre-adult阶段的监测，结合由海洋物理水动力模型驱动的幼虫扩散生物物理粒子追踪模型。例如，挪威采用“交通灯系统”对来自养殖业操作对野生鱼类的风险进行分类，该系统基于多种信息来源，包括来自养殖和野生鱼类的数据、观察哨兵围栏（将鱼类短暂放入围栏中并彻底计数海虱的数据）以及幼虫感染压力的模型预测（Jensen等，2024；Myksvoll等，2018；Sandvik等，2020；van Nes等，2025；Vollset等，2018）。海虱幼虫扩散模型涉及相当高的复杂性，旨在捕捉影响幼虫时空分布的关键过程（Garnier等，2024；Jonsen等，2016, 2014；Moriarty等，2024；Sandvik等，2020）。粒子追踪模拟依赖于详细的水动力模型，这些模型模拟了包括三维水流、盐度和温度在内的多种变量。模拟还结合了海虱生物学的各个方面，特别是幼虫发育时间、死亡率和垂直行为，这些因素可能依赖于所模拟的环境（Garnier等，2025；Moriarty等，2024；Murray等，2022）。这些模型中存在显著的不确定性（Moriarty等，2024；Murray等，2025）。敏感性分析主要研究单个参数或过程，表明预测结果受多种生物和物理参数的影响，如垂直行为（包括游泳速度）、迁徙的环境触发因素和昼夜光周期（Garnier等，2024；Jonsen等，2016）、每个养殖场周围的假设附着半径（Moriarty，2023）、养殖场内的发育情况（Kragesteen等，2023, 2021）以及扩散性（Moriarty，2023）。结合生物和物理过程的不确定性和自然变异性，可以从各种合理的参数化设置中生成多种预测结果。适当优化参数空间是一个巨大的挑战。

集成建模提供了一种替代方案，它考虑了一组合理模型之间的变异性。在集成方法中，一组组成部分模型的结果被用作集成元模型中的预测因子。集成模型在其他领域已被证明非常有效（Araújo和New，2007；Jones和Cheung，2015；Szewczyk等，2025；Wu和Levinson，2021），方法的多样性（Sagi和Rokach，2018；Wu和Levinson，2021）为海虱扩散建模提供了有希望的方向（Moriarty，审稿中）。在最简单的情况下，通过简单平均法将组成部分模型的感染压力进行加权平均（Moriarty，审稿中）。更复杂的技术通过加权输入或估计组成部分与响应变量之间的功能关系来优化预测性能。例如，机器学习方法可以识别预测因子及其相互作用中的复杂非线性模式，从而提高预测准确性（Pichler和Hartig，2023）。值得注意的是，这些模型是根据性能指标为特定目标优化的，因此可以针对特定应用进行定制。此外，多级贝叶斯建模是一种灵活的统计方法，能够自然适应复杂的数据结构和未观测（即潜在）变量，并对不确定性进行概率考虑（B?hn等，2022；Gelman等，2015；Szewczyk等，2021；Szewczyk和McCain，2019）。

在这项“概念验证”研究中，我们旨在评估集成模型是否能够改善海虱感染压力和动态的预测。我们使用20种粒子追踪模型参数化来评估三种集成模型在预测鲑鱼养殖场每条鱼上的平均成年海虱数量方面的性能：EnsAvg使用组成部分的算术平均值；EnsBlend使用具有空间变化权重的组成部分平均值；EnsFcst使用机器学习算法。我们专注于苏格兰西部，使用该地区大西洋鲑鱼养殖场近四年的（2021-2024年）已发布的数据，包括每家养殖场的海虱感染压力以及每种参数化下的模型预测。我们通过重新采样来评估集成性能的普遍性和集成规模的效果。确定最成功的方法将为管理者和决策者提供更可靠的工具集，以实现更可持续的养殖，同时减少对野生种群的影响。

2. **方法**
2.1. **水动力学**
海虱扩散模型使用了WeStCOMS v.2水动力模型（图1）。该模型涵盖了苏格兰西部，包括外赫布里底群岛，模型输出包括2019年4月1日至今的每小时数据（水流、温度、盐度、湍流参数）。WeStCOMS（Aleynik等，2016；Davidson等，2021）基于FVCOM平台（Chen等，2012）构建，该平台使用非结构化三角形网格，能够有效适应沿海地区或复杂地形特征，以及10个随浴深变化的垂直层。沿海海洋模型嵌套在2公里的NEAROMS区域大陆架模型（Dabrowski等，2014）和CEMES-NWS（https://doi.org/10.48670/moi-00054）中，其表面由基于WRF模型（Skamarock等，2019）的SAMS区域高分辨率大气强迫实现（2公里）和1公里全球表面温度（MUR-SST，http://dx.doi.org/10.5067/GHGMR-4FJ04（Armstrong等，2012）驱动。WeStCOMS数据产品可通过https://thredds.sams.ac.uk下载。
图1. 建模范围。蓝色区域表示WeStCOMS2 v.2的模型范围。洋红色开放点表示已许可的养殖场位置。图10中显示了具有示例性垂直分布的地点：1）托里登湖（Loch Torridon）、2）卡伦湖（Loch Carron）、3）阿尔什湖（Loch Alsh）、4）霍恩湖（Loch Hourn）、5）马克岛（Isle of Muck）、6）上林恩湖（Loch Linnhe）、7）斯佩尔夫湖（Loch Spelve）和8）芬湖（Loch Fyne）。

2.2. **鲑鱼和海虱数据**
在苏格兰的鲑鱼养殖场，操作人员至少对每个养殖场的5个围栏中的25条鱼进行每周一次的海虱计数，通常包括chalimus、pre-adult和成年阶段（Moriarty等，2024；SSPO 2015）。苏格兰政府鱼类健康检查员报告并审核了每条鱼上的平均成年雌性L. salmonis数量，并通过在线数据门户公开发布（https://aquaculture.scotland.gov.uk）。2025年4月9日可以获取海水养殖场位置和每月鱼类生物量数据（2002年至今，由苏格兰环境保护局管理），以及每周的海虱处理行动和海虱数量数据（2021年3月29日至今，由海洋局管理），完整数据可追溯到2024年12月。所有位于WeStCOMS v.2时空范围内的养殖场都被纳入分析。当报告的鱼类生物量大于零时，养殖场被视为活跃养殖场，海虱数据位置与已许可的养殖场位置进行了空间关联。当活跃养殖场未报告海虱数据时，数值插值为所有养殖场周中位数中的较小值或总体中位数。每周和每月数据分别与相应的周一或月初对齐，然后海虱数据被线性插值为每日值。活跃养殖场j在第i天的成年雌性海虱数量估计为Aij=Aˉij*Bjw，其中Aˉij*是插值的每条鱼上的平均成年雌性海虱数量，Bj是允许的最大鱼类生物量（Moriarty，2023），w = 240是每吨鱼的估计数量（Aleynik等，2022）。统计建模和验证仅使用了每个生产周期开始直到首次报告抗虱处理期间的实际养殖场海虱计数（图2），以避免处理带来的影响，尽管模型对此选择的敏感性相对较低（图A.7）。

图2. 建模概述。(a) 2021-2024年的模拟使用了WeStCOMS2模型的每小时水动力数据，根据报告的每条鱼上的平均海虱数量和插值到每日值的鱼类生物量来初始化幼虫密度。(b) 在K = 20次模拟中模拟了不同参数化的海虱扩散。对于每个养殖场j在第i天的每个参数化k，使用每日桡足类感染压力IPijk来计算成年等效感染压力AEIPijk，该压力考虑了发育时间和养殖场内的死亡率。(c) 使用AEIP的贝叶斯模型为每种参数化拟合了每条鱼上的平均海虱数量，以与集成模型进行比较。(d) 混合集成（EnsBlend）估计了每种参数化的比例混合权重。平均集成（EnsAvg）使用了每日平均AEIP。预测集成（EnsFcst）使用机器学习方法来提前一周优化每个农场的预测结果。(e) 来自EnsBlend的混合权重后验分布被用来生成整个区域中桡足类分布的集成预测，这些预测基于各自的参数化模型。2.3. 海虱扩散模型海虱幼虫的扩散是使用biotracker v.1.0.0模拟的，该模型是一个拉格朗日粒子跟踪模型，考虑了生命周期发展、死亡率和主动垂直行为。水流、温度、盐度和表面短波辐射数据来自WeStCOMS，其水动力“景观”每小时更新一次。每个粒子的水动力条件是根据其上方和下方最近三个节点或元素质心的距离加权平均值计算得出的，平流过程采用四阶龙格-库塔算法，时间步长为120秒（Adams等人，2021，2015）。水流数据在每次插值时更新到粒子位置，而其他变量和过程（如行为和死亡率）则每120秒更新一次。每个小时从农场j释放的每个海虱幼体初始数量表示为fAij24N，其中f是每只成年雌虫每天产生的幼虫平均数，N是每小时每个农场释放的粒子数量。日繁殖力要么保持恒定，为每只雌虫28.2个幼体（Heuch等人，2002），要么根据文献中的汇总数据（Hamre等人，2019；Kragesteen等人，2021，2023；Samsing等人，2016b；Stien等人，2005）建模为每小时温度的逻辑函数，即f=65.5/(1+e?0.550*(Temp?8.49))+10.5。幼虫发育依赖于积温（°C · d）的累积，两个幼体阶段总共持续γ个积温，而感染性桡足类阶段从γ到150个积温，之后所有海虱都被假设为衰老（Adams等人，2021）。这两个幼体阶段没有区别，因此我们用“幼体”来指代两者。幼体死亡率保持恒定，为m=0.01h?1，或者根据Bricknell等人（2006）的数据建模为盐度的逻辑函数（psu，图A.1b），即m=0.941/(1+e0.594*(Sal.?13.5))+0.0205。粒子根据水流和扩散水平移动，而垂直移动则由水流、扩散和海虱的行为决定。扩散被模拟为Δrd～U(?1,1)6DΔt，其中Δrd是位移，U是从均匀分布中随机抽取的值，D是随模拟变化的扩散系数（Adams等人，2015；Garnier等人，2025；Myksvoll等人，2018）。向下移动是对低盐度的响应，向上游动是对光线的响应（μmol photon s-1 m-2；正向趋光性，与逆向昼夜垂直迁移一致（Hagemann等人，2024）），其中向下刺激具有优先权（Johnsen等人，2016，2014；Myksvoll等人，2018）。每次主动移动时，垂直游速是从标准差设置为平均值1/5的正态分布中抽取的，以模拟个体间和时间上的变化。幼体速度是桡足类速度的一半。对盐度的回避表现为从上阈值处的0%线性增加到下阈值处的100%（Sandvik等人，2020）。桡足类对光刺激比幼体更敏感（Johnsen等人，2016，2014；Myksvoll等人，2018；Sandvik等人，2020）。在表面可见光L0计算为从WeStCOMS提取的全谱辐射的42%的情况下，水柱中光的衰减被设定为k = 0.2，以便在粒子深度z米处的可见光为Lz=L0e?kz。2.4. 构成参数化我们比较了20种粒子跟踪模型参数化（图2b，表1），以捕捉一系列合理的扩散模型。这些参数化模型在2D与3D、扩散率（Dh，Dv）、卵产量（f）、幼体死亡率（m）、发育率（τ）、向下游动的盐度阈值、向上游动的光阈值以及上下游动速度等方面有所不同。我们对每种参数化中的每个参数都使用了独立的随机抽取（图A.2），并施加了三个限制。首先，有四种2D模拟强制所有卵产量和幼体死亡率的组合，以及十六种3D模拟，在这些模拟中这些属性是随机选择的。其次，上盐度阈值被限制为高于下盐度阈值。第三，对游动速度和光阈值的选择应用了0.8的相关系数，这样具有更快上游动速度的模拟倾向于具有更快的下游动速度，而对光敏感度高的幼体具有高桡足类光敏感度的模拟（图A.3）。因此，这些随机参数化但受合理范围限制的扩散模型代表了一系列合理的分布，反映了海虱幼体时空分布的不确定性。表1. 海虱扩散模拟参数化。每个模拟的值是从定义的参数空间中随机选择的。幼体游速设定为桡足类速度的一半。模拟参数化扩散系数（D；m2 s-1）卵产量（f；d-1）幼体死亡率（m；h-1）蜕变为桡足类（γ；°C d）桡足类游速（cm s-1）盐度阈值（psu）光阈值（μmol photon s-1 m-2）水平向上向下幼体桡足类3D.13.76 × 10–3 1.91 × 10–2 f（温度）f（盐度）35.0 0.37 40.58 92 1.2–2 1.84.97 × 10–1 1.76 × 10–4 3D.27.24 × 10° 1.03 × 10–3 f（温度）0.01 37.0 1.61 2.00 20.2–2 0.81.40 × 10–1 7.35 × 10–5 3D.31.21 × 10° 6.57 × 10–2 f（温度）0.01 36.1 0.03 50 0.01 47 21.0–2 4.61.27 × 10–1 8.15 × 10–5 3D.41.32 × 10° 3.41 × 10–5 28.2 0.01 43.9 1.94 1.46 26.5–3 2.05.30 × 10–2 2.20 × 10–6 3D.55.74 × 10–1 6.48 × 10–2 28.2 f（盐度）43.4 1.98 1.76 24.9–2 9.61.33 × 10–1 2.81 × 10–5 3D.61.24 × 10–3 6.17 × 10–4 28.2 0.01 36.8 0.27 10.66 327.7–3 2.01.98 × 10–14.81 × 10–5 3D.73.27 × 10–3 1.42 × 10–2 28.2 0.01 33.6 0.01 53 0.01 83 21.5–2 6.14.58 × 10–1 1.08 × 10–4 3D.82.36 × 10–3 1.10 × 10–4 f（温度）0.01 35.3 0.63 70.96 727.7–3 0.41.62 × 10–1 1.54 × 10–5 3D.93.13 × 10–3 4.54 × 10–2 28.2 0.01 48.0 0.38 50.11 221.1–2 4.41.17 × 10–1 1.11 × 10–5 3D.10 3.91 × 10° 3.22 × 10–3 28.2 f（盐度）30.4 0.69 91.36 23.7–2 4.46.01 × 10–2 6.74 × 10–6 3D.11 8.96 × 10–3 5.96 × 10–5 28.2 f（盐度）47.8 0.24 11.31 26.7–2 7.71.65 × 10–1 1.38 × 10–4 3D.12 1.07 × 10° 2.59 × 10–2 f（温度）f（盐度）49.7 0.27 30.12 427.4–3 2.04.64 × 10–1 1.76 × 10–4 3D.13 7.22 × 10–3 2.01 × 10–5 f（温度）f（盐度）49.4 0.79 40.34 520.7–2 1.51.09 × 10–1 3.33 × 10–5 3D.14 2.55 × 10–1 1.66 × 10–3 f（温度）f（盐度）39.9 0.94 31.39 24.5–2 9.2 2.08 × 10–1 1.42 × 10–4 3D.15 8.98 × 10–12.08 × 10–4 f（温度）0.01 40.7 0.03 60 0.02 2226.6–3 2.04.55 × 10–1 1.60 × 10–4 3D.16 1.43 × 10–2 8.50 × 10–2 28.2 0.01 49.4 0.84 40.66 122.2–2 3.32.49 × 10–16.31 × 10–5 2D.11.16 × 10–1—28.2 0.01 41.8—————2D.21.93 × 10°—f（温度）0.01 41.1—————2D.34.50 × 10–1—28.2 f（盐度）47.1—————2D.41.11 × 10–2—f（温度）f（盐度）39.5—————2.5. 感染压力每个农场的桡足类感染压力（IP）是每天在许可农场位置周围100米半径内的顶部20米水深中每个时间步长的桡足类数量之和。这些总和根据WeStCOMS网格内的100米半径范围（即不包括陆地面积）进行缩放，以估计每日桡足类数量·m-3，并除以24得到每小时率。每个WeStCOMS网格元素的桡足类平均密度·m-3·h-1也以类似的方式按周计算，限制在顶部5米，以对应潜在宿主鱼类的典型深度（Davidsen等人，2008；Eld?y等人，2017；Plantalech Manel-La等人，2009）。为了更好地评估扩散参数化效果，桡足类发育为成虫的过程保持不变。我们假设桡足类附着到成虫阶段之间有固定的五周延迟（Revie等人，2005；Toorians和Adams，2020）。在这个假设下，每日IP用于计算累积的预期成虫等效感染压力（AEIP），考虑了发育时间和死亡率，每个加数都以与桡足类IP相同的单位表示（lice·m-3·h-1）。附着过程没有直接建模，而是通过下面描述的统计模型的回归系数隐式捕获。在农场j的第i天，使用前几天的桡足类感染压力IPij来计算AEIPij：(1)AEIPij=∑t=35T=150IPi?t,js0τ0s1τ1s2τ2s3τ3，其中IPi?t,j是t天前的桡足类感染压力，s是每日存活率（未成熟期：s0=0.970；前成虫期：s1=0.998；成虫期：s2=0.997；怀孕成虫期：s3=0.942），τ是每个生命阶段持续的天数（未成熟期：τ0=15；前成虫期：τ1=20；成虫期：τ2=min?(20,t?τ0?τ1)；怀孕成虫期：τ3=max?(0,t?τ0?τ1?τ2)）。因此，AEIPij提供了从模型预测的35天到150天前桡足类感染压力预期产生的成年海虱密度的相对估计，根据文献中的人口统计率进行了调整（Kragesteen等人，2023；Revie等人，2005；Stien等人，2005；Toorians和Adams，2020）。集成权重对s和τ的变化显示出低敏感性（补充材料）。请注意，遇到和附着过程没有直接表示，IP和AEIP最好解释为暴露的代理指标。对于验证数据集中的每个日-农场组合，都计算了每个构成参数化的AEIP。对于所有统计建模，由于分布严重右偏且界限为零，观察到的每条鱼的平均虱子数量和AEIP进行了四次方根转换，以改善模型诊断。2.6. 构成模型对于每个构成参数化k = 1：K，开发了一个贝叶斯统计模型（图2c），其中第i天在农场j的每条鱼的平均虱子数量y是AEIPijk的四次方根 scale 的函数。因为响应变量是未知数量的鱼的平均海虱计数且没有膨胀，所以使用了障碍对数正态模型，该模型使用逻辑模型来估计Pr(yij = 0)，并对yij > 0的观察值使用对数正态模型，即：(2)yij～hurdlelogNorm(μij,σ,zij)(3)logit(zij)=α0+α1AEIPijk(4)log(μij)=β0+β1AEIPijk，其中μij和σ是每条鱼平均虱子的均值和标准差，zij是yij = 0的概率，α0:1是障碍回归的系数，β0:1是对数正态回归的系数。使用拟合模型的联合后验分布，计算出预测的每条鱼的平均虱子数量为y^ij=(1?zij)×μij。2.7. 集成模型比较了三种集成方法，每种方法都使用了K = 20个构成参数化（图2d）的AEIP。首先，平均集成（EnsAvg）使用了构成AEIP的算术平均值。其次，混合集成（EnsBlend）为每个构成参数估计了空间变化的混合权重。第三，预测集成（EnsFcst）使用机器学习算法并将构成参数作为预测因子。2.7.1. EnsAvg：平均集成平均集成（EnsAvg）展示了一种构建集成的简单方法。集成AEIP是作为每天和每个农场的构成AEIP的算术平均值计算的，并在上述障碍对数正态模型中用作预测因子。2.7.2. EnsBlend：混合集成混合集成（EnsBlend）确定了每个构成参数化的最佳权重。障碍对数正态模型扩展为yij是集成AEIP?的函数，AEIP?是一个潜在的（即未观察到的或中间的）变量，作为构成AEIP的加权平均值计算得出：(5)yij～hurdlelogNorm(μij,σ,zij)(6)logit(zij)=α0+∑k=1KαkAEIPijk(7)log(μij)=β0+β1AEIPij?(8)AEIPij?=∑k=1KpjkAEIPijk(9)pjk=epjk*∑k=1Kepjk*(10)pjk*=sk(latitudej,longitudej)+θjk，其中pjk是农场j的模拟k的权重，满足0≤pjk≤1且∑pj,1:K=1（即每个位置的权重是组成性的：比例之和为1）。通过对pjk的组成约束，使用softmax函数（逻辑函数的多变量泛化）实现了这种约束。未受约束的表示pjk*被建模为使用B样对于每种算法，我们使用了如下所述的10折交叉验证，并对每个超参数实施了最多100个值的网格搜索。调优是分别针对下面详细说明的每个评估指标进行的。2.8. 评估模型性能使用10折交叉验证进行评估，划分时间以确保每个分区代表约4个月。所有预测和观测值都经过四次方根转换。利用交叉验证的预测结果，我们通过计算均方根误差（RMSE）和Spearman’s rho（ρ）来评估每个模型在预测每个农场每周每条鱼观察到的平均虱子数量方面的性能。接收者操作特征曲线下面积（AUCROC）是通过将观测值分类为每条鱼0.5只虱子以上或以下来计算的。指标是根据数据分区计算的：1）按农场划分，以评估给定农场所有周的平均虱子数量预测性能；2）按周划分，以评估给定周内所有农场之间的平均虱子数量的空间模式预测性能。2.9. 重采样以评估泛化能力我们使用重采样来评估集合相对于其组成部分的性能，以及研究对集合大小的敏感性。从20种参数化方法中，我们随机选择了3种、5种、10种或15种不放回的模拟方法，每种方法有100个重复样本（总共400个样本）。对于每种样本大小，都按照上述交叉验证程序拟合EnsBlend和EnsAvg，并计算评估指标。2.10. 软件分析使用R 4.4进行了分析。贝叶斯模型是用Stan（https://mc-stan.org）编写的，并在R中用rstan 2.32.6运行。所有先验分布都提供了少量信息以提高采样效率（表A.1）。模型在Stan中运行，每个链有2000次暖启动和1000个后验样本（总共3000个）。我们确保在采样过程中没有出现发散现象，对于非混合参数，R-hat < 1.01，并且对于所有参数，有效样本量都大于1000。预测集合模型使用tidymodels R包进行参数化。粒子跟踪模拟使用biotracker v.1.0.0（https://doi.org/10.5281/zenodo.14811865）进行，并通过R包biotrackR v.0.1.0（https://doi.org/10.5281/zenodo.14716634）执行和处理。与模拟和分析相关的所有R代码都可以在https://doi.org/10.5281/zenodo.17064093获取。模拟使用了2021年4月1日至2024年12月31日的全部数据跨度来评估组成参数化和集合的性能。还使用相同的参数化方法对2024年3月17日至2024年5月31日的数据进行了额外的模拟，以研究每个农场每小时的桡足类动物的垂直分布。3. 结果3.1. 交叉验证预测组合（EnsFcst）在平均性能上优于所有其他模型（图3）。最佳性能的算法取决于评估指标。在四种情况下选择了随机森林，在另外两种情况下选择了装袋神经网络（农场内：ρ，AUCROC）。PCA仅在原始AEIP值上被选为农场内ρ的指标。下载：下载高分辨率图像（318KB）下载：下载全尺寸图像图3. 交叉验证中的模型性能。性能评估使用了三个指标（均方根误差（RMSE）、Spearman’s ρ和接收者操作特征曲线下面积（AUCROC），这些指标是在每周和每个农场内计算的。更好的性能对应于较低的RMSE，以及较高的Spearman’s ρ和AUCROC。平均值显示了预测组合（黑色：EnsFcst）、混合组合（深棕色：EnsBlend）、平均组合（浅棕色：EnsAvg）和各个单独组成部分参数化（浅蓝色；x：3D；+：2D）。混合组合（EnsBlend）在预测农场之间的每周空间模式方面平均表现优于所有单独的组成部分参数化（图3，“每周平均值”）。在预测农场内的时间模式方面（图3，“农场内平均值”），EnsBlend在RMSE方面表现优于所有组成部分，但在AUCROC和Spearman’s ρ方面排名低于两三个参数化方法。平均组合（EnsAvg）平均而言并未优于最佳组成部分参数化。尽管如此，根据评估指标的不同，它的表现优于75-90%的参数化方法。在各个单独的参数化方法中，3D.7平均表现最好（OptParam = 3D.7），其性能在某些农场和某些周内与3D.3相当。在计算AEIP时（方程（1）的s和τ值范围内，每个参数化方法至少对某些农场和某些周是最优的；图A.5，A.6）。所有模型在捕捉农场之间的每周空间模式方面比捕捉农场内的时间模式更为一致（图3，A.8）。模型在预测哪些农场在给定周的海虱负荷较高或较低方面表现得更好，而在预测给定农场内哪些周的海虱负荷较高或较低方面则表现较差。农场级别的性能存在较大差异，有性能较差的农场尾部，而周级别的性能变化较小（图A.8）。一些地区在所有指标上都表现出一致的高性能（图4a-d，A.9，A.10；深蓝色），例如Loch Torridon（1）、Loch Carron（2）、Upper Loch Linnhe（6）和Loch Fyne（8）。一致的差预测，如Muck岛（5）或Loch Spelve（7）附近，仅限于个别农场而不是区域集群。在某些情况下，指标之间存在矛盾；Loch Hourn（4）和Loch Alsh（3）的农场在ROC-AUC方面表现最差，但在RMSE或ρ方面则不明显。下载：下载高分辨率图像（2MB）下载：下载全尺寸图像图4. 集合交叉验证性能。在交叉验证期间，（a）EnsFcst、（b）EnsBlend、（c）EnsAvg和（d）最佳参数化3D.7的农场级别均方根误差（RMSE）。点表示农场级别的RMSE。背景地图显示了基于农场位置的倒数距离加权表面。颜色越浅表示性能越好，颜色越深表示性能越差。（e）交叉验证期间的每周RMSE与（f）季节趋势。编号区域分别为1）Loch Torridon、2）Loch Carron、3）Loch Alsh、4）Loch Hourn、5）Isle of Muck、6）Upper Loch Linnhe、7）Loch Spelve和8）Loch Fyne。每周性能显示出明显的季节性波动，秋季性能下降（图4e-f，A.9，A.10），这种模式在性能指标和模型之间普遍一致。EnsFcst > EnsBlend > OptParam > EnsAvg的排序在年份和性能指标上平均来说是质量上稳定的（图4f，A.9f，A.10f），其中EnsBlend的表现略优于OptParam和EnsAvg，且变异性较大。EnsBlend在模型中对AEIP的对数正态成分中显示了AEIP与预测的每条鱼成年雌性虱子数量之间更强的关系（图5）。在AEIP = 0时，所有贝叶斯模型预测的每条鱼平均虱子数量可以忽略不计（所有平均值< 0.09；EnsBlend：0.0259 [0.0237–0.0281]，平均值[95%置信区间]，大约每条鱼1只虱子）。在非常低的AEIP以上，预测的每条鱼成年雌性虱子数量趋于稳定在一个AEIP的比例。例如，EnsBlend预测的每条鱼平均虱子数量约为AEIP的20%（95%置信区间：AEIP > 2时为15–24%）。下载：下载高分辨率图像（435KB）下载：下载全尺寸图像图5. 贝叶斯模型的对数正态回归后验。贝叶斯模型的对数正态成分预测了每条鱼成年雌性虱子数量的变化与AEIP的关系，假设每条鱼的平均AF虱数量大于0。（a）后验均值回归及95%可信区间（EnsBlend：深棕色；Ens Avg：浅棕色；组成部分：浅蓝色）作为AEIP的函数。（b）预测的每条鱼AF虱数量占AEIP的比例。3.2. 组成部分参数化之间的变异性组成部分之间的协议强烈依赖于位置（图A.13；图1）。例如，在某些农场，如GORS1，预测值从夏季的高IP到完全没有IP不等。2D和3D参数化之间通常有明显区别，2D模拟在近岸位置预测低IP或没有IP，与3D模型形成鲜明对比。从定性上讲，给定农场的IP要么具有强烈的季节性（例如，APT1每年夏季都有IP峰值），要么与典型的生产周期长度一致（例如，AAC3的峰值间隔约为2年）。3.3. 重采样重采样清楚地表明，与选择单一参数化相比，集合建模是提高海虱扩散模型预测能力的可靠方法。无论评估指标或组成部分的数量如何，EnsBlend在预测农场之间的每周海虱感染空间模式方面几乎普遍优于其组成部分（图6）。对于农场内的时间模式，EnsBlend在RMSE方面同样优于所有组成部分。相比之下，在农场内的Spearman’s ρ和AUCROC方面，EnsBlend很少优于最佳组成部分。随着包含的组成部分数量增加，性能倾向于相对于组成部分提高（图7），表明能够识别表现最好的组成部分。下载：下载高分辨率图像（4MB）下载：下载全尺寸图像图6. 重采样期间的集合排名。20个组成部分分别进行了100次重采样，每种集合大小（3、5、10、15；总共400次重采样），并从中构建了集合。点显示了每个模型在每次重采样中的排名（EnsBlend：深棕色；EnsAvg：浅棕色；组成部分：浅蓝色），其中1表示最佳平均性能。下载：下载高分辨率图像（798KB）下载：下载全尺寸图像图7. 重采样期间的集合性能。20个组成部分分别进行了100次重采样，每种集合大小（3、5、10、15；总共400次重采样），并从中构建了集合。对于每次重采样中的每个指标，使用了组成部分技能分数的平均值和标准差来计算每个模型的z分数，其中0表示组成部分的平均技能，单位是以组成部分技能的标准差表示的。点显示了整体分布，下方的范围显示了中位数和中心80%。无论集合大小或评估方法如何，EnsAvg在几乎所有重采样中的表现都优于中位数组成部分（图6），并且偶尔在较小的集合中优于所有组成部分。相对于组成部分的表现对集合大小不那么敏感（图7）。3.4. 组成部分权重EnsBlend主要由少数几个组成部分主导，大多数组成部分的权重非常小，平均而言只有25%的组成部分在任何农场的权重高于0.1（图8）。参数化3D.7在整个范围内占有较大的权重，3D.3在北部和南部具有较高的权重，3D.11在Skye北部、Harris南部和North Uist对集合的影响更大。参数化3D.2和2D.3也在一些分散的农场中有明显的影响。下载：下载高分辨率图像（1MB）下载：下载全尺寸图像图8. 混合集合权重。图中显示了组成部分参数化的后验平均集合权重地图。仅显示了至少在一个农场上pˉ>0.1的组成部分参数化。点表示个别农场的随机效应，而表面显示了空间样条。集合权重（图8）和整体性能（图3）清楚地表明3D模型优于2D模型。表现最好的两个组成部分（3D.7，3D.3）的参数化具有非常慢的游速（0.01–0.03 cm s-1）和相对较低的盐度阈值（21.5–26.1，21.0–24.6 psu），导致更被动的垂直行为。3.5. EnsBlend桡足类分布EnsBlend估计的桡足类感染压力清楚地遵循了主要海洋 currents的长期平均路径（图9a，c，e）。平均每日IP和95%可信区间的平均每日宽度显示出明显的相关性（图9b，d，f），即高平均值对应于更广泛的可信集合值范围。在Skye南部海岸（图9e，f）等一些地区，平均IP较高且集合内部的一致性较强。在具有复杂水动力学的地区或更受保护的地区（图9d），不确定性也更高，其中垂直动力学的参数化差异对预测的桡足类分布影响更大。下载：下载高分辨率图像（2MB）下载：下载全尺寸图像图9. 2021–2024年的集合虱子密度。图9a，c，e显示了整个模拟期间（a）整个域、（c）Loch Linnhe和（e）Skye附近95%可信区间平均宽度范围内的平均每日IP数量和平均每日桡足类数量（m-3·h-1）。集合感染压力也随季节变化（图10）。虽然罕见，但在夏末期间，中度到高IP在农场中更为常见（IPEns>1的站点：7月至10月平均 = 1.28%，12月至5月平均 = 0.193%），而在春季（0更多的农场在夏季和初秋也没有模型化的IP值（IPEns=0的站点占比：7月至10月平均值为5.77%，12月至5月平均值为1.81%）。下载：下载高分辨率图片（680KB）下载：下载全尺寸图片。图10. 各农场间桡足类感染压力的集合平均值。模型化的感染压力具有明显的季节性，在夏季和秋季压力最大。在2024年4月至5月的详细模拟中，不同农场的集合垂直分布表现出多种模式（图11）。在某些情况下，逆向昼夜垂直迁移占主导地位，白天桡足类浓度较高，夜间扩散到更深的层次。其他地方的桡足类分布较为均匀，浓度随潮汐脉动而变化。春潮和退潮的影响也很明显，表现为大约14天的高浓度和低浓度周期，通常在浓度较高的时期会扩散到更深的水域。下载：下载高分辨率图片（3MB）下载：下载全尺寸图片。图11. 集合垂直分布。集合垂直分布是使用EnsBlend权重计算的，针对2024年4月和5月顶部30米内的每小时粒子密度，分辨率为1米。图中展示了五个农场的分布情况，包括明显的日变化、潮汐变化和春潮-退潮周期。农场位置在图1.4中标出。讨论集合模型明显优于单一的海虱扩散模拟。机器学习方法在特定地点的预测中尤为有效。复杂模拟的参数化不可避免地涉及在多个看似合理的值中进行主观选择。正如重采样所证明的，像EnsBlend这样的模型比单一模拟更能可靠地估计桡足类感染压力的空间分布。当没有合适的数据来拟合复杂的统计模型时，一组扩散模型的平均感染压力通常会比其中任何一个模型表现得更好。因此，在缺乏支持单一参数化的有力证据时，选择集合模型是更为谨慎的选择。建模海虱扩散需要指定许多生物学和物理过程及其复杂的相互作用。实证数据至关重要，持续的实验工作对于加深我们对海虱动态的理解仍然十分必要。然而，这些数据最好被视为定义任何给定过程的合理值范围，而不是确定单一的“真实”值。生物学本质上具有变异性，采样也存在不确定性（Moriarty等人，2024；Shoemaker等人，2020）。与其追求单一的最佳参数集，不如使用模型集合来适应生物学变化、随机性和推断不确定性（Peterson和Bode，2021）。即使有丰富的实证数据来指导模型参数化，集合建模也可能产生更稳健和现实的预测。对于区域模型而言，选择单一参数值会带来额外的复杂性，而集合模型可以克服这些问题。标准常数值（如扩散率等）是否合适取决于当地条件（Visser，1997）。苏格兰西部有多种不同的物理环境，这些环境以复杂的方式相互关联，空间加权的集合模型能够优化在不同条件下的参数化。此外，生物参数在不同种群中无疑会有所不同（Shoemaker等人，2020），包括对最终空间分布有强烈影响的过程（Coates等人，2020）。不同的合理参数选择可能会导致完全不同的结果（图A12）。依赖于环境条件的生物参数还会受到水动力模型准确性的影响，而水动力模型本身也会随空间和时间波动（Aleynik等人，2016）。这里介绍的集合方法旨在优化合理参数值的组合，而不是单一的值。如EnsFcst中使用的机器学习算法明确设计用于优化预测能力，模型超参数经过调整以最大化指定的性能指标（Breiman，2001；Pichler和Hartig，2023；Rubbens等人，2023）。许多机器学习方法的一个关键优势是不假设预测变量和响应变量之间有预定义的关系，允许数据驱动模型的功能形式。这自然地包含了非线性效应和预测变量之间的相互作用（Pichler和Hartig，2023）。选择性能指标还允许模型根据用户的具体目标进行定制（Naidu等人，2023）。贝叶斯模型的几个特点使其非常适合海虱建模的背景。它们非常灵活，能够适应观测值和预测变量之间更为多样的联系，以考虑采样过程和数据特征（Gelman等人，2015；Wisniewski等人，2023）。使用潜在变量也为整合不同数据来源提供了自然途径（Hefley和Hooten，2016；Szewczyk等人，2024a，2021）。在EnsBlend中，贝叶斯方法还为每个时间和空间的集合IP提供了一个后验分布，从而更全面地了解了给定数据下的可能结果范围。例如，利用IP的后验分布估计野生鱼类的暴露情况也可以提供不确定性的估计。虽然可以训练机器学习算法来预测IP，但需要更直接的幼虫IP数据。相比之下，贝叶斯模型可以利用现有数据集推断集合IP，并可以在可用时整合直接数据。EnsBlend和EnsFcst是数据驱动的、实用的方法，可以视为“战术性”模型平均的形式，其目标是实现最佳预测，而不是推断哪个组成模型最有可能是“真实的”。尽管EnsBlend是贝叶斯模型，但它与单纯的贝叶斯模型平均不同（Dormann等人，2018）。虽然组成参数化受到海虱生物学知识的限制，但优化过程明确优先考虑了预测与观测结果的一致性，而不是生物过程的准确表示。在这个框架下，如果结果得到改进，模型参数与现实之间的差异是可以接受的。因此，在解释参数化权重时应保持谨慎。性能与具体的组成参数集以及水动力模型、生物物理模型结构、强制数据和验证数据的误差密切相关。尽管如此，这是基于性能的模型选择或统计程序的一般隐含特征（B?hn等人，2022；Kragesteen等人，2023；例如，Samsing等人，2016a；Sandvik等人，2020）。尽管有这些注意事项，但在组成参数化和EnsBlend权重中仍存在一些明显的趋势。3D模型明显优于2D模型，因为这些模型能够捕捉到与该地区峡湾环流相关的强烈垂直剪切动态。游泳速度较慢的参数化往往表现更好，这与苏格兰Loch Linnhe过去的研究结果相反（Moriarty，2023；Moriarty等人，审稿中），在那里，白天浮游幼体的游泳速度为1.25厘米/秒，桡足类的游泳速度为2.14厘米/秒时，模型与2011年的监控笼数据最为吻合。在这些速度下，海虱大约15分钟内可以从20米深处到达水面，在白天形成一致的聚集。然而，更常见的游泳速度较慢（Garnier等人，2024），包括挪威最近的模型，其中使用的游泳速度为0.05厘米/秒（Johnsen等人，2021；Kragesteen等人，2023；Sandvik等人，2020）。这些较慢的游泳速度会产生更加分散的垂直分布，因为模拟的海虱受到垂直电流和扩散的影响更大。尽管如此，EnsBlend中的权重反映了所包含的具体组成模型之间的相对预测性能，而不是提供直接的生物学洞察。海虱幼虫扩散模型的验证受到数据可用性的限制（Rabe等人，2024）。模型输出的最直接验证来自海虱浮游阶段的样本。然而，海虱幼虫很稀有，需要采集大量水样，并在多样的浮游动物中识别出少量的海虱幼虫。由于海虱分布的斑块性、时间上的轻微错位以及需要更多的颗粒，适当对齐模型预测也并非易事（Murray等人，准备中）。新技术使得这些数据集更易于获取（á Norei等人，2016；Adams等人，2021；Murray等人，准备中；Rabe等人，2024；Zhang等人，2024），未来的创新可能使幼虫监测网络能够直接估计感染压力并生成适合模型验证的区域观测数据集。海虱扩散模型的验证通常依赖于间接数据，即鱼体上的后期阶段计数，正如这里所展示的那样。监控笼通过在一个固定的时间段内将已知数量的鱼放入围栏中并精确计数鱼体上的海虱，可以间接估计桡足类的时间积分相对丰度，前提是关于宿主附着率和鱼体上生存及发育率的假设（Rabe等人，2024）。监控笼为跨地点或时间段的模型预测提供了高度标准化的数据。与直接观察幼虫密度不同，鱼体上的海虱计数与模型输出之间存在时间和多个生态过程的分离，包括附着动态和鱼体上的种群动态。尽管在挪威经常使用（Eliasen等人，2021；Sandvik等人，2021），但在苏格兰的应用仍然有限（Moriarty，2023；Moriarty等人，审稿中；Rabe等人，2024）。然而，由于所需的劳动力和成本，幼虫采样和监控笼在时空范围上受到限制。将模型性能外推到采样发生的区域和时间之外的情况应持谨慎态度。相比之下，苏格兰农场经营者每周进行的海虱计数在地理上覆盖范围较广（Moriarty等人，2024）。农场对成年雌性海虱的计数通常用于确定每个颗粒代表的幼虫数量（Johnsen等人，2016；Moriarty等人，2024；Sandvik等人，2020）。公开可用的成年雌性海虱每周计数构成了大量现有数据（Murray等人，2025）。监管变化和水产养殖业的协调提高了农场监测协议的标准化程度，同时强制性的农场间空间分离减少了空间自相关性，这是机会主义收集的数据集（如公民科学或物种分布模型）的常见限制（Araújo和Guisan，2006；Phillips等人，2009）的常见问题。尽管如此，样本可能受到选择性鱼类处理、农场间和方法随时间的变异或不完整计数的影响。此外，由于只有成年雌性的平均数量被公开报告，将计数与模型化的桡足类IP联系起来需要假设中间阶段的发育和存活情况。由于农场海虱计数用于监管目的，因此可能存在利益冲突，从而进一步影响计数结果（Manheim，2023；Murray等人，2025），尽管业界已经努力标准化收集方法，且监测受到政府监督（Scottish Salmon Producers Organization，2015）。尽管如此，这些数据集涵盖了迄今为止超过四年的时间范围内的大约150个农场，提供了大量的观测数据，使得能够检测到稳健的生态信号，尽管可能存在方法学噪声和偏差。这些数据来源在优势和劣势方面互补。农场计数的时间序列将继续增长；正在进行的工作将通过幼虫采样和监控笼提供更多的感染压力估计。这里介绍的建模仅依赖于农场对成年雌性海虱的监测，未来的工作将探索统计模型，将农场计数与监控笼或直接幼虫采样结合起来，允许根据数据来源进行差异加权（Szewczyk等人，2024a）。未来海虱集合模型的发展有几种明确的方向。传播模型参数化的性能可能具有季节性变化，而不仅仅是空间上的变化。例如，诸如扩散性或盐度等物理状态在一年中的某些时期可能比其他时期表现得更好，或者游动速度可能取决于温度（Crosbie等人，2020年；Moison等人，2012年；Ziarek等人，2011年）。AEIP与每条鱼上预测的虱子数量之间的关系（图5）与之前建模的附着率（约为15-20%）一致（Kragesteen等人，2023年）。然而，通过改进附着过程的模型以及研究抗虱管理措施对鱼的影响（Kragesteen等人，2023年；Tucker等人，2000年），可以直接将其纳入集成模型中，以捕捉这些过程及其环境驱动因素的不确定性，并与鱼体上虱子阶段的观测结果建立更现实的联系。本节中的重采样练习说明了集成方法的一般适用性，但性能评估仅限于20种参数化方案。更全面地探索可能的参数空间可能会产生更稳健的集成模型，从而得出关于最优或集成参数化的更可靠结论。全局敏感性分析通过同时改变所有参数来考虑参数空间内的相互作用，还可以通过识别对模型结果影响最大的参数来为实证研究提供额外信息。

4.1 结论
此处评估的集成模型被认为是一种比寻找单一最优参数化方案更为有力的替代方法。最适合的方法将取决于模型目的和可用数据。多个可行模拟的无权重感染压力平均值可能比选择单一参数化方案更能提高预测准确性，并能提供更稳健的感染压力估计，无论是针对特定地点还是整个区域。在有合理验证数据的情况下，采用经过优化的权重设计更复杂的模型可能会提高预测精度。对于诸如鲑鱼养殖场等特定地点的虱虫感染预测，机器学习方法显得最有前景。集成方法具有灵活性，可以根据具体的建模目标进行定制。集成建模为改进海虱感染压力的预测提供了一种有前途的方法。在复杂的海虱感染动态生物物理系统中，更精细的数据驱动方法为整合多种相互作用不确定性提供了框架，最终有助于实现负责任的管理和监管决策，从而提高养殖鱼类的福利并减少对野生鱼类种群的影响。

附录
附录A. 补充方法和结果

CRediT作者贡献声明
Tim M. Szewczyk：撰写——原始草案、可视化、软件开发、方法论、调查、形式分析、数据管理、概念化。
Dmitry Aleynik：撰写——原始草案、软件开发、方法论、资金获取、概念化。
Helena C. Reinardy：撰写——原始草案、资金获取、概念化。
Kim S. Last：撰写——原始草案、资金获取、概念化。
Andrew C. Dale：撰写——原始草案、方法论、资金获取、概念化。