基于地統(tǒng)計(jì)學(xué)的南海扁舵鰹時(shí)空分布研究

周星星 ，范江濤，于 杰，徐姍楠，蔡研聰，陳作志,

1. 上海海洋大學(xué) 海洋科學(xué)學(xué)院，上海 201306

2. 中國(guó)水產(chǎn)科學(xué)研究院南海水產(chǎn)研究所/農(nóng)業(yè)農(nóng)村部外海漁業(yè)可持續(xù)利用重點(diǎn)實(shí)驗(yàn)室，廣東 廣州 510300

3. 南方海洋科學(xué)與工程廣東省實(shí)驗(yàn)室 (廣州)，廣東 廣州 511458

扁舵鰹 (Auxis thazard) 隸屬于金槍魚科、舵鰹屬，為暖水性中上層集群洄游的小型金槍魚類[1]，廣泛分布于熱帶和亞熱帶海區(qū)，在中國(guó)的東海和南海均有分布。南海蘊(yùn)藏著豐富的金槍魚資源，最新研究表明，以扁舵鰹為主的南海外海的小型金槍魚類資源量為85萬(wàn)噸，但捕撈量很少，屬于輕度開發(fā)的魚類資源，具有較好的漁業(yè)開發(fā)前景[2]。在當(dāng)前近海漁業(yè)資源逐漸枯竭的情況下[3-5]，對(duì)南海扁舵鰹這類小型金槍魚的開發(fā)是轉(zhuǎn)移近海捕撈壓力的方法之一。

魚類資源的時(shí)空分布具有明顯的空間異質(zhì)性特征，可采用地統(tǒng)計(jì)學(xué)方法進(jìn)行研究。相較于經(jīng)典統(tǒng)計(jì)學(xué)方法，地統(tǒng)計(jì)學(xué)能彌補(bǔ)其隨機(jī)獨(dú)立性假設(shè)的局限性，在基于空間自相關(guān)的計(jì)算建模、估值分析方面具有明顯優(yōu)勢(shì)[6]。1985年地統(tǒng)計(jì)學(xué)方法首次應(yīng)用于漁業(yè)領(lǐng)域，用于估計(jì)生物量，而后在魚類資源豐度和空間異質(zhì)性評(píng)估等方面應(yīng)用較為廣泛[7]。目前國(guó)內(nèi)外學(xué)者已對(duì)南海扁舵鰹的形態(tài)鑒定[8]、生物學(xué)特性[1]、食性分析[9]、遺傳結(jié)構(gòu)和遺傳多樣性評(píng)價(jià)[10]、資源評(píng)估及燈光罩網(wǎng)調(diào)查[11-12]等方面做了較多研究，而有關(guān)其時(shí)空分布的研究卻鮮有報(bào)道，目前僅對(duì)閩中、閩東漁場(chǎng)扁舵鰹的時(shí)空分布與溫鹽關(guān)系進(jìn)行了探討。本文利用2016—2017年的南海燈光罩網(wǎng)調(diào)查資料，采用地統(tǒng)計(jì)學(xué)方法進(jìn)行普通克里金插值，對(duì)南海扁舵鰹資源密度指數(shù)空間異質(zhì)性特征進(jìn)行分析，探索其時(shí)空分布變化規(guī)律，為合理開發(fā)和保護(hù)我國(guó)南海扁舵鰹資源提供科學(xué)依據(jù)。

1 材料與方法

1.1 數(shù)據(jù)來源

漁業(yè)數(shù)據(jù)來自2016—2017年4個(gè)航次的燈光罩網(wǎng)調(diào)查數(shù)據(jù) (圖1)，調(diào)查船總噸位421 t，長(zhǎng)43.6 m，配備金屬鹵化物集魚燈 (1 kW)，所用罩網(wǎng)的主尺寸為281.60 m×81.76 m，網(wǎng)衣最大網(wǎng)目35 mm，網(wǎng)囊最小網(wǎng)目20 mm，沉綱配重2 816 kg。調(diào)查期間每晚19點(diǎn)開燈300盞，開燈后約3 h進(jìn)行放網(wǎng)起網(wǎng)作業(yè)，并進(jìn)行漁獲物統(tǒng)計(jì)，記錄內(nèi)容包括調(diào)查時(shí)間、經(jīng)緯度、漁獲物種類、體質(zhì)量等基本信息。4個(gè)航次的調(diào)查時(shí)間依次為2016年4月和10—11月、2017年4—5月和8—9月。

圖1 南海漁業(yè)資源調(diào)查站點(diǎn)圖Fig. 1 Survey station of fishery resources in South China Sea

對(duì)開燈數(shù)量和捕撈量進(jìn)行歸一化處理，作為計(jì)算單位努力捕撈漁獲量 (CPUE) 的基礎(chǔ)。處理公式為：

式中：i為站號(hào)；Xi-nor為第i個(gè)站點(diǎn)的開燈數(shù)量/捕撈量進(jìn)行歸一化；Xi為第i個(gè)站點(diǎn)的捕撈量/開燈數(shù)量；Xmin為該航次的最小開燈數(shù)量/捕撈量；Xmax為該航次的最大開燈數(shù)量/捕撈量。

CPUE作為表征相對(duì)資源豐度的指標(biāo)[13]。計(jì)算公式為：

1.2 方法

結(jié)合經(jīng)典統(tǒng)計(jì)學(xué)和地統(tǒng)計(jì)學(xué)方法分析南海扁舵鰹分布的空間異質(zhì)性，地統(tǒng)計(jì)學(xué)以區(qū)域化變量理論為基礎(chǔ)，要求數(shù)據(jù)必須滿足正態(tài)分布，對(duì)不符合正態(tài)分布的數(shù)據(jù)可通過對(duì)數(shù)、平方根、反正弦平方根、倒數(shù)和Cox-Box等轉(zhuǎn)換，以達(dá)到地統(tǒng)計(jì)學(xué)分析的要求[14-15]。該方法由變異函數(shù)和克里金插值兩部分組成，主要是對(duì)特定空間內(nèi)的變量進(jìn)行變異性建模，并基于此模型進(jìn)行變異性估值[16]。變異函數(shù)是地統(tǒng)計(jì)學(xué)特有的工具，常用的理論模型有球狀模型、高斯模型和指數(shù)模型。其中，球狀模型的一般公式為：

高斯模型的一般公式為：

指數(shù)模型的一般公式為：

采用Garrison[21]的分布重心法，對(duì)各航次站點(diǎn)的經(jīng)緯度以CPUE為權(quán)重進(jìn)行加權(quán)平均，確定南海扁舵鰹資源密度指數(shù)重心變化的軌跡。公式為：

普通克里金插值法既考慮了空間變量的相關(guān)性，又兼顧有效數(shù)據(jù)的數(shù)量，優(yōu)勢(shì)明顯[16]。故本文采用考慮了變量空間相關(guān)性的普通克里金法進(jìn)行插值, 對(duì)不同年份扁舵鰹CPUE進(jìn)行插值時(shí)，全部選取移除一階趨勢(shì)進(jìn)行處理[23]。

采用W-S檢驗(yàn)進(jìn)行正態(tài)性檢驗(yàn)，GS+9.0軟件進(jìn)行變異函數(shù)的計(jì)算和模型的擬合，調(diào)查站點(diǎn)圖和克里金插值圖由Arcgis 10.2軟件繪制，底圖數(shù)據(jù)來源于全國(guó)地理信息資源目錄服務(wù)系統(tǒng)（http://www.webmap.cn/）。

2 結(jié)果

2.1 數(shù)據(jù)檢驗(yàn)與常規(guī)統(tǒng)計(jì)分析

對(duì)各航次南海扁舵鰹資源CPUE數(shù)據(jù)采用W-S擬合優(yōu)度檢驗(yàn)數(shù)據(jù)的正態(tài)性，發(fā)現(xiàn)數(shù)據(jù)不符合正態(tài)分布，進(jìn)行對(duì)數(shù)轉(zhuǎn)換后數(shù)據(jù)均滿足正態(tài)性要求 (表1)。

表1 各航次數(shù)據(jù)W-S正態(tài)性檢驗(yàn)Table 1 W-S normality test of each voyage

對(duì)4個(gè)航次CPUE數(shù)據(jù)進(jìn)行常規(guī)統(tǒng)計(jì)計(jì)算，得到數(shù)據(jù)分布的各項(xiàng)基本特征參數(shù) (表2)，各航次偏度介于1.218 0～3.235 0，均大于0，分布形態(tài)均為右偏；除第2航次外，其余航次峰度值均大于3，呈尖峰分布，低資源密度指數(shù)的海域所占比重較大，高資源密度指數(shù)海域較少；變異系數(shù)CV介于0.879 7～2.034 2，除了第2航次為中等變異情況(0.1＜CV＜1)，其余均屬于強(qiáng)變異程度 (CV＞1)，漁場(chǎng)資源密度指數(shù)差異較大；另外，從各年際和季節(jié)來看，資源密度指數(shù)也存在差異。年際差異程度依次為2017年 (CV均值1.814 1)＞2016年 (CV均值1.236 5)；季節(jié)差異程度依次為春季 (CV均值1.813 8)＞夏季 (CV均值1.594 0)＞秋季(CV均值0.879 7)。

表2 扁舵鰹調(diào)查數(shù)據(jù)基本統(tǒng)計(jì)參數(shù)Table 2 Basic statistical parameters of survey data of A. thazard

2.2 南海扁舵鰹空間異質(zhì)性分析

基于南海扁舵鰹空間異質(zhì)性特征，建立變異函數(shù)模型獲得相關(guān)理論參數(shù) (表3)。結(jié)果表明，球狀模型在各航次空間異質(zhì)性結(jié)構(gòu)中發(fā)生頻率較高；各航次變異函數(shù)塊金值介于0.000 1～0.035 0，差異較小，各航次CPUE空間異質(zhì)性受隨機(jī)因子的影響程度較弱；基臺(tái)值介于0.250 2～1.230 0，變化趨勢(shì)與塊金值相同，均為先下降后增加；一般采用塊金系數(shù)度量樣本空間相關(guān)性，比值表示隨機(jī)因子引起的空間異質(zhì)性中自相關(guān)部分占系統(tǒng)總變異的比例，各航次塊金系數(shù)介于0.000 1～0.077 0，屬于較強(qiáng)的空間自相關(guān)性，在南海扁舵鰹空間分布中，結(jié)構(gòu)性成分起主要作用；變程與資源密度指數(shù)和分布范圍有關(guān)，反映資源密度指數(shù)大小，密度較高的集群變程一般較小，各航次變程值介于1.073 9～2.410 0，無(wú)明顯變化，差異較小。

表3 各航次扁舵鰹資源變異函數(shù)參數(shù)Table 3 Variation function parameters of A. thazard resources in each voyage

2.3 南海扁舵鰹資源密度指數(shù)的時(shí)空變化

基于各航次扁舵鰹空間異質(zhì)性結(jié)構(gòu)分析扁舵鰹洄游的分布特征 (圖2)，各航次漁場(chǎng)分布存在一定差別，但均呈現(xiàn)出明顯的片狀或斑塊狀特征，且其洄游分布路線大致呈西南—東北走向。第1航次有兩個(gè)明顯的資源密度指數(shù)高值區(qū)，主要分布在10°N—12°N海域；第2航次的資源密度指數(shù)高值區(qū)在10°N—14°N的西部海域，由西向東遞減，有兩個(gè)明顯的資源密度指數(shù)低值區(qū)；第3、第4航次的資源密度指數(shù)高值區(qū)均位于東北部海域，由東北向西南遞減，其中第4航次在西南部海域形成了一個(gè)次高值區(qū)，變化梯度較第3航次更加和緩。

圖2 各航次南海扁舵鰹空間異質(zhì)性結(jié)構(gòu)分布圖注： a—d依次對(duì)應(yīng) 1—4 航次。Fig. 2 Distribution of spatial heterogeneity of A. thazard in each voyageNote: a-d correspond to Voyages 1-4.

各航次南海扁舵鰹CPUE重心均分布在10°N—12°N，多靠近島礁附近。以第1航次扁舵鰹CPUE的分布中心為起點(diǎn) (圖3)，第2航次向西北移動(dòng)至12°N附近，第3航次略向西南移動(dòng)，但仍偏北，第4航次又略向西北移動(dòng)，變化幅度較小。同時(shí)，南海扁舵鰹CPUE重心航次間分布差異不顯著 (P＞0.05)，相對(duì)集中在調(diào)查海域的中部偏東南 (表4)。

圖3 各航次南海扁舵鰹CPUE重心移動(dòng)軌跡Fig. 3 Migration trajectory of center of gravity of CPUE of A. thazard

從各航次南海扁舵鰹CPUE分布重心95%的經(jīng)緯度置信區(qū)間 (Bootstrap法，表4) 可知，不同航次的CPUE重心經(jīng)緯度置信區(qū)間多有重疊，且大部分航次重心均位于置信區(qū)間內(nèi)。

表4 各航次南海扁舵鰹CPUE重心的置信區(qū)間 (95%)Table 4 Confidence interval for center of gravity of A. thazard of each voyage (95%)

3 討論

3.1 空間異質(zhì)性特征

本研究采用地統(tǒng)計(jì)方法分析了各航次南海扁舵鰹空間分布的異質(zhì)性結(jié)構(gòu)，發(fā)現(xiàn)其布局總體以低密度海域?yàn)橹?，并具有聚集性特征。空間變異函數(shù)的種類揭示了魚群依賴于當(dāng)前環(huán)境的聚集分布程度[24-25]，南海扁舵鰹空間異質(zhì)性格局以球狀函數(shù)為主，表現(xiàn)出結(jié)構(gòu)性和聚集性較好的空間格局，集聚程度較為明顯。從變程來看，第4航次變程范圍較大，表明其空間自相關(guān)尺度較廣，其余航次的變程范圍均小于4，說明其分布格局呈聚集分布，這與扁舵鰹為中上層集群洄游的魚類屬性有關(guān)，從空間變異函數(shù)類型和變程范圍兩方面均驗(yàn)證了扁舵鰹的聚集分布。此外，從CPUE均值和基臺(tái)值的關(guān)系來看，在一定條件下，資源密度指數(shù)較高的航次其空間異質(zhì)性也高，該結(jié)論也在許多已有研究中得到了驗(yàn)證[26-27]。

3.2 漁場(chǎng)時(shí)空格局

南海扁舵鰹的資源豐度及其分布受到時(shí)間、空間和海洋環(huán)境等多種因素的影響。從4個(gè)航次看，各季節(jié)資源密度指數(shù)差異較大，依次為夏季＞春季＞秋季，這與一般統(tǒng)計(jì)分析結(jié)論相同[28]。3—7月為南海扁舵鰹的繁殖期，種群規(guī)模增長(zhǎng)迅速，故夏季資源密度指數(shù)最大。此外，扁舵鰹還具有明顯的西南-東北洄游特征。結(jié)合其生物學(xué)特性[1,29]發(fā)現(xiàn)，4—6月為南海扁舵鰹的產(chǎn)卵期，洄游至北部海域產(chǎn)卵繁殖，后進(jìn)行索餌洄游至西南部海域[30]。從海洋環(huán)境看，夏季西南季風(fēng)控制南海，季風(fēng)吹動(dòng)表層海水，使得表層洋流由西南向東流動(dòng)，表層豐富的營(yíng)養(yǎng)物質(zhì)和浮游生物隨洋流一起遷移，洋流向大陸架和大陸坡輸送豐富的營(yíng)養(yǎng)物質(zhì)和浮游生物[31]，扁舵鰹向西南方向索餌洄游，資源密度指數(shù)較高。此外，水層深度也會(huì)影響扁舵鰹分布，從圖2知，扁舵鰹多棲息于近岸淺水域[30]，漁場(chǎng)重心靠近島嶼。從海面到海底都有較為充分的陽(yáng)光透射，使浮游生物在光合作用下迅速地進(jìn)行繁殖，給扁舵鰹餌料以豐富的營(yíng)養(yǎng)物質(zhì)，形成高生產(chǎn)力海區(qū)[32-33]。且淺水水域可在風(fēng)浪、潮汐、對(duì)流的影響下，水體混合充分，底層補(bǔ)充到上層，物質(zhì)循環(huán)快，初級(jí)生產(chǎn)力高，整個(gè)水體營(yíng)養(yǎng)好[34]。扁舵鰹多捕獲于近岸淺水海域，不僅與扁舵鰹自身的棲息分布有關(guān)，也與漁業(yè)生產(chǎn)作業(yè)有著直接的關(guān)系[33]。

除自身的生物學(xué)特性外，海洋環(huán)境也會(huì)影響魚類的分布和洄游。由于扁舵鰹具有較強(qiáng)的游泳能力，有能力進(jìn)行長(zhǎng)距離遷移[28]，所以可明顯地觀察到其分布重心在厄爾尼諾-南方濤動(dòng)現(xiàn)象 (El Ni?o-Southern Oscillation, ENSO)、西太平洋暖池、太平洋十年際氣候振動(dòng)等海洋-大氣相互作用系統(tǒng)影響下的變動(dòng)情況[35]。通過對(duì)比各時(shí)期的SSTA3.4指數(shù) (數(shù)據(jù)來自美國(guó)NOAA氣候預(yù)報(bào)中心網(wǎng)站，http://www.cpc.noaa.gov，圖4)，發(fā)現(xiàn)第1航次處于強(qiáng)厄爾尼諾期，第2航次處于弱拉尼娜期，第3、第4航次均處于正常時(shí)期。當(dāng)ENSO事件發(fā)生時(shí)，直接引起熱帶太平洋海域水溫的大規(guī)模變化，對(duì)具有洄游能力的扁舵鰹來說，ENSO引起的水溫結(jié)構(gòu)變化則會(huì)影響其空間分布和洄游。南海扁舵鰹漁場(chǎng)重心在厄爾尼諾期偏東，拉尼娜期偏西，正常時(shí)期位于中間 (圖3)，這與Lehodey等[36]對(duì)中西太平洋鰹魚 (Katsuwonus pelamis) 的研究結(jié)果相似，當(dāng)厄爾尼諾現(xiàn)象發(fā)生時(shí)，鰹魚群體整體向東遷移約4 000 km；當(dāng)拉尼娜現(xiàn)象發(fā)生時(shí)，則反向遷移4 000 km。但扁舵鰹的遷移范圍沒有鰹魚那么大，推測(cè)是因?yàn)轹烎~具有高度洄游能力，扁舵鰹的洄游能力不如鰹魚強(qiáng)。此外，推測(cè)ENSO事件的強(qiáng)度與持續(xù)時(shí)間也會(huì)影響扁舵鰹群體的遷移范圍。

圖4 各航次Nino 3.4 分析Fig. 4 Nino 3.4 index analysis of each voyage

3.3 插值方法選擇

空間插值是基于已知數(shù)據(jù)點(diǎn)對(duì)研究區(qū)域進(jìn)行預(yù)測(cè)的方法，在漁業(yè)資源領(lǐng)域應(yīng)用廣泛[37-38]。通?？杀环譃榇_定性插值和不確定性插值 (地統(tǒng)計(jì)插值)，確定性插值 (如反距離權(quán)重法、自然鄰域法、趨勢(shì)面法和樣條函數(shù)法等) 是以樣本點(diǎn)和插值點(diǎn)之間的距離為權(quán)重進(jìn)行加權(quán)平均，地統(tǒng)計(jì)方法(如克里金插值法) 是以包含空間自相關(guān)的統(tǒng)計(jì)學(xué)方法為基礎(chǔ)，不僅可預(yù)測(cè)表面，還可對(duì)預(yù)測(cè)結(jié)果的準(zhǔn)確性提供某種度量。

不同插值方法對(duì)前提假設(shè)和原始數(shù)據(jù)的要求不同，常用的反距離權(quán)重法適用于均勻分布的點(diǎn)[39]，其優(yōu)點(diǎn)是直觀高效，發(fā)生各向異性時(shí)，會(huì)考慮方向權(quán)重，但研究表明反距離權(quán)重法只考慮插值點(diǎn)與樣本點(diǎn)之間的距離，依賴于前人的經(jīng)驗(yàn)，易受觀測(cè)點(diǎn)數(shù)據(jù)集的影響[40]；普通克里金插值法要求數(shù)據(jù)呈正態(tài)分布，其優(yōu)點(diǎn)為從變量自身特點(diǎn)出發(fā)，考慮觀測(cè)點(diǎn)的整體空間分布情況，可對(duì)插值誤差做出理論估計(jì)，且能給出估計(jì)精度，具有平衡性，結(jié)果更精確，更符合實(shí)際[37]，但區(qū)域化變量離散性太強(qiáng)時(shí)，估值不夠精確，且只限于單變量在空間分布的特征研究[18]。本研究采用的普通克里金插值方法要求區(qū)域變量滿足二階平穩(wěn)或本征假設(shè)，既考慮了平穩(wěn)范圍的大小，又兼顧有效數(shù)據(jù)量，是一種折中方案，與其他插值方法相比，普通克里金插值法適用范圍廣、計(jì)算簡(jiǎn)單，更加符合實(shí)際情況。此外，在進(jìn)行空間插值前，需要對(duì)漁業(yè)數(shù)據(jù)進(jìn)行檢查，降低異常情況的影響，其方法有很多，如劃分小漁區(qū)等，通常越精細(xì)的空間尺度越能體現(xiàn)原始數(shù)據(jù)的空間特點(diǎn)[41]；還可以利用3倍平均值±標(biāo)準(zhǔn)差、箱圖等方法尋找和剔除異常值[42]。

3.4 展望

本研究由于原始調(diào)查站點(diǎn)密度較小，從一些周圍的點(diǎn)到極值點(diǎn)過程中l(wèi)n(CPUE) 漸進(jìn)不顯著，導(dǎo)致內(nèi)插結(jié)果受極值點(diǎn)的影響明顯。今后，在調(diào)查站點(diǎn)設(shè)計(jì)環(huán)節(jié)應(yīng)充分考慮野外實(shí)際情況；考慮積累更長(zhǎng)時(shí)間序列數(shù)據(jù)，適當(dāng)擴(kuò)大采樣范圍和采樣密度，在扁舵鰹棲息范圍內(nèi)多設(shè)置站點(diǎn)，使采樣數(shù)據(jù)更具代表性。