竺山湖是太湖西北部的半封闭湖湾, 南北13 km左右, 东西最长处9.7 km左右, 水深1.9 m左右, 是太湖上游来水的主要湖区^［1］。湖区内渔业资源丰富, 但对该水域的渔业资源状况研究较少。科学的渔业资源评估可以为湖泊的资源利用状况提供数据支撑, 对渔业管理尤为重要^［2-3］。传统的渔业资源评估方法, 如剩余产量、年龄结构等模型均需输入较多的参数, 在实际运用过程中受到诸多限制^［4-6］。因此, 需要一种对参数输入要求简单、数据易获取且能够准确反映资源状况的评估方法。基于长度频率法是在数据缺乏条件下有效评估方法之一^［7］。在较长时间内, 国内外学者大多利用基于长度频率的ELEFAN I方法(electronic length frequency analysis I)^［8］对太湖、长江、南海、太平洋等各大水域^［9-19］鱼类种群参数和资源利用状况进行了研究。2018年, Froese等^［20］提出基于长度频率的贝叶斯生物量估计法(length-based Bayesian biomass, LBB), 该方法仅需要能代表被评估渔业种群的长度频率数据, 数据获取较为简单。利用LBB法, Liang等^［21］和薛国进等^［22］分别提出中国沿海海域14种渔业资源和大长山许氏平鲉(Sebastods schlegelii)处于过度捕捞状态; 王雪辉等^［23］研究发现北部湾二长棘鲷(Parargyrops edita)的开捕体长和渐近体长呈下降趋势, 且当长度频率能够代表资源各个开发阶段的长度组成时, 估算结果能够很好地反映出种群的真实状况。

贝氏䱗(Hemiculter bleekeri)在我国平原地区各大江河、湖泊和池塘广泛分布, 是一种世代交替快、生命周期短的杂食性鱼类^［24］, 在长江下游多个水域中资源量较丰富^［26-27］。以往对贝氏䱗的研究主要集中在繁殖生物学和资源量评估等方面^{［16,25,26,28］}。作为太湖常见小型鱼类, 数量颇多^［29］, 在食物网中处于中间环节, 一方面可利用水生昆虫、浮游动物等多种食物来源; 另一方面又可以作为翘嘴鲌(Culter alburnus)、蒙古鲌(Culter mongolicus)等食鱼性鱼类的饵料生物^［25,30］, 对于维持湖泊生态系统的平衡具有重要意义^［31］, 但对其资源状况目前未见报道。因此, 本研究利用ELEFAN Ⅰ与LBB两种方法对竺山湖贝氏䱗种群参数和资源利用状况进行评估, 以期为太湖的渔业资源科学管理提供科学依据。

竺山湖(119°59′~120°04′E, 31°22′~31°28′N)为太湖西北部湖湾, 在湖区内设置10个采样点逐月进行定点取样(图1)。

图1

图1 　竺山湖采样站点示意图 Fig. 1 Sampling locations of the Zhushan Lake

本研究于2019年10月至2020年9月对竺山湖渔获物按月度进行调查。采样时间为每月下旬, 日落前在采样点放置3条多目复合刺网(网目分别为1.2 cm、2.0 cm、4.0 cm、6.0 cm、8.0 cm、10.0 cm、14.0 cm混合, 长125.0 m、高1.5 m), 次日6:00收网, 收集网内所有渔获物, 并冷冻处理带回实验室, 渔获物分类鉴定参考《太湖鱼类志》^［29］和《江苏鱼类志》^［24］, 随后进行体长(L)和体重(W)的测定, 使用数显游标卡尺测量鱼类体长, 电子天平测量体重, 数据精确到0.01 mm和0.01 g。

采用Microsoft Excel 2019软件拟合贝氏䱗体长与体重的关系式, 其公式为:

式中, L为体长(mm), W为体重(g); a、b为待定参数, a为肥满度, b为幂指数系数, 用以判断贝氏䱗的生长类型。

竺山湖贝氏䱗的体长数据按月以10 mm的间距分组, 一共分为14组, 依次为44.01~54.00 mm, 54.01~64.00 mm……164.01~174.00 mm, 174.01~ 184.00 mm。在FISAT Ⅱ软件中输入贝氏䱗的长度频率数据, 利用ELEFAN I技术^［8］估算其渐近体长(L_∞)和生长系数(K), 渐近体重(W_∞)由公式W_∞=aL_∞^b计算得到。

贝氏䱗理论生长初始年龄(t₀)根据公式^［32］估算。假设其生长遵循von Bertalanffy生长方程^［33］, 因此拟合其体长、体重生长方程公式为:

式中, L_t、W_t分别为t龄鱼的体长(mm)和体重(g); L_∞和W_∞分别为渐近体长(mm)和渐近体重(g); K为生长系数; t₀为理论生长初始年龄, 表示理论上体长和体重等于0时的年龄; b为相关系数。

对生长方程［公式(2)~(4)］中t进行一阶、二阶求导, 具体方法参见《鱼类生态学》^［34］, 可获得其体长、体重的生长速度及加速度曲线方程, 公式分别为:

体长生长速度:

体长生长加速度:

体重生长速度:

体重生长加速度:

体重生长速度达到最大或体重生长加速度为0时的年龄为拐点年龄, 计算公式^［35］为:

采用Pauly^［36］提出的公式估算自然死亡系数(M):

式中, M为自然死亡系数; L_∞为渐近体长; K为生长系数; T为竺山湖采样点年平均温度。

运用体长变换渔获曲线法估算总死亡系数(Z), 公式为^［32］:

式中, Z为总死亡系数, N_t为t龄时贝氏䱗的资源量, N₀为0龄(理论初始生长年龄)时贝氏䱗的资源量, Δt为相应体长组下限生长到上限所需的时间。

捕捞死亡系数(F)根据总死亡系数(Z)和自然死亡系数(M)计算, 公式为:

采用LBB法分析贝氏䱗的长度频率数据, 假设其生长遵循von Bertalanffy生长方程^［33］, 所有参数均采用马尔可夫链-蒙特卡罗方法^［37］计算, 公式如下:

式中, N_L是长度为L的鱼类存活量; ${N_{{L_{{\rm{start}}}}}}$是长度为${L_{{\rm{start}}}}$完全选择的数量, 即鱼类在该长度进入渔具的所有个体均被保留; Z/K是总死亡系数与生长系数的比值。由于长度频率数据不包含绝对丰度的信息, 所以当公式(13)的两边除以它们各自的和时, 公式始终保持不变。这样就可以从右侧分母的和中提取出常数${N_{{L_{{\rm{start}}}}}}$, 然后与分子中的${N_{{L_{{\rm{start}}}}}}$相抵消, 剩下的两个要确定的参数是Z/K和L_∞。

值得注意的是, 在未进行捕捞的自然条件下, Z/K=M/K, ${L_{{\rm{start}}}}$为零, ${N_L}_{_{{\rm{start}}}}$则可以设置为1。因此, 公式(13)可以简化为:

式中, ${P_{L/{L_\infty }}}$是鱼类存活到长度L/L_∞时的概率, 而L/L_∞是M/K的函数, 独立于M、K和L_∞。对于充分开发的种群来说, 达到完全选择体长${L_{{\rm{start}}}}$以外的长度的概率是Z/K的函数。

除了公式(14)中的参数外, 受部分选择影响的捕获数量由公式(16)中的逻辑函数给出的网具选择性函数计算。LBB可以容纳不同类型的选择函数, 因此本研究做了近似于刺网选择的高斯选择^［38］。

式中, S_L是指在长度L处被网具保留的个体, L_c为鱼类与网具相遇时50%被捕获体长, 即开捕体长。α表示网具的陡峭度^［38-39］。

对公式(13)进行拟合以估计L_∞、L_c、α、M/K和F/K:

式中, ${N_{{L_i}}}$是体长组L_i中的个体数, ${N_{{L_{i - 1}}}}$是前一个体长组的个体数。将公式(17)两边同时除以各自的和, 从右侧分母中的和中提取常数F, 得到:

将公式(18)进行拟合到长度频率数据, 得出M/K和F/K的估计值, 可以得到F/M=(F/K)/(M/K)。

LBB估算在Bayesian Gibbs sampler软件JAGS^［40］中实施, 并使用统计语言R执行。将观察到的长度${p_{Li}}$与其期望值$\hat p{\;_{{L_i}}}$进行拟合, 因此模型预测的长度分布${\hat p_{{L_i}}}$为:

式中, $\hat N{\;_{Li}}$是可估计种群动态抑制因子L_∞, M/K和F/K ［公式(18)］以及L_c和α ［公式(16)］的函数。

调查期间在竺山湖共采集到贝氏䱗805尾(表1), 其中9月份的数量最多, 为149尾, 体重为2.17 kg。贝氏䱗体长范围为44.17~181.87 mm, 平均值为(94.24±11.74) mm, 优势体长组在84.01~ 94.00 mm之间, 约占调查总数量的37.76%, 各体长组分布见图2。体重范围为0.88~70.30 g, 平均值为(11.35±5.09) g, 其中15.00 g以下个体约占83.25%, 40.00 g以上的个体仅占0.12%, 各体重组分布见图3。

表1

表1 　竺山湖贝氏䱗的群体结构 Tab. 1 Population structure of Hemiculter bleekeri in Zhushan Lake 采样时间sampling time 尾数number 体长/mm body length 体重/g body weight 范围 range 均值±标准差 $\overline x \pm {\rm{SD}}$ 优势组 dominant 范围 range 均值±标准差 $\overline x \pm {\rm{SD}}$ 优势组 dominant 1月Jan. 48 78.43−122.08 98.89±8.46 84.01−94.00 5.49−22.35 10.14±3.34 0−10.00 2月Feb. 106 73.15−126.90 97.09±11.48 84.01−94.00 4.62−27.38 12.37±5.14 0−10.00 3月Mar. 4 92.46−181.87 117.89±42.76 94.01−104.00 10.40−70.30 26.38±29.30 10.01−20.00 4月Apr. 78 68.07−127.26 90.27±11.30 84.01−94.00 2.10−25.70 11.00±4.73 0−10.00 5月May 52 69.97−106.78 84.34±8.53 74.01−84.00 3.34−19.85 8.53±3.78 0−10.00 6月Jun. 68 50.49−119.37 90.14±11.03 84.01−94.00 1.40−26.04 10.96±4.56 0−10.00 7月Jul. 5 44.17−88.99 73.25±17.00 74.01−84.00 0.88−9.10 5.45±3.00 0−10.00 8月Aug. 36 73.46−122.30 86.59±11.41 74.01−84.00 7.79−22.80 8.57±4.38 0−10.00 9月Sep. 149 85.63−134.46 103.51±9.32 94.01−104.00 7.79−30.18 14.58±4.48 10.01−20.00 10月Oct. 45 75.95−120.70 92.03±7.98 84.01−94.00 6.74−22.92 10.57±3.16 0−10.00 11月Nov. 105 75.88−135.57 94.56±8.17 84.01−94.00 5.20−28.29 11.08±3.17 10.01−20.00 12月Dec. 109 73.17−134.04 92.71±9.58 84.01−94.00 4.30−26.10 10.52±3.87 0−10.00

表1 　竺山湖贝氏䱗的群体结构 Tab. 1 Population structure of Hemiculter bleekeri in Zhushan Lake

图2

图2 　竺山湖贝氏䱗体长分布 Fig. 2 Distributions of the body length of Hemiculter bleekeri in Zhushan Lake

图3

图3 　竺山湖贝氏䱗体重分布 Fig. 3 Distributions of the body weight of Hemiculter bleekeri in Zhushan Lake

根据竺山湖2019年10月至2020年9月的现场调查数据, 贝氏䱗的体长(L)和体重(W)数据拟合得到的关系如图4所示, 关系式为:

W=9.0×10⁻⁵L^3.0707 (R²=0.8775, n=805)

竺山湖贝氏䱗生长为匀速生长, 故本研究运用von Bertalanffy生长方程描述其生长规律, 体长和体重的生长曲线如图5。运用ELEFAN I技术估算其渐近体长L_∞=194.25 mm, 生长系数K= 0.36。根据体长-体重关系式计算出W_∞=95.74 g。再根据Pauly提出的经验公式可估算出理论生长的初始年龄t₀=−0.46 a。因此贝氏䱗的体长和体重生长方程分别为:

图4

图4 　竺山湖贝氏䱗体长-体重关系式 Fig. 4 Length-weight relationship of Hemiculter bleekeri in Zhushan Lake

图5

图5 　竺山湖贝氏䱗生长曲线 Fig. 5 Von Bertalanffy growth curve of Hemiculter bleekeri in Zhushan Lake

L _t =194.25×［1−e^{−0.36(t+0.46)}］

W _t =95.74×［1−e^{−0.36(t+0.46)}］^3.0707

鱼类体长和体重的生长速度和生长加速度用于描述其生长过程变化特征, 分别对体长、体重生长方程进行一次、二次求导得出生长速度和生长加速度方程分别为:

体长生长速度: dL_t/dt=69.93e^{−0.36(t+0.46)}

体长生长加速度: d²L_t/dt²=−25.17e^{−0.36(t+0.46)}

体长生长速度及加速度分析表明(图6), 贝氏䱗的体长生长速度与年龄呈负相关关系, 不具拐点, 即贝氏䱗体长生长速度随着年龄的增加而下降, 其先以较快的趋势降低, 随后降低的速度减缓。体长生长加速度始终为负值, 但随年龄的增长其值一直增加, 最后趋于平缓。

图6

图6 　竺山湖贝氏䱗体长生长速度和加速度变化曲线 Fig. 6 Change of growth rate and acceleration rate in the body length of Hemiculter bleekeri in Zhushan Lake

体重生长速度和加速度曲线方程分别为:

体重生长速度:

dW_t/dt=103.40e^{−0.36(t+0.46)}×［1− e^{−0.36(t+0.46)}］²

体重生长加速度:

d²W_t/dt²=37.22e^{−0.36(t+0.46)}×［1−e^{−0.36(t+0.46)}］×［3e^{−0.36(t+0.46)}−1］

体重生长速度及加速度分析显示(图7), 体重生长速度存在拐点, 拐点年龄t_tp=2.85 a, 即此时贝氏䱗体重生长速度达到最大值, 对应的拐点体长L_t=131.07 mm, 拐点体重W_t=28.61 g。贝氏䱗在2.85 a之前体重生长速度随年龄的增加而增大, 生长加速度为正值; 随后体重生长速度随年龄增长而减小, 生长加速度为负值, 并且随年龄增加逐渐趋近于0。 

图7

图7 　竺山湖贝氏䱗体重生长速度和加速度变化曲线 Fig. 7 Change of growth rate and acceleration rate in the body weight of Hemiculter bleekeri in Zhushan Lake

竺山湖贝氏䱗的体长变换渔获物曲线以“未达到完全补充年龄段和体长接近L_∞的年龄段不能用作回归”为原则, 选择8个数据点进行线性回归分析, 拟合其直线方程为(图8):

ln(N/t)=-2.44t+10.82 (R²=0.8810)

式中, N为种群数量, t为相对年龄。

该拟合的直线方程斜率为−2.44, 其95%的置信区间为−3.33~−1.54, 因此竺山湖贝氏䱗的总死亡系数Z=2.44, 95%的置信区间为1.54~3.33。根据贝氏䱗渐近体长L_∞=194.25 mm, K=0.36和T=16.99 ℃计算得到贝氏䱗的自然死亡系数M= 0.82, 捕捞死亡系数F=Z−M=1.62, 当前开发率E=F/Z=0.66。

图8

图8 　体长变换渔获曲线估算竺山湖贝氏䱗总死亡系数(Z) Fig. 8 Estimation of total mortality parameter (Z) from length-converted catch for Hemiculter bleekeri in Zhushan Lake

基于LBB估算的2019年10月至2020年9月的竺山湖贝氏䱗种群参数中渐近体长L_∞为192.00 mm; 自然死亡系数M与鱼类生长系数K之间的比值M/K=1.88; 捕捞死亡系数F与鱼类生长系数K之间的比值F/K=5.78; 总死亡系数Z与鱼类生长系数K之间的比值Z/K=7.64; 捕捞死亡系数F与自然死亡系数M之比F/M=3.07; 开发率E=0.76 (表2)。

表2

表2 　竺山湖贝氏䱗种群参数估算结果比较 Tab. 2 Comparison of estimated population parameters of Hemiculter bleekeri in Zhushan Lake 方法 method L ∞ /mm M/K F/K Z/K F/M E ELEFAN I 194.25 2.28 4.50 6.78 1.98 0.66 LBB 192.00 1.88 5.78 7.64 3.07 0.76 均值 mean value 193.13 2.08 5.14 7.21 2.53 0.71 注: L∞为渐近体长; K为生长系数; M为自然死亡系数; F为捕捞死亡系数; Z为总死亡系数; E为开发率. Note: L∞ denotes asymptotic length; K denotes growth factor; M denotes natural mortality; F denotes fishing mortality; Z denotes total mortality; E denotes exploitation rate.

表2 　竺山湖贝氏䱗种群参数估算结果比较 Tab. 2 Comparison of estimated population parameters of Hemiculter bleekeri in Zhushan Lake

渔业资源数量变动和利用状况的研究基础在于探究鱼类个体的生长规律, 其中体长和体重为鱼类个体生长的重要指标, 对渔业资源的保护与利用具有重要意义^［35,41］。本研究中贝氏䱗理论生长年龄为负值(-0.46 a), 表明是一种快速生长型的鱼类, 寿命较短^［42］, 体重生长具有拐点, 拐点年龄为2.85 a, 符合小型鱼类的生长特征^［43］。体长-体重关系式与其他水域^［28,44］研究相比, a值较大, 表明竺山湖贝氏䱗群体鱼体肥大^［34］, 这可能与竺山湖水域的地理环境和食物饵料等因素有关^［45-46］。大部分鱼类生长呈现匀速生长类型, 仅少部分鱼类为异速生长类型, b值大多集中在2.5~ 4.0之间^［34］。余文娟^［16］发现长江何王庙故道的贝氏䱗群体生长存在性别差异, 雄性为匀速生长, 而雌性和性别未辨个体的b值显著小于3, 为负异速生长; 刘飞^［28］研究表明赤水河贝氏䱗群体为正异速生长。本研究中贝氏䱗b值接近于3, 表明其在生长过程中具有体形不变、比重不变的特征, 属于匀速生长。因此可见同种鱼类在不同水域中的生长存在一定的差异, 并且同一种群的不同性别也会存在生长差异^［47］, 分析可能与采样时间、所用网具、水体温度和不同水域生态环境等因素有关。

本研究使用多目复合刺网捕获, 且采样周期为12个月, 具有相对较丰富的基础数据, 并利用长度频率数据来估算鱼类种群参数, 避免了年龄鉴定的主观误差, 使得分析结果的精确度更高^{［12,48-49］}。渔业资源评估的准确性又与鱼类生长系数K值、自然死亡系数M值等种群参数密切相关。本研究对参数的估算采用ELEFAN I方法和LBB比较分析, 操作较为简便, 尤其适用于生长迅速、短寿命的小型鱼类^［9］。

本研究中两种方法估算的贝氏䱗渐近体长相差不大, ELEFAN I技术估算的生长速率K=0.36, 表明其属于快速生长的鱼类, 与FishBase数据库和其他水域贝氏䱗的研究结果相似^［16,28］, 并且e^−K=0.70<1, 符合理论要求^［49-50］。一般鱼类的K值与M值存在正比例关系, 且受鱼类生长状况、栖息环境等多种因素的影响^［51］, M/K值范围在1.5~2.5之间时M值较为合理^［35］。由表2可看出, 基于ELEFAN I技术估算的M/K值(2.28)高于LBB估算的结果(1.88), 而Z/K、F/K、和F/M值均略低于LBB的估算值。相似的现象在王雪辉等^［23］利用LBB估算北部湾二长棘鲷参数的研究中也有报道, 其研究发现M/K、Z/K和F/K值等参数与其他独立评估方法评估的值也存在一定的差异。推测存在偏差的原因可能与长度频率数据不能代表资源开发范围的组成有关^［20,52］。

长期以来, 太湖渔业捕捞投入使用的网具主要有高踏网、飞机网、网簖和丝网等, 渔业资源承受着较大的捕捞压力^［53］。研究表明，在高强度的捕捞条件下，鱼类容易出现资源衰退迹象^［54］, 渔业资源整体呈现出低龄化、小型化现象^［55-56］。本研究中贝氏䱗生长拐点年龄为2.85 a, 文献记载其产卵群体以1~2 a为主^［16］, 表明贝氏䱗在产卵之后仍能继续生长; 值得注意的是, 贝氏䱗繁殖期集中在4—7月^［16,24］, 而2020年以前太湖在9月开捕, 开捕时贝氏䱗已产卵, 保证了群体资源的补充。因此，推测竺山湖水域中高强度捕捞可能使得繁殖群体中低龄鱼比例较高, 进而会导致贝氏䱗资源量较大但规格普遍偏小, 消耗水域内大量的浮游动物资源。优化鱼类群落结构、可持续利用小型鱼类资源是湖泊渔业管理的一个重要目标^［30］。因此, 在太湖退捕政策实施后, 应通过增加食鱼性鱼类(如翘嘴鲌、蒙古鲌等)来对其资源量进行控制, 减小其对浮游动物的牧食, 促进上行效应^［57］, 改善湖区的渔业生态环境、提高水产品的经济效益。

本研究利用ELEFAN I和LBB两种方法, 对竺山湖贝氏䱗的体长-体重特性、生长-死亡参数以及利用状况进行了分析。贝氏䱗体长-体重关系式为W=9.0×10⁻⁵L^3.070, 呈匀速生长; 体长和体重方程分别为L_t=194.25×［1−e^{−0.36(t+0.46)}］、W_t=95.74× ［1−e^{−0.36(t+0.46)}］^3.0707; ELEFAN I和LBB估算的相对总死亡系数分别为6.78和7.84, 种群的实际开发率大于0.5。尽管本研究利用基于体长频率的两种方法对贝氏䱗种群参数和利用现状进行了评估, 但仍存在不足: (1)认为年龄差异是导致体长不同的原因，未考虑相同年龄的个体体长差异; (2)对于体长频率法，样本是否具有代表性至关重要，渔获物样本体长分布频率的偏差可能影响评估结果。因此，后续研究会继续调查、采集、补充样本量，同时结合年龄分析，以期更准确地了解贝氏䱗的生长特性和资源利用状况。