Forest fire progress monitoring using dual-polarisation Synthetic Aperture Radar (SAR) images combined with multi-scale segmentation and unsupervised classification

Study area

As shown in Fig. 2, the Thomas wildfire occurred on the coast of southern California, USA, a region with a Mediterranean climate and a mosaic of hills, mountains, plains, and terraces, with a high degree of development of shrubs and evergreen forests. The fire began on 5 December 2017, and was fully contained on 12 January 2018. The total burned area was approximately 114 078 ha, making it the largest fire in recent modern California history at that time (Zhang et al. 2019).

Fig. 2.

Map of the study area.

Sentinel-1 dataset

The Sentinel-1 satellite consists of two satellites (Sentinel-1A and Sentinel-1B), which as successfully launched on 3 April 2014 and 25 April 2016. They have carried the same type of C-band polarised SAR system, which has the characteristics of a short revisit period, wide coverage area, high resolution, and all-weather/all-weather imaging. We chose a 5-scene Sentinel-1 IW-mode SLC image with VV + VH dual-polarisation mode. For information on Sentinel-1 image parameter, see Supplementary Table S1.

Reference data

The spatial resolution of Sentinel-2 multi-spectral remote sensing data is 10, 20, and 60 m, while the temporal resolution is 5 days. Cloud-free Sentinel-2A remote sensing images before and after the Thomas wildfire were selected as reference data to validate the feasibility and reliability of the burned area extraction method in this paper. The normalised burn ratio (NBR) is a widely used index for estimating burn severity and monitoring burned areas. It is calculated using the near-infrared band (NIR, wavelengths 0.85–0.88 μm) and the short-wave infrared band (SWIR2, wavelengths 2.11–2.29 μm). The calculation formula is as (Sismanis et al. 2023):

Calculating the difference in NBR (dNBR) before and after a wildfire can better estimate the area burned. The dNBR is calculated as:

where NBR_pre and NBR_post are NBR images before and after the wildfire, respectively.

Supplementary Fig. S1. shows a map of pre-fire vegetation types in the Thomas Wildfire and the normalized difference vegetation index (NDVI) map. Reference ranges for the Thomas wildfire were obtained by combining California Fire Perimeters vector data and visual interpretation of 12 January 2018, Sentinel-2 imagery. Land cover maps are from the US National Land Cover Database (NLCD-2011), based on 30-m Landsat data (Homer et al. 2015). As can be seen from the land cover map, the primary vegetation type in the burn area is shrubs, followed by evergreen forests and mixed forests. NDVI was calculated using post-fire Sentinel-2 imagery, and the NDVI map indicates that the lower right portion of the burned area has much less vegetation cover than the upper left portion.

Constructing feature datasets

The radar signal before the burn was mainly derived from the vegetation canopy, and the forest fire burned the canopy, the understorey, and the surface vegetation successively. The radar signal reaches the ground, interacts with the soil surface, and returns the radar signal, and the main factor affecting the radar backscatter value is soil moisture. Therefore, this paper uses the backscatter parameters of different polarisation methods to construct the Backscatter Burn Ratio (BBR), which is calculated by:

where γ_pre represents the Sigma value of the backward scattering coefficient or the Gamma value of the normalised backward scattering coefficient before the fire; and γ_post represents the Sigma or Gamma value after the fire.

Forest fires destroy pre-existing vegetation structures, altering the backscattering mechanisms of the burned area. Before the fire, the target scattering echo mainly originates from the canopy of the vegetation, characterised by volume scattering. After the fire, the vegetation canopy is significantly reduced, and the target scattered echoes mainly originate from the soil surface, characterised by surface scattering. Polarisation target decomposition, the primary method for extracting polarisation scattering features, can identify different types of scattering mechanisms. First, the H-A-Alpha polarisation decomposition method is used to obtain the polarisation parameters such as entropy, scattering angle, and anisotropy. Second, the polarisation burn difference is constructed based on the pre-fire and post-fire polarisation parameters. The formula for its calculation is:

where D_δ denotes the difference between the pre-fire and post-fire polarisation parameters. δ_pre and δ_post are the entropy, scattering angle, or anti-entropy of the pre-fire and post-fire, respectively. δ̅_pre and δ̅_post are the average values of the pre-fire and post-fire polarisation parameters in the burned area, respectively.

From the imaging characteristics of SAR images, it is known that the grey scale of SAR images reflects the backward scattering characteristics of ground objects to radar waves. Different ground objects with similar backscattering coefficients will show the same or similar grey values in SAR images, thus presenting a more complex grayscale distribution. Therefore, using only image grey-scale-based features to extract the burned area does not obtain satisfactory classification results. We use the texture information of the target’s interior to compensate for the lack of grayscale features and thus improve the accuracy. Texture analysis is used to study the spatial distribution relationship between pixels. In this paper, eight standard features such as mean, contrast, variance, dissimilarity, homogeneity, correlation, information entropy, and angular second moments are extracted using the Grey Level Co-occurrence Matrix (GLCM) method for Sentinel-1 remote sensing images before and after the fire. We used a window size of 9 × 9, an offset distance of 1, and a grey quantisation level of 64 as the parameters for GLCM calculation. At the same time, four directions (0°, 45°, 90°, and 135°) were selected to calculate the average values of the texture feature parameters. We construct the texture feature variation parameter based on the texture features of different polarisation methods with the expression:

where P_c is the difference between the texture feature parameters before and after the fire. c̅_pre and c̅_post denote the average values of the texture features before and after the fire in four directions: 0°, 45°, 90°, and 135°, respectively.

In this paper, 13 feature parameters were calculated using dual-polarised Sentinel-1A images. For detailed feature names and feature meanings, see Table S2.

Feature selection

We have constructed too many feature factors from the perspective of backscatter intensity, scattering mechanism, and SAR texture features. These features are more or less correlated with each other, and some of them have mutual constraints. Therefore, using all these features for burned area extraction will generate data redundancy and affect work efficiency. At the same time, the classification results are often uncertain due to the problem of ‘dimensional catastrophe.’ Therefore, before performing the burning area extraction, we reduce the data redundancy between features by calculating the correlation coefficients between the feature factors.

The Pearson correlation coefficient (Dufera et al. 2023) is a statistic used to reflect the degree of linear correlation between two features, reducing the potential co-linearity of the features. Because of this, this paper uses the Pearson correlation coefficient to determine the correlation between features, aiming to determine whether there is a high degree of correlation between features. The formula calculates the Pearson correlation coefficient:

where COR is the Pearson correlation coefficient; x and y are the feature factors; cov(x, y) is the covariance of the feature factors x and y; var(x) and var(y) are the variances of the feature factors.

Multi-scale segmentation and optimal segmentation scale selection

Multi-scale segmentation (Du et al. 2016; Fu et al. 2022) is a top-down region merging technique. The merging process follows a minimal heterogeneity criterion, which considers the spectral and shape characteristics of the images, whose shape characteristics include smoothness and compactness. The heterogeneity criterion is calculated as:

where f is the heterogeneity criterion; h_colour is the spectral feature; h_shape is the shape feature; h_smooth is the smoothness; and h_compact is the compactness. ω₁ and ω₂ are scaling factors; λ_i denotes the weight occupied by the i band; and σ_i denotes the standard deviation of the i band. E denotes the number of image pixels contained in the boundary contour of the image object; N denotes the number of image pixels contained in the image object; and L denotes the length of the boundary that contains the smallest outer rectangle within the range of the image object.

From Eqns 7 to 11, it can be seen that spectral features and shape features, smoothness, and compactness are complementary, which means that an increase in the ratio coefficient of one side will correspondingly decrease the ratio coefficient of the other side. Consequently, we used the iterative preference method to determine the shape factor and compactness factor of the study object. We chose the optimal segmentation scale of the study area with the help of the estimation of scale parameters (ESP) method. ESP takes the local variance (LV) of the homogeneity of the image object as the average standard deviation of the segmented object layer, and when the rate of change (ROC) of the LV peaks, the segmentation scale corresponding to that point is the optimal segmentation scale. ROC is calculated as:

where L_i and L_i−1 denote the average standard deviation of the i object layer and i − 1 object layer of the target layer, respectively.

Unsupervised classification algorithms and parameter selection

K-means clustering (Sinaga and Yang 2020) is an unsupervised classification method that is easy to implement and efficient and is the most widely used of all clustering algorithms. This method uses distance as a similarity criterion and assumes a higher degree of similarity the closer two objects are to each other (Abid et al. 2021). The k-means clustering algorithm is a commonly used cluster analysis algorithm with an iterative solution, and the specific steps are as:.

Assume that the dataset is X = {x₁, x₂, …, x_i, …, x_n}, the number of clusters is K, and the cluster centre is C = {c₁, c₂, …, c_j, …, c_k}.

One of the main problems of the K-means clustering algorithm is determining the optimal number of clusters (the K parameter). Silhouette coefficients are a way of assessing how good or bad the clustering results are. For a sample point in the data set, the Silhouette coefficients are calculated as:

where a represents the average distance of the sample x_i to other samples in the same cluster, and the smaller a is indicates that the sample x_i is the more it should be clustered into that cluster. b represents the sample x_i to the minimum of the average distance from the sample points in the other clusters, assuming that the number of clusters is K, i.e. b = min{b₁,b₂, …, b_k}, where 1 ≤ k ≤ K. The larger the value of S(x_i), the tighter the cluster in which the sample point x_i is located. The size of the silhouette coefficient is between [−1, 1], and the average silhouette coefficient for each cluster is calculated as:

s_k represents the mean value of the silhouette coefficients of the sample points of the kth cluster, where n_k shows the number of sample points of the kth cluster, and a larger value of s_k indicates a better clustering effect, and vice versa.

Data set preparation

First, we normalised all the feature factors to the range of [0−1] to ensure the comparability and consistency of the data to improve the accuracy and reliability of the research results in this paper. Then, using Pearson’s correlation coefficient, a random selection of 100 000 sample points within the study area was used to calculate the correlation between the featured factors. Fig. 3 shows the correlation coefficients for each feature factor.

Fig. 3.

Correlation coefficients of the featured factors.

Feature importance is a measure used to evaluate the usefulness of a feature in the model classification process. This paper uses feature importance ranking to reveal the influence of different feature factors on the extraction of burned areas. Fig. 4 shows the importance ranking of feature factors.

Fig. 4.

Importance ranking of feature factors.

As shown in Fig. 3, the correlation coefficients between the feature factors of the same feature type are relatively high. For example, the Pearson correlation coefficient between P_variance, P_contrast, P_{dissimilarity} feature factors, and the other seven texture feature factors is the lowest at 0.52, and the Pearson correlation coefficient of the three polarisation decomposition feature factors is the lowest at 0.70. The feature factors are ranked in importance to select more representative feature factors. Combining Pearson’s correlation coefficient and the importance ranking of the features, the texture parameters P_mean and P_variance, the polarisation decomposition parameters D_alpha and D_entropy, and the backscattering intensity parameters BBR_gamma and BBR_sigma are finally selected in this paper, and the six selected features are used for the unsupervised, object-oriented extraction of burned area.

Fig. 5 illustrates the spatial distribution of the six feature factors at different times for the Thomas wildfire. The changes within the burned area are not uniform due to the complexity of the terrain, among other reasons, and show different degrees of variability on different parameters. On some parameters, the change in the burned areas is so pronounced that it can be differentiated from the unburned areas, while on some parameters, the change in the burned areas is less pronounced and can easily be confused with the unburned areas. As seen in Fig. 5, P_mean and P_variance can recognise the range of the combustion zone very well. In addition, other feature factors selected by feature selection also see an approximate range of the combustion zone.

Fig. 5.

Characteristic factor maps of the Thomas fire.

Multi-scale segmentation experiment

The optimal segmentation scale of the study area was determined based on the estimation of scale parameters (ESP) analysis. This experiment uses 20 as the initial segmentation scale, the step size is 1, and the weights of each band are all 1 for multi-scale segmentation to obtain the ROC-LV graph (Fig. 6).

Fig. 6.

Multi-scale segmentation results.

As shown in Fig. 6, the ROC-LV curves have multiple peaks, such as 62, 74, 88, 97, etc. We chose localised images of the study area for the test to avoid the subjectivity of segmentation scale selection. The ground object category samples are outlined through visual interpretation as the standard for accuracy evaluation, and the over-segmented and under-segmented objects in the local images are analysed. The optimal segmentation scale that meets the study area is finally selected as 97. The segmentation results obtained under the optimal parameters are superimposed on the optical image for display (Fig. S2), showing the multi-scale segmentation algorithm can maximise the retention of the authenticity of the segmented image, and the segmentation boundary better reflects the difference between the combustion and non-combustion regions.

The silhouette coefficient algorithm determines the value of K

Using 100 000 randomly selected sample points from 15 January 2018, Thomas Fire SAR dataset, the optimal number of clusters (K-value) was determined using the silhouette coefficient method. Fig. S3a illustrates the trend of the silhouette coefficient scores calculated for values of K taken from 2 to 20. The highest values of silhouette coefficients are generated by lower k values (Fig. S3a). The next highest value, silhouette score of 0.248, was found when K = 6. When the number of clusters is two, the extraction results of combustion area are shown in Fig. S3b. From Fig. S3b, it can be seen that a large number of non-burned areas are wrongly extracted as burned areas when K = 2, and the extraction results are not good. Therefore, the clustering number used in this study was six.

Burned area extraction and accuracy evaluation

Visual comparison of the pixel-based and object-based classifications (Fig. 7) showed severe fragmentation in the extraction of pixel-based classification. When performing pixel-level classification, the classification algorithm assigns each pixel to its corresponding category, and the salt-and-pepper noise causes the brightness values of some pixels to be inconsistent with their actual categories, leading to misclassification. Therefore, pixel-based classification is affected by the salt-and-pepper noise, which leads to poor fire trace ground extraction. Object-oriented burned area extraction suppressed the salt-and-pepper noise using a multi-scale segmentation method. As a result, the burned area extraction results of this paper’s method have a more continuous distribution and better integrity without apparent fragmentation. The experimental results show that the method in this paper can capture the shape and boundary of the burned area more accurately, which provides a reliable basis for subsequent fire monitoring and management.

Fig. 7.

Classification results.

This paper uses a multi-scale segmentation approach to combine neighbouring pixels into larger objects and then performs burned area extraction on these objects. As can be seen from a1 to c3 of Fig. 8, the classification results of this paper’s method are complete in shape, without the phenomenon of many fragmented patches inside the burned area, and the distribution of the burned area is more accurate. As shown from d1 to d3 in Fig. 8, roads, some cultivated land, and bare land have similar backscattering characteristics with the burned area, resulting in their misclassification as burned area. In contrast, the classification method in this paper is processed by multi-scale segmentation, and there is no situation in which many geomorphic features are misclassified.

Fig. 8.

Detailed map of the burning area extraction results. a1–d1 are optical images, a2–d2 are pixel-based classification results, and a3–d3 are this paper’s method’s burned area extraction results.

Agreement (burned) in Fig. 9 indicates the area where the burned area extracted based on unsupervised classification of SAR images agrees with the burned area extracted based on optical images. S1 only shows the area where SAR images are commission extracted compared to the burned area based on optical image recognition. S2 only indicates where SAR images are omitted to be removed compared to the burned area based on optical image recognition.

Fig. 9.

A comparison of burned area extraction results based on SAR and optical images for the Thomas wildfire, 15 January 2018. (a) Pixel-based classification, and (b) object-orientated classification. S1 refers to Sentinel-1 and S2 to Sentinel-2.

The method of accuracy assessment using the confusion matrix depends on the selected samples, which vary from person to person when selected. Therefore, to evaluate the reliability of the experimental results of this paper more objectively, this paper uses the burned area extracted by Sentinel-2 as the reference data to quantitatively assess the results of the burned area removed on 15 January 2018, and the results of the quantitative evaluation of the Thomas wildfire are in Table 1.

In this paper, we fully used the various feature change parameters caused by forest fire and combine its scattering intensity, polarisation decomposition, and texture feature change parameters to realise the unsupervised extraction of the burning area. As seen in Table 1, the overall accuracy of the polarised SAR image forest fire progress monitoring method combined with multi-scale segmentation and unsupervised classification proposed in this paper is higher than that of the pixel-based classification method. The extraction accuracy of this our method is 91.67%, and the omission errors and commission errors are 8.31 and 7.61%, respectively. Compared with the pixel-based unsupervised classification method, the overall accuracy of this paper’s method is improved by 33.70%, and omission errors and commission errors are reduced by 33.72 and 1.15%, respectively. These results indicate that the method in this paper exhibits higher accuracy and stability in combustion area extraction and achieves significant improvement compared to the traditional pixel-based method.

The difference in feature extraction performance between single-temporal and multi-temporal data

Burning area extraction models were constructed for the two cases using multi-temporal phase and post-fire single-temporal phase polarisation SAR images, respectively. The experimental results are in Fig. 10. As seen in Fig. 10, the use of single-temporal phase feature parameters can only describe the scattering mechanism of a certain aspect due to the influence of complex burning conditions, which has a limited ability to explain forest fires. The forest fire monitoring model constructed using multi-temporal polarisation SAR images can extract the burned area accurately.

Fig. 10.

Burned area extraction results for different feature combinations.

The statistical analysis of the burned areas extraction results based on multi-temporal phase features and post-fire single-temporal phase features, respectively (Table 2). The results show that the burning area extraction accuracy based on multi-temporal phase polarisation data is improved by about 30% compared to post-fire single-temporal phase polarisation data, and the commission extraction error and omission extraction error are reduced by about 10 and 32%, respectively. The scattering characteristics of some areas after forest fire burning are relatively similar to those of some unburned ground objects, which are difficult to be distinguished from the single-temporal data after forest fire. Compared with using only post-fire data, the use of multi-temporal phase data not only expresses the difference between the burned area itself and other ground objects but also makes full use of the changes in the burned area before and after the forest fire, providing reliable features for the subsequent extraction of the fire trails.

Optimal feature combination for burned area extraction

The existing remote sensing extraction of burned areas mainly focuses on a single SAR feature type, and the mining of multiple types of features, such as backscatter intensity, polarisation decomposition, texture, etc. of SAR images is not yet sufficient (Shama et al. 2023). In this paper, multiple types of feature parameters caused by forest fires are analysed in depth, and the monitoring capability of forest fires with different polarisation feature combinations is compared. The results of burned area extraction for different feature combinations are in Fig. 10. The forest fire changed parameters such as soil moisture, biomass, roughness, and dielectric constant, which led to the change of SAR backscattering intensity (Belenguer-Plomer et al. 2019, 2021). At the same time, the forest fire burned down the forest canopy destroying the original vegetation structure, thus changing the scattering mechanism of the burned area, which is manifested by a large reduction of the volume-scattered component and an increase of the surface-scattered component from the soil surface. Therefore, as can be seen from the line comparison in Fig. 10, there are limitations in extracting the burned area with only a single type of feature. By combining three types of effective parameters (scattering intensity variation parameter, scattering mechanism variation parameter, and texture feature variation parameter), the burned area extraction results have significantly improved.

Results of the burned area extracted based on different feature combinations were analysed statistically, which indicated: (1) that due to the complex texture of the study area, it is difficult to distinguish the burned area from the non-burned area by using only a single feature type; (2) the extraction accuracy is improved by about 10% after fully considering the parameters of multiple types of feature changes caused by forest fires; and (3) the error of commission extraction and omission of extraction is reduced by about 7 and 10%, respectively. Using different feature combinations can improve the ability to obtain surface information, thus improving the extraction accuracy of the burned area.