基于MaxEnt模型和ArcGIS精准预测湖南省杉木良种在湖北省同一适宜引种生态区

胡 超，于 静

（1.湖北省林业局林木种苗管理总站，湖北武汉 430079；2.岭南生态文旅股份有限公司，湖北武汉 430062）

林木良种需种植在适宜生态区域，如自然条件不适宜，良种也很难丰产和稳产。林木良种不适应引种区自然条件会造成巨大损失。20世纪70年代，在油茶（Camellia oleifera）种植中，较多地方在引种前缺乏科学预判，盲目引种，导致幼林生长不良、成林产量低，人力、物力等方面损失较大[1]。

传统的林木良种引种适宜生态区判断仅凭主观经验，如适宜生态区内，水平气候带和垂直气候带都分布着特有类型的森林植被。良种由北向南或由西向东适宜种植的范围大于相反方向的范围，适宜的海拔高度不超过300～500 m[1]。湖北成功引种广东和福建的马尾松（Pinus massoniana）种子，打破“马尾松南种北移的幅度不能超过2° ～ 3°”的定论[2]。研究表明马尾松能成功北移6° ～ 7°育苗[2]。温度、降水和土壤等主要环境因子相似，即为林木良种同一适宜引种生态区。

杉木（Cunninghamia lanceolata）是湖北省主要造林树种之一。杉木良种数量较多，且生长周期长，对所有杉木良种都进行引种试验的可行性不大。使用MaxEnt 模型和ArcGIS 软件精准分析杉木良种同一适宜引种生态区十分必要[3]。基于MaxEnt模型的同一适宜生态区的研究中，环境因子数据常来自世界气候-全球气候数据库网站，仅有19 个环境因子，空间分辨率仅为5arc-min[4-6]，效果较差。为了获得更精准的引种效果，本研究应用MaxEnt模型[4-5]，以1 hm2为1 单元，选取34 个重要环境因子，采用中国2 160 个基准地面气象观测站的观测数据，推算模拟无测站区域的环境资源分布情况，精准预测湖南省杉木良种在湖北省的同一适宜引种生态区，旨在为湖北省杉木良种造林工作和经营管理提供理论依据。

1 数据与方法

1.1 基本信息来源

湖南省审定杉木良种为靖州县排牙山国有林场马园种子园的杉木优良种源种子、杉木1 代种子园种子和杉木2 代种子园种子；会同县林业科学研究所的杉木优良种源种子、杉木1 代种子园种子和杉木2 代种子园种子；攸县林业科学研究所的杉木初级种子园种子和杉木2 代种子园种子；资兴市天鹅山国有林场的杉木2 代种子园种子；江华瑶族自治县江华国有林场的杉木2 代种子园种子；城步县林木良种场的杉木初级种子园种子；湖南省林业科学院的杉木靖全（01 ～ 27 个无性系）、杉木靖半（01～07 个无性系）、杉木会无（01～25 个无性系）、杉木会半（01～08 个无性系）、杉木会全（01～10 个无性系）、杉木江全（01 ～ 03 个无性系）、杉木攸半（01 ～02个无性系）、杉木攸全（01～18个无性系）、湘林所杉木无性系（Y1 ～ Y8）和湘林所杉木家系（A ～ D、A1、A3、A4）。湖南省杉木良种审定信息来源于湖北省林业局林木种苗管理总站。

累年年日照时数（Bio1）、累年年均气温（Bio3）、累年年极端最低气温（Bio4）、累年年极端最高气温（Bio5）、累年最寒冷月（1月）均气温（Bio6）、累年最热月（7月）均气温（Bio7）、累年年日最低气温≤0 ℃平均日数（Bio8）、累年年日最高气温≥35.0 ℃平均日数（Bio9）、累年年均气温日较差（Bio10）、累年年平均5 cm 地温（Bio13）、累年年均降水量（Bio14）、累年年最大降水量（Bio15）、累年年最小降水量（Bio16）、累年月最长连续无降水日数（Bio17）、累年月最长连续降水日数（Bio18）、累年月最长连续降水量（Bio19）、累年年最大日降水量（Bio20）、累年年均相对湿度（Bio21）、累年年均风速（Bio24）、累年年极大风速（Bio25）、累年年日最大风速≥5.0 m∕s日数（Bio26）和 累年年日 最大风速≥10.0 m∕s日 数（Bio27）数据获取于中国气象科学数据共享服务网。≥0 ℃积温（Bio11）、≥10 ℃积温（Bio12）、湿润指数（Bio22）和干燥度（Bio23）数据获取于中国科学院资源环境科学数据中心。太阳辐射日均值（Bio2）数据获取于国家青藏高原科学数据中心。土壤类型（Bio28）、土壤酸碱度（Bio29）、土壤阳离子交换能力（Bio30）、土壤深度（Bio31）、土壤有机碳含量（Bio32）、土壤沙含量（Bio33）和土壤有效水含量（Bio34）数据获取于中国西部环境与生态科学数据中心。

中国行政区划数据、中国海拔高度（DEM）数据获取于中国科学院资源环境科学数据中心和湖北省林业调查规划院。

1.2 信息数据处理

1.2.1 分布数据处理

为避免样点数据在某个地理空间上过度聚集，在湖南省杉木适生范围内，采用ArcGIS 10 的Create fishnet 工具生成空间为30行×30列的网格数据，以1个网格为1个采样单元，对杉木良种的分布数据进行采样（图1）[6]。根据选育单位确定杉木良种适宜的自然地理环境条件，如适宜海拔范围为1 500 m以下。在Excel 表中，剔除高程＜0 m、高程＞1 500 m、土壤厚度＜30 cm 和异常值的采样点，分布记录共496 条。按照MaxEnt 软件的“Samples”要求整理数据，将分布点以“物种+ 经度+ 纬度（西经、南纬的值为负，经纬度为十进制小数格式）”另存为csv格式文件。

图1 湖南省杉木产区采样点分布Fig.1 Distribution of sampling points in C.lanceolata pro⁃duction area of Hunan province

1.2.2 环境因子处理

地形因子（经度、纬度和高度）与环境因子有较好的回归关系，利用中国2 160 个基本、基准地面气象观测站的观测数据，推算模拟无测站区域的环境资源分布情况。建立Bio1～Bio10、Bio13～Bio21 和Bio24～Bio27 等23 个环境因子的空间分布模型，表达式为：

式中，Y为环境因子要素；λ为经度；φ为纬度；h为海拔高度（m）；函数f（λ,φ,h）为气候学方程；ε为残差项，可视为小地形因子（坡度、坡向等）和下垫面对环境的影响。将f（λ,φ,h）展成三维二次趋势面方程[7]：

式中，b0～b9为待定系数，利用SAS 9.4建立逐步回归优化模型，模拟23 个环境因子的宏观趋势项，分别建立23个环境因子的小网格推算模型，各模型均通过可信度α=0.01的显著性检验（表1）。

表1 环境因子的小网格推算模型Tab.1 Small grids reckoning models of environmental factors

在中国海拔高度（DEM）数据支持下，在ArcGIS 10 里，用23 个环境因子的小网格推算模型，将环境因子Bio1～Bio10、Bio13～Bio21 和Bio24～Bio27 分别插值为100 m × 100 m 网格的基础数据[8-9]。用IDW 法分别插值其残差项为100 m ×100 m 网格的修正数据。采用Spatial Analyst 工具→数学→逻辑→加，将每个环境因子的基础数据和修正数据叠加，即为环境因子栅格数据。将23个环境因子栅格数据用投影栅格工具统一为地理坐标系D_WGS_1984。以湖北省和湖南省矢量边界为掩膜，裁剪出这23个环境因子的栅格数据图层。最后，用栅格转ASCⅡ工具将这23 个环境因子的栅格数据转换保存为MaxEnt模型所需要的ASCⅡ格式文件。

在ArcGIS 10 里，将下载的Bio11、Bio12、Bio22、Bio23 和Bio28～Bio34 等11 个环境因子数据通过重采样工具使其像元大小与Bio1 ～ Bio10、Bio13 ～Bio21、Bio24～Bio27等23个环境因子一致[10]。11个环境因子数据统一为地理坐标系D_WGS_1984。以湖北省和湖南省矢量边界为掩膜，裁剪出这11个环境因子的栅格数据图层。最后，用栅格转ASCⅡ工具将这11 个环境因子的栅格数据转换保存为Max⁃Ent模型所需要的ASCⅡ格式文件。

续表1 Continued

1.3 模型构建

1.3.1 MaxEnt软件建模

将之前导出的杉木良种分布数据（csv 格式）的文件，通过Browse 加载到MaxEnt 软件“Samples”模块。把ASCⅡ格式文件的34 个环境数据加载到MaxEnt 软件“Environmental layers”模块。使用au⁃to features 选项，根据自动特征规则进行计算，所有的要素类型都将用到。结果以comulative 类型和ASCⅡ格式输出，并定义其输出位置。将界面选择settings 里的“Random test percentage”设置为25，随机选取75%的样本点数据作为训练数据[11]；设置settings 中的replicates，本试验选择3 次重复作为平行试验，最大迭代次数设为500 次，收敛阀值设为0.000 01，取值范围设为0 ～ 100[12]。选择“Do jack⁃knife to measure variable importance”衡量所有变量的重要性，采用MaxEnt软件分别对每个环境影响因子进行刀切图绘出。

进入新世纪以来，随着人们环保意识的不断提升，我国人民对环境保护的重视程度越来越高，在市政工程建设中，都市生态环境保护也成为一项十分重要的工作。因此在一些发达城市的市政工程建设中，会有意识的采取一些都市生态保护措施，发挥出了不俗的作用。但与此同时，不重视都市生态保护的现象依旧十分严重，如决策过于随意、超预算投资、资金不合理浪费等。针对这些情况，我国必须进一步优化市政工程管理，为都市生态健康发展提供良好的保障。

1.3.2 ROC曲线绘制

绘制响应曲线评价模型精度。ROC 曲线以真阳性率为纵坐标轴（敏感性，实际存在且被预测为存在的比率），以假阳性率（特异性，实际不存在但被预测为存在的比率）为横坐标；AUC 值为ROC 曲线与横坐标轴围成的面积值，为0 ～ 1。AUC 值越大，表示与随机分布相距越远，环境因子变量与预测的杉木良种同一适宜引种生态区间的相关性越大，即模型预测效果越好，反之说明模型预测效果越差。AUC 值在0.5～0.6、0.6～0.7、0.7～0.8、0.8～0.9和0.9～1分别表示模拟效果失败、较差、一般、好和非常好[4,8]。34 个环境因子预测模型的训练和测试样本的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 模型优化

在使用MaxEnt 模型进行较大空间范围的杉木良种同一适宜引种生态区预测时，环境因子变量过多、变量空间共线性过强将导致模型的复杂性增加，随机误差增大。过多低贡献率的环境因子变量会导致模型运行结果的准确性降低。因此，需要对环境因子进行筛选或降维[13]。

1.4.1 筛选贡献率高的环境因子变量

在34个环境因子中，对适宜引种生态区分布贡献较大的环境因子变量有Bio1、Bio2、Bio7、Bio8、Bio16～Bio18、Bio20、Bio21和Bio32等10个因子，累计 贡献 率 为91.4%（表2）。Bio3 ～ Bio6、Bio9 ～Bio15、Bio19、Bio22 ～ Bio31、Bio33 和Bio34 等24 个环境因子变量的贡献率都小于1%，对分布影响有限，对这24个环境因子变量进行剔除[14]。

表2 各环境因子变量的贡献率Tab.2 Contribution rates of environmental factor variables

1.4.2 筛选正规化训练增益高的环境因子变量

采用刀切法测定各环境因子变量权重。刀切法即每次都忽略一个环境因子变量，基于剩下的环境因子变量对适宜引种生态区进行预测，以MaxEnt模型自带程序画出的柱形图为依据评估环境因子变量的重要性。红色条带代表所有变量的贡献；深蓝色条带越长，说明该变量越重要；浅蓝色条带代表除该变量外，其他所有变量组合的贡献。10 个贡献较大的环境因子变量中，Bio7 对应的深蓝色条带很短，说明它的增益值接近于0，表明它对预测适宜引种生态区并不重要，所以，剔除环境因子变量Bio7（图3）。

图3 刀切法的环境因子变量重要性分析Fig.3 Analysis on importances of environmental factors variables by Jackknife method

1.4.3 筛选多重共线的环境因子变量

用ArcGIS 软件的值提取至点工具，提取有效分布点的环境因子变量数值，用SPSS软件对贡献较大的9 个主导环境因子变量进行Spearman 相关分析，检验环境因子变量间的多重共线性，9 个主导环境因子变量间的相关系数|r| ＜0.8（表3）。保留Bio1、Bio2、Bio8、Bio16 ～ Bio18、Bio20、Bio21 和Bio32 等9个主导环境因子变量重新建模[14]。

表3 关键环境因子变量的相关系数Tab.3 Correlation coefficients of key environmental factor variables

2 结果与分析

2.1 模型的有效性

用选出的9 个主导环境因子变量重新建模，模型的训练和测试样本的AUC 值分别为0.753 和0.723（图4），均在0.7～0.8，表明重建模型的适用性及模拟精度均达到合格水平，与主导环境因子变量间的相关性较大，预测同一适宜引种生态区的结果合格，可以据此进行引种推广。

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

2.2 杉木良种同一适宜引种生态区预测和适生等级划分

运用MaxEnt 模型进行3 次重复试验，选取重复试验中AUC值最高的图层导入ArcGIS软件，进行适宜等级划分和可视化表达（图5）。MaxEnt模型输出的数据为ASCⅡ格式，用ArcGIS的“ASCⅡto Raster”功能，输出数据类型选FLOAT，使该结果可在ArcGIS 软件中显示[14]。利用“Reclassify”功能,划分分布值等级及相应分布范围，并使用不同颜色表示。划分标准为存在概率＜0.05 为不适生区；0.05 ≤存在概率＜0.33 为低适生区；0.33 ≤存在概率＜0.66 为中适生区；存在概率≥0.66 为高适生区[4,8]。整体来看，湖南省杉木良种在湖北省的中适生区面积为622 618 hm2，主要分布在鄂东的黄梅县和通城县；鄂中的监利县、石首市和松滋市；鄂西的宜都市、长阳县、五峰县、巴东县、建始县、鹤峰县、恩施市、宣恩县、咸丰县和来凤县。湖南省杉木良种在湖北省的低适生区面积为5 942 618 hm2，主要分布在鄂东的蕲春县、浠水县、武穴市、黄州区、鄂州市、阳新县、大冶市、江夏区、蔡甸区、汉南区、咸安区、通山县、崇阳县、赤壁市和嘉鱼县；鄂中的仙桃市、洪湖市、潜江市、荆州区、沙市区、江陵县、公安县、枝江市和松滋市；鄂西的枝江市、当阳市、夷陵区、点军区、秭归县、利川市和来凤县。在低适生区引种杉木良种时，需选择适宜的小生境。无高适生区。

图5 湖南省杉木良种在湖北省同一适宜引种生态区分布Fig.5 Suitable ecological introduction distributions of C.lanceolata superior varieties from Hunan province in Hubei province

2.3 杉木良种地理分布与生物气候变量的关系

用刀切法检测9个主导环境因子变量对于分布增益的贡献。结果表明，累年年最小降水量（Bio16）对杉木分布的增益最大，当年最小降水量为600 ～1 174.35 mm 时，分布值随年最小降水量的升高而增大（表4，图6）；其次是累年年日照时数（Bio1），当年日照时数为1 172.509～1 300 h时，分布值随年日照时数增加而增大；当年日照时数为1 300～2 761.266 h时，分布值随年日照时数减少而减小（图7）。

图6 累年年最小降水量（Bio16）反馈曲线Fig.6 Annual minimum precipitation(Bio16)feedback curve

图7 累年年日照时数（Bio1）反馈曲线Fig.7 Annual sunshine hours(Bio1)feedback curve

表4 主导环境因子变量的贡献率Tab.4 Contribution rates of dominant environmental factor variables

3 讨论与结论

林木良种如果没有在适宜地区种植，可能导致育苗和造林失败。传统的杉木良种引种同一适宜生态区都是以乡（镇）、县、市和省等行政单位为单元。影响杉木生长的光、热、水和气等环境因子受太阳辐射和大气环流的影响，在地面上呈地带性分布。由于山体起伏，垂直森林地带由断断续续呈孤岛状分布的地块组成。为获得精准的引种效果，本研究以1 hm2为单元，运用MaxEnt 模型精准预测同一适宜引种生态区。

传统的林木引种是以单个树种划出同一适宜生态区。随着自然条件演变和科学技术的发展，转基因抗性育种、种间和远缘杂交育种等遗传改良工作被广泛开展，每年都有新的林木良种通过审定。在相同的立地条件下，同一树种不同良种间的生长差异显著。为获得精准的引种效果，本研究以单个良种划出同一适宜生态区，运用MaxEnt模型进行精准预测。

在大尺度空间上，气候因子是限制杉木良种同一适宜引种生态区范围的关键因素。MaxEnt 模型运行结果表明，累年年最小降水量和累年年日照时数是影响杉木良种分布的最重要气候因子，反映了杉木喜阴好湿的生物特性。本研究利用主导气候因子与分布概率间的反馈曲线获得各主导气候因子的阈值，当年最小降水量大于920 mm、年日照时数为1 300 ～ 1 550 h 时，杉木良种存活概率较高。此结果有利于了解杉木良种与气候的关系，对于精准预测适宜引种生态区具有重要意义。

通过运用MaxEnt 模型对湖南省杉木良种在湖北省同一适宜引种生态区进行分析，证明MaxEnt模型在林木良种引种应用方面的可行性及可信度，同时结合刀切法探讨对杉木良种生长影响最显著的环境因子，为良种适生性分析提供更进一步的技术支撑。基于MaxEnt和Arc GIS 精准预测适宜引种生态区，能避免引种的盲目性，获得良好的引种效果，该方法值得推广应用。