基于Bayes 狀態(tài)空間剩余產(chǎn)量模型JABBA 的南極48 區(qū)犬牙魚(yú)資源評(píng)估

符武月 程 前 劉崇煥 魏藝?guó)Q 劉漢超

（1、大連海洋大學(xué) 漁業(yè)資源，遼寧 大連116023 2、大連市現(xiàn)代海洋牧場(chǎng)研究院，遼寧 大連116023）

南極犬牙魚(yú)主要有小鱗犬牙魚(yú)（Dissostichus Eleginoides）和莫氏犬牙魚(yú)（Dissostichus Mawsoni）兩個(gè)經(jīng)濟(jì)種，主要分布在南極及附近海域。犬牙魚(yú)經(jīng)濟(jì)價(jià)值極高，是南極海域最重要的經(jīng)濟(jì)種類(lèi)。該漁業(yè)開(kāi)始于1980 年代,1990 年代后初迅速擴(kuò)大,并存在大量的非法、未報(bào)和不受管制的(IUU)漁獲量，為此南極海洋生物資源養(yǎng)護(hù)委員會(huì)（CCAMLR）加強(qiáng)了對(duì)南極犬牙魚(yú)的研究與管理，每年對(duì)主產(chǎn)犬牙魚(yú)的南極48、58 和88 區(qū)進(jìn)行年漁獲量限制。南極海域48 區(qū)犬牙魚(yú)產(chǎn)量主要來(lái)自小鱗犬牙魚(yú)，莫氏犬牙魚(yú)產(chǎn)量很少，1998 年以來(lái)該海區(qū)犬牙魚(yú)產(chǎn)量基本穩(wěn)定在2 000-4 000t 左右，CCAMLR 每年對(duì)48 區(qū)犬牙魚(yú)漁獲量均進(jìn)行了限制[1]，而進(jìn)行漁獲量數(shù)量限制的科學(xué)基礎(chǔ)是漁業(yè)資源評(píng)估。

本文將利用CCAMLR 提供的漁業(yè)統(tǒng)計(jì)公報(bào)中48 區(qū)的2004-2018 年歷年漁獲量和捕撈努力量數(shù)據(jù)[2]，采用基于Bayes狀態(tài)空間剩余產(chǎn)量模型JABBA 對(duì)48 區(qū)犬牙魚(yú)進(jìn)行資源評(píng)估，旨在為科學(xué)合理利用該漁業(yè)資源提供科學(xué)依據(jù)。

1 材料與方法

1.1 材料

本文所用數(shù)據(jù)來(lái)自CCAMLR 提供的漁業(yè)統(tǒng)計(jì)公報(bào)數(shù)據(jù)庫(kù)中的48 區(qū)2004-2018 年歷年漁獲量和捕撈努力量數(shù)據(jù)[2]。具體的歷年漁獲量和捕撈努力量數(shù)據(jù)取自按漁具、魚(yú)種和分區(qū)報(bào)表，該區(qū)漁獲量主要來(lái)自小鱗犬牙魚(yú)，而莫氏犬牙魚(yú)產(chǎn)量占極少。捕撈方法主要有4 種漁具，分別是近底層延繩釣、籠壺、底拖網(wǎng)和中層拖網(wǎng)，而該區(qū)產(chǎn)量?jī)H來(lái)自延繩釣和籠壺，籠壺產(chǎn)量占比不足1%，關(guān)于籠壺的捕撈努力量以延繩釣為標(biāo)準(zhǔn)進(jìn)行了標(biāo)準(zhǔn)化處理。歷年漁獲量數(shù)據(jù)以小鱗犬牙魚(yú)和莫氏犬牙魚(yú)兩種產(chǎn)量之和進(jìn)行統(tǒng)計(jì)，歷年總捕撈努力量以各漁具標(biāo)準(zhǔn)化后的數(shù)據(jù)進(jìn)行統(tǒng)計(jì)。歷年漁獲量（t）、捕撈努力量（Thousand hooks）和CPUE 數(shù)據(jù)如表1 所示。

1.2 研究方法

1.2.1 JABBA 程序包模型評(píng)估方法

表1 南極48 區(qū)犬牙魚(yú)歷年產(chǎn)量（t）和CPUE（t / thousand hooks）

本文利用基于Bayes 狀態(tài)空間的剩余產(chǎn)量模型R 程序包JABBA v1.2[3]對(duì)南極48 區(qū)犬牙魚(yú)進(jìn)行資源評(píng)估。JABBA 程序包主要內(nèi)置選項(xiàng)包括：（1）擬合多種CPUE 的時(shí)間序列及其相關(guān)的標(biāo)準(zhǔn)差；（2）估計(jì)模型的過(guò)程方差；（3）估計(jì)可選的單個(gè)或成組CPUE 的時(shí)間序列觀測(cè)方差；（4）通過(guò)設(shè)置BMSY/K 及其形狀參數(shù)m 確定魚(yú)類(lèi)生產(chǎn)量函數(shù)，如Fox, Schaefer or Pella-Tomlinson 模型等。關(guān)于程序包JABBA v1.2 的詳細(xì)描述，包括模型、模型參數(shù)的先驗(yàn)及其誤差診斷工具等見(jiàn)Winker et al.(2018a)[3]。

1.2.2 漁業(yè)資源重建時(shí)間的確定

對(duì)于能夠正常補(bǔ)充的漁業(yè)資源，從當(dāng)前資源量到達(dá)能夠提供MSY 的資源量（Bmsy）所需的平均時(shí)間是B/Bmsy、內(nèi)稟增長(zhǎng)率r和實(shí)際使用的捕撈死亡系數(shù)F 的函數(shù)，資源重建時(shí)間可用如下模型估算[10]：

2 結(jié)果與討論

本文Bayes 模型參數(shù)估計(jì)采用兩條鏈，迭代次數(shù)為10 000次。從參數(shù)初始值開(kāi)始迭代到參數(shù)馬氏鏈的平穩(wěn)分布，通常要經(jīng)過(guò)一定的迭代次數(shù)n，從n 開(kāi)始進(jìn)入平穩(wěn)分布，則1-n 的樣本不能作為該參數(shù)的后驗(yàn)分布樣本，這部分樣本要丟棄且稱(chēng)為Burn in。另外馬氏鏈的樣本之間并非相互獨(dú)立，可能存在大量的自相關(guān)問(wèn)題。所以達(dá)到平穩(wěn)分布后的樣本通常要間隔5-10個(gè)樣本取一個(gè)稱(chēng)為T(mén)hin，使Thin 后的樣本進(jìn)入該參數(shù)的后驗(yàn)分布統(tǒng)計(jì)，以盡可能地消除自相關(guān)問(wèn)題。模型參數(shù)收斂診斷采用可視化的MCMC 迭代蹤跡圖。

模型各參數(shù)MCMC 取樣的軌跡圖（圖1）顯示三條鏈混合得很好，并且可以看出各參數(shù)取樣經(jīng)過(guò)一定迭代次數(shù)以后收斂到的平穩(wěn)分布。 模型的對(duì)數(shù)殘差值（圖2）顯示CPUE 的時(shí)間序列與模型計(jì)算的資源量指數(shù)差異很小，均方根誤差（RMSE）僅為15.8%；另外觀測(cè)的CPUE 時(shí)間序列（圓圈與標(biāo)準(zhǔn)差范圍條）與模型計(jì)算的資源量指數(shù)（實(shí)線(xiàn)）（圖3）擬合的很好。

圖1 模型參數(shù)MCMC 取樣蹤跡圖Fig 1 Trace plots for the model parameter drawn from MCMC samples.

圖2 CPUE 指數(shù)殘差診斷圖Fig 2 Residual diagnostic plots of CPUE indices.

各參數(shù)的后驗(yàn)分布與先驗(yàn)分布比較顯示于圖4，模型各主要參數(shù)與漁業(yè)資源管理重要估計(jì)量匯總于表2 中。內(nèi)稟增長(zhǎng)率r的邊際后驗(yàn)分布的中位數(shù)為0.14，其95%置信區(qū)間的下限與上限分別為0.068 和0.283，與其先驗(yàn)分布密度差異較小；環(huán)境容納量K 的邊際后驗(yàn)分布的中位數(shù)為92 521 t，其95%置信區(qū)間的下限與上限分別為39 693 t 和318 457 t，比其先驗(yàn)分布密度窄；時(shí)間序列的起始年2004 的資源量與環(huán)境容納量的比值B2004/K 的邊際后驗(yàn)分布的中位數(shù)為0.782，其95%置信區(qū)間的下限與上限分別為0.53 和1.005，比其先驗(yàn)分布密度更窄。

最大持續(xù)產(chǎn)量MSY 的邊際后驗(yàn)分布的中位數(shù)為4 247t，其95%置信區(qū)間的下限與上限分別為2 114t 和14 891t，支持獲得MSY 的資源生物量Bmsy的邊際后驗(yàn)分布的中位數(shù)為37 007t，其95%置信區(qū)間的下限與上限分別為15 876t 和127 377t，最優(yōu)捕撈死亡系數(shù)Fmsy中位數(shù)為0.118,而歷年實(shí)際捕撈死亡系數(shù)的平均值與Fmsy比值僅為0.3-0.4，捕撈死亡系數(shù)最大的2008 年該比值也僅為0.5-0.6，這說(shuō)明南極犬牙魚(yú)捕撈強(qiáng)度相對(duì)較小（圖5）；2004-2018 年資源量基本穩(wěn)定在60 000t-72 000t 左右，最高年為2004 年為71 429 t，最低年為2010 年也達(dá)到了61 570t，歷年資源量的平均值與支持獲得MSY 的資源生物量Bmsy之比高達(dá)1.8-1.9，該比值最低的2010 年也達(dá)到了1.678，這說(shuō)明資源量水平較好，需加大捕撈強(qiáng)度以使資源量進(jìn)一步下降到資源量的 最 有 水 平 Bmsy（圖5-6）。 從2004-2018 年 的F/Fmsy、B/Bmsy和 資 源量的軌跡圖可以看出，歷年的捕撈強(qiáng)度呈緩慢下降趨勢(shì)，但均大大小于最優(yōu)捕撈死亡系數(shù)Fmsy，歷年的資源生物量基本穩(wěn)定在60 000 t-70 000 t左右，歷年的資源量也都大大地超過(guò)最優(yōu)資源量Bmsy。

2018 年南極48 區(qū)犬牙魚(yú)資源量為原始資源量的75.2%，據(jù)CCAMLR估計(jì)，2015 年羅斯海的犬牙魚(yú)資源量為原始資源量的75%，這與本次評(píng)估結(jié)果基本吻合。關(guān)于南極48 區(qū)犬牙魚(yú)資源重建，由于2018 年資源量與最優(yōu)資源量之比B/Bmsy已達(dá)1.879，捕撈不足導(dǎo)致資源存在一定浪費(fèi)，需適當(dāng)加大捕撈強(qiáng)度使資源量降到最有水平Bmsy。假如2019 年適當(dāng)加大捕撈死亡系數(shù)F=0.3，相當(dāng)于捕撈努力量釣鉤，按漁業(yè)資源重建時(shí)間模型可得=2.81 年，即若從2019 年開(kāi)始投入90 079 000 釣鉤，2.81 年后的2022 年資源量即可回到最佳資源量Bmsy水平，此后持續(xù)使用Fmsy可連續(xù)獲得最大持續(xù)產(chǎn)量MSY。

從捕撈死亡系數(shù)和資源生物量與其對(duì)應(yīng)的最優(yōu)值之比兩個(gè)指標(biāo)對(duì)南極海域48 區(qū)的南極犬牙魚(yú)資源利用進(jìn)行評(píng)價(jià)（圖7），南極48 區(qū)犬牙魚(yú)的資源利用情況基本呈順時(shí)針的趨勢(shì)，即2004-2005 年資源量和捕撈死亡系數(shù)均稍有下降，2005-2009 年捕撈死亡系數(shù)連續(xù)稍稍增大而資源量連續(xù)稍稍下降，2009-2011年捕撈死亡系數(shù)呈連續(xù)緩慢下降且資源量基本穩(wěn)定，2011-2018年捕撈死亡系數(shù)基本穩(wěn)定且資源量連續(xù)緩慢提升。從資源利用是否捕撈過(guò)度來(lái)看，該資源利用有95.3%的可能不存在捕撈過(guò)度，而存在捕撈過(guò)度概率僅為1.3%，介于二者之間的可能性為3.4%。可以得出該資源利用不存在捕撈過(guò)度情況。

圖3 觀測(cè)CPUE 與模型計(jì)算CPUE 指數(shù)擬合圖Fig 3. Time-series of observed (circle and SE error bars) and predicted (solid line) CPUE.

圖4 模型各參數(shù)的先驗(yàn)與后驗(yàn)分布比較圖Fig 4. Prior and posterior distributions of various model and management parameters.

圖5 歷年F/FMSY 和B/BMSY 變化動(dòng)態(tài)（灰色陰影區(qū)為95%置信區(qū)間）Fig 5 Trends in harvest rate relative to FMSY and biomass relative to BMSY. (Shaded grey area indicates 95%credibility intervals)

總之，2004-2018 年南極48 區(qū)犬牙魚(yú)歷年資源生物量水平總體較高，基本穩(wěn)定在60 000-70 000t，且遠(yuǎn)遠(yuǎn)高于最優(yōu)資源量37 000-40 000t；捕撈死亡系數(shù)相對(duì)較小，基本穩(wěn)在0.04-0.07之間，也遠(yuǎn)遠(yuǎn)低于最優(yōu)捕撈死亡系數(shù)的0.118；2004-2018 年資源利用不存在過(guò)度捕撈的情況。建議2019-2021 適當(dāng)加大捕撈死亡系數(shù)到F=0.3，這樣2022 年南極48 區(qū)犬牙魚(yú)資源量可以回到最優(yōu)資源量Bmsy，為了持續(xù)獲得最大持續(xù)產(chǎn)量MSY，從2022 年開(kāi)始可以一直使用最優(yōu)捕撈死亡系數(shù)Fmsy=0.118，這樣2022 年以后即可持續(xù)獲得最大持續(xù)產(chǎn)量MSY。當(dāng)然由于環(huán)境和魚(yú)類(lèi)生物學(xué)參數(shù)的波動(dòng)，需要每年進(jìn)行漁業(yè)資源評(píng)估以提供更切實(shí)可行的漁業(yè)資源管理方案。

表2 模型各參數(shù)的估計(jì)中位數(shù)及其95%置信區(qū)間Table 2. Summary of posterior quantiles presented in the form of marginal posterior medians and associated the 95%credibility intervals of parameters.

圖6 模型估計(jì)的歷年資源量動(dòng)態(tài)（噸）（灰色陰影區(qū)為95%置信區(qū)間）Fig 6 Retrospective analysis for stock biomass (t).

圖7 歷年資源利用狀態(tài)與評(píng)價(jià)圖Fig 7. Chart for the resource utilization status and evaluation over the years.