Remote sensing technology for rapid extraction of burned areas and ecosystem environmental assessment

Study area

Muli County is located on the southeastern edge of the Qinghai-Tibet Plateau (longitude 100°03′–101°41′ east, latitude 27°40′–29°10′ north), at the boundary between the Qinghai-Tibet Plateau and the Yunnan-Guizhou Plateau. Its main landform features are mountains, canyons and mountain plains. The terrain is complex, and the altitude difference a are enormous. It belongs to a typical alpine and canyon area. Muli County has a subtropical monsoon climate, with rainy summers and less rain in winter. Muli is rich in forest resources, with a forest coverage rate of 67.3%, including Picea asperata, Quercus semicarpifolia, Pinus densata, etc., and more than 10 kinds of national first and second-class protected tree species. The volume of standing trees is 117 million m², accounting for 10% of Sichuan Province and 1% of the whole country (Liu, 2014), and the vegetation is vertically distributed.

The annual average temperature is about 14.0 °C, the temperature difference between day and night is significant, and solar radiation is powerful. The four seasons are not pronounced, and the dry and wet seasons are distinct. The dry season is from December to March of the following year., the climate is dry, with little precipitation, long sunshine hours, and uneven spatial and temporal distribution of precipitation, often suffering from natural disasters such as droughts, debris flows and floods (Li et al., 2021d, 2021c). The wet season is from June to September, the precipitation is concentrated, the rainfall is abundant, the daily temperature difference is slight, and the climate is suitable. The forest in the study area is widely distributed, the terrain is complex and changeable, and the forest fire prevention situation is complicated. Forest fires have occurred in the region for three consecutive years, causing massive casualties, economic losses and loss of forest resources. The fire areas are shown in Fig. 2:

Data

Reference data

The reference data are the fire survey data of the Sichuan Forestry Survey and Design Institute. It was vector data that was visually interpreted from sentinel-2 remote sensing images and corrected using field verification. The reference data mainly includes the forest loss area (heavily burned area), the forest where fire passed area (forest loss is not apparent, mildly burned area), and various survey data details are shown in Table 3.

dNBR computation

The NBR (Normalized Burn Ratio) is a remote sensing index obtained by replacing the red band in the NDVI (Normalized Difference Vegetation Index) equation with the SWIR band. Since the reflectance change in the NIR in the burned areas is more significant than that in the SWIR compared with the non-burned areas, it can be used as an essential indicator factor to extract burned areas. The NBR values range from 1 to −1 and negatively correlate with forest fire intensity. However, the NBR index has some advantages in the extraction of burned areas, but using the mono-temporal NBR is still affected by low reflectivity areas such as water bodies, which is prone to significant errors. The dNBR (delta Normalized Bun Ratio) (Chen, Loboda & Hall, 2020) utilized a multi-temporal that is the difference of NBR index in the post-fire images and pre-fire images, and it can effectively reduce the influence of factors such as clouds and aerosols compared with the NBR index.

Moreover, it can be seen that the NBR reflectance of burned areas showed an increase after the occurrence of forest fires. In contrast, the SWIR reflectance of burned areas is extremely sensitive to water bodies, and the vegetation humidity in the forest decreases after burning, so the higher the intensity of the forest fire, the more significant the decrease in humidity within the forest and the decrease in humidity then leads to the decrease the SWIR reflectance. Finally, the NBR index shows a decrease. Therefore, the dNBR index of the burned area images is significantly higher than the low reflectance areas such as water bodies. In summary, the dNBR index of multi-temporal can be used to distinguish the burned areas. The NBR’s and the dNBR’s mathematical formulas are as Eqs. (1) and (2):

(1) $$NBR={\displaystyle \frac{\rho {}_{nir}-{\rho }_{swir}}{\rho {}_{nir}+{\rho }_{swir}}}$$

(2) $$dNBR=NB{R}_{pre}-NB{R}_{post}$$where ${\rho }_{nir}$ represents the NIR band reflectance of sentinel-2 (Band 8); ${\rho }_{swir}$ represents the SWIR band reflectance of sentinel-2 (Band 11); $NB{R}_{pre}$ represents the pre-fire NBR index; $NB{R}_{post}$ represents the post-fire NBR index.

OTSU thresholding for image segmentation

Threshold segmentation is one of the most commonly used image segmentation methods; it is simple and effective. Moreover, scholars have made researched the methods of automatically obtaining thresholds in image segmentation work. The OTSU algorithm, also known as the maximum interclass variance method, is an automatic algorithm for the global binarization of clustered images proposed by Japanese scholar Zhenyuki Otsu in 1979 (Otsu, 1979). It can obtain more ideal segmentation results, has the characteristics of high segmentation efficiency and broad application scenarios, the results are stable, and the segmentation quality can be guaranteed, so it has been widely used (Fan et al., 2016). It works as a discrete simulation of the one-dimensional Fisher discriminant. In computer vision and image processing work, this algorithm assumes that an image contains only two classes of pixels (background pixels and target pixels) after bimodal histogram processing and then minimizes their intra-class variance and maximizes inter-class variance by calculating the best threshold that can distinguish between the two classes of pixels.

The OTSU algorithm is computed as follows: Assuming that an image has several pixels, the probability of occurrence in each gray-level pixel $P_i$ is Eq. (3):

(3) $$P_i=n_i/N,(i=0,1,2,3\dots L-1),{\sum }_{i=0}^{L-1}{P}_i=1$$where the grayscale $i$ takes a range $\left[0,L-1\right]$, the number of pixels in the grayscale $i$ level is $n_i$. The image was divided into background pixels $C_0$ and target pixels $C_1$ using thresholding, $C_0$ consisting of pixels with grayscale values $\left[0,k\right]$ and $C_1$ pixels with grayscale values $\left[t+1,L-1\right]$. Then, the probability of occurrence of each gray ${\mu }_t$ is Eq. (4):

(4) $${\mu }_t={\sum }_{i=0}^{L-1}i{P}_i$$

The occurrence probability of $C_0$ and $C_1$ pixels are Eqs. (5) and (6):

(5) $${\omega }_0={\sum }_{i=0}^{L-1}{P}_i$$

(6) $${\omega }_1={\sum }_{i=t+1}^{L-1}{P}_i=1-{\omega }_0$$

The average gray ${\mu }_0$ and ${\mu }_i$ are Eqs. (7)–(9):

(7) $${\mu }_0={\sum }_{i=0}^{t}i{P}_i/{\omega }_0$$

(8) $${\mu }_1={\sum }_{i=t+1}^{L-1}i{P}_i/{\omega }_1$$

(9) $${\mu }_t={\omega }_0{\mu }_0+{\omega }_1{\mu }_1$$

The interclass variance ${\delta }_t^2$ is Eq. (10):

(10) $${\delta }_t^2={\omega }_0({\mu }_0-{\mu }_t{)}^2+{\omega }_1({\mu }_1-{\mu }_t{)}^2={\omega }_0{\omega }_1({\mu }_0-{\mu }_1{)}^2$$

The optimal threshold to distinguish between the two classes of pixels is the one that takes the value in the interval when the maximum value is reached.

After a forest fire, there is a sharp increase in the maximum interclass variance within a reasonable image element window, and the OTSU algorithm was applied to extract the dNBR threshold automatically. It can be used as an essential parameter in active fire detection.

Building decision trees

Decision tree learning is an essential method in the field of statistics, data mining and machine learning, mainly used as a prediction model to analyze test and target variables in a round-robin manner and finally to predict the category of the samples. The training process is to obtain rules for classification by summarizing expert experience, simple statistical analysis of mathematical models, induction, etc. A training set is divided into several subsets, and the recursion is repeated in the resulting subsets. This top-down decision tree induction method is one of the greedy algorithms and is one of the most commonly used training methods, which can be applied in predicting continuous variables and remote sensing classification.

In this work, we constructed a decision tree classification model with the dNBR based on the sensitive feature analysis of burned areas. First, we divided burned areas into the mildly burned areas (fire traced area) and the heavily burned area (forest loss area) and calculated the first dNBR threshold based on OTSU, which was used to distinguish the heavily burned area (forest loss area). Then a second dNBR threshold is calculated in the rest area based on OTSU, which is used to distinguish between the mildly burned areas (fire traced area) as well as the rest area, and the decision tree classification model is shown in Fig. 3.

Accuracy verification

The accuracy verification of burned areas identification used the following calculation formula (Eq. (11)):

(11) $$\mathrm{P}=(1)-\left|{\displaystyle \frac{\mathrm{a}-\mathrm{b}}{\mathrm{b}}}\right|\times 100\mathrm{\%}$$

P represents the accuracy of the method, a represents the number of image elements of the identified burning area, and b represents the number of image elements visually interpreted after the field survey respectively.

To better verify the experimental results, we used the kappa coefficient to verify the accuracy. Its equation is as follows.

(12) $$\mathrm{K}\mathrm{a}\mathrm{p}\mathrm{p}\mathrm{a}={\displaystyle \frac{{\mathrm{p}}_{\mathrm{o}}-{\mathrm{p}}_{\mathrm{e}}}{1-{\mathrm{P}}_{\mathrm{e}}}}$$where ${\mathrm{P}}_{\mathrm{o}}$ is the sum of the diagonals of the confusion matrix divided by the total number of samples, and ${\mathrm{P}}_{\mathrm{e}}$ is the sum of the “product of actual and predicted numbers” for each of the categories, divided by the “square of the total number of samples”.

Construction of a comprehensive ecological index

Therefore, the remote sensing ecological index (RSEI) is rapid monitoring and evaluation of the ecological method based on natural factors, which integrates moisture, greenness, heat and dryness information and can also be used to evaluate the ecological and environment of the burned area.

Greenness

In this work, the normalized vegetation index (NDVI) was selected to characterize the greenness of the study area, which can reflect the plant biomass, leaf area index and vegetation cover in the region, calculated as Eq. (13):

(13) $$NDVI={\displaystyle \frac{{\rho }_{nir}-{\rho }_{red}}{{\rho }_{nir}+{\rho }_{red}}}$$where, ${\rho }_{nir}$ and ${\rho }_{red}$ represent the reflectance of Landsat images in the NIR and red bands, respectively.

Wetness

The wetness index is characterized based on the moisture component obtained from the tassel cap transformation, and the wetness can reflect the surface water body conditions, especially the moisture of the soil and the calculation method for wetness extraction based on Landsat8/OLI data is shown in Eq. (13) (Baig et al., 2014). Before extraction, the water bodies were masked using the Modified Normalized Difference Water Index (MNDWI) (Sun, Wu & Tan, 2012), so the wetness reflects the actual land surface moisture conditions.

(14) $$\mathrm{W}\mathrm{e}\mathrm{t}=0.1511{\mathrm{b}}_2+0.1973{\mathrm{b}}_3+0.3283{\mathrm{b}}_4+0.3407{\mathrm{b}}_5-0.7117{\mathrm{b}}_6-0.4559{\mathrm{b}}_7$$where: ${\mathrm{b}}_2,{\mathrm{b}}_3,{\mathrm{b}}_4,{\mathrm{b}}_5,{\mathrm{b}}_6,{\mathrm{b}}_7$ represents the reflectance of Landsat8/OLI images in bands 2, 3, 4, 5, 6, and 7, respectively; Wet represents the humidity component of Landsat8/OLI images.

Heat

This work used Landsat 8 surface temperature to invert the land surface temperature (LST) to indicate the hotness indicator. The dual-band feature of Landsat 8 thermal infrared data offers the possibility of using a split-window algorithm to invert the surface temperature. Therefore, single-channel surface temperature inversion methods are mainly used in surface temperature inversion work to estimate surface temperature. This article chose the single window algorithm (SMW) proposed by Ermida et al. (2020). The SMW algorithm is based on an empirical relationship between TOA brightness temperatures in a single TIR channel and LST and utilizes simple linear regression. The model consists of a linearization of the radiative transfer equation that maintains an explicit dependence on surface emissivity, and it is calculated as Eq. (15).

(15) $$LST=A_i{\displaystyle \frac{T_b}{\epsilon }+B_i{\displaystyle \frac{1}{\epsilon }+C_i}}$$where $T_b$ is the TOA brightness temperature in the TIR channel, and $\epsilon$ is the surface emissivity for the same channel. The algorithm coefficients $A_i$; $B_i$; and $C_i$ are determined from linear regressions of radiative transfer simulations performed for 10 classes of Total Column Water Vapor (TCWV) that values from NCEP/NCAR reanalysis data are available on GEE (I = 1, …, 10), ranging from 0 to 6 cm in steps of 0.6 cm, with values of TCWV above 6 cm being assigned to the last class.

(4) Dryness

The dryness index represents an indicator of land dryness, which also indicates the degree of soil dryness and the land is affected by weathering and sanding. The dryness index consists of the normalized building bare soil index (NDBSI) incorporating the bare soil index (SI) and the building index (IBI), and the burned areas are mainly affected by the SI, and most burned areas are far from the towns and not easily affected by the IBI. Rikimaru’s model was used to calculate the following Eq. (16).

(16) $$SI={\displaystyle \frac{({\rho }_{swir1}+{\rho }_{red})-({\rho }_{blue}+{\rho }_{NIR})}{({\rho }_{swir1}+{\rho }_{red})+({\rho }_{blue}+{\rho }_{NIR})}}$$where, ${\rho }_{blue},{\rho }_{red},{\rho }_{swir1},{\rho }_{NIR}$ denotes the reflectance of the blue band, red band, near-infrared band, and short-wave infrared band1 corresponding to the OLI image, respectively.

The proposed ecological index should be able to appear both as a single indicator and to combine the information from the above four indicators. Therefore, we used principal component analysis (PCA) (Gewers et al., 2021), a multivariate statistical method, to construct the ecological index. It is a multi-dimensional data compression technique that can select a few critical variables by an orthogonal linear transformation. Its most significant advantage is that the weights of the integrated indicators are not determined artificially; instead, they are determined automatically and objectively according to the contribution of each indicator to each principal component, thus avoiding the bias of the results caused by the different weight settings from a person.

Since the four indices are not uniform, it is necessary to normalize the original values of the four indices to a range between [0, 1] to avoid the imbalance of each index weight in the principal component analysis, and the normalization formula is as Eq. (17):

(17) $$BI=(I_i-I_{min})/(I_{max}-I_{min})$$where, BI is the normalized image element value of a factor, $I{}_i$ is the image element value of a factor, $I{}_{max}$ and $I{}_{min}$ the maximum and minimum values of the factor, respectively, the normalized indexes have a scale between 0 and 1. The normalized metrics are combined into a new image, and a principal component analysis is performed on the new image to obtain the first principal component (PC1) and related statistics. The obtained PCl is normalized to obtain the RSEI.

The extraction results of burned areas

According to the OTSU calculation, the heavily burned threshold and mildly burned threshold for each fire was obtained and then can get the heavily burned areas and mildly burned areas, and the two areas were added to get the total burned area. The results are shown in the following Table 4.

From Table 4, we can see that the accuracy of burned area extracted by combining the dNBR and OTSU thresholding image segmentation method for three fires was 97.69%, 89.32%, and 96.37%, respectively. The method performs better in the accuracy of small burned areas, but the performance is average in the extraction of the large burned area; the accuracy of heavily burned area extraction is higher than that of mildly burned.

Combined with the actual visual interpretation areas, the burned areas map derived from GEE using the OTSU threshold image segmentation method was compared, as shown in Fig. 4.

As can be seen, the heavily burned areas that were extracted matched the visual interpretation area, and the proposed method’s accuracy is 98.4%, 91.89% and 98.96%, respectively. In contrast, the mildly burned areas are mainly distributed near the heavily burned areas, and the extracted burned areas are complete with clear boundaries and close to the visual interpretation results. However, it is difficult to intuitively extract the mildly burned areas and discriminate them from the remote sensing images.

Ecological environment evaluation results

As can be seen from Table 5, NDVI and wetness both played a positive role in promoting the ecological environment of the study area during the period, while SI and LST played a negative role together, which is consistent with the actual situation. The PC1’s contribution in the calculation of the RSEI was more outstanding than 90%, indicating that the PC1 had concentrated most of the four indicators’ information in the three fires; meanwhile, the four component indicators’ contribution to PC1 was relatively stable, so PC1 was chosen to create the RSEI instead of the four component indicators.

The RESI’s mean values were calculated from the results of PCA in before and after the fire as shown in the following Table 6.

We have concluded that the RSEI decreases significantly after three fires. To finely analyze the changes of RSEI, we calculated the difference of RESI before and after the fire for each pixel. Where −0.2 to −0.05 is a slight increase, −0.05 to 0.05 is no significant change, 0.05–0.20 is a slight decrease, and greater than 0.20 is a significant decrease. As shown in Fig. 5 below:

As we can see the Fig. 5, the following phenomena were present in all three fires: The post-RSEI of mildly burned and heavily burned was significantly lower than that of pre-RSEI, indicating that the fires brought different degrees of damage to the ecological environment; meanwhile, the difference in RSEI of heavy burned was higher than that of mildly burned, indicating that the damage of RSEI in the area of heavy burned was more serious, and the vegetation restoration measures in the area should be considered to reduce the impact of burned on the ecological environment as early as possible.

The monitoring analysis of the ecological environment

To further investigate the changes in the local ecological conditions, the normalized RSEI ecological index was divided into five levels according to the 0.2 value interval, from low to high: excellent: 0.8–1.0; good: 0.6–0.8; moderate: 0.4–0.6; fair: 0.2–0.4; poor: 0–0.2. We calculated the area and proportion of ecological classes before and after the three fires. And the transfer matrix was used to analyze the ecological changes as follows.

In 2019, we can see (Table 7) that the ecological environment of the burned areas changed: after the burn, the ecological environment had no excellent, but the poor ecological environment had an area of 0.09 hm²; In the mildly burned area, the ecological environment changed from good to medium with an area of 34.26 hm², the changing area accounted for 49.8%, and all the excellent ecological environment changed to good. The ecological environment of the heavily burned area had changed more, and it has a poor with an area of 2.59 hm², the excellent ecological environment (1.08 hm²) changed to a medium (0.99 hm²) ecological environment with accounted for 95.12%, which means that the ecological environment of the heavily burned area has changed more prominently.

In 2020, there was no poor ecological environment before the fire (Table 8), but poor ecological conditions appeared after the fire with an area of 13.88 hm². The areas of the excellent ecological environment also decreased from 1,370.5 to 57.06 hm². In the mildly burned areas, the ecological environment was concentrated in the good ecological environment before the fire, but after the fire, its ecological environment changed to medium grade, accounting for 49.73% (3,962.3 hm²), and the ecological environment changed from excellent to good grade, and it accounted for 92.30%. In the heavily burned area, the ecological environment evaluation was concentrated in the good grade, but it was concentrated in the medium grade after the fire. Changes in the ecological environment were concentrated in the pre-fire good and excellent ecological environment, and most of the excellent ecological environment was changed to the medium or poor ecological environment, with accounted for 98.90%, and most of the good ecological environment was changed to the medium ecological environment with accounted for 72.90%.

Also, we can see (Table 9) that the ecological environment of the burned area changed in 2021: before the fire, there was no poor ecological environment, but after the fire, there was a poor ecological environment with an area of 27.38 hm², while the excellent ecological environment changed from 32.68 to 0.63 hm² with accounted for 61.51%, and the excellent ecological environment changed to good or medium with accounted for 97.93%; In the heavily burned area: medium ecological environment all changed to the worse ecological environment with accounted for 100%, and good ecological environment changed to the worse or medium ecological environment with accounted for 99.02%, excellent ecological environment mainly changed to the medium and worse ecological environment with accounted for 100%.

In summary, there was no poor ecological environment in the three pre-fires, but the ecology showed a poor level after being burned in high fire intensity in 2020 and no poor level after the 2019 and 2021 fires. So in terms of area burned and ecological changes, the fires in 2020 were more severe than those in 2019 and 2020. On the other hand, regardless of the heavily burned area and mildly burned area, the pre-fire ecology had a good area, but its post-fire ecology change area was also the largest, and the post-fire ecology had the most moderate area. It shows that the ecological environment changed from good to moderate after the fire.