基于GAM 和GWR 模型分析環(huán)境因子對魚類分布的影響-以山東近海多鱗鱚為例

簡盈 ，張?jiān)评?，宋業(yè)暉 ，張崇良，紀(jì)毓鵬，任一平*

(1.中國海洋大學(xué) 水產(chǎn)學(xué)院，山東 青島 266003；2.青島海洋科學(xué)與技術(shù)試點(diǎn)國家實(shí)驗(yàn)室 海洋漁業(yè)科學(xué)與食物產(chǎn)出過程功能實(shí)驗(yàn)室，山東 青島 266237；3.海州灣漁業(yè)生態(tài)系統(tǒng)教育部野外科學(xué)觀測研究站，山東 青島 266003)

1 引言

多鱗鱚（Sillago sihama），隸屬于鱸形目（Perciformes）、鱚科（Sillaginidae）、鱚屬（Sillago），為體長約20 cm 的小型海洋魚類，多分布在河口和近岸水域[1]。多鱗鱚屬于沿海暖水性底層魚類，多棲息于沙泥和沙礫的海域，其肉質(zhì)鮮美、商業(yè)價值高[2]。多鱗鱚是底棲生物食性魚類，主要攝食多毛類和端足類，在生態(tài)系統(tǒng)物質(zhì)循環(huán)和能量流動中有重要意義[3]。目前對多鱗鱚的研究主要集中在基因、攝食、生長和死亡參數(shù)估算、早期發(fā)育形態(tài)等方面[1,3-5]。有關(guān)多鱗鱚棲息地的研究罕見報道。棲息地為魚類的生長提供適宜的生物及非生物環(huán)境，直接或間接影響魚類的行為、豐度、滲透壓調(diào)節(jié)、新陳代謝等[6-7]。此外，棲息地的喪失還可能導(dǎo)致生態(tài)系統(tǒng)的失調(diào)。

廣義可加模型（Generalized Additive Model,GAM）通常用于研究物種分布和環(huán)境因素之間的關(guān)系，如李迎冬等[8]使用GAM 探究不同季節(jié)山東半島南部海域小黃魚（Larimichthys polyactis）分布與環(huán)境因子的關(guān)系；朱文斌等[9]基于GAM 分析了浙江沿岸日本鳀（Engraulis japonicus）幼魚資源量受環(huán)境因子影響發(fā)生的變動。GAM 在分析數(shù)據(jù)之前，通常假定環(huán)境變量對物種分布的影響是恒定的，即假設(shè)變量間具有空間平穩(wěn)性或者同質(zhì)性[10]。然而，在一個大的空間尺度上，空間平穩(wěn)性并不總是存在，特別是在海洋生態(tài)系統(tǒng)，水生生物的高度流動性和生物與環(huán)境變量之間的空間變化關(guān)系，使空間非平穩(wěn)性比空間平穩(wěn)性可能性更高[11]。

地理加權(quán)回歸（Geographically Weighted Regression,GWR）模型本質(zhì)上是改進(jìn)的全局回歸模型，是空間數(shù)據(jù)分析方法之一。根據(jù)采樣點(diǎn)的坐標(biāo)，GWR 模型在每個空間位置擬合出一個局部回歸方程，從而得出各解釋變量與響應(yīng)變量關(guān)系的局部回歸系數(shù)[12]。與GAM 相比，GWR 模型更適合用于研究物種分布與環(huán)境因子間的空間異質(zhì)性，而不是在整個研究區(qū)域內(nèi)只得到一個全局平均回歸參數(shù)。GWR 模型被廣泛應(yīng)用于陸地上相關(guān)的研究，如楊逍遙等[13]使用GWR模型研究西安市共享單車分布情況與各影響因素的關(guān)系；丁亞鵬等[14]基于GWR 模型分析了伊河流域土壤有機(jī)碳空間分布與影響因素間的關(guān)系，并探討影響因素的空間異質(zhì)性；而其在海洋生物上的研究較少。Windle 等[15]最早應(yīng)用GWR 模型探討了海洋物種大西洋鱈魚（Atlantic cod）分布與環(huán)境因子和生物因子的關(guān)系，并討論了GWR 模型在漁業(yè)生態(tài)學(xué)研究上的優(yōu)勢。

本研究利用2016 年秋季在山東近海開展的漁業(yè)資源底拖網(wǎng)調(diào)查數(shù)據(jù)，使用GAM 和GWR 模型探討環(huán)境因子對多鱗鱚分布的影響，旨在了解多鱗鱚數(shù)量空間分布特征，比較GAM 與GWR 模型的分析效果，以期為開展?jié)O業(yè)生物空間分布研究提供新思路、新方法，為多鱗鱚資源的科學(xué)養(yǎng)護(hù)和管理提供基礎(chǔ)。

2 材料與方法

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

多鱗鱚樣品來自2016 年秋季（10 月）在山東近海進(jìn)行的漁業(yè)資源底拖網(wǎng)調(diào)查，調(diào)查海域?yàn)?5°～38.5°N，118°～124°E。將調(diào)查海域分為萊州灣海域（A 區(qū)）、煙威漁場及其鄰近海域（B 區(qū)）、山東近海南部及其鄰近海域（C 區(qū)），3 個海域同時開展調(diào)查（圖1）。采用定點(diǎn)設(shè)站的方法來進(jìn)行漁業(yè)資源調(diào)查，A 區(qū)設(shè)置48 個站位，B 區(qū)設(shè)置51 個站位，C 區(qū)設(shè)置63 個站位。調(diào)查船為220 kW 單拖漁船，拖速為3.0 kn，每站拖曳時間約為1 h，拖網(wǎng)調(diào)查均安排在白天進(jìn)行。拖網(wǎng)網(wǎng)具網(wǎng)口高度約為7.53 m，寬約為15 m，囊網(wǎng)網(wǎng)目大小為17 mm。同時，采用RBR-XR-420 型CTD 采集每個站位的水深（Depth）、底層海水溫度（Sea Bottom Temperature，SBT）、底層海水鹽度（Sea Bottom Salinity，SBS）。本次調(diào)查均按《海洋調(diào)查規(guī)范 第6 部分：海洋生物調(diào)查》（GB/T 12763.6-2007）[16]進(jìn)行樣品鑒定，獲取多鱗鱚的密度和生物量等數(shù)據(jù)。依據(jù)拖網(wǎng)時間1.0 h和拖速3.0 kn 進(jìn)行調(diào)查數(shù)據(jù)標(biāo)準(zhǔn)化處理，獲得各調(diào)查站位多鱗鱚相對資源量（單位：g/h）[17]。

圖1 山東近海漁業(yè)資源底拖網(wǎng)調(diào)查區(qū)域Fig.1 Bottom trawl survey areas for fishery resources in the Shandong coastal waters

2.2 分析方法

為了避免環(huán)境變量之間的共線性，通過方差膨脹因子（Variance Inflation Factor,VIF）函數(shù)檢驗(yàn)環(huán)境變量之間的關(guān)系，篩選出適合加入模型的因子。一般認(rèn)為，當(dāng)VIF＞4 時，變量存在共線性[18]，因此，在建模前去除導(dǎo)致共線性的解釋變量。

本研究應(yīng)用GAM 對多鱗鱚相對資源量和環(huán)境因子進(jìn)行分析[19]。GAM 利用非參數(shù)的方法靈活處理數(shù)據(jù)間的復(fù)雜聯(lián)系，通過樣條平滑函數(shù)反映響應(yīng)變量與解釋變量之間的關(guān)系。GAM 的一般表達(dá)式[20]為

式中，Y是多鱗鱚相對資源量（單位：g/h）；g（Y）為聯(lián)系函數(shù)，本文采用自然對數(shù)；xj為解釋變量，即各站位的環(huán)境因子；α為函數(shù)截距；fi（xj）為平滑函數(shù)；ε表示隨機(jī)誤差。

GWR 模型既能描述響應(yīng)變量與解釋變量之間的關(guān)系，又能反映空間異質(zhì)性。它將數(shù)據(jù)的地理位置嵌入到回歸參數(shù)中，通過每個空間位置分配的權(quán)重不同，從而在每個觀測點(diǎn)都進(jìn)行局部回歸。GWR 模型的一般表達(dá)式[21]為

式中，ui是第i個采樣點(diǎn)的緯度數(shù)據(jù)；vi是第i個采樣點(diǎn)的經(jīng)度數(shù)據(jù)，即（ui,vi）是第i個采樣點(diǎn)的經(jīng)緯度坐標(biāo)；β0(ui,vi)是第i個采樣點(diǎn)回歸常數(shù)；βk是第i個采樣點(diǎn)的第k個待估回歸參數(shù)；xki為第i個點(diǎn)上第k個變量的觀測值；t為解釋變量的個數(shù)；εi是第i個采樣點(diǎn)的隨機(jī)誤差。

GWR 模型的空間權(quán)重函數(shù)決定著局部回歸系數(shù)等參數(shù)值，本文采用Gauss 函數(shù)法來確定權(quán)重函數(shù)和帶寬，其公式[22]為

式中，i為預(yù)測點(diǎn)；j為觀測點(diǎn)；wij為觀測點(diǎn)j在預(yù)測點(diǎn)i的權(quán)重；dij為i、j兩點(diǎn)間的距離；θ為帶寬，是描述權(quán)重與距離之間函數(shù)關(guān)系的非負(fù)數(shù)衰減參數(shù)。帶寬的大小影響著GWR 模型的分析結(jié)果，如局部回歸系數(shù)值等。本研究基于赤池信息準(zhǔn)則（Akaike Information Criterion，AIC）的自動黃金分割搜索方法選取最優(yōu)帶寬，即帶寬最好的模型AIC 值最低。

在模型擬合過程中，依照AIC[23]、決定系數(shù)（R2）[15]和偏差解釋率（Deviance Explained,DE）[15]評價模型效果，根據(jù)逐步回歸的方法依次將因子加入模型中。AIC 值越低，R2和DE 越高，表明模型的擬合效果越好。

AIC 公式為

式中，D為偏差（殘差平方和）；φ為離差參數(shù)（方差）；df為有效自由度。

以上模型構(gòu)建過程均在R4.0.3 軟件中實(shí)現(xiàn)，其中GAM 由mgcv 包構(gòu)建，GWR 模型由spgwr 包構(gòu)建，各回歸系數(shù)空間變化的顯著性檢驗(yàn)使用GWmodel 包實(shí)現(xiàn)。

3 結(jié)果

3.1 模型構(gòu)建

根據(jù)VIF 檢驗(yàn)，水深、底層海水鹽度和底層海水溫度的VIF 值分別為1.57、2.17、1.51，均小于4，故本研究選取水深、底層海水鹽度和底層海水溫度3 個環(huán)境變量參與構(gòu)建模型。

根據(jù)AIC、DE 和R2，采用逐步回歸分析，最終選擇水深（Depth）、底層海水鹽度（SBS）和底層海水溫度（SBT）作為解釋變量，多鱗鱚相對資源量作為響應(yīng)變量共同構(gòu)建最優(yōu)GAM。表達(dá)式為

式中，s()為自然樣條平滑函數(shù)。其迭代過程如表1所示，最優(yōu)GAM 的AIC 值為168.113，偏差解釋率為31.87%，R2值為0.285。

表1 GAM 變量篩選及影響因子的參數(shù)分析Table 1 The variable screening process for GAM and parameters analysis

GWR 模型篩選過程如表2 所示。通過比較各個GWR 模型的AIC 值和決定系數(shù)，Depth+SBT 模型的AIC 值最小，為126.296，R2值最大，為0.449，最優(yōu)帶寬為0.113。最優(yōu)GWR 模型的表達(dá)式為

表2 多鱗鱚空間分布的環(huán)境影響因子的GWR 模型篩選過程Table 2 Forward-selection procedure of GWR model for environmental influencing factors on spatial distribution of Sillago sihama

3.2 多鱗鱚空間分布

山東近海A 區(qū)、B 區(qū)、C 區(qū)共調(diào)查了162 個站位，其中A 區(qū)、B 區(qū)、C 區(qū)分別捕獲到多鱗鱚的站位有41 個、10 個、11 個。多鱗鱚主要分布在山東近岸海域，在122.5°E 以東的遠(yuǎn)岸海域未捕獲到多鱗鱚。大部分站位的多鱗鱚相對資源量在0～50 g/h，在150～215 g/h 的站位有2 個，分布在A 區(qū)，即37.5°～38.5°N，119°～120.5°E 的海域范圍內(nèi)（圖2）。

圖2 山東近海多鱗鱚的空間分布Fig.2 Spatial distribution of Sillago sihama in the Shandong coastal waters

3.3 GAM 和GWR 模型對比

利用AIC 和R2指標(biāo)對模型表現(xiàn)性能進(jìn)行評估。結(jié)果顯示，由相同解釋變量構(gòu)建的GWR 模型與GAM在不同解釋變量組合中模型性能表現(xiàn)不同（表1，表2）。此外，分別對比具有相同解釋變量組合的GWR模型和GAM 的AIC 和R2指標(biāo)發(fā)現(xiàn)，GWR 模型具有較低的AIC 值和較高的R2值?？梢姡珿WR 模型比GAM具有更優(yōu)的擬合效果。

GAM 結(jié)果顯示，在底層海水鹽度、水深和底層海水溫度3 個環(huán)境因子中，水深對多鱗鱚相對資源量的影響最大，其偏差解釋率達(dá)到23.50%（表1）。多鱗鱚相對資源量隨著深度的增加而降低（圖3）。底層海水溫度對多鱗鱚相對資源量的影響小于水深的影響，隨著底層海水溫度的升高，多鱗鱚相對資源量整體呈現(xiàn)先增加后降低的趨勢，在19.8℃左右相對資源量達(dá)到頂峰。而底層海水鹽度的偏差解釋率最低，為11.00%（表1）。隨著底層海水鹽度的升高，多鱗鱚相對資源量逐漸升高（圖3）。

圖3 GAM 中不同環(huán)境因子對山東近海多鱗鱚相對資源量的影響Fig.3 Effects of different environmental variables on relative abundance of Sillago sihama in the Shandong coastal waters in GAM

GWR 模型分析結(jié)果顯示，多鱗鱚的相對資源量與水深、底層海水溫度這兩個環(huán)境因子之間的關(guān)系在空間上是非平穩(wěn)的。從各參數(shù)上看，底層海水溫度和水深對多鱗鱚相對資源量在空間上變化影響顯著（p＜0.01）。截距局部回歸系數(shù)的變化范圍為-13.722～7.805，67.2%為負(fù)影響，32.8%為正影響。底層海水溫度局部回歸系數(shù)的范圍為-0.349～0.755，對多鱗鱚相對資源量主要是正影響，正相關(guān)占68.8%，負(fù)相關(guān)占31.2%。水深局部回歸系數(shù)的范圍是-0.028～0.029，負(fù)相關(guān)占76.6%，正相關(guān)占23.4%，對多鱗鱚相對資源量的負(fù)影響多于正影響（表3）。

表3 最優(yōu)GWR 模型局部參數(shù)估計(jì)匯總統(tǒng)計(jì)Table 3 Summary statistics of the local parameter estimates for the optimal GWR model

根據(jù)GWR 模型R2可以看出最優(yōu)GWR 模型在近岸水域的表現(xiàn)更佳，近岸水域水深和底層海水溫度的決定系數(shù)值達(dá)到0.30～0.40，遠(yuǎn)岸決定系數(shù)值大部分為0.10～0.20，少部分為0～0.10（圖4）?？傮w上，近岸海域的模型擬合度較高，遠(yuǎn)岸模型擬合度較低。

圖4 最優(yōu)GWR 模型決定系數(shù)（R2）的空間分布Fig.4 Spatial distribution of R2 values in the optimal GWR model

為了進(jìn)一步分析GWR 模型中水深和底層海水溫度對多鱗鱚空間分布的影響，分別繪制了水深和底層海水溫度的局部回歸系數(shù)空間分布（圖5）。由圖可以看出，局部回歸系數(shù)隨空間位置改變而改變，同一因子對不同地點(diǎn)多鱗鱚相對資源量的影響程度不同。在調(diào)查區(qū)域中，水深主要與多鱗鱚相對資源量形成負(fù)相關(guān)關(guān)系，底層海水溫度主要與多鱗鱚相對資源量呈正相關(guān)關(guān)系。對比水深和底層海水溫度的局部回歸系數(shù)看出，水深和底層海水溫度局部回歸系數(shù)高值區(qū)均在A 區(qū)，并且底層海水溫度局部影響程度高于水深，最高可達(dá)0.90，而水深局部回歸系數(shù)最高僅為0.04。

圖5 GWR 模型中水深（a）、底層海水溫度（b）局部回歸系數(shù)值的空間分布Fig.5 Spatial distribution of regression coefficient values for depth (a) and sea bottom temperature (b) in the GWR model

4 討論

4.1 多鱗鱚的空間分布

本研究中發(fā)現(xiàn)多鱗鱚聚集區(qū)域雖然存在空間上的變化，但集中分布于近岸海域，這與以往的研究相一致[24-25]。值得注意的是，多鱗鱚主要分布在35°～38.5°N，122.5°E 以西的近岸海域，這可能與多鱗鱚的棲息習(xí)性有關(guān)。有調(diào)查研究發(fā)現(xiàn)，山東海域近岸底質(zhì)類型主要為砂泥和砂礫，該特點(diǎn)與多鱗鱚喜棲息底質(zhì)類型相一致[2,26-27]。此外，山東近岸海域具備水深較淺、營養(yǎng)鹽豐富和餌料充足等優(yōu)勢，可以為多鱗鱚的生長提供有利的條件[26,28]。夏季是多鱗鱚繁殖期，產(chǎn)卵后開始分散索餌，之后隨著溫度降低逐漸向南洄游且遷移至較深水域[25,29]?？梢?，多鱗鱚的洄游習(xí)性對其相對資源量分布的時空變化也有影響。

多鱗鱚相對資源量較大的站位多分布在37.5°～38.5°N，119°～120.5°E 范圍內(nèi)，其原因可能是萊州灣海岸線較長且黃河、小清河等10 余條入海河流帶來了豐富的營養(yǎng)鹽，提高了該海域的初級生產(chǎn)力[30]。此外，多鱗鱚喜食底棲生物，主要攝食多毛類[24]。結(jié)合史書杰[26]對渤海底棲生物豐度和分布的研究，本研究發(fā)現(xiàn)多鱗鱚相對資源量高值區(qū)與底棲生物（其中多毛類豐度最大）資源分布的高值區(qū)重疊?？梢?，萊州灣可以為多鱗鱚提供較佳的攝食環(huán)境。

4.2 環(huán)境因子對多鱗鱚空間分布的影響

GAM 和GWR 模型結(jié)果顯示，山東近海多鱗鱚分布主要受到水深、底層海水溫度這兩個環(huán)境因子的影響。底層海水鹽度與多鱗鱚資源分布的相關(guān)性較低，這可能與多鱗鱚屬于廣鹽性魚類有關(guān)。沿岸、紅樹林、河口、沙洲等鹽度變化較大海區(qū)均有多鱗鱚分布，其原因可能是在高低鹽度變化交界的海域，其適宜鹽度梯度有利于刺激多鱗鱚成魚的性腺成熟[31-32]。

GAM 結(jié)果顯示，水深與山東近海多鱗鱚資源分布的相關(guān)性最高，且隨著水深的增加，多鱗鱚相對資源量呈現(xiàn)降低的趨勢。多鱗鱚群體主要分布在0～60 m海域[4,33]。黃海水深10～30 m 處水溫變化較大且其初級生產(chǎn)力高值區(qū)主要集中于中上層水域[8,34]。可見，水深差異會直接影響溫度、海洋初級生產(chǎn)力等其他因素直接或間接對魚類的空間分布造成影響。GWR模型結(jié)果也顯示，多鱗鱚相對資源量與水深呈顯著負(fù)相關(guān)關(guān)系。萊州灣（A 區(qū)）大部分水深在10 m 以內(nèi)，最大水深為18 m，而黃海海域中部（B 區(qū)、C 區(qū)）水深在30～70 m 上下，由西向東逐漸加深[35]?？傊?，萊州灣水深小于黃海中部，更有利于多鱗鱚棲息。

底層海水溫度與秋季山東近海多鱗鱚資源分布的相關(guān)性也較高。GWR 模型結(jié)果顯示，萊州灣多鱗鱚分布與底層海水溫度主要呈現(xiàn)正相關(guān)關(guān)系，而在黃海中部（B 區(qū)、C 區(qū)）主要呈現(xiàn)負(fù)相關(guān)關(guān)系，這表明多鱗鱚分布與底層海水溫度存在空間非平穩(wěn)性關(guān)系。GAM 根據(jù)所有調(diào)查數(shù)據(jù)得出的結(jié)論在某些區(qū)域不適用的原因可能是由于區(qū)域間的底層海水溫度差別較大，也有可能在某些區(qū)域底層海水溫度不是影響多鱗鱚分布的主要影響因子[12]。多鱗鱚屬于暖水性底層魚類，秋季B 區(qū)和C 區(qū)北部受到黃海冷水團(tuán)的影響，使其周圍海域溫度降低，致使多鱗鱚呈零星分布[28,36]。

4.3 模型分析

本研究比較分析了GAM 和GWR 模型在擬合多鱗鱚相對資源量與環(huán)境因子關(guān)系方面的效果。研究表明，GAM 的表現(xiàn)劣于GWR 模型。值得注意的是，兩個模型的最優(yōu)模型的解釋變量組成并不一致，如最優(yōu)GWR 模型中缺少底層海水鹽度這一解釋變量?？梢?，模型的解釋變量組合對GAM 和GWR 模型的擬合效果很重要。GAM 使用平滑函數(shù)來擬合非線性關(guān)系時，容易過度擬合[27]，而GWR 模型是基于線性相關(guān)的，過度擬合風(fēng)險較小。此外，GAM 通常用于分析響應(yīng)變量與解釋變量間的非線性關(guān)系，得出全局統(tǒng)一規(guī)律，在處理復(fù)雜數(shù)據(jù)方面具有靈活方便的特點(diǎn)[12]。正如本研究中GAM 的結(jié)果可以直觀地分析出適宜水深和溫度范圍。但在某些情況下，變量之間的關(guān)系不是通用的，而是取決于位置[37]。故GAM 在處理類似水深等可能存在空間非平穩(wěn)性的變量時，不能全面準(zhǔn)確反映自身的變異規(guī)律和可能的響應(yīng)機(jī)制。

對比本研究中GAM 和GWR 模型結(jié)果，可以看出GWR 模型更有利于反映出空間上變量間復(fù)雜的關(guān)系，考慮空間關(guān)系的局部特征[38]。它能夠根據(jù)不同的帶寬在不同的觀測點(diǎn)納入不同的近鄰點(diǎn)進(jìn)行數(shù)據(jù)分析得出各解釋變量與響應(yīng)變量關(guān)系的局部回歸系數(shù)，側(cè)重于闡釋空間異質(zhì)性及其規(guī)律性特征，這是傳統(tǒng)分析方法難以實(shí)現(xiàn)的[39]。同時，GWR 模型除考慮響應(yīng)變量與各解釋變量在局部范圍內(nèi)的特定關(guān)系外，還能考慮到該關(guān)系在全局尺度上的變化，有利于增強(qiáng)參數(shù)估計(jì)值對空間變化的靈敏性[22]。可見，GWR 模型在處理空間非平穩(wěn)性數(shù)據(jù)時具有巨大的優(yōu)越性。此外，生態(tài)數(shù)據(jù)并不都符合正態(tài)分布和線性關(guān)系的假設(shè)，GWR 模型代表使用相對較新的工具來研究空間非平穩(wěn)性，在漁業(yè)數(shù)據(jù)空間動態(tài)分析中表現(xiàn)出了發(fā)展?jié)摿15]。

綜上所述，本研究比較了GAM 和GWR 模型對山東近海多鱗鱚資源分布的擬合效果，當(dāng)物種分布與環(huán)境變量存在空間非平穩(wěn)性時，建議使用GWR 模型對物種分布進(jìn)行研究，但在實(shí)際應(yīng)用中，還應(yīng)根據(jù)數(shù)據(jù)可獲得性和實(shí)驗(yàn)?zāi)康牟煌x擇合適的模型。此外，本研究中多鱗鱚主要分布于萊州灣，其相對資源量與水深、底層海水溫度呈非線性關(guān)系和存在空間異質(zhì)性。進(jìn)一步的研究，還可以從食物網(wǎng)角度分析不同物種間相互作用對多鱗鱚分布的影響，進(jìn)而根據(jù)山東近海多鱗鱚的分布特征和生物學(xué)特征規(guī)劃保護(hù)區(qū)，這對多鱗鱚的資源管理、保護(hù)以及可持續(xù)利用具有重要意義。