基于卫星遥感的森林火灾风险预警研究

李哲全，张 贵，谭三清，吴 鑫

（中南林业科技大学 林学院，湖南 长沙 410004）

卫星遥感技术提供了地球表面的全球观测信息，使森林火灾风险预警不再受地面气象站数据的限制，能够更快速的获得森林火灾风险因子信息，且卫星遥感数据具有高实时性，可满足森林火灾风险预警的时效性。

Cipriani 等[1]将高程、坡向、坡度、土地利用类型、道路和建筑物等因子进行危险等级划分，通过火灾结构指数（SFI）进行叠加，生成研究区的森林火灾风险图，进而确定森林火灾发生概率最高的区域；刘琳[2]使用了NDVI、TVDI、植被覆盖类型3 个因子建立了火险综合预报指数模型对重庆市进行森林火灾监测预警；王双等[3]采用Logistic 回归分析方法，地表温度数据、植被指数数据、气象数据、人文数据、植被类型数据、森林火灾监测数据等，构建森林火灾风险概率模型，并对湖南省森林火灾风险等级进行划分；贾旭[4]将各个火险因子依据是否在短时间内发生变化，分为可变因子与稳定因子两大类，并使用层次分析法来分开计算因子的权重，后根据权重进行计算生成内蒙古的火险区划；Adab 等[5]利用地形因子、植被湿度因子、人类活动因子建立了火灾结构指数、火灾风险指数与混合火灾指数对伊朗南部地区进行了森林火灾监测预警；Babu 等[6]提出了一种新的静态火险指数（SFDI），该指数利用土地覆盖类型数据（MCD12Q1）、数字高程模型数据（DEM）和热点数据（MOD14）对印度某一地区进行统计分析计算SFDI。马文苑等[7]使用2010—2017年卫星监测热点数据，基于逻辑斯蒂模型和随机森林模型分析气象、地形、植被和人类活动对山西省林火发生的影响，选取山西省主要林火驱动因子，建立林火发生概率模型。

综上所述，目前大部分国内外学者将气象因子、地形因子、植被状况划分为若干等级，对应不同等级的森林火灾风险区划，依次建立了森林火险预警模型，但由此得到的预警结果是一个静态的结果，而本文利用卫星遥感数据与DEM 数据提取森林火灾风险预警因子并构建森林火灾风险预警指标体系与模型，可进行森林火灾风险的动态预警以满足森林防火工作的时效性。

1 材料与方法

1.1 数据源

研究中使用的数据为MOD02 数据、DEM 数字高程数据、2009—2019年广东省森林火灾卫星热点数据。

MOD02 是MODIS 经校正后的原始数据，分为3 种类型：MOD02QKM 分辨率为250 m，包含了1、2 波段的数据；MOD02HKM 分辨率为500 m，包含了1～7 波段的数据，其中1、2 波段的数据被重采样至500 m 的分辨率；MOD021KM 分辨率为1 km，包含了MODIS 36 个波段的信息，其中1～7 波段的数据都被重采样至1 km 分辨率。

DEM 数据包含了高程在内的各种地形因子如坡度、坡向、坡度变化率等因子在内的线性和非线性组合的空间分布。

1.2 研究方法

1.2.1 森林火灾风险预警指标体系的构建

国内外研究成果表明，森林火灾的发生主要与地表温度、地表植被状况、地形地貌、降水量、可燃物含水率等因素有密切关系。例如NDVI 能有效地表征地面可燃物状态；TVDI 可以表征地表温度与土壤湿度两个森林火灾风险因子；NDII7 反映可燃物含水率；森林火灾的发生与高程、坡度、坡向等地形因子存在不同程度的关系；发生森林火灾的历史状况也与森林火灾风险存在密切关系。利用遥感技术对森林火灾进行风险预警具有高时效性，森林火灾风险程度是多因子耦合的结果。考虑到遥感数据提取森林火灾风险预警因子的可获得性，选择NDVI、NDII7、TVDI、高程、坡度、坡向、森林火灾历史数据等作为森林火灾风险预警因子，构建森林火灾风险预警模型，实现森林火灾的实时预警。

NDVI 是反映地表植被覆盖信息的一种植被指数，能够反映植被的绿度、光合作用强度，现主要被用于植被的动态监测。其值的计算通过近红外波段（841～876 nm）与红光波段（620～670 nm）两个通道反射率之差除以它们的和，计算公式如下：

式中：Rnir和Rred是指经过大气校正的MODIS 影像数据中的近红外波段和红光波段的地物反射率。

湿度，即可燃物含水率。覃先林[8]的研究工作己经论证了短波红外波段与可燃物含水率呈现负相关关系，并且与其他植被指数相比，第7 波段的归一化红外指数（Normalized difference infrared index7，NDII7）对可燃物含水率更为敏感（图1），从图1可看出植被的NDVI、NDII7 值和NDWI 均大于裸地，且植被和裸地的NDII7 值相差大于两者的NDVI 值之差。因此，将NDII7作为可燃物含水率的指标之一，该指数由MODIS的两个波段计算，其公式如下：

图1 不同地物类型在不同波段和不同含水率指数的变化曲线Fig.1 Variation curves of different surface feature types in different wave bands and different water content indexes

式中：ρi为第i波段的反射率。

TVDI 是以地表温度（Ts）和植被指数（NDVI）构建的二维空间的Ts-NDVI 到的，它主要用于旱情的监测、土壤湿度的反演和反应植被缺水状况。在文中用于表示温度与森林火灾之间的关系，因为TVDI 值越高，表示区域内温度与植被指数之比越大，越容易发生森林火灾，其计算公式如下：

式中：Ts为地表温度（℃）；Tmin为相应NDVI 下最低陆地表面温度（℃）；Tmax为相应NDVI 下最高陆地表面温度（℃），即所对应的湿边和干边方程；a、c、b、d为干湿边方程的回归系数，分别为干湿边方程的截距和斜率[9-10]。

高程、坡度、坡向可基于DEM 数据进行提取；历史森林火灾因子则通过2009～2019年森林火灾卫星热点数据计算获得。

1.2.2 森林火灾风险预警模型

自然灾害风险指数。根据自然灾害的形成与发展机理，自然灾害风险是由致灾因子的危险性与受灾体的暴露性和脆弱性综合作用形成的，可通过自然灾害风险指数来表现风险程度。

自然灾害风险指数的计算公式如下：

式中：xi为第i风险项的取值；wi为第i风险项的权重值。

森林火灾动态危险性指数。为了建立一种简单易操作的森林火灾动态危险性评估模型，可提出如下假设:1）可燃物含水率可预测森林火灾发生的可能性；2）植被覆盖率可作用一个辅助指标提供相应的参考；3）地表温度会影响可燃物含水率的大小，因此地表温度可作为其中的一个动态因素来预测森林火灾的发生；4）本模型不考虑风所带来的影响。使用基于李晓恋[11]研究提出的森林火灾动态危险性指数（式7）进行森林火灾风险预警，公式如下：

式中：FDDI为森林火灾动态危险性指数，NDII7为MODIS 第7 波段归一化红外指数，用于表示可燃物含水率，T为地表温度，Pv为植被覆盖率。

1.2.3 地表温度反演

劈窗算法请参考文献[12]。

2 结果与分析

2.1 地表温度反演

本研究利用大气水汽含量来估算地表温度反演中所需的大气透过率。通过MODIS 影像的第2与19 波段的反射率来估算得到大气水汽含量，公式如下：

式中：W为大气水汽含量，ρi为第i波段的反射率，α、β为常量，依据地物类型来取值α=0.020，β=0.651。

将W代入公式分别计算第31 与32 波段的大气透过率。

依据上述关系，对大气透过率进行估算，结果如图2所示。

图2 研究区第31、32 波段的大气透过率（2021年2月14日MODIS 数据）Fig.2 Atmospheric transmittance of the 31st and 32nd wavebands in the study area(MODIS data of February 14,2021)

地表比辐射率是进行地表温度遥感反演的基本参数，主要取决于地表覆盖物质的结构。计算公式如下：

式中：εi为第i波段的地表比辐射率，Pv为植被覆盖度，可通过NDVI 来计算，Rv、Rs为植被与裸土的辐射比率。εiv、εis为植被和裸土在第i波段的地表比辐射率，赋值为常数ε31v=0.986 72、ε32v=0.989 90、ε31s=0.967 67、ε32s=0.977 90。

Rv与Rs不仅取决于温度变化，而且取决于植被覆盖度，且后者的影响更大。因此，可建立它们与植被覆盖度之间的关系来进行估计。计算公式如下：

植被覆盖度Pv主要是通过NDVI 来估计：

式中：NDVIv、NDVIs为高植被覆盖像元与裸土像元的NDVI 值，通常取常数NDVIv=0.9、NDVIs=0.15。

依据上述关系，对地表比辐射率进行估算，结果如图3显示。

图3 研究区第31、32 波段的地表比辐射率（2021年2月14日MODIS 数据）Fig.3 Surface emissivity of the 31st and 32st wavebands in the study area(MODIS data of February 14,2021)

dε为估计相互作用校正项，由于热辐射相互作用在植被与裸土各占一半时达到最大，所以根据Sobrino[13]的研究提出的经验公式来计算，公式为：

依据劈窗算法得到广东省的地表温度反演，见图4。

图4 2021年2月14日广东省地表温度反演结果Fig.4 Retrieval results of land surface temperature in Guangdong province on February 14,2021

2.2 森林火灾风险预警指标计算

基于式（1），提取MOD02 中对应波段，计算得到区域NDVI 值。

基于式（2），提取MOD02 中对应波段，计算得到区域NDII7 值。

图5 2021年2月14日广东省NDVI 值分布Fig.5 Distribution of NDVI in Guangdong province on February 14,2021

图6 2021年2月14日广东省NDII7 值分布Fig.6 Distribution of NDII7 in Guangdong province on February 14,2021

基于式（3）～（5），使用反演的地表温度与NDVI 计算得到区域TVDI 值。

图7 2021年2月14日广东省TVDI 值分布Fig.7 Distribution of TVDI in Guangdong province on February 14,2021

高程。高程与森林火灾的发生和发展有明显的关系。对于高程的划分，中科院1987年规范了基于全国区域的高程划分，其范围为＜1 000 m，1 000～3 500 m，3 500～5 000 m，＞5 000 m，该划分对于广东省来说间距过大，不能很好的反应广东省的地形特征。经统计，广东省实际的高程分布区间在[-57～1 874 m]间，因此参考周成虎等[14-15]、哈凯等[16]的研究成果，以200 m为一级进行划分，将广东省的高程分为[-57～200 m)，[200～400 m)，[400～600 m)，[600～800 m)，[800～1 000 m)，[1 000～1 200 m)，[1 200～1 400 m)，[1 400～1 874 m]8 个划分区，并依据分区对广东省的森林火灾热点数量进行划分统计。表1为基于高程的森林火灾热点数分布表。

表1 基于高程的森林火灾热点数分布Table 1 Distribution table of forest fire hot spot number based on elevation

坡度。坡度使用DEM 数据在ArcMap10.4 生成坡度图，并根据中国科学院地理研究所编制的中国1:1 000 000 地貌图制图规范，将研究区坡度划分为[0°～5°)、[5°～10°)、[10°～15°)、[15°～20°)、[20°～25°)、[25°～30°)、[30°～35°)、[35°～40°)、[40°～90°)9 个类别，再根据2009年至2019年的森林火灾热点数据进行统计分析，得出不同坡度的森林火灾热点分布频率。经统计得，坡度在[0°～25°]之间，广东省森林火灾热点数量很高，这是由于坡度小的地方人类的活动越多，失火的概率越大，反之亦然。但是当森林火灾发生时，火势会向坡度高的地方蔓延。表2为基于坡度的森林火灾热点数分布表。

表2 基于坡度的森林火灾热点数分布Table 2 Distribution table of number of forest fire hot spots based on slope

坡向。在森林火灾预警研究中，坡向可作为日照时长的指标。表3为基于坡向的森林火灾热点数分布表。

表3 基于坡向的森林火灾热点数分布Table 3 Distribution table of number of forest fire hot spots based on aspect

历史森林火灾因子。基于2009—2019 森林火灾卫星热点数据，生成区域的历史森林火灾发生区域图，将发生过森林火灾的区域赋值为1，未发生森林火灾的区域赋值为0（图8）。

图8 广东省2009—2019年森林火灾发生区域Fig.8 Forest fire occurrence area in Guangdong province from 2009 to 2019

2.3 森林火灾风险预警

选择TVDI、NDVI、NDII7、高程、坡度、坡向与森林火灾历史数据作为预警因子。

传统森林火灾预警模型是由于早期对各预警因子间关系的研究较少而进行的一种简单模型计算。随着对森林火灾的研究加深，各预警因子间的关系逐渐清晰，并得到量化，且由于传统模型对因子的直接叠加过于粗糙，该模型的预警精度已不能保证。

传统森林火灾风险预警模型基于已分类并进行归一化处理的预警因子进行直接叠加得到森林火灾风险分布图（图9）。

将森林火灾动态危险性指数公式进行适宜于本研究预警方法的改进（式17），公式如下：

利用SPSS19.0 软件对高程、坡度、坡向进行显著性检验，得到高程、坡度、坡向的权重（表4）。最终构建出森林火灾风险预警模型。公式如下：

表4 地形因子权重Table 4 Weight of terrain factors

基于改进后的森林火灾动态危险性指数模型（式18），叠加预警因子的栅格图层，在ArcMap10.4 软件中进行栅格计算，得到森林火灾风险分布图（图10）。

由图10可知，在2021年2月15日当天，广东省各市均存在高森林火灾风险区，其中广州市、佛山市、中山市、东莞市、汕头市、深圳市的高森林火灾风险区面积较大。

对比图9与图10。在基于传统森林火灾预警模型得到的预警图（图9）中，中高森林火灾风险区域占比面积偏大，对于具体森林火灾风险程度的危险等级变化较小，且对比基于改进后的森林火灾动态危险性指数模型计算得到的预警图（图10），对中高风险区的描绘较粗略，这是由于传统预警模型是将各因子进行了归一化处理后再进行因子叠加运算的结果。

图9 2021年2月15日广东省森林火灾风险分布（基于传统森林火灾风险预警模型）Fig.9 Distribution of forest fire risk in Guangdong province on February 15,2021 (Based on the traditional forest fire risk early warning model)

对改进后的森林火灾动态危险性指数模型进行ROC 检验。利用模型计算得出的预测概率，在SPSS19.0 软件中绘制ROC 曲线并对AUC 值进行计算，得到AUC 值为1.0，大于0.7，证明该模型能够有效区分森林火灾高风险区与低风险区，模型拟合效果较好（图10）。

图10 2021年2月15日广东省森林火灾风险分布（基于2021年2月14日MODIS 数据）Fig.10 Distribution of forest frie risk in Guangdong province on February 15,2021 (Based on MODIS data of February 14,2021)

3 结论与讨论

3.1 结 论

1）从森林火灾热点数量与地形因子之间关系看出，2009—2019年广东省森林火灾热点多集中于坡度低于25°，海拔低于300 m 的区域；坡向因子对区域内森林火灾热点个数影响较弱，在各坡向上，森林火灾热点分布较均匀。

2）使用TVDI、NDVI、NDII7、高程、坡度、坡向与森林火灾历史数据因子拟合的森林火灾动态危险性指数FDDI 来反映森林火灾发生的风险概率，结果值越大，表明NDVI 值越小、TVDI 值越小、NDII7 值越小，森林火灾风险程度大；结果值越小，表明NDVI 值越大、TVDI 值越大、NDII7 值越大，森林火灾风险程度越小。

图11 ROC 曲线Fig.11 ROC Curve

3.2 讨 论

1）传统森林火灾预警方法多使用地面气象数据进行监测预警，而地面气象数据受气象站位置的制约，当地面气象站距离越小时，精度越高；反之亦然。且这些预警方法在操作层面上是以点为对象，在处理时需进行空间插值，会导致不同的插值结果，泛用性不高。本文使用了遥感数据可实时分析监测预警森林火灾风险，使林业管理部门可以高效的做出预防措施。

2）使用MODIS 数据得到的最终预警图为1 km 分辨率，对于小范围的森林火灾可能预警效果不佳，会影响预警精度。因此在后续的研究中可结合使用高空间分辨率的影像进行预警研究，以提升风险预警的准确性。

3）该预警模型的优点是简单易执行，能够实现实时的森林火灾风险预警，为防火工作提供直接的信息，并结合了地形因子，解决了原森林火灾动态危险性指数模型的局限性，但仅仅加入地形因子还是不够的，森林火灾的发生还涉及其他预警因子，例如人为活动、雷击等火源因子、风速风向等因子。下一步研究方向为进行诸如人为活动、风速风向等预警因子的遥感提取研究，以满足森林火灾风向预警的时效性，并以此建立更加完善的基于卫星遥感的森林火灾风险预警模型，提高森林火灾动态预警的预警精度，合理安排森林防火工作，降低损失。