贝类(主要为滤食性贝类，如牡蛎、蛤、扇贝等)是我国海水养殖的重要品种，据统计2020年我国贝类养殖总产量高达1 480.0万t，占海水养殖总产量的69.3%以上(农业农村部渔业渔政管理局等, 2021)。滤食性贝类具有较高的经济价值，然而，近些年由于经济效益的驱使，一些海域养殖密度不断增加，养殖规模持续扩大，这些粗犷式的养殖方式会对养殖生物本身及海洋生态系统带来一系列的负面影响，如死亡率增加、生长速率减慢、浮游植物衰减、浮游植物群落结构及底层环境的改变等(Prins et al, 1997; Velasco et al, 2009)。为更好地保障贝类养殖业的绿色可持续发展，摸清养殖贝类与环境的相互作用关系、明确贝类的养殖容量显得尤为重要。

Inglis等(2000)将贝类养殖容量分为4类，即物理容量、养殖容量、生态容量和社会容量。在生态系统水平上，贝类养殖的容量被定义为不会对生态系统的结构和功能造成显著影响(即不可接受压力)的最大养殖密度(McKindsey et al, 2006)。因此，在评估养殖贝类的生态容量时，需要一种方法从整体水平上评估养殖贝类对生态系统结构、功能以及其他功能群的影响。根据营养动力学原理，Ecopath模型从物质和能量平衡角度出发，充分考虑了物种间的捕食关系、竞争关系与生态转化效率，是评估养殖贝类生态容量的有效工具(Jiang et al, 2005; Byron et al, 2011; 林群等, 2013)。目前，利用Ecopath模型法评估贝类生态容量在世界范围已经得到广泛应用。Byron等(2011)利用Ecopath模型确定了Narragansett湾养殖牡蛎的生态容量为297 t/km²。Gao等(2020)利用Ecopath模型得出桑沟湾长牡蛎(Crassostrea gigas)的生态容量为目前生物量的1.8倍，为976 t/km²。Jiang等(2005)利用Ecopath模型评估了Tasman和Golden湾新西兰翡翠贻贝(Perna canaliculus)的生态容量，得出该海域贻贝的生态容量为65 t/km²，远低于其养殖容量(310 t/km²)，一旦超过此生态容量，该海域生态系统结构与功能将发生显著性的改变。虽然已有学者基于Ecopath模型评估了贝类养殖活动对胶州湾海域生态系统的影响(Han et al, 2017、2018; 马孟磊等, 2018)，也有学者基于饵料收支平衡关系的生态模型、Dame指标和Herman模型等方法估算了胶州湾菲律宾蛤仔的养殖容量(即产量容量而非生态容量)(刘学海等, 2015; 董世鹏等, 2020)，但至今未见该海域菲律宾蛤仔生态容量的相关报道。

除了具有重要的经济价值外，滤食性贝类还通过直接摄食水体颗粒物或间接改变营养元素流动，在碳汇扩增、净化水质、缓解富营养化等方面发挥着重要的生态服务功能(Dowd, 2003; 方磊等, 2011; Filgueira et al, 2015; 唐启升等, 2016; Alonso et al, 2021)。在近海养殖生态系统中，浮游植物形成的初级生产力非常高，滤食性贝类一方面通过对浮游植物大量摄食，将颗粒有机碳转化为粪便和假粪等形式，加速了有机碳向海底的输送过程；另一方面通过吸收海水无机碳形成碳酸钙贝壳，在实现碳的长期封存中发挥着重要功能(任黎华, 2014)。

菲律宾蛤仔(Ruditapes philippinarum)是一种广温、广盐的海洋双壳贝类，在世界分布广泛，是一种非常重要的海水贝类渔业开发品种。胶州湾位于山东半岛南岸(总面积约为350 km²)，是我国北方重要的菲律宾蛤仔养殖基地。研究发现，菲律宾蛤仔的养殖活动已对胶州湾底栖生态系统造成了一定程度的干扰(麻骜等, 2014; 丁敬坤等, 2020)。本研究以胶州湾菲律宾蛤仔为对象，构建胶州湾海域的Ecopath模型，评估湾内菲律宾蛤仔的生态容量，并分析菲律宾蛤仔的碳收支情况，研究结果可为菲律宾蛤仔的增养殖可持续发展及其生态服务功能的发挥提供科学指导。

1 材料与方法 1.1 Ecopath模型构建

Ecopath模型包括一组线性方程，描述了胶州湾生态系统各功能群的生产量：

$ \begin{gathered} {B_i} \times {(P{\text{/}}B)_i} \times E{E_i} = \hfill \\ \;\;\;\;\;\;\;\;\; \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\sum\limits_{j = 1}^n {{B_j}} \times {(Q{\text{/}}B)_j} \times D{C_{ij}} + {Y_i} + {E_i} + BA \hfill \\ \end{gathered} $

(1)

式中，B_i和B_j分别为被捕食功能群i和捕食功能群j的生物量；(P/B)_i为功能群i的生产量/生物量或者生态系统稳定情况下的总死亡率；EE_i为功能群i的生态转化效率，被定义为生态系统中该功能群的利用比例；(Q/B)_j表示捕食功能群j的食物消耗量/生物量；DC_ij是捕食功能群j胃含物中被捕食者i的比例；Y_i是功能群i的捕捞产量；E_i为净迁移率；BA_i是功能群i的生物量积累率。对于每个功能群i，需要输入DC_ij和4个基本参数(B、P/B、Q/B和EE)中的至少3个以构建模型。

参照Han等(2017、2018)和马孟磊等(2018)建立的模型，采用Ecopath with Ecosim (EwE) 6.5软件构建胶州湾的营养通道模型。胶州湾Ecopath模型由23个功能群组成，包括碎屑、浮游植物、浮游动物、棘皮动物、多毛类、菲律宾蛤仔、头足类、虾类、蟹类和不同类型鱼类等，各功能群B (t/km²)、P/B (/yr)、Q/B (/yr)的值及其组成的食物矩阵参照Han等(2017、2018)和马孟磊等(2018)的研究。各功能群的年捕捞量参照《山东渔业统计年鉴》以及生态特征相似的黄渤海海域的相关报道(山东省海洋与渔业厅, 2019; Gao et al, 2020)。在当前胶州湾Ecopath模型的基础上，不断提高模型中菲律宾蛤仔的生物量，观察生态系统中包括碎屑、浮游植物等在内的其他功能群的变化，当系统中任意其他功能群的生态营养效率大于1时，推算得出胶州湾菲律宾蛤仔的生态容量(Jiang et al, 2005; Byron et al, 2011; 林群等, 2013)。

1.2 Ecosim动态模拟

在Ecopath静态模型的基础上，Ecosim加入时间动态模块，评估不同菲律宾蛤仔养殖密度对胶州湾生态系统造成直接或者间接的潜在影响，基本方程如下：

$ d{B_i}{\text{/}}dt = {g_i}\sum\limits_j {{Q_{ji}}} - \sum\limits_j {{Q_{ji}}} + {I_i} - ({M_i} + {F_i} + {e_i}) \times {B_i} $

(2)

式中，dB_i/dt为功能群i在单位时间内的生物量变化；g_i为净生长效率；Q_ji为消耗率，其中B_i分为易捕食部分和不易捕食部分(Walters et al, 1997)；I_i是迁入率；M_i是自然死亡率；F_i是捕捞死亡率，e_i是迁出率。

Ecosim中每个功能群的易捕食部分(v)为捕食者和被捕食者间的营养控制关系(即上行控制关系，下行控制关系和混合关系)，根据Cheung等(2002)、Christensen等(2004)和Kluger等(2016)的研究做相应调整(表1)：

$v_{i}=0.1515 \times \mathrm{TL}_{i}+0.0485 $

(3)

式中，TL_i是通过Ecopath模型中得到的每个功能群的营养级，v是Ecosim中的易捕食部分，介于0到1之间，其中0代表上行控制关系，1代表下行控制关系。

在计算时真正所需的易捕食部分(v_new)是通过如下方程获得的(Kluger et al, 2016)：

$ \mathrm{ln} (v_{new}) = 2.301 985 × v_{i} + 0.001 051 $

(4)

1.3 不同菲律宾蛤仔生物量密度情景设置

将模型中菲律宾蛤仔的生物量分别设置为当前(情景1)、达到生态容量(情景2)、生态容量10倍(情景3)以及无养殖情况(情景4)，探究不同情景下胶州湾生态系统结构和功能特征的动态变化。生态系统总体特征采用如下指标进行表征：1)系统总流量(total system throughput)，即生态系统能量流动的总和，可反映生态系统的规模大小；2)容量(capacity)，为生态系统总流量和优势度的综合产物，代表了优势度的上限；3)优势度(ascendency)，可综合反映生态系统总流量和平均交互信息；4)平均交互信息(average mutual information)，给出了生态系统内物质交换网络的相关信息；5)循环指数(Finn’s cycling index)，量化了生态系统的再循环利用程度，是衡量系统成熟度和稳定性的重要指标；6)熵(entropy)，给出了生态系统内流量的总数和多样性的相关信息；7)能量转换效率(transfer efficiency)，即输出与被摄食之和与总流量的比值，体现了该营养级被生态系统的利用效率；8) Kempton's Q指数，为衡量生态系统生物量多样性的相对指数。

1.4 胶州湾菲律宾蛤仔个体及种群水平的碳收支

基于吴亚林(2018)有关不同温度下不同规格菲律宾蛤仔的滤水及呼吸情况，牟新悦(2018)有关胶州湾颗粒有机碳的分布情况，郭永禄(2005)有关2龄菲律宾蛤仔的平均湿重情况(约为5.95 g)以及周毅等(2002)有关菲律宾蛤仔贝壳及软体部组织的含碳量，在个体水平(1个生长周期，约730 d)和种群水平(当前生物量和生态容量)上评估菲律宾蛤仔的碳收支情况，以此分析胶州湾菲律宾蛤仔的碳汇能力。

2 结果与分析 2.1 胶州湾菲律宾蛤仔的生态容量分析

利用构建的胶州湾Ecopath模型，评估胶州湾养殖菲律宾蛤仔的生态容量。当菲律宾蛤仔的生物量达到239.9 t/km²时，胶州湾生态系统处于平衡状态，菲律宾蛤仔一旦超过此生物量，碎屑的生态营养效率首先超过1，此时Ecopath模型失去平衡，由此得出胶州湾菲律宾蛤仔的生态容量为239.9 t/km²。

胶州湾菲律宾蛤仔当前、达到生态容量、超过生态容量及无养殖情境下生态系统总体特征参数的变化见表2。由表2可以看出，当胶州湾菲律宾蛤仔生物量达到239.9 t/km² (即生态容量)时，系统总消费量、总呼吸量、总生产量和总生物量有所升高，而总输出量、流向碎屑总量、系统净生产量、总初级生产量/总生物量有所降低，其余生态系统能流指数差别基本不大；当菲律宾蛤仔生物量为此生态容量10倍时，生态系能流特征参数均发生了显著变化，其中，总输出量、流向碎屑总量和系统净生产量甚至出现了负值；而当胶州湾生态系统无菲律宾蛤仔养殖时，流向碎屑总量、总初级生产量/总呼吸量、能量流入碎屑比例等显著提高，浮游植物能流占比显著降低。

从表3可以看出，随着菲律宾蛤仔生物量的增加，胶州湾生态系统有机碎屑转换效率、初级生产者转换效率以及总流量转换效率均明显提高，其中，有机碎屑转换效率提升最大，当达到菲律宾蛤仔生态容量时，营养级Ⅰ~Ⅱ之间的有机碎屑转换效率、初级生产者转换效率和总流量转换效率是未养殖状态下的3.0、1.4和2.1倍。

2.2 不同菲律宾蛤仔生物量密度下胶州湾生态系统结构和功能的变化

基于Ecosim的动态模拟结果可知，随着菲律宾蛤仔生物量密度的增加，胶州湾生态系统的规模不断扩大，当达到生态容量时，系统的总流量、容量和优势度较当前分别增加1 672.7 t/(km²·yr)、1 190.8 lowbits和4 363.0 flowbits；同时在情景2中，生态系统循环指数(12.7%)和平均交互网络信息参数(1.1 bits)较当前分别提高了103.0%和26.8% (图1)。然而，菲律宾蛤仔生物量密度的增加会对胶州湾生态系统多样性造成一定的影响，具体表现为生态系统熵由2.9降低至2.6，Kempton's Q则由0.7降低至0.6。在情景1和2中，系统能量转换效率变化不大，在Ecosim模拟的前2年系统能量转换效率呈波动性下降，随后基本保持稳定。值得注意的是，当菲律宾蛤仔的生物量密度设为生态容量10倍时，胶州湾生态系统失去平衡，除Kempton's Q指数外，其余参数均表现为负值或接近于0；而在无菲律宾蛤仔养殖情况下，生态系统各特征参数随时间未发生明显改变，故均未在图中展示。

菲律宾蛤仔生物量密度的增加对其他功能类群的生物量也会造成一定程度的影响，包括浮游植物、小型底层鱼类、虾类、其他底层软体动物等功能群的生物量出现一定程度的下降，浮游动物等功能群的生物量变化基本不大(图2)。

2.3 胶州湾菲律宾蛤仔个体及种群水平的碳收支

在1个生长周期(约730 d)内，菲律宾蛤仔约摄取3 310.1 mg碳，其中，约有46.2%的碳沉降至海底，约有40.6%的碳用于呼吸和钙化重新释放进入水体，约有13.2%的碳通过收获移出(图3)。

若以胶州湾菲律宾蛤仔当前生物量计算，胶州湾菲律宾蛤仔每年将有0.6万t碳被沉积，有0.3万t碳因收获被移出；而当菲律宾蛤仔生物量达到生态容量时，胶州湾每年将有1.5万t碳被沉积，有0.6万t碳被移出(图4)。

3 讨论

本研究利用Ecopath模型评估了胶州湾菲律宾蛤仔的生态容量，并进一步利用Ecosim模块动态分析了菲律宾蛤仔生物量扩大对胶州湾生态系统的影响。同样在先前的研究中，Jiang等(2005)和Gao等(2020)也基于Ecopath模型法及增加养殖贝类生物量使模型失衡得到生态容量的假设，分别评估了不同养殖海域贻贝、牡蛎等贝类养殖的生态容量。本研究评估出胶州湾菲律宾蛤仔的生态容量为239.9 t/km² (表2)。基于已有文献资料(Han et al, 2017; 马孟磊等, 2018; 丁敬坤等, 2020)，当前胶州湾菲律宾蛤仔的平均生物量密度整体尚未达到其生态容量，但在胶州湾部分养殖区，菲律宾蛤仔生物量可高达3 444.8 t/km² (范颖, 2015)，已远超出了菲律宾蛤仔的生态容量。

由系统总流量、容量、优势度、循环指数、能量转换效率和Kempton's Q等生态系统特征参数的变化可知，菲律宾蛤仔生物量的持续扩大会对胶州湾生态系统造成一定的影响。当胶州湾菲律宾蛤仔生物量从当前增加至生态容量时，生态系统总流量、容量和优势度明显增加(图1)，说明此时胶州湾生态系统的规模及其抵抗外界环境变动的能力增强。根据Odum (1969)和Finn (1976)的报道，一个成熟的生态系统应为系统内再循环比例较高而熵值较低。本研究中，当胶州湾菲律宾蛤仔生物量增加至239.9 t/km²时，生态系统循环指数较当前提高了103.0%，而熵值较当前降低了10.4%，意味着此时胶州湾生态系统具有更高的成熟度与稳定性。然而随着菲律宾蛤仔生物量密度的增加，平均交互网络信息的表征参数有所提高(由0.9提升至1.1 bits)，说明生态系统此时受到更多限制，能量沿着更加固定的路径流动(例如通过菲律宾蛤仔的路径) (Ulanowicz et al, 1997)。与此同时，生态系统的Kempton's Q指数出现一定程度的下降，表明较高的菲律宾蛤仔生物量密度会对生态系统物种多样性产生一定的负面影响。当胶州湾菲律宾蛤仔生物量持续增加至生态容量10倍时，随着模拟时间的延长，生态系统各特征参数接近于0或出现负值(数据未展示)，生态系统发生崩溃，说明系统无法承载如此高密度的菲律宾蛤仔养殖。

浮游植物是滤食性贝类的主要饵料来源，菲律宾蛤仔生物量的增加直接导致了浮游植物类群生物量的减少，并因食物竞争关系间接影响了其他底层软体动物等竞争者种群的生物量(图2)。鉴于上述因菲律宾蛤仔生物量增加带来的直接或间接影响，胶州湾生态系统无法容纳高于菲律宾蛤仔生态容量的生物量。浮游动物是滤食性贝类的主要食物竞争者(Lam-Hoai et al, 2001; Gao et al, 2020)，在本研究中，菲律宾蛤仔生物量的增加对浮游动物生物量的影响较小。在胶州湾生态系统中，虾类是浮游动物的主要捕食者，同时作为菲律宾蛤仔的饵料竞争者，菲律宾蛤仔生物量的增加间接导致了虾类种群生物量的减少，从而可能降低了其对浮游动物的摄食压力，导致浮游动物种群生物量变化较小，这与李颖虹(2010)的实测数据一致。

本研究还发现，在无菲律宾蛤仔养殖活动时，浮游植物能流占比显著降低，说明浮游植物的生长可能因贝类摄食活动的减少而受到一定程度的限制(如周转率降低等)。同时发现，当无菲律宾蛤仔养殖活动时，生态系统总初级生产量/总呼吸量和能量流入碎屑比例显著提升，意味着此时大量初级生产力并没有沿食物链传递而是沉积形成了碎屑；胶州湾生态系统来源于有机碎屑的转换效率，特别是营养级Ⅰ~Ⅱ之间的碎屑转换效率也显著降低(由15.4降低至5.5)，说明此时有机碎屑在营养级之间的传递不通畅，而胶州湾适量的菲律宾蛤仔养殖可有效提升有机碎屑转换效率及生态系统总转换效率，减少生态系统的能量浪费。

滤食性贝类具有较高的滤水、钙化与生物沉积速率，这些生理特征使贝类具有重要的生物碳汇功能(唐启升等, 2013、2016)。尽管有学者指出，滤食性贝类并不是纯粹的固碳者，在表现出生物碳汇的同时，还因呼吸和钙化过程释放碳而表现出部分碳源的特征(Chauvaud et al, 2003)。但越来越多的监测结果显示，养殖贝类呼吸产生的CO₂并未直接从海洋释放进入大气，由贝类呼吸产生的CO₂可被生态系统中浮游植物或大型藻类吸收，一部分碳还可被贝类重新利用形成贝壳(刘启珍等, 2010; 蒋增杰等, 2012; Jiang et al, 2015; Zhao et al, 2019)。在本研究中，胶州湾菲律宾蛤仔通过呼吸和钙化每年约释放0.4和0.2万t碳，通过收获每年约移出0.3万t碳，通过生物沉积每年约有0.6万t碳沉入海底，这里菲律宾蛤仔收获移出碳的计量基本采用与海洋行业标准《养殖大型藻类和双壳贝类碳汇计量方法碳储量变化法》相同的方法(该行业标准未考虑其释放以及沉积的碳)。从整体上来看，相比菲律宾蛤仔释放产生的碳，有更多的碳被移出和沉积，养殖菲律宾蛤仔发挥了重要的渔业碳汇功能。然而受水动力等因素的影响，养殖海域沉积物–水界面发生着剧烈的有机碳矿化作用，至于菲律宾蛤仔生物沉积产生的有机碳中，有多少比例被永久埋藏还需要进一步的量化与分析。