Satellite-Observed Soil Moisture as an Indicator of Wildfire Risk

Wildfires are a concerning issue in Canada due to their immediate impact on people’s lives, local economy, climate, and environment. Studies have shown that the number of wildfires and affected areas in Canada has increased during recent decades and is a result of a warming and drying climate. Therefore, identifying potential wildfire risk areas is increasingly an important aspect of wildfire management. The purpose of this study is to investigate if remotely sensed soil moisture products from the Soil Moisture and Ocean Salinity (SMOS) satellite can be used to identify potential wildfire risk areas for better wildfire management. We used the National Fire Database (NFDB) fire points and polygons to group the wildfires according to ecozone classifications, as well as to analyze the SMOS soil moisture data over the wildfire areas, between 2010–2017, across fourteen ecozones in Canada. Timeseries of 3-day, 5-day, and 7-day soil moisture anomalies prior to the onset of each wildfire occurrence were examined over the ecozones individually. Overall, the results suggest, despite the coarse-resolution, SMOS soil moisture products are potentially useful in identifying soil moisture anomalies where wildfire hot-spots may occur.


Introduction
Wildfires are considered to be a natural disturbance to terrestrial ecosystems, as well as an integral part of ecosystems' natural evolution. However, large wildfires, particularly when due to human activity, can be considered hazardous to ecosystems as they interrupt an ecosystem's natural evolution. Nevertheless, whether it is natural or human induced, large wildfires significantly affect regional climate, the land surface hydrology, biogeochemical processes, and wildlife habitats, as well as the livelihoods of people. Recent studies have attributed the increase of large wildfire occurrences to a warming climate due to increased greenhouse effect [1][2][3][4]. The warming trend is expected to continue in the coming decades, potentially increasing the vulnerability of ecosystems to wildfires and related climate disturbances [5].
The increase in vulnerability of ecosystems to wildfires in a warming climate is of great concern to Canada since the country holds some of the largest forest ecozones in the world. For example, the Boreal Shield, being the largest ecozone in Canada, covers almost 20% of the total landmass, accounting for the majority of recorded wildfire area burned annually across Canada. On average, approximately two million hectares of forests are subject to fire annually across Canada [6,7]. The total carbon emission from these wildfires is estimated at 27 Tg carbon per year [6,[8][9][10][11][12]. Within Canada, managed forests, which are subjected to more intensive fire suppression measures than unmanaged forests, are primarily located in the southern extent of the country, while the majority of unmanaged forests are located away from population centers [13]. Given the risks involved in managing large wildfires, identifying potential wildfire risk areas or "hot-spots" in advance is one of the most important aspect of wildfire

Data and Methods
In this study, we used fire data from the National Fire Database [56] for an 8-year period, from 2010-2017. The NFDB, maintained by Natural Resources Canada (NRC) via the Canadian Wild-Land Fire Information system, is composed of annual fire points and polygons derived from various regional fire management agencies [7,57]. The fire points in the NFDB indicate the geographical locations of ignition points, and the polygons represent the perimeters of burned areas. The spatial distribution of all fires between 2010-2017, over the ecozones analyzed in this study, are illustrated in Figure 1. Using this spatial data set, the fire point data (date and area) were used to extract corresponding SMOS SM data corresponding to time and location of each fire.
Remote Sens. 2020, 12, x FOR PEER REVIEW 4 of 14 we tested the hypothesis that a negative shift (i.e., negative anomaly of SM) had occurred when SM observations were drawn from time periods prior to fire initiation. We evaluated the statistical significance of shifts to the histogram using a one sample Kolmogorov-Smirnov (KS) normality test. The total number of all fires, as well as the number of SMOS grid cells within each of the ecozones, is illustrated in Figure 2a. Based on Figure 1 and 2a, it is clear that the majority of the fires occurred within the Boreal ecozones (i.e., Boreal plain and Boreal shield) and in the Montane cordillera. In some ecozones, such as the Boreal plain and Montane cordillera, there are few SMOS grid pixels, which implies a relatively large number of fires occurred within each SMOS grid pixel in that ecozone (see Figure 2a). Figure 2b illustrates the ratio of the number of fires exceeding a specific size (200 ha and 3600ha) to the total number of fires in each ecozone.  Processing of SMOS brightness temperature observations provided a soil water content estimate for the top 5 cm (here, we used Level 2 Soil Moisture V650, released on 15 November 2017. SMOS data is freely available at https://smos-diss.eo.esa.int/oads/access/). For the development of the SMOS product, radiometric observations made from SMOS instrument were processed to Level 1 brightness temperature (Tb) with~40-km spatial resolution and reported on a 15-km hexagonal grid called Discrete Global Grid (DGG; ISEA 4H9 grid). The brightness temperature was then processed to a SM product and reported on the same grid (15-km ISEA grid). Therefore, the associated grid of L2 SM product was 15 km; however, the radiometric resolution of the instrument was~40 km. The accuracy of SMOS L2 SM products, compared to several in-situ data sets, approached an unbiased root mean square error of 0.04 m 3 m −3 [44,45,58]. We also obtained soil moisture estimates from SMAP Soil Moisture [59,60] retrieval, based on radiometer data scaled to a 9 km EASE-2 grid (here we have used version R16010. SMAP data is freely available at https://nsidc.org/data/smap/smap-data.html.). Given the longer time series of SMOS, most of our analyses (as described below) focused on the SMOS time series. The accuracy of SM from SMOS is reduced over mountainous regions due to the influence Remote Sens. 2020, 12, 1543 4 of 14 of complex topography on local incidence angles [61,62]. Therefore, the cordillera ecozones were excluded from the overall results.
Prior to statistical analysis, the fire points, as well as the SMOS SM data, were grouped according to the ecozones in Canada. The ecozones of Canada are characterized regions grouped by various abiotic and biotic factors [63]. To aggregate the wildfires and SMOS SM data, we used the geographic location information (latitude and longitude) of the wildfires from the NFDB database, matched-up to that from SMOS SM metadata (based on area of SMOS grid), and then grouped according to their respective ecozones. The SMOS SM data was quality controlled prior to the analysis. Anomalously high and low soil moisture values were excluded from the analysis (i.e., SM values below the sensitivity of the satellite sensor. In general, for silty sand and silty loam soils texture types, moisture values above 0.6 m 3 /m 3 indicate saturated or near saturated soil. Therefore, we excluded SM values > 0.6 m 3 /m 3 .). Temporal gaps in the SM data filled using MATLAB's 'pchip' function, which stands for Piecewise Cubic Hermite Interpolating Polynomial, a shape-preserving piecewise cubic interpolation method. The pchip interpolation method is local, which means the interpolation within a subinterval is determined by only the nearby four data points, the two data points on either side of that interval, and is not influenced by the data points farther away. Another advantage of using pchip is that it has no overshoots, less oscillation if the data are not smooth, and is computationally less expensive (All analysis was performed using MATLAB R2017a ver. 9.2). Over the SMOS time series (2010-2017) and for each pixel in each of the ecozones, we calculated daily standardized anomalies using the interannual daily mean (µ) and standard deviation (σ) (i.e., the daily mean and standard deviation were computed across 8 years for each pixel). For example, we calculated the standardized anomaly (SA) of the SMOS observed SM for each day (i) as: SA i = (SM i − µ). For our analysis, we examined the soil moisture anomaly conditions for three time periods prior to the fire initiation: 3, 5, and 7 days prior to initiation. Although these intervals were chosen arbitrarily, it maximized the SMOS SM data sample size since the average SMOS satellite repeat cycle is less than 3 days. Over these time periods, we constructed histograms of soil moisture anomalies observed over the three time periods. Histograms of soil moisture anomalies for each fire, within the three time periods, grouped over each ecozone, were prepared for our statistical analysis. Given that we had drawn the SM anomaly samples from a normally distributed sample (given the standardization shown above), we tested the hypothesis that a negative shift (i.e., negative anomaly of SM) had occurred when SM observations were drawn from time periods prior to fire initiation. We evaluated the statistical significance of shifts to the histogram using a one sample Kolmogorov-Smirnov (KS) normality test.
The total number of all fires, as well as the number of SMOS grid cells within each of the ecozones, is illustrated in Figure 2a. Based on Figures 1 and 2a, it is clear that the majority of the fires occurred within the Boreal ecozones (i.e., Boreal plain and Boreal shield) and in the Montane cordillera. In some ecozones, such as the Boreal plain and Montane cordillera, there are few SMOS grid pixels, which implies a relatively large number of fires occurred within each SMOS grid pixel in that ecozone (see Figure 2a). Figure 2b illustrates the ratio of the number of fires exceeding a specific size (200 ha and 3600 ha) to the total number of fires in each ecozone.

Results
To illustrate the potential sensitivity of SMOS soil moisture products prior to a fire initiation, in Figure 3, we show some selected soil moisture time series prior to the onset of fire. In Figure 3, we included a SMAP time series for reference only and to potentially highlight the consistency between two passive microwave L-band sensors shown to be sensitive to soil moisture. We did not complete the histogram analysis of the SMAP data set as the time series is relatively short for the SMAP sensor (since 2015). Four randomly selected fires (unique fire identifier is indicated in square brackets on the figure subtitles) of various burned area sizes in four major ecozones were selected. The first fire located in the Boreal Plain Ecozone (Figure 3a 13), also known as the "Puntzi Lake fire", took place in mid-summer 2015, burned approximately 8089ha, and was located in the Montane Cordillera ecozone (Figure 3d; here, we include a Montane Cordillera fire as the SMOS quality flags associated with topographic impacts on data quality were not present for the specific pixel containing the fire). It is clear from all four cases (Figure 3a-d), that both SMAP and SMOS reveal a drying trend, at least one week prior to fire initiation, which indicates the two satellite systems' sensitivity to the drying of the surface prior to the onset of a fire. However, it should be noted that both sensors vary in terms of magnitude of the

Results
To illustrate the potential sensitivity of SMOS soil moisture products prior to a fire initiation, in Figure 3, we show some selected soil moisture time series prior to the onset of fire. In Figure 3, we included a SMAP time series for reference only and to potentially highlight the consistency between two passive microwave L-band sensors shown to be sensitive to soil moisture. We did not complete the histogram analysis of the SMAP data set as the time series is relatively short for the SMAP sensor (since 2015). Four randomly selected fires (unique fire identifier is indicated in square brackets on the figure subtitles) of various burned area sizes in four major ecozones were selected. The first fire located in the Boreal Plain Ecozone (Figure 3a , also known as the "Puntzi Lake fire", took place in mid-summer 2015, burned approximately 8089 ha, and was located in the Montane Cordillera ecozone (Figure 3d; here, we include a Montane Cordillera fire as the SMOS quality flags associated with topographic impacts on data quality were not present for the specific pixel containing the fire). It is clear from all four cases (Figure 3a-d), that both SMAP and SMOS reveal a drying trend, at least one week prior to fire initiation, which indicates the two satellite systems' sensitivity to the drying of the surface prior to the onset of a fire. However, it should be noted that both sensors vary in terms of magnitude of the dryness. The sudden "jumps" in soil moisture, particularly in the Hudson Plain ( Figure 3c) and in the Montane Cordillera fires (Figure 3d), are due to small local rainfall events and the uncertainty of the products ( [45,59]; an unbiased RMSE of approximately ± 0.04 m 3 m −3 ). Although illustrative of soil moisture time series prior to fire and potentially for the application of these products for identifying potential wildfire hotspots, Figure 3 is anecdotal. Therefore, we expand our analysis to all 48,754 of the fires recorded in the National Fire Database.
Remote Sens. 2020, 12, x FOR PEER REVIEW 6 of 14 dryness. The sudden "jumps" in soil moisture, particularly in the Hudson Plain ( Figure 3c) and in the Montane Cordillera fires (Figure 3d), are due to small local rainfall events and the uncertainty of the products ( [45,59]; an unbiased RMSE of approximately +/-0.04m 3 m -3 ). Although illustrative of soil moisture time series prior to fire and potentially for the application of these products for identifying potential wildfire hotspots, Figure 3 is anecdotal. Therefore, we expand our analysis to all 48,754 of the fires recorded in the National Fire Database.  Figure 4 illustrates histograms of the soil moisture anomalies 3-days, 5-days, and 7-days (column-wise, left to right) prior to the occurrence of fire, for Atlantic Maritime, Boreal Plain, Pacific Maritime, Taiga Plain, Boreal Shield, and Taiga Shield (row-wise, top to bottom) ecozones, respectively. These seven ecozones are selected based on relatively large numbers of fires (as shown in Figure 2a). It is clear from Figure 4 that the histogram distributions are quite different in each ecozone, although the shape of the distributions are somewhat consistent at 3-days, 5-days, and 7days prior to the occurrence of fire. The KS normality test of the distribution indicates that, for all ecozones except the Taiga Shield, at 3-days prior to the fire (Figure 4p), the test statistic (ks) is greater than the critical value (cv) rejecting the null-hypothesis at 5%. In other words, the KS test suggests, that, for all ecozones (except that previously mentioned), the soil-moisture anomalies are not normally distributed (hv =1). All the KS test statistics are displayed on the top-left of each sub-figure in Figure 4.
For all the ecozones, except the Atlantic Maritime (Figure 4a-c), the histograms indicate a negative shift suggesting drier soil conditions. The skewness (Sk) of the histogram are computed using  Figure 4 illustrates histograms of the soil moisture anomalies 3-days, 5-days, and 7-days (column-wise, left to right) prior to the occurrence of fire, for Atlantic Maritime, Boreal Plain, Pacific Maritime, Taiga Plain, Boreal Shield, and Taiga Shield (row-wise, top to bottom) ecozones, respectively. These seven ecozones are selected based on relatively large numbers of fires (as shown in Figure 2a). It is clear from Figure 4 that the histogram distributions are quite different in each ecozone, although the shape of the distributions are somewhat consistent at 3-days, 5-days, and 7-days prior to the occurrence of fire. The KS normality test of the distribution indicates that, for all ecozones except the Taiga Shield, at 3-days prior to the fire (Figure 4p), the test statistic (ks) is greater than the critical value (cv) rejecting the null-hypothesis at 5%. In other words, the KS test suggests, that, for all ecozones (except that previously mentioned), the soil-moisture anomalies are not normally distributed (hv = 1). All the KS test statistics are displayed on the top-left of each sub-figure in Figure 4.
MATLAB's skewness function. The skewness (Sk) parameter displayed on the top-right of each histogram indicates that almost all histograms, except the Atlantic Maritime (Figure 4a-c) and Boreal Shield (Figure 4m-o), show a positive skewness, suggesting a left leaning distribution, i.e., negative soil moisture anomaly. In Atlantic Maritime (Figure 4a-c) and Boreal Shield (Figure 4m-o), the histogram indicates a slight positive shift (negative skewness), suggesting relatively wetter conditions prior to fire initiation.  In Atlantic Maritime (Figure 4a-c) and Boreal Shield (Figure 4m-o), the histogram indicates a slight positive shift (negative skewness), suggesting relatively wetter conditions prior to fire initiation.
Previous research demonstrates that, although there is a much larger number of small fires, many are subjected to fire suppression measures and thus represent much less area burned than observed in large fires [7,64,65]; therefore, we focused on larger fires (>200 ha) for the subsequent analysis. The contrast between total number of fires vs the number of large fires is evident in Figure 2 and is consistent with previously published studies (see also Figure 4 in Reference [7]). For example, Figure 2a indicates that, although there were approximately 3000 fires that occurred in Atlantic Maritime ecozone during the study period, only a very small fraction were large fires. Similarly, less than 5% of thẽ 12,000 fires in the Boreal ecozones (Boreal Plain and Boreal Shield) were considered large fires. To examine the importance of SM anomaly on large fires, we recomputed the histogram ( Figure 5) after excluding all the fires with burned area < 200 ha. The histogram ( Figure 5) shows the normalized SM anomaly at 3-days, 5-days, and 7-days (column-wise, left to right) prior to the occurrence of fire, for the Hudson Plain, Taiga Plain, Boreal Shield, and Taiga Shield ecozones (row-wise, top to bottom). Some ecozones were excluded as the number of large fires (>200 ha) were very small over the time period resulting in insufficient data to form a histogram. Figure 5 demonstrates that for all ecozones except the Taiga shield (Figure 5j-l), the soil moisture anomalies show a shift towards a negative anomaly (dry soil moisture anomaly as also indicated a positive skewness values) observable for 3-days, 5-days, and 7-days (columns left-to-right) prior to the fire initiation. The SMOS observed SM anomaly does not show similar sensitivity over Taiga shield, which may be related to numerous open water surfaces over the area, which will have an influence on the SMOS retrieval accuracy [44].
Remote Sens. 2020, 12, x FOR PEER REVIEW 8 of 14 Previous research demonstrates that, although there is a much larger number of small fires, many are subjected to fire suppression measures and thus represent much less area burned than observed in large fires [7,64,65]; therefore, we focused on larger fires (>200ha) for the subsequent analysis. The contrast between total number of fires vs the number of large fires is evident in Figure  2 and is consistent with previously published studies (see also Figure 4 in Reference [7]). For example, Figure 2a indicates that, although there were approximately 3000 fires that occurred in Atlantic Maritime ecozone during the study period, only a very small fraction were large fires. Similarly, less than 5% of the ~12,000 fires in the Boreal ecozones (Boreal Plain and Boreal Shield) were considered large fires. To examine the importance of SM anomaly on large fires, we recomputed the histogram ( Figure 5) after excluding all the fires with burned area < 200ha. The histogram ( Figure 5) shows the normalized SM anomaly at 3-days, 5-days, and 7-days (column-wise, left to right) prior to the occurrence of fire, for the Hudson Plain, Taiga Plain, Boreal Shield, and Taiga Shield ecozones (rowwise, top to bottom). Some ecozones were excluded as the number of large fires (>200 ha) were very small over the time period resulting in insufficient data to form a histogram. Figure 5 demonstrates that for all ecozones except the Taiga shield (Figure 5j-l), the soil moisture anomalies show a shift towards a negative anomaly (dry soil moisture anomaly as also indicated a positive skewness values) observable for 3-days, 5-days, and 7-days (columns left-to-right) prior to the fire initiation. The SMOS observed SM anomaly does not show similar sensitivity over Taiga shield, which may be related to numerous open water surfaces over the area, which will have an influence on the SMOS retrieval accuracy [44].  Distributions of fires with burned area greater than 3600 ha are illustrated in Figure 6. Histograms could only be created for three ecozones: the Taiga Plain, the Boreal Shield, and the Taiga Shield. Overall, the figure suggests that SMOS SM anomaly observable 3, 5, and 7 days prior to the fire particularly over the Boreal Shield captures the drying surface conditions prior to the occurrence of fire, as indicated by the shift in the distribution towards negative anomalies. Capturing the SM anomaly over the Boreal Shield ecozone is particularly important since this ecozone is one of the largest and most diverse in-terms of species and has a wide range of climatic and ecosystem conditions [7,13]. The anomaly distributions illustrated in Figures 6d-f and 5g-i suggest the SMOS SM product was sensitive to dry surface conditions in advance of fire in this ecozone.
Remote Sens. 2020, 12, x FOR PEER REVIEW 9 of 14 Distributions of fires with burned area greater than 3600ha are illustrated in Figure 6. Histograms could only be created for three ecozones: the Taiga Plain, the Boreal Shield, and the Taiga Shield. Overall, the figure suggests that SMOS SM anomaly observable 3, 5, and 7 days prior to the fire particularly over the Boreal Shield captures the drying surface conditions prior to the occurrence of fire, as indicated by the shift in the distribution towards negative anomalies. Capturing the SM anomaly over the Boreal Shield ecozone is particularly important since this ecozone is one of the largest and most diverse in-terms of species and has a wide range of climatic and ecosystem conditions [7,13]. The anomaly distributions illustrated in Figure 6d-f,5g-i suggest the SMOS SM product was sensitive to dry surface conditions in advance of fire in this ecozone.

4.Discussion
The results of this study demonstrated that there is skewness in the SMOS observed soil moisture distributions towards negative anomalies prior to the fire ignition for the majority of ecoregions. By comparing the histograms of respective ecozones in Figure 4,5, we can see that Figure  5 shows a relatively greater positive skewness than Figure 4, which implies SMOS SM anomalies are sensitive to the drier surface conditions prior to larger wildfires. Within the regions analyzed, the SMOS product includes quality flags on the retrieval. The quality flags denote potential issues associated with the product. These flags may be related to regions of high topographic variation, denser forest canopies, or rainfall during the retrieval; when these flags are present, it is recommended that the data be interpreted with caution. While the present study does not explicitly separate the influence of various science flags on SMOS SM distributions over the boreal forest, the majority of the regions analyzed included science flags that suggest that the retrieval is impacted by vegetation. Although there are some studies that demonstrate that SMOS has some correlation with ground observations of soil moisture in the boreal forest in Canada, these studies are limited to only one or two pixels (e.g., 50,54). Therefore, it was somewhat of an open question regarding the value of product in this environment for the purposes of identifying surface drying prior to wildfire initiation.

Discussion
The results of this study demonstrated that there is skewness in the SMOS observed soil moisture distributions towards negative anomalies prior to the fire ignition for the majority of ecoregions. By comparing the histograms of respective ecozones in Figures 4 and 5, we can see that Figure 5 shows a relatively greater positive skewness than Figure 4, which implies SMOS SM anomalies are sensitive to the drier surface conditions prior to larger wildfires. Within the regions analyzed, the SMOS product includes quality flags on the retrieval. The quality flags denote potential issues associated with the product. These flags may be related to regions of high topographic variation, denser forest canopies, or rainfall during the retrieval; when these flags are present, it is recommended that the data be interpreted with caution. While the present study does not explicitly separate the influence of various science flags on SMOS SM distributions over the boreal forest, the majority of the regions analyzed included science flags that suggest that the retrieval is impacted by vegetation. Although there are some studies that demonstrate that SMOS has some correlation with ground observations of soil moisture in the boreal forest in Canada, these studies are limited to only one or two pixels (e.g., 50,54). Therefore, it was somewhat of an open question regarding the value of product in this environment for the purposes of identifying surface drying prior to wildfire initiation.
We recognize that there are several limitations to this study overall. As illustrated in our results (Figures 4-6), there are a significant number of fires that also initiated under wetter soil moisture anomalies. Soil moisture is among a multitude of factors necessary for wildfire; therefore, further analysis of vegetation conditions, fine fuel abundance, and connectivity would be necessary to further outline the role of satellite observed soil moisture anomalies on individual fires. Furthermore, it is very likely a soil moisture threshold may exist in several regions; below this critical moisture threshold, the likelihood of fire ignition is possible [66]; in many cases, this threshold may not be associated with a negative soil moisture anomaly. This study focused on retrieved SMOS SM products (e.g., SMOS SM Level-2 used in this study), in regions of high vegetation where the product may be impacted by significant vegetation water content. Therefore, using a forward model better adapted to each ecozone in conjunction with a land surface model, which assimilates the satellite brightness temperature (along with site specific parameters and uncertainties), may provide better soil moisture estimates. Finally, the spatial and temporal resolution of data could be improved by using advanced downscaling techniques, which may result in better accuracy pin-pointing the hot-spots when also combined with vegetation density.

Conclusions
It is well known that dry land surface conditions precede wildfire and that methods to identify potential wildfire hot-spots in advance are useful for predicting wildfire occurrence. However, it is not known the degree to which satellites sensitive to SM are useful for this purpose, particularly over the various boreal ecozones within Canada. In this study, surface SM conditions preceding wildfire occurrences across various Canadian ecozones were analyzed for the period 2010-2017, using SMOS satellite-derived SM data. Results show that, over the majority of the ecozones, histograms created from SM anomalies at 3, 5, and 7 days prior to the fire show a negative shift suggested that the satellite-based SM products were sensitive to surface drying prior to the occurrence of the wildfire. We would recommend further evaluation of these products over Canada, particularly if they are integrated into operational models.