Satellite Observations of Fire Activity in Relation to Biophysical Forcing Effect of Land Surface Temperature in Mediterranean Climate

: The present work is aimed at gaining more knowledge on the nature of the relation between land surface temperature (LST) as a biophysical parameter, which is related to the coupled effect of the energy and water cycles, and ﬁre activity over Bulgaria, in the Eastern Mediterranean. In the ecosystems of this area, prolonged droughts and heat waves create preconditions in the land surface state that increase the frequency and intensity of landscape ﬁres. The relationships between the spatial–temporal variability of LST and ﬁre activity modulated by land cover types and Soil Moisture Availability (SMA) are quantiﬁed. Long-term (2007–2018) datasets derived from geostationary MSG satellite observations are used: LST retrieved by the LSASAF LST product; ﬁre activity assessed by the LSASAF FRP-Pixel product. All ﬁres in the period of July–September occur in days associated with positive LST anomalies. Exponential regression models ﬁt the link between LST monthly means, LST positive anomalies, LST-T2 (as a ﬁrst proxy of sensible heat exchange with atmosphere), and FRP ﬁre characteristics (number of detections; released energy FRP, MW) at high correlations. The values of biophysical drivers, at which the maximum FRP (MW) might be expected at the corresponding probability level, are identiﬁed. Results suggest that the biophysical index LST is sensitive to the changes in the dynamics of vegetation ﬁre occurrence and severity. Dependences are found for forest, shrubs, and cultivated LCs, which indicate that satellite IR retrievals of radiative temperature is a reliable source of information for vegetation dryness and ﬁre activity. MSG satellite effect vegetation a basis at a regional scale. is on two key aspects of wild ﬁre problems that are not systematically studied in the an evaluation the importance of LST mean and difference precursor/s of vegetated land surface dryness the related pre-ﬁre ﬁre exploration of variability of these of vulnerability in across the delineated into main land quantiﬁcation of ﬁre activity using long-term records of satellite information from Meteosat, ground observations, and SVAT model data of Soil Moisture Availability (SMA) to the land cover as reference data the Eastern Mediterranean region The speciﬁc objectives be


Introduction
Fire is a global phenomenon closely linked to climate variability and human practices with critical regional implications [1,2]. It is an important process in the modulation of the Earth system through the links between weather, climate, and vegetation and has the potential to impact on the global climate system by changing the ability of the surface to absorb and emit energy [3,4].
Climatic variability is a major driver of fire in many terrestrial ecosystems, as reflected in Bradstock's conceptual model of four climatic "switches" that influence fire regimes by controlling fuel amount, fuel moisture, and fire weather at contrasting temporal scales [1]. Fire regimes are also affected by other controls such as landscape-scale patterns of vegetation, topography, and human activities [5]. Climate is connected to fires at two distinct temporal scales [6]. Short-term climatic anomalies (from months to years) affect fires by modifying vegetation growth and fuel moisture before the fire and by influencing weather during the period of fire spread. In addition, climate has more indirect, long-term (decadal or longer) effects on the distribution of major vegetation types, which in turn constrain the landscape-scale mosaics of fuels and vegetation.
Soil moisture (SM) has been identified as a key variable in understanding and predicting wildfire hazards [7][8][9][10][11]. Soil moisture, defined as the water contained in the unsaturated soil zone [12], not only influences vegetation growth conditions and consequently the accumulation of wildfire fuel, but also determines the vegetation moisture content and hence the flammability of the vegetation [7,9,13].
Land surface temperature (LST) is one of the most important parameters controlling physical processes responsible for the land surface balance of water, energy, and CO 2 [14,15]. In the context of wildfire studies, the LST consideration has two main aspects: first, as a pre-fire indicator of land surface state (SM deficit in root zone and vegetation dryness, evapotranspiration, leaf temperature), and second, as a post-fire characteristic of fireinduced environmental changes. The variations in the spatial distribution of LST as a result of fires are usually the focus of research works that characterize fire intensity and burn severity. Higher post-fire LST of burned areas has been reported in remotely-sensed images [3,4,[16][17][18], that is accounted mainly due to a decrease in transpiration and an increase in the Bowen ratio (β = sensible heating/latent heating) [19].
Although relationships between drought and fire seem quite interrelated, only a few studies have explored LST as a pre-fire indicator [20,21]. Short-term analyses based on LST daily anomalies have been performed to predict fire occurrence [20].
The vegetation cover partitions the incoming solar radiation into sensible and latent heat fluxes depending on its structure as well as affects the surface roughness, which can in turn alter heat and moisture transport; the type of vegetation and its seasonal dynamics affect the land's vulnerability to fire ignition and spread. That is why spatial-temporal patterns of LST can contribute to the monitoring of processes that structure ecosystem development and may be associated with fire occurrences to help fire management. Satellite measurements are a source of information for the accurate estimation of LST at the global and regional level, thus helping to evaluate the land surface-atmosphere exchange processes and can serve as a valuable metric of surface state [4,14].
In this study, we consider data for skin surface temperature retrieved by satellite measurements to explore the significance of LST as a biophysical driver, which (in combination with SM, air temperature, and humidity) controlled long-term wildfire activity during the study period of 2007-2018 over the Eastern Mediterranean region, as seen in coherent satellite active fire observations. Studies have used satellite remote sensing for fire monitoring from the early 1990s. With the rapid development of remote sensing technology, satellites are increasingly used for fire monitoring [22][23][24][25][26] and applied at large scales. In recent years, many studies have been committed to identifying the spatial and temporal characteristics of fires from the global and regional perspectives with the use of satellite data [27][28][29][30][31]. All these studies were performed by using observations by low earth orbit (LEO) satellite sensors with the Very High Resolution Radiometer (AVHRR), Moderate Resolution Imaging Spectroradiometer (MODIS), and Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER).
Geostationary (GEO) fire products were first generated over the Americas using the Geostationary Operational Environmental Satellite Visible Infrared Spin Scan Radiometer Atmospheric Sounder (GOES-VAS), e.g., in [32], and this led to the development of the long-standing GOES WildFire Automated Biomass Burning Algorithm (GOES WFABBA) product [33]. Roberts et al. [34] first demonstrated the retrieval of Fire Radiative Power (FRP) from GEO data and developed a full "fire thermal anomaly" active fire detection and FRP retrieval algorithm for GEO systems. This was initially applied to data from Meteosat Second Generation [35], and an operational version is now used to generate a series of geostationary active fire detection and FRP retrieval products that span much of the globe, including Meteosat over Africa and Europe [36,37], GOES-East and West over the Americas [38], and Himawari over Asia [39]. Wooster et al. [40] summarized the history of achievements in the field of active fire remote sensing, reviewed the physical basis of the approaches used, the nature of the active fire detection and characterization techniques deployed, and highlighted some of the key current capabilities and applications.
However, since every observation by the satellite sensors lasts only an instant, what they actually measured was the rate of release of FRE per unit of time or Fire Radiative Power (FRP) in MW [24,41,42]. It was shown that the FRE released by a fire is directly proportional to the biomass consumed as well as the smoke emitted [43]. Measurement of fire radiative energy (FRE) release rates or power (FRP) from space offers opportunities for the characterization of biomass burning and emissions in a quantitative manner. In essence, it has enabled the rating of fires as well as an estimation of regional FRP fluxes, which reflect the relative concentrations of biomass burning activities.
Mediterranean regions are some of the most affected by wildfires, and remote-sensed information about fire activity as provided by the SEVIRI instrument on board Meteosat-8 is especially valuable for forest and civil protection activities [44]. On the other hand, the Mediterranean is one of the most responsive regions to climate change, as evidenced by "pronounced warming" and significant decreases in spring and summer precipitation, which have led to regime shifts toward more arid climates [45,46]. Triggered by large-scale atmospheric forcing, Mediterranean regional heat waves are often amplified by surface preconditioning such as negative soil moisture anomalies and vegetation stress [47].
The research works of Sifakis et al. [48], Amraoui et al. [44], and Di Biase and Laneve [49] were aimed at understanding the spatial and temporal patterns of landscape fires in the Mediterranean, applying different remote sensing data from GEO sensors. Compared to LEO systems, GEO products offer higher temporal resolutions but coarser spatial resolutions, and each sensor only provides data over a specific region of the Earth. The observations from the geostationary meteorological satellite MSG provide a valuable source of information about fire occurrences in the Mediterranean region. Sifakis et al. [48] reported that during the summer period of 2007, MSG-SEVIRI data successfully detected 82% of the fire events in Greek territory, with less than 1% false alarms. Using data from MSG-1 satellite, Amraoui et al. [44] performed an analysis of the spatial distribution of fire events in the months of July and August during the period of 2007-2009. Around half of the fire pixels were detected in croplands, and the remaining half was evenly distributed between forest and shrub. Based on the analysis of the low, mid, and upper atmospheric fields of geopotential, temperature, relative humidity, and wind in two extreme events of fire activity that struck Greece and Italy on 24-25 July and 22-27 August 2007, they suggested a conceptual model for meteorological conditions that favor the occurrence of severe wildfire episodes in Italy and the Balkan Peninsula.
The current study is the first one to use the FRP product retrieved by geostationary satellite observations to investigate fire activity in a recent long-term period for the Eastern Mediterranean.
In Mediterranean ecosystems, prolonged droughts and heat waves have created the preconditions for the increasing frequency and intensity of forest fires [50], the underlying mechanism being the reduction of live and dead fuel moisture content as a response of the soil-plant system to the increased vapor pressure deficit [51]. Vegetation response varies with species as well as with forest structure and soil/terrain characteristics, and it is determined by evapotranspiration. The moisture of dead fuels is affected by weather variations as well, and it is regulated through evaporation [52]. The climate and weather patterns of the Mediterranean region highlight the value of having an improved understanding of the relationships between drought and wildfire; more specifically, an understanding of how drought is related to fire danger outputs. A number of studies have examined the link between drought indicators and wildfire occurrence [53][54][55]. Many drought indices (see [56]) are driven by stand and climate variables of precipitation and/or temperature, but more recent developments include variables that express conditions at the land surfaceatmosphere interface such as vegetation health [57], soil moisture deficit [58][59][60][61][62], actual evapotranspiration [63], and evaporative demand [56,64,65]. These last indexes are important because they can reflect the outcome of the coupled physical processes of energy-water cycles on the land surface.
In our previous work [66], long-term LST anomalies retrieved by MSG SIVIRI were used to identify drought-prone areas in the Eastern Mediterranean. Song [21] showed that there is no linear relationship between LST anomaly and fire occurrences, but the time series of LST anomaly before fire events show a significant trend. Resolving LST-fire interactions facilitate fire predictions and the issuance of early warnings.
We use long-term spatially and temporally consistent satellite observations from the SEVIRI sensor of geostationary MSG satellite to assess the biophysical forcing effect of vegetation fires on a climatic basis at a regional scale. The study is focused on two key aspects of wild fire problems that are not systematically studied in the literature: (i) an evaluation of the importance of LST (monthly mean values, anomalies, and difference with air temperature, T2m) as precursor/s of vegetated land surface dryness and the related pre-fire signals of vegetation stress and fire occurrence; (ii) an exploration of the variability of these biophysical drivers of vulnerability in biomass burning across the climatic gradient, delineated into main land cover types, and the quantification of relations with fire activity using long-term records of satellite information from Meteosat, ground observations, and SVAT model data of Soil Moisture Availability (SMA) to the land cover as reference data over the Eastern Mediterranean region (Bulgaria). The specific objectives to be met are the following:

•
To characterize spatial-temporal patterns of fire activity (July-September) using longterm satellite data records from SEVIRI observations (2004-2019) in terms of the number of biomass burning detections and the severity of burning (FRP, MW) according to the LSASAF FRP-Pixel product; • To use LSASAF LST product data to statistically investigate and evaluate the relationship of the biophysical parameters of LST, LST anomaly, and the difference between skin and air temperatures (LST-T2m) to the occurrence and severity of wild fires on a short-term climatic basis (2007-2018); • To characterize the wild fire vulnerability of the main vegetation types (forest, shrublands, cultivated) in relation to the LST and SMA warm and dry anomalies.

Materials and Methods
This work is based on the methodology developed in our work [66], adapted to assess environmental control on fire activity and its spatial-temporal variability by using satellite observations from geostationary Meteosat. As a biophysical index related to the forcing of fire activity, the satellite-derived IR skin temperature from SEVIRI observations with the EUMETSAT LSASAF LST product is proposed. Based on long-term data records, LST mean summer monthly values, and their anomalies and deviation from the air temperature, T2m were applied as indicators of the coupled energy and water cycles. Fire-promoting SM anomaly patterns derived by the SVAT modeling approach were used as a reference for the evaluation of the LST-fire intertwining. The pattern of fire activity over Bulgaria was characterized by the evaluation of SEVIRI-based satellite detections of biomass burning according to the LSASAF FRP-Pixel product. Benefiting from the geostationary satellites' high temporal resolution, the relationship between LST and fire regimes was studied and the underlying biophysical processes were explored.
Statistical comparative analyses of long-term data records of the LSASAF product, SVAT meteorological model outputs, and ground observations (meteorological parameters and actual fires) were applied following the methodology in [66], which was further extended (in Section 2.5). Overlaid graphical analyses of the summer season (June-September 2007-2018) dynamics of biophysical indexes in relation to fire characteristics (number of detections and severity) were performed to study the sensitivity of the approach. All datasets taken from information sources with different spatial resolutions were placed over the grid of the ECMWF global NWP model, IFS version O1280 (about 9 km spatial resolution). Only values with 1400 points that cover the region of Bulgaria were considered.

The LSASAF FRP-PIXEL Product
In order to characterize fire activity in this study, radiative energy released by fire events, as provided by the LSA SAF Fire Radiative Power-Pixel (FRP) product [36,67], was used. The LSASAF FRP-PIXEL product provided information on the location, timing, and fire radiative power (FRP, in MW) output of landscape fires ("wildfires") detected every 15 min across the full Meteosat disk at the native spatial resolution of the SEVIRI sensor. The FRP provided information on the measured radiant heat output of detected fires; the heat output was a direct result of the combustion process whereby carbon-based fuel is oxidized to CO 2 with the release of a certain "heat yield" [24]. FRP calculation relies on the Middle Infrared (MIR) method and assumes FRP to be proportional to the difference between the observed fire pixel radiance measured from the MSG SEVIRI radiometer and the "background" radiance that would have been observed at the same location in the absence of fire. The LSASAF algorithm for deriving FRP from SEVIRI radiance measurements is fully described in [37]. The type of fire detection algorithm that is proposed for use with SEVIRI is based on the principles applied to generate active fire detections within the MODIS Fire Products [24].
While this FRP methodology offers many advantages, among the uncertainness reported, two factors are likely to reduce the satellite-measured FRP that is emitted by the fire:

•
The first concerns surface fires in forests, where an unknown amount of radiant energy may be intercepted (scattered and absorbed) by the forest canopy. • Second, although atmospheric effects perturb MIR wavelength observation far less than those in shorter wavelengths, allowing FRP retrieval through even dense smoke and plumes, the radiative impact of the absorptive black carbon released during combustion may result in some underestimation of the FRP by the satellite measurement.
The FRP algorithm is subject to other sources of uncertainties, essentially due to the characteristics of SEVIRI, as described in [37]. Among these limitations, one can state that the FRP-derived value is quite sensitive to the fire location within the pixel. A fire located at the center of the instrument's instantaneous field-of-view will elevate the pixel temperature much more than a fire located far away from its center. Additionally, abnormally low radiances might be observed surrounding a fire pixel due to the negative lobes of the point spread function. As a result, the background temperature might be colder, thereby increasing the estimated FRP value. The rather coarse spatial resolution of SEVIRI may cause the non-detection of smaller/less intense fires that MODIS can detect, for instance, but SEVIRI misses. The non-linearity of the SEVIRI 3.9 channel above 310 K and the saturation of that band above 335 K are responsible for the error in the FRP estimation. As a result, FRP values derived from SEVIRI underestimate those derived from MODIS by about 40% [34].
A significant remote sensing constraint for fires is the presence of cloudiness in fire weather situations, as all satellite methods used for monitoring fires have limitations that tend to cause important biases in the final product. The algorithms cannot detect existing fires due to the elimination of cloudy pixels or the assumption of some detection in partially cloudy pixels as "false" fire detections based on "contextual" tests.

Biophysical Indexes
In our work [66], the spatial-temporal variability of land surface state dry anomalies was characterized by the land surface temperature with the use of LST retrievals from the MSG satellite measurements and SMA. Using long-term data records, the same biophysical drivers of drought occurrence have been further explored here in terms of their forcing effects on fire activity through the related biophysical processes. In addition to land surface radiometric temperature characterized by its monthly values and calculated anomalies, the blended parameter, which is the difference between skin surface temperature and air temperature at 2 m above the earth's surface (LST-T2m), is proposed to be used as a first proxy of sensible heat flux exchange. The patterns of these biophysical indexes were defined as the mean monthly values for June, July, August, and September. The SMA anomalies during the summer months were considered in parallel to explore the fire forcing effect of soil moisture deficit.
To quantify the relations, statistical analyses and modeling were applied. Since the considered biophysical parameters and the procedure of their application have been explained in detail in [66], only the main points are marked here, and the new elements are further explained in Section 2.5.

The LSASAF LST Product
For the purposes of the current study, the LSASAF LST product (MLST, LSA-001) has been applied [68]. The LST retrieval was based on clear-sky measurements from the MSG system in the thermal infrared window (MSG/SEVIRI channels IR10.8 and IR12.0), every 15 min within the area covered by the MSG disk. Theoretically, LST values could be determined from MSG 96 times per day, but in practice, fewer observations were available due to cloud cover.
The 12-year time series for the region of Bulgaria for the period of June-September (2007-2018), which included only days/locations with satellite fire detections, was constructed. Datasets for LSASAF LST at 0900 UTC and at 1200 UTC (averaged over 5 measurements within ±30 min of these two time slots to avoid the limitations of cloudy pixels) were inferred at MSG pixel bases. For some climatic assessments, the averaged 3 × 3 pixels LST values around the grid points were also processed. The difference between LST and the air temperature, T2m (observations at the SYNOP network at 0900 UTC and 1200 UTC) was also constructed for the same studied period.

SVAT Model and SMAI
For the identification of fire-promoting SM anomaly patterns, the "SVAT_bg" model developed at the National Institute of Meteorology and Hydrology (NIMH) of Bulgaria was used [61,69,70]. This is a simple 1D site-scale meteorological model that exploits the concept of one layer of vegetated land surface and two levels of moisture availability along the root zone depth. One of the main "SVAT_bg" model output parameters is soil moisture (SM) for different root zone soil depths (20, 50, 100 cm). For assessing land surface state, the SMA concept was adopted to serve as an information source for "warnings" of environmental constraints. Based on the "SVAT_bg" model-derived SM, a quantitative SMA index, SMAI was developed and operationally calculated at a site scale for the region of each NIMH synoptic station. At each point of the used grid, the SMAI values at the nearest synoptic station were considered for the calculation of the anomalous water content in the unsaturated root zone. SMAI was designed as a 6-level threshold scheme to account for moistening conditions [66]. The site-scale assessment of mean monthly SMAI anomalies at synoptic stations was compared with the overlaid FRP Pixel detections from corresponding MSG pixels covering the region of the station applying the statistical analyses.

Ground Observations of Forest Fires
Ground observations of the number of actual forest fires were considered according to the database of the State Forest Agency of Bulgaria (SFA), the responsible national institution, and were used for comparison with satellite detections.

Target Region and Land Cover
The study area includes the whole territory of Bulgaria, Southeastern Europe, located within latitude circles of 40.25 N and 45.0 N and meridians of 20.5 N and 29.5 E. The region falls under Mediterranean climate influences, characterized by dry summers and mild, wet winters, both with irregular precipitation distribution. Information about land covers (LC) provided by the ESA-CCI Land Cover Map product at a spatial resolution of 300 m [71] was used. Maps were updated on an annual basis. The LC classification system from ESA-CCI includes 24 vegetation types. For the purposes of this study, we reclassified these types into the following three main vegetation groups ( Figure 1): forests (classes 50,60,61,62,71,72,80,81,82,90), scrubland (11,12,40,100,110,120,121,122), and cultivated areas (10,20,30,130). Other land cover types, including permanent snow and ice (220) and barren and water bodies (200, 210), where there are limited fires, were not analyzed in our study. Urban areas (190) were also excluded.

Numerical Analyses
To characterize the fire distribution pattern, a dataset of FRP detections (location, timing, and energy) for the period of June-September (2004-2019) was constructed [72,73]. All pixels with at least one fire detection were considered. Fire regime was characterized by the frequency of detections and severity according to the released energy from biomass burning (FRP fire radiative energy, MW). Each fire pixel was associated with the coordinates of the nearest point on the grid, where LST, LST anomalies, and (LST-T2m) had been derived. The procedure was performed for monthly (June, July, August, and September) accumulated energy released from fires in the corresponding grid point.
In order to validate the adopted methodology, monthly mean anomalies of LST and SMA that had been inferred in [66] and visualized in color-coded maps were evaluated for consistency with FRP fire pixel spatial distribution over Bulgaria by using qualitative comparative analysis. Examples for July 2007 are presented in Figure 2. Figure 2a shows that the higher positive LST anomalies (in red) correspond to a higher number of FRPdetected fire pixels. The stronger negative SMA anomalies (indicated by reddish colors) are related to a higher number of fire detections (Figure 2b). The consistency of the data allows quantitative studies to be performed further.
Based on long-term records (June-September 2007-2018), stochastic graphical analysis was performed. The consistency in the behavior of the fire activity and the biophysical drivers LST, LST anomalies, (LST-T2m), and SMA anomalies was analyzed in terms of their mean, spatiotemporal variability on a monthly and annual basis, as well as their anomalous distribution and relations. The summer seasonal dynamics of biophysical conditions for different LC vegetation types were studied and their relation to vulnerability of fire ignition/spread was evaluated.
A graphical description of the locality, the spread and skewness groups of numerical data of FRP detections, released energy from biomass burning, LST, LST anomalies, and the (LST-T2m) temperature difference was performed through their boxplot quartiles. Two types of regression analyses were performed; the first one made use of estimates by the method of least squares of the conditional mean of the response variable across values of the predictor variables. The second method used Quantile Regression (QR) estimates of the conditional median of the response variable, which has received increasing attention in recent years and is applied in many areas, including data analysis in the natural sciences. To identify regions most vulnerable to fires in Bulgaria, the spatial pattern of the accumulated fire numbers detected by the LSASAF FRP-Pixel product for the period of June-September 2004-2019 was considered. A statistical evaluation of the spatial distribution of fire pixels was performed by applying the Mann-Kendall statistic test. To perform all these comparisons, the R language for statistical computing was used [74].
All numerical evaluations were performed for the three main vegetation groups: forest, shrubland, and cultivated, as well as considered all together in a sample "all land cover". Such a reclassification into three LC categories has been used in other studies of fire activity over the Mediterranean using satellite data, e.g., in [44].
A summary of all data used for comparative analyses in the current study are presented in Table 1.

Results
Fire regimes were described through statistical distributions of frequency and severity over the studied area in the summer season during a 16-year time period. Using 12 years of data for LST and related biophysical parameters, the environmental determinants of fire regimes were assessed by exploring how environmental drivers operating over a range of scales affected the spatial and temporal patterns of these fires.

Active Fire Monitoring from Space
The results from the study of fire activity derived from satellite observations over Bulgaria are presented in maps of accumulated fire detections for each month of the July-September period and each year of the whole 2004-2019 period. Examples for selected years are shown in Figure 3 (see also [72,73]). The fire severity was assessed following the color-indicated released      Figure 5a shows the trend lines for LST and the energy released from biomass burning for a sample of all land cover types without classifying them ("all LC"); the linear regression line fits the relation between the parameters at high R 2 values of 0.68. The values of R 2 for forest-shrubs-cultivated types vary between 0.61 and 0.68, and the validity of these conclusions passed the significant test from 0.001 up to 0.003 levels (Figure 5b-d). The FRP (blue lines) decreases with the decrease in LST (red lines) within individual diapason for biomes (e.g., 32-29 • C for forest and 36-34 • C for shrubs/cultivated, as seen in Figure 5c,d). A high level of fit between the trends in the number of fire pixels and the LST as a biophysical driver is also experienced (Supplementary Material). In this case R 2 varies within 0.67-0.62 for shrubs and forest vegetation types, at high significant levels, indicating a descending trend for the fire accuracy along with the LST decrease.
To clarify the nature of this linear structure, it is examined by linear regression on the individual components in each data pair. Again, the results are similar (values in the bottom left panel in Figure 5), although at the limit of reliability. For at least one group in each pair, a linear model with a similar slope of about −0.25 is obtained, and its significance level is slightly above 0.10. Therefore, for at least one pair, a significant confirmation of the estimated downward trend line was obtained.

Statistical Analyses
Monthly wildfire activity (2007-2018) characterized by the number of fire pixel according to the LSASAF FRP-Pixel product, and the ratio between the total energy released (MW) and the total number of fire pixels detected for a specific month and LC type are shown in Table 2. The number of fire pixels in forest and shrubs in June represents 6.5-7.6% percent of the fire pixels detected in July. Fires in cultivated LC in June seem lower by about 5% than their amount in July due to the specific vegetation in the managed LC areas still being in the growing phase. For September, the wildfires (forest, shrubs) from the forest fires in August are at 37% and 50% from the shrub fires in the same month. However, in the cultivated LC, the relatively large number of fire pixels are detected in September (75% from their number in August), which is likely due to some controlled agriculture burning activities. The ratio between the total FRP energy released and the accumulated number of fire pixels at a specific LC type was introduced to serve as a measure of mean Detected Fire Intensity (DFI) observed by the satellite. The DFI from the forest fire pixels is the lowest one, ranging from 35 MW/pixel in June to 96 MW/pixel in August. For shrubs and cultivated LC types, the DFI range is 72-101 MW/pixel and 70-107 MW/pixel, respectively.
The satellite-derived FRP are likely to estimate a specific portion of the actual fireemitted FRP (see Section 2.1). Data in Table 2 show that shrubs LC is characterized by the largest amount of total energy released by fires for the summer season. In order to assess how the portion of satellite-detected fire energy depends on the vegetation type during the fire season, the DFI values of forest and cultivated LC are divided by the DFI obtained for shrubs in each month (Table 3). It can be seen that the ratio between the DFI parameters of forest and shrubs is 0.48 in June, then increases gradually to 0.95 in August and becomes 0.89 in September. At the same time, the ratio between the DFI values for cultivated-shrubs LC is close to 1.0, ranging from 0.96 in June to 1.06 in August-September. This result comes from the different ability of satellite observations (FRP product) to estimate actual emitted energy from wildfires at different LC types, and this also depends on the biophysical status of the vegetation along the fire season.

Box Plots Analyses
Boxplots of the distribution of fire radiative energy (MW) according to the FRP-Pixel product and biophysical indexes for forest and shrubs LCs are shown in Figures 6 and 7.  The plots in Figure 6 show different aggregation of data by months. Through a visual comparison of the corresponding notches, it can be seen that there is sufficient reason to accept the statistical hypothesis that the medians are significantly different. The red lines show the median and first and fourth quartile averages of the same statistics for all individual boxplots presented and are used to visually assess the degree of difference depending on the month.
The boxplots in Figure 6a,b demonstrate the variability of the monthly mean values of FRP (MW) for forest and shrubs LC during the summer. The statistical analyses confirm that for July, the ratio between the median of total FRP box plots for shrubs and forest fires is higher than the corresponding ratio for the DFI parameter in Table 3: shrub fires were detected with about 50% higher released FRP (MW), while the DFI parameter was only 35% larger than this for forest fires. This effect seems to be due to higher energy emitted by some forest fires in July.
On the contrary, for June, the medians of total FRP are not much different for these two LC types (Figure 6a,b), but the DFI parameter is twice higher for shrubs than the corresponding DFI value obtained for forest fires (Table 3). This suggests the occurrence of more forest fires under the canopy of dry dead forest fuel with much less energy measured by the satellite and that a certain amount of radiant energy may have been intercepted (scattered and absorbed) by the forest canopy.
The summer course of the FRP median corresponds well to the LST statistical behavior for the forest LC (Figure 6a,c). In June, optimal SMA is present in the deep forest root zone that leads to fully open stomata conductance and unlimited evapotranspiration, leading to decreasing canopy temperature and minimum LST. For shrubs (Figure 6b,d), there is a disagreement: the LST median maximum appears in August, but the maximum of the total FRP median is in July. The reason could be that in shrub lands, the SMA easily reaches the wilting point, in which the stomata conductance strongly decreases in reaction to the increasing water stress conditions through the transpiration regulation mechanism. These lead to an increasing shrub canopy LST in August.
All fires in the period of July-September occur in days associated with positive values of LST anomalies (averaged in the area of 3 × 3 pixels) around the same location (Figure 7a,b). All 25% of the fire detections in the days associated with negative anomalies of LST occur in June in forest and shrubs LCs. This result suggests that in June, some fires may mostly occur in dead fuel where the canopy temperature is below the historical average for today's day-of-year.
The trend of increasing (LST-T2m) from June to September (Figure 7c,d) corresponds to the gradual decrease of SMA to vegetation cover on the one hand, and to the increasing LST on the other. This is consistent with the result on the disagreement in the behavior of LST versus FRP in August for shrubs LC (Figure 6b,d).

Quantile Regression
Results from the quantile regression analysis of the relationship between log FRP (MW) and biophysical parameters (LST, LST anomalies, (LST-T2m), and SMA anomalies) are presented on Figure 8a-d.
These graphs can quite comprehensively characterize the dependences between each pair of quantities. QR examples of the "all LC" type are here presented. Figure 8a shows at which LST values the maximum FRP (MW) might be expected and the corresponding probability level the FRP might reach the maximum. Accordingly, the QR in Figure 8a shows how the log(FRP, MW) maximum depends on the predictor LST, with a confidence of 5%, 25%, and other levels of significance. The QR lines indicate that the high values of the log(FRP, MW) > 5 can be expected to occur at LST = 35 • C with a 25% probability level; extreme values of log(FRP, MW) > 7 are possible with 5% reliability at LST above 40 • C, and with 10% reliability above 45 • C.
On the plot of Figure 8b, the probability level of FRP maxima at a specific temperature difference (LST-T2m) is illustrated. Very high values of log(FRP, MW) > 6 can be expected to occur at (LST-T2m) = 5 • C with 10% reliability; extreme values of log(FRP, MW) > 7 are possible at a 5% significance level for (LST-T2m) above 12 • C. This illustrates how the log(FRP, MW) maximum depends on the predictor (LST-T2m), with confidence levels of 5%, 10%, and other levels of significance.
The QR between FRP and LST anomalies (Figure 8c) shows that very high values of log(FRP, MW) > 6 can be expected to occur at negative LST anomalies, with a probability of less than 10%. Extreme values of log(FRP, MW) > 7 are not possible even at a 5% significance level, with negative LST anomalies. In parallel, higher FRP energy occurs at negative SMAI anomalies (Figure 8d). The decreasing trend of energy released, along with the decreasing values of negative or low positive anomalies, i.e., decreasing drought (agri-cultural and ecological), at corresponding probability levels are shown. Extreme values of log(FRP, MW) > 7 are possible at 10% significance with negative SMAI anomalies. Such QR graphs have been developed for forest, shrubs, and cultivated vegetation types and confirm similar relations for the fire activity over the whole country. The results obtained for some smaller samples over specific regions of Bulgaria show an intersection of the probability lines that is an indication of insufficient data.

Correlation Analyses
The results of the least-squares method applied for the quantitative evaluation of the relation between biogeophysical indexes (LST, LST anomaly, (LST−T2m) temperature difference, including the SMAI anomaly at 20, 50, 100 cm depth) and the severity of wildfire events are presented in Table 4. Monthly means LST at around 0900 UTC, which is 1200 local time (selected as a less cloudy period), are considered. For all land cover types (all biomes included in a sample), exponential regression lines fit the relations at high R 2 values and high significance levels ( Figure 9).  The SMA index is especially efficient as a measure of fuel dryness, taken at 50 cm or 100 cm soil depth, for all cover types in the period of July-August. With the progressive establishment of dry conditions during August and September, the correlation decreases when the SMA at 20 cm and 50 cm becomes fully exhausted, and the SMAI is equivalent to a constant value of zero.
Exponential regression models better fit the FRP (MW)-LST relations from Table 4. Examples for cultivated LC are shown in Figure 9a-c, indicating high significance levels (p < 0.007) of the dependences for July, August, and September. Similar exponential models for the regression of FRP (MW) versus the LST anomalies, R 2 at about 0.6 (p < 0.003), are presented in (Figure 9d-f). Results confirm a statistically significant high correlation between biophysical drivers (LST, LST anomalies) and the wildfire severity.
All these results suggest that the biophysical index LST and the related parameters LST anomalies and (LST-T2m) are sensitive to the dynamics of the occurrence and severity of vegetation fires in all studied LC types.

Spatial Pattern and Trends
To characterize fire activity, plots of monthly sums of FRP detections are considered. Figure 10 shows a color-coded map of accumulated fire pixels in the period of June-September 2004-2019. The fire numbers are categorized into the following classes: 1-30 (low, yellow color), 30-90 (moderate, orange color), 90-150 (high, red color), and >150 (extreme, dark red). Thus, the spatial distribution of the spots with higher fire activity is revealed. The "+" sign on the map indicates that the trend passed the Mann-Kendall significance test at the 5% level; thus, regions with a statistically significant positive trend in fire activity over the studied period are localized. The long-term trends of fire activity (number of detections) are different when it is delineated into different land cover types (superposition of Figures 1b and 10). These considerations show that there is no positive trend of fire activity in regions with forest land cover types. There are limited hot spots for fires, located in the northern and southeastern parts of Bulgaria, indicating an increasing trend in fire activity that is mostly scattered in the shrublands. Positive trends are also seen in limited areas of cultivated land cover types, but in this case, it is rather a result of human activities than a climate trend.

Discussion
The results of our study confirm the driving role of land surface temperature and the related biophysical parameters in the spatial-temporal evolution of fire activity in Eastern Mediterranean as well as in forcing the vulnerability of forest, shrubs, and cultivated land cover types to wildfire. Based on the 2007-2018 summer time series, the trend analyses (Section 3.2.1) show the synchronized behavior of annual LST values with the released FRP energy (MW) from biomass burning (respectively to the number of fire pixels). The analysis of satellite retrievals from geostationary satellite data show a statistically significant declining trend for the wildfire characteristics, consistent with the course of LST for all considered vegetation types (Figure 3).
The box-plot statistics (Figure 6a,c) show that the median of the FRP (MW) energy released from forest fires in the summer course corresponds well to the behavior of the LST median, which confirms the relationship between the biophysical index LST and fire activity. All fires in the period of July-September occurred in days associated with positive values of LST anomalies, while the other 25% of the fire detections in days associated with negative anomalies of LST occurred in June at forest and shrubs LCs (Figure 7a,b).
A high correlation between the monthly mean values of the biophysical indexes (LST, LST anomalies, LST-T2m, and SMA) and FRP-detected fire pixels was obtained. Exponential regression models fit the relations for all LC types at high significance levels (Section 3.3.3). This result is a step forward from the findings of Song [21], who stressed on the spatial-temporal information provided by LST anomalies (2001-2019 time series over Australia) on fire activity, confirming a significant trend before the fire events but not establishing a linear relationship with fires.
The values of LST, LST anomalies, and LST-T2m, for which the maximum FRP (MW) values might be expected at the corresponding probability level, were statistically estimated (Section 3.3.2) as being able to contribute to fire preparedness activities. For example, extreme values of log(FRP, MW) > 7 were possible at a 5% significance level for (LST-T2m) above 12 • C. A previous study reported that forest fires in August over Bulgaria might occur under land surface drought conditions defined by middy (MSG LST-T2m difference) >12 • C prior to the fire [62].
These dependences obtained for the Eastern Mediterranean confirm and extend the reported results in the literature that land surface state can significantly contribute to the enhancement of weather impacting on the fire environment such as drought (long-term and/or short-term), terrain, and fuel conditions [62,75,76]. Plants constitute the main ignition material in the landscape, and their moisture content plays an important role because it may serve to retard ignition or mitigate the propagation of a fire [77][78][79][80]. For this reason, the fuel moisture content is a common component of fire danger and related fire regime assessments; e.g., in [81,82]. Skin temperature is an important disclosure of fuel moisture content, and the results reported confirm the role of LST as a biophysical index of fire activity. The physical mechanism of LST-SMA relations and related fire occurrence covers a range of coupled biophysical processes: as a moist soil surface dries out, more of the incoming solar energy is reflected, and a larger fraction of the absorbed energy is used to heat the air and soil [83]. The heat flow into the soil increases at first, then decreases as the soil becomes very dry [84]. This results in increasing land surface temperatures, which also influences the rate of drying and evapotranspiration. As a result, high LST is a cause and product of dry periods; in other words, drought begets drought and extreme fire activity.
Land surface temperature explains 80% of the variance in air temperature, and vegetation density also plays an important role in explaining the air temperature variance [85]. Since only data from days with detected fires in the summer months were considered, most of the results of the statistical analyses are related to conditions of limited SMA for the vegetation cover.
In such conditions, the (LST-T2m) parameter is a measure of the evapotranspiration rate, and relations are specified by atmosphere/land weather and climate conditions as well as by physiological and structural features of the LC [86,87] Hence, the trend of increasing (LST-T2m) from June to September (see Figure 7c,d) is related to gradually decreasing evapotranspiration during the summer months due to the gradual decrease of SMA to vegetation cover on the one hand, and/or the increasing LST on the other.
All these results suggest that the biophysical index LST and the related parameters of LST anomalies and (LST-T2m) are sensitive to the dynamics of vegetation fire occurrence and severity. The dependences are valid for forest, shrubs, and cultivated LCs at high significant levels and indicate that satellite IR retrievals of radiative temperature is a reliable source of information for vegetation dryness and fire occurrence. Some existing mismatches are due to highly inhomogeneous land cover from one side and the complicated nature of their biogeophysical relations. LST is directly linked to precipitation, cloudiness, and solar irradiance, while SMAI is influenced by the cumulative effect of these meteorological parameters and with functional links to vegetation. In the case of when the period of moisture depletion is long, the LST-SMA relation is weak because the LST continues to increase; on the other hand, after full SM depletion, the "SVAT_bg" model does not assess any further SM changes. In terms of the model simulations, full SM depletion denotes an interruption of liquid and water vapor flow because SM has reached its constant minimum value (at specific soil/climate), which is set to be the maximum hygroscopic value [69].
In addition, there are some limitations in the remote sensing approach for assessing the relation of fire activity to the distribution of its biophysical drivers. Such cases may be due to different abilities of the satellite observations and the derived FRP product to estimate actual emitted energy from wildfire in different LC types, and this also depends on the biophysical status of the vegetation during the fire season. This ability of FRP is significantly different as compared to wildfires in the forest and shrubs LCs in June (reported here, Table 2 and Figure 6a,b) due to the specific physical properties of the forest canopy, which leads to larger amounts of radiant energy to be intercepted. This in turn makes some of the low-energy forest fires impossible to be detected by the satellite. In August and September, with the decrease in SMA and the fraction of vegetation cover, the canopy biophysical properties changed so that a larger amount of energy released by forest fires could be evaluated by satellite measurements, and most of the forest fires could be detected as well.
Although landscape fires in Mediterranean countries are mainly caused by human activities related to agricultural practices, this study confirms that fire regime components (occurrence and severity) exhibit characteristic spatial and temporal patterns that reflect differences in the relative importance of various environmental drivers. In this regard, an identification of the regions most vulnerable to biomass burning is of importance for some prevention activities. Different land covers have widely differing flammability, which depends on species composition, stand age and density, microclimate, and soil conditions [31,88]. For this reason, a detailed discrimination of the land cover types is needed in studying short-term weather and land surface state influences on the fire ignition and spread. In our study, the LC typology is simplified to consider three types: forests, scrubland, and cultivated; this is because relationships between monthly mean values of the physical indexes and fire activity are explored over long-term 12-year periods. Moreover, the forest LC is considered to be spatially stable throughout all 12 years while more land cover changes may occur in the category scrubland during the observation period, and more detailed classification would be difficult to apply. Our qualitative comparative analysis of the ESA CCI LC maps for 2006-2015 does not show any significant changes in the land cover over Bulgaria for the test period.
Recent concerns about the potential increases of forest fires under climate change underline the importance of fire-climate feedbacks [89][90][91]. In this context, an important result of our study based on a dataset of satellite FRP detections from the last 16 years is that no positive trend in fire activity for forest land cover type has been observed over Bulgaria.
There are significant increasing trends in fire activity for single spots in the low land areas of shrubs LC types in the northern and southeastern part of Bulgaria. A reasonable explanation for these results is that the forests are a more stable environment over a longterm period, which also confirms the validity of our approach for the simplification of the LC variety into three main types for the purposes of the current climate study.

Conclusions
This study quantified the relationship between the spatial-temporal variability of land surface temperature (LST) and fire activity on a short-term climatic scale over the Eastern Mediterranean (Bulgaria), accounting for physical properties such as land cover and soil moisture that, combined with LST, provide a valuable metric for surface state. Based on satellite observations from the geostationary Meteosat LSASAF-FRP-Pixel product (summer season 2007-2018), the distribution of fire activity on a monthly/annual basis in relation to the biophysical forcing effect of LST (assessed by IR MSG satellite measurements) was observed in land cover types of different fire vulnerability (forest, shrubland, cultivated).
A synchronized annual behavior between LST and FRP with a declining trend (2007-2018) in wildfire characteristics was seen. Exponential models fit the relationships of LST monthly means, LST anomalies, and LST-T2, as a first proxy of sensible heat exchange with atmosphere, as well as SMA with FRP fire characteristics (MW). All fires in the period of July-September occurred in days associated with positive LST anomalies.
The values of LST, LST anomalies, and LST-T2m, for which the maximum FRP energy (MW) might be expected at corresponding probability levels, were estimated as being capable of contributing to fire preparedness activities. For example, monthly mean extremes of log(FRP, MW) > 7 could be observed with a 5% reliability at an LST monthly mean above 40 • C and (LST-T2m) above 12 • C, and with a 10% reliability above 45 • C. The reported biophysical forcing effects of LST on vegetation fires provided more understanding of the relationships between drought and wildfire; more specifically, of how drought is related to fire danger outputs. Since forest fires have the potential to affect regional climate through changes in the energy budget [91][92][93], this knowledge is especially important for the Mediterranean region, where land-atmosphere coupling has become one of the important aspects of global environmental change.
To advance the added value of this study by using LST as a biophysical index of drought for fire management, further evidence in support of this approach is observed: First, the use of related FRP products from the MODIS sensor to validate the usefulness of this approach for short-term applications; Second, the evaluation and adaptation of the LST applications in the operational mode for the analyses of real fire situations in the scope of early warnings of fire risk assessment, thus contributing as the final step in the fire management practice.
Acknowledgments: Data processing and preparation of illustration material in Sections 2.5 and 3.1 were performed by Andrey Kulishev within the frame of the LSASAF project implementation. LSASAF is acknowledged for providing the requested LST and FRP data. Data in Section 2.3 have been kindly provided by the State Forest Agency of Bulgaria. Thanks are due to the anonymous reviewers for the valuable suggestions and comments.