胡 超 于 靜

（1.湖北省林業(yè)局林木種苗管理總站，湖北 武漢 430079； 2.嶺南生態(tài)文旅股份有限公司，湖北 武漢 430062）

林木良種是有適宜生態(tài)區(qū)域要求的，如果自然條件不適宜，再好的良種也達(dá)不到豐產(chǎn)、穩(wěn)產(chǎn)。因林木良種不適應(yīng)引種區(qū)自然條件而造成巨大損失的教訓(xùn)是深刻的：20 世紀(jì)70 年代，各地在油茶Camellia oleifera生產(chǎn)發(fā)展過(guò)程中調(diào)購(gòu)種子比較隨意，較多地方因?yàn)槌乩韰^(qū)域引種造林，引種前沒(méi)有進(jìn)行科學(xué)預(yù)判，盲目性的引種，導(dǎo)致幼林生長(zhǎng)不良、成林產(chǎn)量很低，在人力、物力等方面都造成了較大的損失[1]。

傳統(tǒng)的林木良種引種適宜生態(tài)區(qū)憑主觀經(jīng)驗(yàn)判斷較多，如，宜林范圍內(nèi)每個(gè)按水平分布的氣候帶和垂直氣候帶都分布著特有類(lèi)型的森林植被。經(jīng)緯度由北向南，由西向東調(diào)運(yùn)范圍大于相反方向的范圍，海拔不超過(guò)300～500 m，但是，1958年，湖北引種廣東、福建馬尾松Pinus massoniana種子成功，用事實(shí)改變了過(guò)去專(zhuān)家認(rèn)為“馬尾松南種北移的幅度不能超過(guò)2°～3°”的定論[2]。1979年李傳志[2]論證馬尾松一次北移6°～7°育苗是可以成功的。所以，溫度、降水、土壤等主要環(huán)境因子相似，即為林木良種同一適宜引種生態(tài)區(qū)。

杉木Cunninghamia lanceolata是湖北省主要造林樹(shù)種之一。傳統(tǒng)的杉木良種引種適宜生態(tài)區(qū)也是憑主觀經(jīng)驗(yàn)判斷較多。杉木良種數(shù)量較多，且生長(zhǎng)周期長(zhǎng)，像農(nóng)作物良種一樣，對(duì)所有杉木良種都進(jìn)行引種試驗(yàn)的可行性不大，所以，采用較為可靠的預(yù)測(cè)方法對(duì)引種的可行性提供數(shù)據(jù)參考尤為必要。目前，可用方法包括了MaxEnt 和ArcGIS 分析法等[3]。鑒于此，本研究擬廣東杉木良種為引種對(duì)象，應(yīng)用MaxEnt 生態(tài)學(xué)模型[4-5]，以100×100=10 000 m2，即1 hm2為單元，用34個(gè)環(huán)境因子預(yù)測(cè)廣東杉木良種在湖北省同一適宜引種生態(tài)區(qū)，以期為引種應(yīng)用提供參考。

1 數(shù)據(jù)與方法

1.1 基本信息來(lái)源

廣東省審定杉木良種信息來(lái)源于湖北省林業(yè)局林木種苗管理總站（表1）。

表1 廣東省杉木產(chǎn)區(qū)審定杉木良種Tab. 1 Cunninghamia lanceolata superior varieties approved by Cunninghamia lanceolata production area in Guangdong province

34 個(gè)環(huán)境因子數(shù)據(jù)獲取于中國(guó)氣象科學(xué)數(shù)據(jù)共享服務(wù)網(wǎng)、中國(guó)科學(xué)院資源環(huán)境科學(xué)數(shù)據(jù)中心、國(guó)家青藏高原科學(xué)數(shù)據(jù)中心、中國(guó)西部環(huán)境與生態(tài)科學(xué)數(shù)據(jù)中心（表2）。

表2 廣東省杉木良種在湖北省同一適宜引種生態(tài)區(qū)環(huán)境因子Tab. 2 Environmental factors of the Hubei province identical suitable introduction ecological distribution of superior varieties of Cunninghamia lanceolata in Guangdong province

中國(guó)行政區(qū)劃數(shù)據(jù)、中國(guó)海拔（DEM）數(shù)據(jù)

獲取于中國(guó)科學(xué)院資源環(huán)境科學(xué)數(shù)據(jù)中心和湖北省林業(yè)調(diào)查規(guī)劃院。

1.2 信息數(shù)據(jù)處理

1.2.1 分布數(shù)據(jù)處理 為避免樣點(diǎn)數(shù)據(jù)在某個(gè)地理空間上過(guò)度聚集，在廣東省杉木適生范圍內(nèi)，用ArcGIS 的Create fishnet 工具生成空間為40 行×40 列的格網(wǎng)數(shù)據(jù)，以1 個(gè)格網(wǎng)作為1 個(gè)采樣單元對(duì)杉木良種的分布數(shù)據(jù)進(jìn)行采樣（圖1）[6]。根據(jù)選育單位確定的杉木良種適宜的自然地理環(huán)境條件范圍，如，適宜海拔范圍為1 500 m 以下，在Excel 表中，剔除高程小于0 m、高程大于1 500 m、土壤厚度小于30 cm 和異常值的采樣點(diǎn)，全部采樣分布記錄共320 條。按照MaxEnt 軟件的“Samples”的要求整理數(shù)據(jù)，將分布點(diǎn)以“物種+經(jīng)度+緯度（西經(jīng)、南緯的值為負(fù)，經(jīng)緯度為十進(jìn)制小數(shù)格式。）”另存為CSV 格式文件。

圖1 廣東省杉木產(chǎn)區(qū)采樣點(diǎn)分布Fig. 1 Distribution of sampling points in Cunninghamia lanceolata production area of Guangdong province

1.2.2 環(huán)境因子處理 地形因子（經(jīng)度、緯度、高度）與環(huán)境因子有較好的回歸關(guān)系，利用中國(guó)2 160 個(gè)基本、基準(zhǔn)地面氣象觀測(cè)站的觀測(cè)數(shù)據(jù)（1981—2010 年），推算模擬無(wú)測(cè)站區(qū)域的環(huán)境資源分布情況。建立Bio1～Bio10、Bio13～Bio21、Bio24～Bio27 等23 個(gè)環(huán)境因子的空間分布模型，其表達(dá)式為：

式中，Y為環(huán)境因子要素；λ為經(jīng)度；φ為緯度；h為海拔(m)；函數(shù)f(λ,φ,h)為氣候?qū)W方程；ε為殘差項(xiàng)，可視為小地形因子(坡度、坡向等)及下墊面對(duì)環(huán)境的影響。將f(λ,φ,h)展成三維二次趨勢(shì)面方程[7]。

式中b0～b9為待定系數(shù)，利用SAS9.4 建立逐步回歸優(yōu)化回歸模型，模擬23 個(gè)環(huán)境因子的宏觀趨勢(shì)項(xiàng)，分別建立23 個(gè)環(huán)境因子的小網(wǎng)格推算模型（表3）。

表3 環(huán)境因子的小網(wǎng)格推算模型Tab. 3 Small grids reckoning models of regionalization indexes of environmental factors

在中國(guó)海拔（DEM）數(shù)據(jù)支持下，在ArcGIS里，用23 個(gè)環(huán)境因子的小網(wǎng)格推算模型，將環(huán)境因 子Bio1～Bio10、Bio13～Bio21、Bio24～Bio27 分別插值為100 m × 100 m 網(wǎng)格的基礎(chǔ)數(shù)據(jù)[8-9]。用IDW 法分別插值其殘差項(xiàng)為100 m × 100 m 網(wǎng)格的修正數(shù)據(jù)。用Spatial Analyst 工具→數(shù)學(xué)→邏輯→加，將每個(gè)環(huán)境因子的基礎(chǔ)數(shù)據(jù)和修正數(shù)據(jù)疊加相加為環(huán)境因子?xùn)鸥駭?shù)據(jù)。23 個(gè)環(huán)境因子?xùn)鸥駭?shù)據(jù)用投影柵格工具統(tǒng)一為地理坐標(biāo)系D_WGS_1984。以湖北省和廣東省矢量邊界為掩膜，裁剪出這23 個(gè)環(huán)境因子?xùn)鸥駭?shù)據(jù)圖層。最后，用柵格轉(zhuǎn)ASCII 工具將這23 個(gè)環(huán)境因子?xùn)鸥駭?shù)據(jù)轉(zhuǎn)換保存為MaxEnt 所需要的ASCII 格式文件。

在ArcGIS 里，將下載的Bio11、Bio12、Bio22、Bio23、Bio28～Bio34 等11 個(gè)環(huán)境因子數(shù)據(jù)通過(guò)重采樣工具使其像元大小與Bio1～Bio10、Bio13～Bio21、Bio24～Bio27 等23 個(gè)環(huán)境因子一致[10]。11 個(gè)環(huán)境因子數(shù)據(jù)統(tǒng)一為地理坐標(biāo)系D_WGS_1984。以湖北省和廣東省矢量邊界為掩膜，裁剪出這11 個(gè)環(huán)境因子?xùn)鸥駭?shù)據(jù)圖層。最后，用柵格轉(zhuǎn)ASCII 工具將這11個(gè)環(huán)境因子?xùn)鸥駭?shù)據(jù)轉(zhuǎn)換保存為MaxEnt 所需要的ASCII 格式文件。

1.3 模型構(gòu)建

1.3.1 MaxEnt 軟件建模 （1）物種數(shù)據(jù)：將之前導(dǎo)出的杉木良種分布數(shù)據(jù)(csv 格式)的文件，通過(guò)Browse 加載到MaxEnt 軟件“Samples”模塊。

（2）環(huán)境數(shù)據(jù)：把ASCII 格式文件的34 個(gè)環(huán)境數(shù)據(jù)加載到MaxEnt 軟件“Environmental layers”模塊。

（3）參數(shù)設(shè)置：使用auto features 選項(xiàng)，根據(jù)自動(dòng)特征規(guī)則進(jìn)行計(jì)算，所有的要素類(lèi)型都將用到。結(jié)果以comulative 類(lèi)型和ASCII 格式輸出，并定義其輸出位置。設(shè)置界面的選擇 settings里‘Random test percentage’ 設(shè) 置 為25，隨 機(jī) 選取75%的樣本點(diǎn)數(shù)據(jù)作為訓(xùn)練數(shù)據(jù)[11]，‘replicates’選擇3 次重復(fù)作為平行試驗(yàn)，最大迭代次數(shù)設(shè)為500 次，收斂閥值設(shè)為0.00 001，取值范圍0～100[12]。 選 擇‘Do jackknife to measure variable importance’衡量所有變量的重要性，MaxEnt 軟件分別對(duì)每一個(gè)環(huán)境影響因子進(jìn)行刀切圖繪出。

1.3.2 ROC 曲線(xiàn)繪制 繪制響應(yīng)曲線(xiàn)(response curves)評(píng)價(jià)模型精度。ROC 曲線(xiàn)以真陽(yáng)性率為縱坐標(biāo)(敏感性，實(shí)際存在且被預(yù)測(cè)為存在的比率)，以假陽(yáng)性率(1-特異性，實(shí)際不存在但被預(yù)測(cè)為存在的比率)為橫坐標(biāo)，AUC 值指 ROC 曲線(xiàn)與橫坐標(biāo)圍成的面積值，值域?yàn)?～1。AUC 值越大表示與隨機(jī)分布相距越遠(yuǎn)，環(huán)境因子變量與預(yù)測(cè)的物種地理分布之間的相關(guān)性越大，即模型預(yù)測(cè)效果越好，反之說(shuō)明模型預(yù)測(cè)效果越差。AUC 值在 0.5～0.6，0.6～0.7，0.7～0.8，0.8～0.9，0.9～1 分別表示模擬效果失敗、較差、一般、好、非常好[4,8]。34 個(gè)環(huán)境因子預(yù)測(cè)模型的訓(xùn)練樣本和測(cè)試樣本的AUC 值 達(dá) 到0.802 和0.739( 圖2)，AUC 均 值 在0.7～0.8 之間，說(shuō)明模型預(yù)測(cè)效果一般。

圖2 初始模型的ROC 曲線(xiàn)分析及AUC 值Fig. 2 ROC curve and AUC value for the initial model

1.4 模型優(yōu)化

在使用MaxEnt 模型進(jìn)行較大空間范圍的杉木良種同一適宜引種生態(tài)區(qū)預(yù)測(cè)時(shí)，如果環(huán)境因子變量過(guò)多、變量空間共線(xiàn)性過(guò)強(qiáng)，將導(dǎo)致模型的復(fù)雜性增加，隨機(jī)誤差增大。所以，過(guò)多低貢獻(xiàn)率的環(huán)境因子變量會(huì)導(dǎo)致模型運(yùn)行結(jié)果的準(zhǔn)確性降低。因此，需要對(duì)環(huán)境因子進(jìn)行篩選或降維[13]。

1.4.1 篩選貢獻(xiàn)率高的環(huán)境因子變量 在34 個(gè)環(huán)境因子中，對(duì)于杉木良種同一適宜引種生態(tài)區(qū)分布貢獻(xiàn)較大的環(huán)境因子變量有：Bio3、Bio6、Bio8、Bio10、Bio11、Bio17～Bio19、Bio31，累計(jì)貢獻(xiàn)率為93.8%。Bio1、Bio 2、Bio 4、Bio 5、Bio 7、Bio 9、Bio12～Bio16、Bio20～Bio30、Bio32～Bio34 等25 個(gè)環(huán)境變量的貢獻(xiàn)率都小于1%（表4），對(duì)杉木的種植分布影響有限，對(duì)這25 個(gè)環(huán)境因子變量進(jìn)行剔除[14]。

表4 各環(huán)境因子變量的貢獻(xiàn)率Tab. 4 Contribution rate of each environmental factor variable

1.4.2 篩選正規(guī)化訓(xùn)練增益高的環(huán)境因子變量刀切法(jackknife test)測(cè)定各環(huán)境因子變量權(quán)重。刀切法就是每次都忽略一個(gè)環(huán)境因子變量，然后基于剩下的環(huán)境因子變量來(lái)對(duì)杉木良種同一適宜引種生態(tài)區(qū)進(jìn)行預(yù)測(cè)，然后MaxEnt 軟件自帶程序畫(huà)出柱形圖作為依據(jù)評(píng)估環(huán)境因子變量的重要性。紅色條帶代表所有變量的貢獻(xiàn)；深藍(lán)色的條帶越長(zhǎng)，說(shuō)明該變量越重要；淺藍(lán)色的條帶長(zhǎng)度代表除該變量以外，其他所有變量組合的貢獻(xiàn)。Bio3、Bio6、Bio8、Bio10、Bio11、Bio17～Bio19、Bio31 等6 個(gè)貢獻(xiàn)較大的環(huán)境因子變量中，Bio31 對(duì)應(yīng)的深藍(lán)色條帶很短（圖3），說(shuō)明它本身的增益值幾乎接近于0，表明它對(duì)預(yù)測(cè)杉木良種同一適宜引種生態(tài)區(qū)并不是重要環(huán)境因子變量，所以，剔除環(huán)境因子變量Bio31。

圖3 刀切法的環(huán)境因子變量重要性分析Fig. 3 Analysis of the importance of environmental factors variables in the Jackknife method

1.4.3 篩選多重共線(xiàn)的環(huán)境因子變量 用GIS 軟件的值提取至點(diǎn)工具提取有效分布點(diǎn)的環(huán)境因子變量數(shù)值，用SPSS 軟件對(duì)貢獻(xiàn)較大的Bio3、Bio6、Bio8、Bio10、Bio11、Bio17～Bio19 等9 個(gè)主導(dǎo)環(huán)境因子進(jìn)行Spearman 相關(guān)分析（表5），檢驗(yàn)環(huán)境因子變量之間的多重共線(xiàn)性，Bio8 分別與Bio3、Bio6、Bio10 和Bio11 的相關(guān)系數(shù)|r|≥0.8，分別對(duì)比初始模型中二者的貢獻(xiàn)率，Bio3、Bio6、Bio10 和Bio11 貢獻(xiàn)率較小，所以，剔除貢獻(xiàn)率較小的變量Bio3、Bio6、Bio10 和Bio11，提高模型模擬的精度[14]。

表5 關(guān)鍵環(huán)境因子變量的相關(guān)系數(shù)Tab. 5 Correlation coefficient of key environmental factor variables

2 結(jié)果與分析

2.1 模型的有效性

用剩余的Bio8、Bio17～Bio19 等4 個(gè)主導(dǎo)環(huán)境因子變量重新建模，重建模型的訓(xùn)練樣本和測(cè)試樣本的AUC 值達(dá)到0.7563 和0.750(圖4)，AUC均值在0.7～0.8 之間，表明重建模型適用性及模擬精度均達(dá)到合格水平，與主導(dǎo)環(huán)境因子變量之間的相關(guān)性較大，預(yù)測(cè)同一適宜引種生態(tài)區(qū)的結(jié)果合格，可以據(jù)此進(jìn)行引種推廣。

圖4 重建模型的ROC 曲線(xiàn)分析及AUC 值Fig. 4 ROC curve and AUC value for the reconstruction model

2.2 杉木良種同一適宜引種生態(tài)區(qū)預(yù)測(cè)及適生等級(jí)劃分

MaxEnt 進(jìn)行3 次重復(fù)試驗(yàn)，選取重復(fù)試驗(yàn)中，AUC 值最高的圖層導(dǎo)入ArcGIS 軟件進(jìn)行適宜等級(jí)劃分和可視化表達(dá)（圖5）。MaxEnt 模型輸出的數(shù)據(jù)為ASC Ⅱ格式，用ArcGIS 的ASCII to Raster 功能，輸出數(shù)據(jù)類(lèi)型選FLOAT，使該結(jié)果可在 ArcGIS 中顯示[14]。利用“Reclassify”功能，劃分分布值等級(jí)及相應(yīng)分布范圍，并使用不同顏色表示，劃分標(biāo)準(zhǔn)為： 存在概率＜0.05 為不適生區(qū)；0.05 ≤存在概率＜0.33 為低適生區(qū)； 0.33 ≤存在概率＜0.66 為中適生區(qū)； 存在概率≥0.66 為高適生區(qū)[4,8]。整體來(lái)看，廣東省杉木產(chǎn)區(qū)的杉木良種在湖北省的低適生區(qū)面積為1 294 849 hm2，主要分布在：鄂東南的黃梅縣、武穴市、蘄春縣、陽(yáng)新縣、大冶市、通山縣、崇陽(yáng)縣、通城縣、咸安區(qū)和赤壁市；鄂西南的長(zhǎng)陽(yáng)縣、秭歸縣、興山縣、巴東縣、建始縣、鶴峰縣、恩施市、宣恩縣和來(lái)鳳縣。低適宜區(qū)域在引種杉木良種時(shí)，需要選擇適宜的小生境。

圖5 廣東省杉木良種在湖北省同一適宜引種生態(tài)區(qū)分布Fig. 5 Distribution of the Hubei province identical suitable introduction ecological distribution of superior varieties of Cunninghamia lanceolata in Guangdong province

2.3 杉木良種地理分布與生物氣候變量的關(guān)系

用刀切法(Jackknife Test)檢測(cè)4 個(gè)主導(dǎo)環(huán)境因子變量對(duì)于分布增益的貢獻(xiàn)，結(jié)果(表6)表明：年日最低氣溫≤0.0℃日數(shù)(Bio8)對(duì)杉木分布的增益最大，分布值隨年日最低氣溫≤0.0℃日數(shù)的升高而減小(圖6)。月最長(zhǎng)連續(xù)降水日數(shù)(Bio18)也對(duì)杉木分布的影響較大，當(dāng)月最長(zhǎng)連續(xù)降水日數(shù)為14～21 日，分布值隨月最長(zhǎng)連續(xù)降水日數(shù)的升高而增大 (圖7)。月最長(zhǎng)連續(xù)降水量(Bio19)也對(duì)杉木分布的影響較大，當(dāng)月最長(zhǎng)連續(xù)降水量為118.971～280 mm，分布值隨月最長(zhǎng)連續(xù)降水量的升高而增大；當(dāng)月最長(zhǎng)連續(xù)降水量為280～883.572 mm，分布值隨月最長(zhǎng)連續(xù)降水量的升高而減小(圖8)。

圖6 年日最低氣溫≤0.0℃日數(shù)(Bio8)反饋曲線(xiàn)Fig. 6 Days with annual daily minimum temperature below 0℃ (Bio8) feedback curve

圖7 月最長(zhǎng)連續(xù)降水日數(shù)(Bio18)反饋曲線(xiàn)Fig. 7 Monthly longest continued precipitation days(Bio18) feedback curve

圖8 月最長(zhǎng)連續(xù)降水量(Bio19)反饋曲線(xiàn)Fig. 8 Monthly maximum continued precipitation (Bio19)feedback curve

表6 主導(dǎo)環(huán)境因子變量的貢獻(xiàn)率Tab. 6 Contribution rate of dominant environmental factor variable

3 討論與結(jié)論

如果沒(méi)有做到“適地區(qū)適良種”，則可能導(dǎo)致林木良種育苗和林木良種造林的失敗。主要原因是不同地理、氣候和土壤等環(huán)境因子的質(zhì)或量不同，因而對(duì)該林木良種所要求的生態(tài)條件的滿(mǎn)足程度不同。所以，只有在適生地域內(nèi)，在適宜的立地條件下，選擇適合的林木良種，才能發(fā)揮林木良種造林的優(yōu)良特性，實(shí)現(xiàn)速生豐產(chǎn)。否則，林木良種表現(xiàn)不好，甚至不如一般的當(dāng)?shù)仄贩N。

基于MaxEnt 生態(tài)位模型的同一適宜生態(tài)區(qū)研究中，環(huán)境因子數(shù)據(jù)常來(lái)自于世界氣候—全球氣候數(shù)據(jù)庫(kù)網(wǎng)站，僅有19 個(gè)環(huán)境因子，空間分辨率僅為5arc-min[5-6,10,14]。本研究選取34 個(gè)重要環(huán)境因子，用中國(guó)2 160 個(gè)基準(zhǔn)地面氣象觀測(cè)站的觀測(cè)數(shù)據(jù)，推算模擬無(wú)測(cè)站區(qū)域的環(huán)境資源分布情況，相較于全球氣候數(shù)據(jù)庫(kù)環(huán)境因子數(shù)量更多，數(shù)據(jù)點(diǎn)更密集。

傳統(tǒng)的杉木良種引種同一適宜生態(tài)區(qū)都是以鄉(xiāng)鎮(zhèn)、縣、市、省等行政單位為單元。然而，影響杉木成活生長(zhǎng)的光、熱、水、氣等環(huán)境因子，受太陽(yáng)輻射、大氣環(huán)流的影響，在地面上呈地帶性的分布。由于山體起伏，垂直森林地帶在實(shí)際上并不都是連續(xù)的，而是由斷斷續(xù)續(xù)地呈孤島狀分布的地塊組成。所以，本研究是以100×100=10 000 m2，即1 hm2為單元。

傳統(tǒng)的林木引種是以單個(gè)樹(shù)種劃出同一適宜生態(tài)區(qū)[5-6,10,14-16]。然而，隨著自然條件演變和科學(xué)技術(shù)的發(fā)展，轉(zhuǎn)抗性基因育種、種間和遠(yuǎn)緣雜交育種等遺傳改良工作在廣泛開(kāi)展，每年都有新的林木良種通過(guò)審定。在相同的立地條件下，同一樹(shù)種，不同良種之間的生長(zhǎng)好壞是有顯著差異的。為了獲得精準(zhǔn)的引種效果，本研究以單個(gè)良種劃出同一適宜生態(tài)區(qū)精準(zhǔn)預(yù)測(cè)廣東省杉木良種在湖北省同一適宜生態(tài)區(qū)。

通過(guò)運(yùn)用MaxEnt 生態(tài)位模型對(duì)廣東省杉木良種在湖北省同一適宜引種生態(tài)區(qū)進(jìn)行分析，證明了MaxEnt 模型在林木良種引種應(yīng)用方面的可行性以及可信度，同時(shí)結(jié)合刀切法探討對(duì)杉木良種生長(zhǎng)影響最顯著的環(huán)境因子，這對(duì)廣東省杉木良種適生性分析提供了更進(jìn)一步的技術(shù)支撐。如同諸葛亮需要精準(zhǔn)的預(yù)測(cè)天氣，才能草船借到箭一樣。還如同哥倫布需要精準(zhǔn)的指南針，才能航海發(fā)現(xiàn)新大陸一樣。需要基于MaxEnt 和ArcGIS 精準(zhǔn)預(yù)測(cè)廣東省杉木良種在湖北省同一適宜引種生態(tài)區(qū)，才能避免引種廣東省杉木良種的盲目性，從而獲得良好的引種效果，值得推廣應(yīng)用。