首页 > 分享 > 基于GIS空间技术和MaxEnt模型预测川西松材线虫病入侵风险

基于GIS空间技术和MaxEnt模型预测川西松材线虫病入侵风险

花匠小妙招
2024-09-15 16:57

松材线虫病又称松树萎蔫病（pine wilt disease），是一种由松材线虫（Bursaphelenchus xylophilus）引起的松树毁灭性病害[1]。松材线虫病最早发现于北美洲，20世纪初期传入亚洲（日本），我国在1982年首次发现该病害[2]。经过近40年的扩散传播，松材线虫逐渐适应我国生态环境与松树树种多样性，近年危害程度有加深加快的趋势[3]。根据国家林业和草原局2022年第6号公告：截止2021年，我国松材线虫病疫情已波及19个省（区、市）742县级行政区，发生面积超171万hm2（http://www.forestry.gov.cn）。松材线虫病主要通过松墨天牛（Monochamus alternatus）和云杉花墨天牛（M. saltuarius）等媒介昆虫在松树体内传播[4]。天牛羽化后需要飞往健康松树补充营养，在该过程中即有可能将附着在媒介昆虫体表的松材线虫传播到健康松树上[2]。我国所有的松树种类均有可能成为松材线虫的寄主，近40年来，松材线虫病已在我国累计致死松树10亿株，造成千亿元级巨大的直接和间接经济损失，严重威胁森林生态系统结构稳定性与功能发挥[5]，已成为我国头号林业检疫性有害生物[6]。

已有研究表明，在全球范围内松材线虫病传播受气温和降水格局的影响[7]，尤其是夏季均温高于20 ℃的地区，其松材线虫病传播感染风险显著增加[8]。当前全球气候变暖以及降水格局（极端干旱）变化等潜在全球气候变化将进一步增加松材线虫病传播的风险[9]。在我国，松材线虫病疫情有北扩西进的趋势，这与我国温带和青藏高原地区是全球气候变化敏感区可能存在关联性，但相关区域或局域尺度的影响因子尚不清楚。松材线虫的定殖与传播是一个复杂的病害系统，除受大尺度气候格局影响外，还受许多局域尺度环境和和生物因子的影响[10−11]。例如，地形因子已被证明通过改变水热和土壤条件直接影响松树的生长与分布，松树健康与否会间接影响媒介昆虫天牛的生存与定殖，最终影响松材线虫的传播[12−14]。另外，不同植被类型及其群落植物物种组成也会在很大程度上决定松材线虫病的发生概率[15−16]。居民点聚集、离林区道路距离以及人为利用枯死松木等也被证明是导致松材线虫病发生和传播的关键因素[12−13]。

当前，对松材线虫及媒介昆虫天牛（Monochamus spp.）的适宜生境和风险评估主要采用基于生态位理论的模型模拟，包括生物气候分析系统、基于规则集的遗传算法模型、生态位因子分析模型以及最大熵模型等[17−19]。其中最大熵（maximum entropy，MaxEnt）模型以其操作简易、运算时间短、精准性高、适合小样本的优点备受科研人员青睐，是预测和评价物种生境空间分布模式的有效工具[20−21]。近几年我国有学者利用松材线虫和中间宿主松墨天牛的地理分布数据，基于MaxEnt模型成功预测了松材线虫病的潜在发生区[22−24]，为有关部门开展松材线虫病害的早期监测、预警及防控提供科技支撑。

我国西南高山林区是针叶林的主要分布区，川西理县米亚罗林区则是西南林区的重要组成部分[25]，其位于长江上游，当前针叶林正处于快速生长阶段，保护该林区高山森林对于保障长江中下游生态安全、实现“3060”双碳目标具有重要作用[26−27]。2018年国家林业和草原局印发的《松材线虫病防治技术方案》指出要发挥和加强和科技对松材线虫病的支撑与保障作用，通过监测与预测松材线虫适生区与高发区，预防和扭转松材线虫病疫情快速扩散蔓延态势。然而，关于西南青藏高原东缘川西林区松材线虫病的研究还未开展，据川西地区森防站（理县森林病虫害防治检疫站）监测发现川西林区已有云杉花墨天牛分布，且存在松树枯死多点发生的情况，具备松材线虫传播所需条件。随着川西经济发展与基础设施建设，松木产品相关物流往来急剧增加，可能削弱川西高山峡谷区对松材线虫及其传播媒介天牛的天然屏障作用，为松材线虫病向西扩张提供可能，因此开展高山林区松材线虫病入侵风险评估与预测的相关研究工作迫在眉睫[28]。

综上所述，考虑到松材线虫自然入侵与传播至少需要同时具备媒介昆虫（云杉花墨天牛）和寄主（松属乔木）两个条件，本文以理县米亚罗林区为研究对象，基于24个云杉花墨天牛和55个枯死松树地理分布位点监测数据，结合GIS空间技术提取相应的气候环境数据，利用MaxEnt模型构建与预测该林区松材线虫病传播与发生的潜在分布区，解析其受局域环境条件与人为活动的影响程度，补充松材线虫病广域研究未能充分考虑局域环境与生物因子的不足，为川西林业部门精确制定松材线虫病防疫政策与措施提供科学依据。

1. 研究区概况与研究方法

1.1 研究区概况

川西理县（30°54′43″ ~ 31°12′12″N，102°33′46″ ~ 103°30′30″E）地处青藏高原东南缘，位于川西高原岷江上游、邛崃山系东侧，是四川盆地向高原过渡的典型区域，境内多高山峡谷区，海拔落差达4 500 m（最低海拔为1 442 m，最高海拔5 922 m）。气候属季风性气候，县城多年平均气温为10.7 ℃，降雨主要集中在夏季，年平均降水量约为600 mm。由于地貌复杂多样，塑造形成了丰富的植被类型，拥有较高的生物多样性，是我国生物多样性热点地区之一。历史上该地区形成以原始暗针叶林为特征的典型植被，主要针叶优势树种有岷江冷杉（Abies fargesii var. faxoniana）、云杉（Picea asperata）、铁杉（Tsuga chinensis）、岷江柏木（Cupressus chenginana）等[25]。在20世纪50—90年代，大面积的原始森林被国家采伐利用，同期开展以云杉、川西云杉（Picea likiangensis var. rubescens）、紫果云杉（Picea purpurea）和落叶松（Larix gmelinii）为主的人工造林，当前采伐迹地植被正处于全面恢复阶段[29]。

1.2 数据收集与整理

四川理县林业和草原局承担国家级林业有害生物中心站点测报工作，自2005年开始林业病虫害监测、建立数据库、采集病虫害标本等工作。本文从2005—2019年理县森林病虫害防治检疫数据库中提取到云杉花墨天牛24个分布点以及55个枯死松树分布点空间信息，并根据这些分布点的地理位置设置半径为100 m的缓冲区，当同一类型多个地理分布点缓冲区存在相互重叠的情况，则随机保留一个，避免后续模型分析过度拟合（图1）。

图 1 四川理县云杉花墨天牛和枯死松树分布情况

Figure 1. Distributing pattern of Monochamus saltuarius and dead pine trees in Lixian County of Sichuan Province

气候变量是物种生态位模型研究中常用的环境变量。从世界气候数据中心（https://www.worldclim.org/）下载当前国内通用的1990—2019年平均气候数据（共包括19个气候因子），利用ArcGIS软件提取理县气候环境变量的栅格数据。考虑到气候变量对云杉花墨天牛分布的制约以及变量之间的空间自相关性，利用ArcGIS中的SDMtoolbox插件剔除空间自相关性绝对值 ≥ 0.80的气候因子[30]，最终保留13个气候变量。这些变量分别为年均温、平均气温日较差、昼夜温差与年温差比值、温度季节变化、年温度变化范围、年降水量、最湿月降水量、最干月降水量、降水量季节性变异系数、最湿季降水量、最干季降水量、最暖季平均降水量、最冷季平均降水量。此外，地形地貌因子也是影响物种分布的关键因素，利用地理空间数据云数据库（http://www.gscloud.cn/）中提取理县地形因子（高程、坡度、坡向）信息；从全球地理信息公共产品数据库（http://www.globeland30.org/）中提取理县地表覆盖类型信息（包括森林、草地、灌木地、苔原、冰川和永久积雪、耕地、水体和人造地表）。

除了以上自然因子外，人为活动和干扰也是导致云杉花墨天牛和松材线虫传播的主要途径。本文将云杉花墨天牛和枯死松树地理分布点距水源、村庄（居民点）和道路（包括林区道路）的距离作为人为活动强弱的量化指标，相关数据来源于全国地理信息资源目录服务系统（http://www.webmap.cn/）。

1.3 数据处理、模型构建与风险评估

以理县的高程数据（digital elevation model，DEM）为底图利用ArcGIS将所有环境变量的栅格统一边界、坐标系统和像元，将其转换成ASCII格式[31]。将筛选出的24个云杉花墨天牛自然分布点、55个枯死松树的自然分布点和20个生物与非生物变量导入MaxEnt软件中，随机选取80%的云杉花墨天牛和松树自然分布点作为模型训练集，剩余的20%作为测试集，以Logistic格式输出结果，得到每一栅格的发生指数[30]。MaxEnt模型中调控倍率（regularization multiplier）、最大迭代次数（maximum iterations）和收敛阈值（convergence threshold）等都保持默认，对模型进行20次自举法（bootstrap）重复，以平均指数作为最终模型结果[4]。运用Jackknife刀切法确定20个生物与非生物因子的贡献率，模型将自动计算并分别生成了云杉花墨天牛和枯死松树的工作特征曲线（Receiver operation characteristic curve，ROC）以及曲线下面积值（Areas under curves，AUC），用于评估模型预测结果的精准度。当AUC值介于0.9 ~ 1则为优、0.8 ~ 0.9则为良、0.7 ~ 0.8则为中、0.6 ~ 0.7则为差、0.5 ~ 0.0则为模拟失败。

将MaxEnt模型分析生成数据通过ArcGIS转换后，采用自然断点分级法生成阈值后对模型结果进行重分类，以确定云杉花墨天牛适生区和枯死松树易发区[32−33]，将适生区分为高度适生区（HSI ≥ 0.394）、中度适生区（0.105 ≤ HSI < 0.394）和低度适生区（HSI < 0.105）[34]；同时将枯死松树发生区划分为高发区（HSI ≥ 0.410）、易发区（0.121 ≤ HSI < 0.410）和低发区（HSI < 0.121）[35]。最后，耦合云杉花墨天牛适生区和枯死松树发生区数据，利用加权求和方法划分松材线虫病潜在高风险区、中风险区和低风险区[36]。

2. 结果与分析

2.1 模型模拟结果

MaxEnt模型应用ROC分析法模拟发现云杉花墨天牛适生区的AUC均值为0.993，表明模型拟合程度高，模型的预测结果为优，能较准确地预测云杉花墨天牛的适生区（图2a）。枯死松树的自举法（bootstrap）重复得到AUC均值为0.969（图2b），模型拟合程度高，模型的预测结果为优，能较准确地预测枯死松树的发生区。

图 2 MaxEnt模型应用ROC分析法检验云杉花墨天牛适生区和枯死松树发生区的模拟结果

Figure 2. Simulating results of potential target areas for M. saltuarius and occurred areas for dead pines examined using ROC curve in the MaxEnt models

2.2 云杉花墨天牛适生区与枯死松树发生区的影响因子

MaxEnt模型分析发现距居民点距离、年均温、最湿月降雨量和最湿季度降雨量是云杉花墨天牛分布区的主要影响因子，它们对模型的贡献率达到79%；而降水量季节性变异系数、最干季降水量和最冷季平均降水量对云杉花墨天牛适生区模型模拟的贡献率为0%。枯死松树发生区的主要影响因子有人为干扰因子（距居民点距离、距最近道路距离）、最湿季降水量以及年均温，综合贡献率为76.3%，是枯死松树分布的主要影响因子。平均日较差、温度季节变化和最干季降雨量的贡献率为0，对枯死松树发生区无影响（表1）。

表 1 四川理县生物与非生物因子对云杉花墨天牛适生区和枯死松树发生区MaxEnt模拟的贡献率

Table 1. Contribution rates of biotic and abiotic factors to MaxEnt model simulations of suitable areas for M. saltuarius and occurred areas for dead pine trees in Lixian County of Sichuan Province

影响因子
Impact factor贡献率 Contribution rate/%影响因子
Impact factor贡献率 Contribution rate/% 云杉花墨天牛
M. saltuarius枯死松树
Dead pine tree云杉花墨天牛
M. saltuarius枯死松树
Dead pine tree 年均温
Annual mean temperature 20.40 5.50 最干季降水量
Precipitation of the driest month 0.00 0.00 平均气温日较差
Mean temperature diurnal range 0.10 0.00 最暖季平均降水量
Mean precipitation of the warmest season 0.40 2.90 昼夜温差与年温差比值
Isothermality 0.50 1.90 最冷季平均降水量
Mean precipitation of the coldest season 0.00 0.10 温度季节变化 Temperature seasonality 1.00 0.00 地表覆盖情况
Land surface coverage condition 1.20 1.30 年温度变化范围
Annual temperature variation range 1.00 1.10 海拔 Elevation 7.00 2.10 年降水量 Annual precipitation 0.30 5.30 坡度 Slope 1.40 1.40 最湿月降水量
Precipitation of the wettest month 19.00 4.20 坡向 Aspect 3.90 0.50 最干月降水量
Precipitation of the driest month 0.10 2.10 距最近居民点距离
Distance to the nearest residential point 28.60 29.00 降水量季节性变异系数
Variation coefficient of precipitation seasonality 0.00 0.20 距最近道路距离
Distance to the nearest road 2.00 24.90 最湿季降水量
Precipitation of the wettest month 11.00 16.90 距水源距离
Distance to the water source 2.00 0.70

从地理与空间分布看，云杉花墨天牛适生区主要集中在海拔1 900 ~ 2 300 m、坡度为34° ~ 65°的区域（图3a、c）。当年均温为7.8 ~ 10.1 ℃、最湿月降水量介于125 ~ 130 mm、最湿季降水量介于345 ~ 358 mm且附近1.5 km内有居民点时，云杉花墨天牛危害最容易发生（图4a、c、e；图3e）。枯死松树发生区主要集中在海拔1 900 ~ 3 000 m、坡度40° ~ 65°的山地，当年均温为5.7 ~ 10.0 ℃、最湿季降水量介于340 ~ 360 mm、且距居民点2 km内、距道路1 km内时极易发生松树枯死（图3b、d、f；图4b、d、f）。

图 3 四川理县云杉花墨天牛与枯死松树的预测概率对海拔、坡度与居民点的响应曲线

红线表明MaxEnt模型中云杉花墨天牛或枯死松树发生概率随各因子的变化趋势，蓝色部分为95%置信区间。下同。The red line indicates variation trend of occurrence probability of M. saltuarius or dead pine trees in the MaxEnt models as response to abiotic or biotic factors; the blue part represents 95% confidence intervals. The same below.

Figure 3. Response curves of predicted probability for M. saltuarius and dead pine trees to elevation, slope and residential variables in Lixian County, Sichuan Province

图 4 四川理县云杉花墨天牛与枯死松树的预测概率对气候因子与人为活动的响应曲线

Figure 4. Response curves of predicted probability for M. saltuarius and dead pine trees to climatic and anthropogenic variables in Lixian County, Sichuan Province

2.3 云杉花墨天牛与枯死松树分布区预测及松材线虫发生风险评估

在理县县域内利用MaxEnt模型模拟预测云杉花墨天牛适生区，发现高适生区面积约为5 859 hm2，占针叶林总面积3.92%，主要以下孟乡、通化乡与杂谷脑镇的风险尤高。云杉花墨天牛中适生区面积约为14 506 hm2，占针叶林总面积的9.71%，主要集中在上下孟乡、通化、蒲溪两乡以及甘堡乡到朴头乡公路沿线。整体看理县县域多为云杉花墨天牛的低适生区（图5）。

图 5 四川理县云杉花墨天牛适生区的MaxEnt模型预测

Figure 5. Prediction of MaxEnt model for potential target areas of M. saltuarius in Lixian County, Sichuan Province

MaxEnt模型模拟发现枯死松树高发区和易发区面积分别为23 721 hm2和39 985 hm2，各占理县针叶林总面积的15.87%和26.75%，且在理县全县11个乡镇均有分布，多集中在公路沿线两侧（图6）。

图 6 四川理县枯死松树发生区的MaxEnt模型预测

Figure 6. Prediction of MaxEnt model for potential occurred areas of dead pines in Lixian County, Sichuan Province

叠加云杉花墨天牛适生区与枯死松树发生区数据，赋予其相同权重的MaxEnt模型模拟即可预测川西理县松材线虫病发生的风险等级与规模。结果发现川西理县松材线虫潜在高风险区面积为10 616 hm2，呈带状分布于县域内的下孟乡、上孟乡、蒲溪乡、甘堡乡和杂谷脑镇等10个乡镇，约占理县针叶林总面积的7.1%。松材线虫病中风险区面积为55 445 hm2，主要集中在朴头乡到米亚罗镇公路（道路）沿线地区，约占针叶林总面积的37.10%（图7）。

图 7 四川理县耦合云杉花墨天牛与枯死松树的松材线虫病潜在风险MaxEnt模型预测

Figure 7. Prediction of MaxEnt model for potential risk of pine wilt disease based on integration of M. saltuarius and dead pines in Lixian County, Sichuan Province

3. 讨　　论

自然条件下松材线虫扩散与传播受媒介昆虫与寄主和气候环境等共同驱动[37]。Hirata等[9]对全球松材线虫病在未来气候变化背景下的潜在分布进行预测，指出四川是松材线虫病入侵的适生区。但是该研究仅从省级报道有松材线虫病例而将该省作为一个松材线虫病发生的最小单元，这在全球尺度分析可行，但从局域尺度则显然不甚合理，并未考虑四川省域内巨大的海拔落差以及丰富的高山峡谷植被类型，同时未考虑其他局域环境因子的影响[25]。因此，从县域尺度开展松材线虫病潜在发生规模的风险评估研究对于疫情精准防控更加有参考意义[38]。本文同时耦合云杉花墨天牛及其潜在寄主枯死松树的多年监测数据，首次将气候变量、地表覆盖类型、地形变量、人为因素等用于MaxEnt模型构建，研究和评估松材线虫病在川西米亚罗林区发生的潜在风险及其驱动因子，从县域尺度预测松材线虫病高发区，这将有助于提前做好疫情防控[36]。

松材线虫病媒介昆虫天牛的天然传播能力有限（通常仅局限于数公里内），大部分疫情传播还是以人为传播为主，表现在人为调运带疫的苗木、松材、松木包装箱、电缆盘、光缆盘及松木制品等导致松材线虫远距离传播[1]。有证据表明，我国超过75%的松材线虫病疫区是人为传播造成的[2]。因此，松材线虫病防控的第一要务是监测与检疫，在天然林区或人工难以达到的地区，主要以监测为主，一旦发现病木应第一时间上报并做好采样检测工作；在人为活动频繁或市、县、乡镇等物流流通主干道，做好检验检疫工作，做到精准防控[39]。可能得益于青藏高原高寒气候与川西高山峡谷区的生态屏障作用，在一定程度上阻挡了松材线虫病的扩散与传播，目前未发现松材线虫入侵。但是，四川盆地及其周边43县已然沦为松材线虫病疫区，理县紧邻四川盆地，随着该地区旅游开发和基础设施建设（例如高速公路、水电站等）以及乡村居民楼房建设等人为活动，很大程度上增加了川西林区松材线虫病传播与扩散的风险[40]。

本研究通过MaxEnt模型模拟发现距居民点距离远近、年均温和最湿季降水量是云杉花墨天牛适生区和枯死松树发生区的关键影响因子。具体地讲，距居民点约1.5 km内，年均温为7.8 ~ 10.1℃，最湿季降水量为348 ~ 358 mm时云杉花墨天牛或枯死松树发生的概率最高（图5、6）。本文相关结果与云南省云杉花墨天牛适生性空间模拟、松材线虫病诱导的松树枯死等预测结论相一致[24,41]。不同的是，本文还发现云杉花墨天牛适生区和枯死松树分布区多位于坡度34° ~ 65°的阴坡半阴坡（图3），这可能与川西地区多高山峡谷导致多陡坡地形地貌有关，陡坡地段土壤水分与养分含量通常低于平缓地段的相应含量[25]，生长在陡坡地段的松树长势常弱于其在缓坡地段的长势，因而具有更高的死亡风险，这是导致川西地区林木枯死主要集中在陡坡（35°以上）的原因。此外，松树衰弱木或濒死木是云杉花墨天牛产卵到羽化前生长发育阶段的关键宿主（羽化后飞往健康大树补充营养）[2]，所以该地区也是云杉花墨天牛的适生区。将云杉花墨天牛适生区和枯死松树分布区数据进行加权求和，模拟预测评估松材线虫病发生的潜在风险，发现川西理县松材线虫入侵高风险区面积为10 616 hm2（约占全县针叶林总面积7.1%，图7），并标识出这些高风险地区的具体位置，为当地林业主管部门制定相应的措施开展精确的松材线虫病防范提供科学依据。

松材线虫病常具有潜伏期，若等到林木表现出明显的征兆，则松材线虫已经具有一定的传播规模[2]，本研究发现松材线虫病潜在的高发区主要集中在公路和居民区附近，最科学的防范办法是从源头（外源输入）阻止松材线虫的媒介昆虫和松科病木流入，因此要严把检疫关，对经过该县的涉木（木材、木制品或外包装）物流进行严格的检验检疫，同时对居民点开展科普宣传，对发现枯死松木应及时上报林草部门并彻底清理枯死木，避免私自利用（作为薪材等）引起的人为传播[42]。川西理县亚高山林区当前具有大面积的针叶人工纯林，同时已监测发现有松材线虫媒介昆虫（云杉花墨天牛和小灰长角天牛），一方面建议加强针叶林病虫害防治，提高林木健康程度以及森林稳定性，另一方面建议开展天牛诱捕监测和防治工作。此外有研究指出，随着气候变化与人类活动的加剧，松材线虫病的遗传多样性逐渐提高，传播风险进一步提升，防控难度增大[43−44]。

值得注意，本研究也有一定的局限性。本文仅使用云杉花墨天牛分布点与枯死松树发生点数据，缺乏将这些发生点所在森林群落结构（如树种组成、林龄、种群密度等）等相关群落生态学指标纳入分析[29]；另外本研究未能明确哪些枯死松树点位是由云杉花墨天牛引起的，这会对加权预测松材线虫病发生风险产生一定的误差。虽然当前模型模拟精度达到统计要求，但未能从生物学机理上解释松材线虫发生潜在风险的驱动机制[45]。因此，下一步建议将生态系统结构（物种、遗传与生态系统多样性）与功能（生物量与生产力、稳定性与抗性）量化并纳入MaxEnt模型当中，从机理和过程上提高松材线虫病风险评估与预测的准确性，为松材线虫防治提供科学依据[46]。

4. 结　　论

本文利用四川理县云杉花墨天牛和枯死松树的地理分布点、气候环境数据、地形地貌数据以及人为活动数据，基于最大熵（MaxEnt）模型对云杉花墨天牛和枯死松树的空间分布进行了预测，以期耦合松材线虫入侵和传播所需中间媒介和寄主来预测评估川西理县林区松材线虫病发生的潜在风险。模型预测发现沿道路两侧以及居民点（区）附近1.5 km内是云杉花墨天牛与枯死松树发生的集中区，年均温和最湿季降水量也将显著影响云杉花墨天牛和枯死松树发生概率，这将在一定程度上决定松材线虫病潜在传播风险。模型模拟预测理县松材线虫病感染高风险区面积为10 616 hm2，约占全县针叶林总面积7.1%。鉴于人为活动与物流运输等对理县松材线虫病入侵与传播的较大潜在风险，建议以林草部门为主体，协调公安、交通、卫生等部门严把道路沿线检验检疫关口，同时加强林区病虫害与气候环境变化监测，科学精准防控松材线虫病入侵川西林区。