胡 超，于 靜

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

林木良種需種植在適宜生態(tài)區(qū)域，如自然條件不適宜，良種也很難豐產和穩(wěn)產。林木良種不適應引種區(qū)自然條件會造成巨大損失。20世紀70年代，在油茶（Camellia oleifera）種植中，較多地方在引種前缺乏科學預判，盲目引種，導致幼林生長不良、成林產量低，人力、物力等方面損失較大[1]。

傳統(tǒng)的林木良種引種適宜生態(tài)區(qū)判斷僅憑主觀經驗，如適宜生態(tài)區(qū)內，水平氣候帶和垂直氣候帶都分布著特有類型的森林植被。良種由北向南或由西向東適宜種植的范圍大于相反方向的范圍，適宜的海拔高度不超過300～500 m[1]。湖北成功引種廣東和福建的馬尾松（Pinus massoniana）種子，打破“馬尾松南種北移的幅度不能超過2° ～ 3°”的定論[2]。研究表明馬尾松能成功北移6° ～ 7°育苗[2]。溫度、降水和土壤等主要環(huán)境因子相似，即為林木良種同一適宜引種生態(tài)區(qū)。

杉木（Cunninghamia lanceolata）是湖北省主要造林樹種之一。杉木良種數(shù)量較多，且生長周期長，對所有杉木良種都進行引種試驗的可行性不大。使用MaxEnt 模型和ArcGIS 軟件精準分析杉木良種同一適宜引種生態(tài)區(qū)十分必要[3]?；贛axEnt模型的同一適宜生態(tài)區(qū)的研究中，環(huán)境因子數(shù)據(jù)常來自世界氣候-全球氣候數(shù)據(jù)庫網站，僅有19 個環(huán)境因子，空間分辨率僅為5arc-min[4-6]，效果較差。為了獲得更精準的引種效果，本研究應用MaxEnt模型[4-5]，以1 hm2為1 單元，選取34 個重要環(huán)境因子，采用中國2 160 個基準地面氣象觀測站的觀測數(shù)據(jù)，推算模擬無測站區(qū)域的環(huán)境資源分布情況，精準預測湖南省杉木良種在湖北省的同一適宜引種生態(tài)區(qū)，旨在為湖北省杉木良種造林工作和經營管理提供理論依據(jù)。

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

1.1 基本信息來源

湖南省審定杉木良種為靖州縣排牙山國有林場馬園種子園的杉木優(yōu)良種源種子、杉木1 代種子園種子和杉木2 代種子園種子；會同縣林業(yè)科學研究所的杉木優(yōu)良種源種子、杉木1 代種子園種子和杉木2 代種子園種子；攸縣林業(yè)科學研究所的杉木初級種子園種子和杉木2 代種子園種子；資興市天鵝山國有林場的杉木2 代種子園種子；江華瑤族自治縣江華國有林場的杉木2 代種子園種子；城步縣林木良種場的杉木初級種子園種子；湖南省林業(yè)科學院的杉木靖全（01 ～ 27 個無性系）、杉木靖半（01～07 個無性系）、杉木會無（01～25 個無性系）、杉木會半（01～08 個無性系）、杉木會全（01～10 個無性系）、杉木江全（01 ～ 03 個無性系）、杉木攸半（01 ～02個無性系）、杉木攸全（01～18個無性系）、湘林所杉木無性系（Y1 ～ Y8）和湘林所杉木家系（A ～ D、A1、A3、A4）。湖南省杉木良種審定信息來源于湖北省林業(yè)局林木種苗管理總站。

累年年日照時數(shù)（Bio1）、累年年均氣溫（Bio3）、累年年極端最低氣溫（Bio4）、累年年極端最高氣溫（Bio5）、累年最寒冷月（1月）均氣溫（Bio6）、累年最熱月（7月）均氣溫（Bio7）、累年年日最低氣溫≤0 ℃平均日數(shù)（Bio8）、累年年日最高氣溫≥35.0 ℃平均日數(shù)（Bio9）、累年年均氣溫日較差（Bio10）、累年年平均5 cm 地溫（Bio13）、累年年均降水量（Bio14）、累年年最大降水量（Bio15）、累年年最小降水量（Bio16）、累年月最長連續(xù)無降水日數(shù)（Bio17）、累年月最長連續(xù)降水日數(shù)（Bio18）、累年月最長連續(xù)降水量（Bio19）、累年年最大日降水量（Bio20）、累年年均相對濕度（Bio21）、累年年均風速（Bio24）、累年年極大風速（Bio25）、累年年日最大風速≥5.0 m∕s日數(shù)（Bio26）和 累年年日 最大風速≥10.0 m∕s日 數(shù)（Bio27）數(shù)據(jù)獲取于中國氣象科學數(shù)據(jù)共享服務網?！? ℃積溫（Bio11）、≥10 ℃積溫（Bio12）、濕潤指數(shù)（Bio22）和干燥度（Bio23）數(shù)據(jù)獲取于中國科學院資源環(huán)境科學數(shù)據(jù)中心。太陽輻射日均值（Bio2）數(shù)據(jù)獲取于國家青藏高原科學數(shù)據(jù)中心。土壤類型（Bio28）、土壤酸堿度（Bio29）、土壤陽離子交換能力（Bio30）、土壤深度（Bio31）、土壤有機碳含量（Bio32）、土壤沙含量（Bio33）和土壤有效水含量（Bio34）數(shù)據(jù)獲取于中國西部環(huán)境與生態(tài)科學數(shù)據(jù)中心。

中國行政區(qū)劃數(shù)據(jù)、中國海拔高度（DEM）數(shù)據(jù)獲取于中國科學院資源環(huán)境科學數(shù)據(jù)中心和湖北省林業(yè)調查規(guī)劃院。

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

1.2.1 分布數(shù)據(jù)處理

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

圖1 湖南省杉木產區(qū)采樣點分布Fig.1 Distribution of sampling points in C.lanceolata pro?duction area of Hunan province

1.2.2 環(huán)境因子處理

地形因子（經度、緯度和高度）與環(huán)境因子有較好的回歸關系，利用中國2 160 個基本、基準地面氣象觀測站的觀測數(shù)據(jù)，推算模擬無測站區(qū)域的環(huán)境資源分布情況。建立Bio1～Bio10、Bio13～Bio21 和Bio24～Bio27 等23 個環(huán)境因子的空間分布模型，表達式為：

式中，Y為環(huán)境因子要素；λ為經度；φ為緯度；h為海拔高度（m）；函數(shù)f（λ,φ,h）為氣候學方程；ε為殘差項，可視為小地形因子（坡度、坡向等）和下墊面對環(huán)境的影響。將f（λ,φ,h）展成三維二次趨勢面方程[7]：

式中，b0～b9為待定系數(shù)，利用SAS 9.4建立逐步回歸優(yōu)化模型，模擬23 個環(huán)境因子的宏觀趨勢項，分別建立23個環(huán)境因子的小網格推算模型，各模型均通過可信度α=0.01的顯著性檢驗（表1）。

表1 環(huán)境因子的小網格推算模型Tab.1 Small grids reckoning models of environmental factors

在中國海拔高度（DEM）數(shù)據(jù)支持下，在ArcGIS 10 里，用23 個環(huán)境因子的小網格推算模型，將環(huán)境因子Bio1～Bio10、Bio13～Bio21 和Bio24～Bio27 分別插值為100 m × 100 m 網格的基礎數(shù)據(jù)[8-9]。用IDW 法分別插值其殘差項為100 m ×100 m 網格的修正數(shù)據(jù)。采用Spatial Analyst 工具→數(shù)學→邏輯→加，將每個環(huán)境因子的基礎數(shù)據(jù)和修正數(shù)據(jù)疊加，即為環(huán)境因子柵格數(shù)據(jù)。將23個環(huán)境因子柵格數(shù)據(jù)用投影柵格工具統(tǒng)一為地理坐標系D_WGS_1984。以湖北省和湖南省矢量邊界為掩膜，裁剪出這23個環(huán)境因子的柵格數(shù)據(jù)圖層。最后，用柵格轉ASCⅡ工具將這23 個環(huán)境因子的柵格數(shù)據(jù)轉換保存為MaxEnt模型所需要的ASCⅡ格式文件。

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

續(xù)表1 Continued

1.3 模型構建

1.3.1 MaxEnt軟件建模

將之前導出的杉木良種分布數(shù)據(jù)（csv 格式）的文件，通過Browse 加載到MaxEnt 軟件“Samples”模塊。把ASCⅡ格式文件的34 個環(huán)境數(shù)據(jù)加載到MaxEnt 軟件“Environmental layers”模塊。使用au?to features 選項，根據(jù)自動特征規(guī)則進行計算，所有的要素類型都將用到。結果以comulative 類型和ASCⅡ格式輸出，并定義其輸出位置。將界面選擇settings 里的“Random test percentage”設置為25，隨機選取75%的樣本點數(shù)據(jù)作為訓練數(shù)據(jù)[11]；設置settings 中的replicates，本試驗選擇3 次重復作為平行試驗，最大迭代次數(shù)設為500 次，收斂閥值設為0.000 01，取值范圍設為0 ～ 100[12]。選擇“Do jack?knife to measure variable importance”衡量所有變量的重要性，采用MaxEnt軟件分別對每個環(huán)境影響因子進行刀切圖繪出。

進入新世紀以來，隨著人們環(huán)保意識的不斷提升，我國人民對環(huán)境保護的重視程度越來越高，在市政工程建設中，都市生態(tài)環(huán)境保護也成為一項十分重要的工作。因此在一些發(fā)達城市的市政工程建設中，會有意識的采取一些都市生態(tài)保護措施，發(fā)揮出了不俗的作用。但與此同時，不重視都市生態(tài)保護的現(xiàn)象依舊十分嚴重，如決策過于隨意、超預算投資、資金不合理浪費等。針對這些情況，我國必須進一步優(yōu)化市政工程管理，為都市生態(tài)健康發(fā)展提供良好的保障。

1.3.2 ROC曲線繪制

繪制響應曲線評價模型精度。ROC 曲線以真陽性率為縱坐標軸（敏感性，實際存在且被預測為存在的比率），以假陽性率（特異性，實際不存在但被預測為存在的比率）為橫坐標；AUC 值為ROC 曲線與橫坐標軸圍成的面積值，為0 ～ 1。AUC 值越大，表示與隨機分布相距越遠，環(huán)境因子變量與預測的杉木良種同一適宜引種生態(tài)區(qū)間的相關性越大，即模型預測效果越好，反之說明模型預測效果越差。AUC 值在0.5～0.6、0.6～0.7、0.7～0.8、0.8～0.9和0.9～1分別表示模擬效果失敗、較差、一般、好和非常好[4,8]。34 個環(huán)境因子預測模型的訓練和測試樣本的AUC 值分別為0.767 和0.700，均在0.7 ～0.8，說明模型預測效果一般（圖2）。

圖2 初始模型的ROC曲線及AUC值Fig.2 ROC curves and AUC values of initial model

1.4 模型優(yōu)化

在使用MaxEnt 模型進行較大空間范圍的杉木良種同一適宜引種生態(tài)區(qū)預測時，環(huán)境因子變量過多、變量空間共線性過強將導致模型的復雜性增加，隨機誤差增大。過多低貢獻率的環(huán)境因子變量會導致模型運行結果的準確性降低。因此，需要對環(huán)境因子進行篩選或降維[13]。

1.4.1 篩選貢獻率高的環(huán)境因子變量

在34個環(huán)境因子中，對適宜引種生態(tài)區(qū)分布貢獻較大的環(huán)境因子變量有Bio1、Bio2、Bio7、Bio8、Bio16～Bio18、Bio20、Bio21和Bio32等10個因子，累計 貢獻 率 為91.4%（表2）。Bio3 ～ Bio6、Bio9 ～Bio15、Bio19、Bio22 ～ Bio31、Bio33 和Bio34 等24 個環(huán)境因子變量的貢獻率都小于1%，對分布影響有限，對這24個環(huán)境因子變量進行剔除[14]。

表2 各環(huán)境因子變量的貢獻率Tab.2 Contribution rates of environmental factor variables

1.4.2 篩選正規(guī)化訓練增益高的環(huán)境因子變量

采用刀切法測定各環(huán)境因子變量權重。刀切法即每次都忽略一個環(huán)境因子變量，基于剩下的環(huán)境因子變量對適宜引種生態(tài)區(qū)進行預測，以MaxEnt模型自帶程序畫出的柱形圖為依據(jù)評估環(huán)境因子變量的重要性。紅色條帶代表所有變量的貢獻；深藍色條帶越長，說明該變量越重要；淺藍色條帶代表除該變量外，其他所有變量組合的貢獻。10 個貢獻較大的環(huán)境因子變量中，Bio7 對應的深藍色條帶很短，說明它的增益值接近于0，表明它對預測適宜引種生態(tài)區(qū)并不重要，所以，剔除環(huán)境因子變量Bio7（圖3）。

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

1.4.3 篩選多重共線的環(huán)境因子變量

用ArcGIS 軟件的值提取至點工具，提取有效分布點的環(huán)境因子變量數(shù)值，用SPSS軟件對貢獻較大的9 個主導環(huán)境因子變量進行Spearman 相關分析，檢驗環(huán)境因子變量間的多重共線性，9 個主導環(huán)境因子變量間的相關系數(shù)|r| ＜0.8（表3）。保留Bio1、Bio2、Bio8、Bio16 ～ Bio18、Bio20、Bio21 和Bio32 等9個主導環(huán)境因子變量重新建模[14]。

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

2 結果與分析

2.1 模型的有效性

用選出的9 個主導環(huán)境因子變量重新建模，模型的訓練和測試樣本的AUC 值分別為0.753 和0.723（圖4），均在0.7～0.8，表明重建模型的適用性及模擬精度均達到合格水平，與主導環(huán)境因子變量間的相關性較大，預測同一適宜引種生態(tài)區(qū)的結果合格，可以據(jù)此進行引種推廣。

圖4 重建模型的ROC曲線及AUC值Fig.4 ROC curves and AUC values of reconstruction model

2.2 杉木良種同一適宜引種生態(tài)區(qū)預測和適生等級劃分

運用MaxEnt 模型進行3 次重復試驗，選取重復試驗中AUC值最高的圖層導入ArcGIS軟件，進行適宜等級劃分和可視化表達（圖5）。MaxEnt模型輸出的數(shù)據(jù)為ASCⅡ格式，用ArcGIS的“ASCⅡto Raster”功能，輸出數(shù)據(jù)類型選FLOAT，使該結果可在ArcGIS 軟件中顯示[14]。利用“Reclassify”功能,劃分分布值等級及相應分布范圍，并使用不同顏色表示。劃分標準為存在概率＜0.05 為不適生區(qū)；0.05 ≤存在概率＜0.33 為低適生區(qū)；0.33 ≤存在概率＜0.66 為中適生區(qū)；存在概率≥0.66 為高適生區(qū)[4,8]。整體來看，湖南省杉木良種在湖北省的中適生區(qū)面積為622 618 hm2，主要分布在鄂東的黃梅縣和通城縣；鄂中的監(jiān)利縣、石首市和松滋市；鄂西的宜都市、長陽縣、五峰縣、巴東縣、建始縣、鶴峰縣、恩施市、宣恩縣、咸豐縣和來鳳縣。湖南省杉木良種在湖北省的低適生區(qū)面積為5 942 618 hm2，主要分布在鄂東的蘄春縣、浠水縣、武穴市、黃州區(qū)、鄂州市、陽新縣、大冶市、江夏區(qū)、蔡甸區(qū)、漢南區(qū)、咸安區(qū)、通山縣、崇陽縣、赤壁市和嘉魚縣；鄂中的仙桃市、洪湖市、潛江市、荊州區(qū)、沙市區(qū)、江陵縣、公安縣、枝江市和松滋市；鄂西的枝江市、當陽市、夷陵區(qū)、點軍區(qū)、秭歸縣、利川市和來鳳縣。在低適生區(qū)引種杉木良種時，需選擇適宜的小生境。無高適生區(qū)。

圖5 湖南省杉木良種在湖北省同一適宜引種生態(tài)區(qū)分布Fig.5 Suitable ecological introduction distributions of C.lanceolata superior varieties from Hunan province in Hubei province

2.3 杉木良種地理分布與生物氣候變量的關系

用刀切法檢測9個主導環(huán)境因子變量對于分布增益的貢獻。結果表明，累年年最小降水量（Bio16）對杉木分布的增益最大，當年最小降水量為600 ～1 174.35 mm 時，分布值隨年最小降水量的升高而增大（表4，圖6）；其次是累年年日照時數(shù)（Bio1），當年日照時數(shù)為1 172.509～1 300 h時，分布值隨年日照時數(shù)增加而增大；當年日照時數(shù)為1 300～2 761.266 h時，分布值隨年日照時數(shù)減少而減?。▓D7）。

圖6 累年年最小降水量（Bio16）反饋曲線Fig.6 Annual minimum precipitation(Bio16)feedback curve

圖7 累年年日照時數(shù)（Bio1）反饋曲線Fig.7 Annual sunshine hours(Bio1)feedback curve

表4 主導環(huán)境因子變量的貢獻率Tab.4 Contribution rates of dominant environmental factor variables

3 討論與結論

林木良種如果沒有在適宜地區(qū)種植，可能導致育苗和造林失敗。傳統(tǒng)的杉木良種引種同一適宜生態(tài)區(qū)都是以鄉(xiāng)（鎮(zhèn)）、縣、市和省等行政單位為單元。影響杉木生長的光、熱、水和氣等環(huán)境因子受太陽輻射和大氣環(huán)流的影響，在地面上呈地帶性分布。由于山體起伏，垂直森林地帶由斷斷續(xù)續(xù)呈孤島狀分布的地塊組成。為獲得精準的引種效果，本研究以1 hm2為單元，運用MaxEnt 模型精準預測同一適宜引種生態(tài)區(qū)。

傳統(tǒng)的林木引種是以單個樹種劃出同一適宜生態(tài)區(qū)。隨著自然條件演變和科學技術的發(fā)展，轉基因抗性育種、種間和遠緣雜交育種等遺傳改良工作被廣泛開展，每年都有新的林木良種通過審定。在相同的立地條件下，同一樹種不同良種間的生長差異顯著。為獲得精準的引種效果，本研究以單個良種劃出同一適宜生態(tài)區(qū)，運用MaxEnt模型進行精準預測。

在大尺度空間上，氣候因子是限制杉木良種同一適宜引種生態(tài)區(qū)范圍的關鍵因素。MaxEnt 模型運行結果表明，累年年最小降水量和累年年日照時數(shù)是影響杉木良種分布的最重要氣候因子，反映了杉木喜陰好濕的生物特性。本研究利用主導氣候因子與分布概率間的反饋曲線獲得各主導氣候因子的閾值，當年最小降水量大于920 mm、年日照時數(shù)為1 300 ～ 1 550 h 時，杉木良種存活概率較高。此結果有利于了解杉木良種與氣候的關系，對于精準預測適宜引種生態(tài)區(qū)具有重要意義。

通過運用MaxEnt 模型對湖南省杉木良種在湖北省同一適宜引種生態(tài)區(qū)進行分析，證明MaxEnt模型在林木良種引種應用方面的可行性及可信度，同時結合刀切法探討對杉木良種生長影響最顯著的環(huán)境因子，為良種適生性分析提供更進一步的技術支撐?；贛axEnt和Arc GIS 精準預測適宜引種生態(tài)區(qū)，能避免引種的盲目性，獲得良好的引種效果，該方法值得推廣應用。