Forty-Year Fire History Reconstruction from Landsat Data in Mediterranean Ecosystems of Algeria following International Standards

: Algeria, the main fire hotspot on the southern rim of the Mediterranean Basin, lacks a complete fire dataset with official fire perimeters, and the existing one contains inconsistencies. Preprocessed global and regional burned area (BA) products provide valuable insights into fire patterns, characteristics, and dynamics over time and space, and into their impact on climate change. Nevertheless, they exhibit certain limitations linked with their inherent spatio-temporal resolutions as well as temporal and geographical coverage. To address the need for reliable BA information in Algeria, we systematically reconstructed, validated, and analyzed a 40-year (1984–2023) BA product (NEALGEBA; North Eastern ALGeria Burned Area) at 30 m spatial resolution in the typical Mediterranean ecosystems of this region, following international standards. We used Landsat data and the BA Mapping Tools (BAMTs) in the Google Earth Engine (GEE) to map BAs. The spatial validation of NEALGEBA, performed for 2017 and 2021 using independent 10 m spatial resolution Sentinel-2 reference data, showed overall accuracies > 98.10%; commission and omission errors < 8.20%; Dice coefficients > 91.90%; and relative biases < 3.44%. The temporal validation, however, using MODIS and VIIRS active fire hotspots, emphasized the limitation of Landsat-based BA products in temporal fire reporting accuracy terms. The intercomparison with five readily available BA products for 2017, by using the same validation process, demonstrated the overall outperformance of NEALGEBA. Furthermore, our BA product exhibited the highest correspondence with the ground-based BA estimates. NEALGEBA currently represents the most continuous and reliable time series of BA history at fine spatial resolution for NE Algeria, offering a significant contribution to further national and international fire hazard and impact assessments and acts as a reference dataset for contextualizing future weather extremes, such as the 2023 exceptional heat wave, which we show not to have led to the most extreme fire year over the last four decades.


Introduction
Fire is an intrinsic process of the Earth's system [1], which is projected to increase in burned areas (BAs) owing to climate change and escalating anthropogenic effects [2].Among all forest disturbances, fire is the major forest destructive agent in the Mediterranean Basin [3], where heat wave-related fires may increase in the future [4].When looking across the African side of the Mediterranean Basin, Algeria is the main fire hotspot [5].Throughout history, this country has witnessed an unprecedented series of large extreme fires [6], with a record-breaking heatwave in the summer of 2023, which affected ecological and socioeconomic assets [7,8], even resulting in human casualties.These fires may seriously degrade forest habitats in this country, a large part of which may not be restored [9].However, accurate and spatially explicit BA data about this region, as in several Middle Eastern and North African countries, are scarce or even lacking [10].
Remote sensing has become the most efficient tool for addressing all fire management aspects, including the generation of BA products [11].Unlike ground-based fire datasets that are often biased, incomplete or exhibit inconsistencies [12][13][14], satellite-derived BA products provide spatially and temporally consistent and reliable information about fires on regional and global scales [15].In practice, several pixel-level BA products have been employed in a wide range of research works, including global fire trends [16], characterization of fire regimes [17,18], climate impact on fire patterns [19,20], fire emission modeling [21] and fire model benchmarking [22], and to derive global databases of single fire events [23].
Early global BA products were based on the coarse resolution data from the SPOT-VGT, ERS2-ATSR2, ENVISAT-AATSR, NOAA-AVHRR, PROBA-V, and MODIS sensors [15].In the last few years, major efforts have been made to develop comprehensive global and regional BA products, according to two major programs: the ESA Fire Disturbance Climate Change Initiative (FireCCI) and the NASA MODIS Land Science Team.The current global BA products from the ESA FireCCI project include FireCCI51 (2001-2020; 250 m), which derives from the MODIS surface reflectance imagery coupled with thermal anomaly observations [24]; FireCCILT11 (1982-2018; 5 km) from AVHRR LTDR [25], including new developments; and sensors as in FireCCIS310 (2019; 300 m) from Sentinel-3 SYN coupled with VIIRS active fire hotspots [26].Perspectives with newly delivered medium-resolution sensors have been evaluated regionally, such as FireCCISFD11 and FireCCISFD20 (2016 and 2019, respectively; 20 m) from the Sentinel-2 MSI coupled with the active fire data for sub-Saharan Africa [27,28].On the other hand, MCD64A1 collection 6.1 (2000-present; 500 m) is NASA's current standard global BA product, which derives from MODIS daily surface reflectance imagery combined with MODIS active fire data [29].
Other available coarse-resolution products from different agencies include the Copernicus Climate Change Service Burned Area product, version 1.1 (C3SBA11) (2017-2022; 300 m) from Sentinel-3 OLCI data [30], and the European Forest Fire Information System (EFFIS) BA product (2000-present; 250 m and 20 m) from the MODIS and Sentinel-2 imagery [31].Additional efforts have been made to provide finer-resolution BA products-a major end-user request [32].Landsat-based BA mapping includes the novel 30 m resolution Global Annual Burned Area Maps (GABAM 1985-2019; 30 m), which derived from the Landsat dense time-series data by means of a global automated BA mapping approach [33] in the Google Earth Engine (GEE) [34].In the same context, albeit with limited spatial coverage, other products include the Monitoring Trends in Burn Severity (MTBS 1984-2022; 30 m) across the whole of the U.S. [35], and the Landsat Collections 1 and 2 BA products for CONUS (1984-2022; 30 m) [36].
Although freely accessible to the scientific community with widespread use on different spatio-temporal scales, the above-mentioned BA products exhibit certain limitations.These are mainly caused by the inherent coarse spatial resolution that results in very high omission rates of small burned patches [27,28,[37][38][39], particularly in the Mediterranean Basin where smaller fires happen [10,40], poor temporal fire reporting accuracy to prevent a fire seasonality analysis [33], and limited spatio-temporal coverage [35,36], which restrict their usage in other areas of the globe.
In Algeria, the available ground-based fire dataset provides invaluable information that can hardly be obtained by satellite-based systems.This includes, among others, the exact date and time of ignition/intervention/extinction, burned vegetation type and species, origin, and causes of fires.Nevertheless, this dataset is acknowledged to be incomplete, lacks fire perimeters, and displays discrepancies in fire extent terms.This is attributed to the visual estimation of fire-affected areas, which is often conservative, especially in inaccessible areas.Furthermore, a standardized BA estimation protocol across local forest services in the country is lacking.
Considering these limitations in both national statistics and the performance of existing BA datasets, the development of a reliable and long-term BA product for such an insufficiently investigated part of the Mediterranean Basin would strengthen future research into forest management plans and for understanding Mediterranean fire hazards in this southern part of the Mediterranean Basin [5,13].This would allow for accurate in-depth analyses of the fire regime over lengthy periods, and for us to learn the factors that underlie fire occurrence and propagation in this region with a Mediterranean climate, but with substantial socio-economic and political differences compared to the more-studied Euro-Mediterranean side.
In recent years, several BA-mapping approaches have been developed for different study regions using medium-resolution data [41,42], including the BA Mapping Tools (BAMTs) [43].By leveraging the powerful capabilities of the GEE's cloud computing platform [34], the BAMTs constitute not only a significant stride as innovative, time-efficient, and resource-conserving tools for accurate multi-year BA mapping [10,39] but also the creation of independent reference data for validation exercises [44][45][46].
Based on these premises, we aim to exploit these efficient tools for systematically reconstructing the fire history in NE Algeria.Specifically, we (1) generate a BA product from the Landsat Collection-2 Surface Reflectance (LC2SR) product covering the 1984-2023 period; (2) assess its spatio-temporal consistency; and (3) provide pieces of evidence for a significantly revised BA estimate compared to existing BA products (GABAM, FireCCI51, C3SBA11, MCD64A1, and EFFIS) and a ground-based fire dataset.This work constitutes the mandatory initial step for creating a spatially explicit BA database following international standards for the whole of Algeria to further open a major research field for fire hazard, impacts, and vulnerability assessments that lead to firefighting and fire management policies [47].

Study Area
The study area is in NE Algeria, on the southern rim of the Mediterranean Basin, and spans longitudes from 3.71 • E to 8.68 • E and latitudes from 36.21 • N to 37.08 • N. The total study area covers 17,036 km 2 .It includes the coastal "Wilayas" (administrative provinces) of Annaba, Béjaïa, El-Tarf, Jijel, Skikda, and Tizi Ouzou in a narrow strip that stretches 440 km in length and up to 96 km in depth, with altitudes ranging from m.s.l. to 2291 m (Figure 1a,b).Landscapes are typically Mediterranean, dominated by forests, shrubs, and grasslands, where the mean annual rainfall varies from 725 to 924 mm (Figure 1c and Table 1).The main forest types are composed of oaks, such as cork-oak (Quercus suber L.), Algerian oak (Quercus canariensis Willd.),African oak (Quercus afares Pomel.), and holm oak (Quercus ilex L.), and of conifers, such as maritime pine (Pinus pinaster Aiton.),Aleppo pine (Pinus halepensis Miller.), and Atlas cedar (Cedrus atlantica (Endl.)Carrière.).Our study area is acknowledged as being one of the regional plant biodiversity hotspots in the Mediterranean Basin [48].It harbors the densest, richest, and most pristine forest ecosystems, which are reflected in four National Parks (Djurdjura, Gouraya, Taza, and El-Kala) that form part of the UNESCO World Network of Biosphere Reserves (WNBR) [49].However, this region experiences the highest incidence of fires in the country [50,51] owing to high fuel availability, extreme climate conditions in summer (prolonged drought and Saharan heat waves) and increasing anthropogenic activities [5].
this region experiences the highest incidence of fires in the country [50,51] owing to high fuel availability, extreme climate conditions in summer (prolonged drought and Saharan heat waves) and increasing anthropogenic activities [5].

Landsat and Sentinel-2 Imagery
The NASA/USGS Landsat program has launched eight Earth Observation (EO) satellites to date.These satellites have been continuously collecting data on the Earth's land  The NASA/USGS Landsat program has launched eight Earth Observation (EO) satellites to date.These satellites have been continuously collecting data on the Earth's land surface for over 50 years, the world's longest EO data series (1972-present day) [54].Different multispectral and thermal sensors have been developed, including Landsat 1-5 MSS , Landsat-4 TM (1982-1993), Landsat-5 TM (1984-2012), Landsat-7 ETM + (1999-2021), and the currently operational Landsat-8 OLI/TIRS (2013-present day) and Landsat-9 OLI-2/TIRS-2 (2021-present day) [55].Landsat instruments from TM to OLI-2/TIRS-2 eras have been generating scenes over 185 km swath width using the Worldwide Reference System-2 (WRS-2) with seven, eight and 11 spectral bands at 15 m (panchromatic), 30 m (visible, NIR, SWIR, coastal/aerosol, and cirrus), and 60-120 m (thermal bands) spatial resolutions every 16 days, reduced to 8 days when two satellites are combined (currently, Landsat 8 and 9).The ESA/Copernicus Sentinel-2 Mission comprises a constellation of two identical EO satellites, Sentinel-2A and 2B, which have been simultaneously orbiting the Earth since July 2015 and July 2017, respectively.The Sentinel-2 MSI sensor captures scenes of 110 × 110 km 2 using the Military Grid Reference System (MGRS) with 13 spectral bands at 10 m (visible and NIR), 20 m (red edge and SWIR), and 60 m (atmospheric bands) every 10 days, or every 5 days with both satellites combined [56].In this work, the LC2SR images from the last five Landsat satellites  and Sentinel-2 Level-1C Top of Atmosphere (TOA) reflectance images (2017 and 2021) were used, both available as image collections in the Earth Engine Data Catalogue (EEDC) "https://developers.google.com/earth-engine/datasets(accessed on 20 October 2022)".

Global Burned Area Products
To evaluate our generated BA product, a comparison to different existing BA products was conducted.For this purpose, five readily available BA products with different spatial resolutions were used: 1.
The GABAM is the first and only medium-resolution BA product at 30 m resolution with a global scale and long-term BA data to date.A novel automatic pipeline [33] was applied to generate the annual global BA maps from the Landsat time-series on the GEE platform [34] from 1985 to 2019.Nevertheless, 1986Nevertheless, , 1988Nevertheless, , 1990Nevertheless, , 1991Nevertheless, , 1993Nevertheless, , 1994Nevertheless, , 1997 and 1999 are unavailable.Yearly GABAM composites were downloaded as 10 The FireCCI51 provides monthly global BA maps at 250 m spatial resolution based on a hybrid algorithm coupling daily surface reflectance imagery and the active fire data from MODIS for 2001-2020 [24].The FireCCI51 pixel product in GeoTIFF format was downloaded at https://doi.org/10.5285/58f00d8814064b79a0c49662ad3af537,accessed on 18 April 2023.

3.
The C3SBA11 delivers monthly global BA maps at 300 m spatial resolution from 2017 onwards, using an adaptation of the FireCCI51 BA algorithm to Sentinel-3 OLCI data [26].

Active Fire Products
Active fire observations from the MODIS and VIIRS sensors were leveraged to examine the temporal fire reporting accuracy of our BA product by comparing detected burn dates to active fire dates.In line with this, the standard MCD14ML collection 6.1 product from MODIS onboard Terra and Aqua provides the systematic detection of active fires worldwide at a 1 km resolution from November 2000 to the present day [57].In contrast, the 375 m standard VIIRS VNP14IMGML collection 1 product offers better accuracy in detecting active fires from January 2012 onwards [58].These datasets were provided as feature point data in the ESRI Shapefile format by the Fire Information for Resource Management System (FIRMS) "https://firms.modaps.eosdis.nasa.gov/download/(accessed on 4 April 2023)".

Ground-Based Fire Dataset
Ground-based annual BA estimates over the 1985-2023 period for the six wilayas under study were provided by the Directorate General of Forests (DGF), which operates under the authority of the Algerian Ministry of Agriculture and Rural Development (MADR), Algeria.The DGF oversees fire prevention and control and develops national strategies and action plans for fire management across the country.Fire statistics are collected from 40 local forest services of the different wilayas that are regularly affected by fires [50].However, detailed single-fire datasets have been recently compiled only for the 2012-2023 period, where detailed information is provided: smallest administrative division, forest name or locality, latitudes/longitudes, date and time of ignition/intervention/extinction, burned vegetation type, and more.Here, we note that the minimum recorded fire size is not standardized across years and wilayas-a common feature that has been previously reported in Spain [59], France [60] or Tunisia [13], for which a standard threshold of 0.1 ha is usually assumed.

Land Cover Map
The ESA WorldCover project produced 10 m resolution global land cover maps for 2020 and 2021 from Sentinel-1 and Sentinel-2 data with 11 land cover classes conforming to the UN-FAO's Land Cover Classification System (LCCS) [52].The land cover map of 2021 for our study area was obtained in GeoTIFF format from ESA at https://worldcover2021. esa.int/downloader, accessed on 6 May 2023.

Burned Area Generation
In 2014, Bastarrika et al. [61] developed a Landsat interactive BA mapping procedure, named the BA Mapping Software (BAMS).This software allowed us to supervise the BA mapping between single pre-and post-fire images by adopting a two-phase supervised strategy [62].This approach has been widely employed to create fire reference perimeters to validate global and regional BA products [63][64][65], and to operationally reconstruct historical fires in Greece (1984-1991 and 1999-2011) [66], Argentina (1999-2011) [67] and Tunisia (1984-2010) [13].However, it heavily depended on the system's computational capacity and had several limitations, including the inability to work with higher-resolution Sentinel-2 data and maintenance difficulties due to changes in the Landsat data format and metadata, which led to frequent malfunctioning and uncontrolled commissions over agricultural and cloud-shadowed areas [43].To address these limitations, Roteta et al. [43] proposed the BA Mapping Tools (BAMTs).The current version (BAMTs v1.7) presents a set of six tools implemented in the GEE environment [34], which include BA mapping (BA Cartography tool), the selection of validation areas (VA tool) and image dates (VA Dates tool), reference data creation (RP tool), validation (Assessment tool), and image data visualization (Image Viewer).The JavaScript API scripts for these tools with a user guide are publicly available at https://github.com/ekhiroteta/BAMT,accessed 20 October 2022).
To generate our BA product, yearly pre-and post-fire consecutive multitemporal composites covering the study area were derived using 10,406 surface reflectance images from the LC2SR product for the 1984-2023 period using the BA Cartography tool.Clouds and cloud shadows were masked using Landsat pixel_qa band.For each year, the pre-fire period was defined from 1 January to 30 April, and the post-fire period from 1 May to 31 December , following Roteta et al. [43].This was done to employ the maximum number of cloud-free scenes and thus reduce BA omissions and capture all the fires with detection delays caused by the Landsat revisit time (8-16 days).The BA Cartography tool is a semiautomatic BA algorithm that performs a supervised classification of Landsat or Sentinel-2 data based on a Random Forest (RF) machine-learning classifier [68].The RF model is widely applied for BA extractions from medium-resolution data [33,43,69].
In addition to six reflectance bands over the visible near-infrared (VNIR) and shortwavelength infrared (SWIR), three spectral indices are computed for each Landsat image during the pre-and post-fire periods: the normalized difference vegetation index (NDVI) [70]; the normalized burn ratio (NBR) [71]; and the normalized burn ratio 2 (NBR2) [72].The date with the lowest NBR value (the date of the strongest burn signal) was used to generate image mosaics, meaning that each pixel has its specific mosaicking date.Then, a difference image was generated by subtracting the post-and pre-fire image composites.Subsequently, the delineation of burned patches was achieved through an iterative training process of the RF classifier, starting with an initial set of burned and unburned training sample polygons, defined based on both temporal composites and the difference image (post-fire-pre-fire).The classification output was visually inspected with each iteration, and new polygons were incrementally added to the training dataset of the model to reduce commission and omission errors.This strategy continued until a satisfactory delineation accuracy of the burned patches was obtained.A shapefile vector layer was used as a mask to limit the processing extent to our study area where an average of 38,894 burned and 45,723 unburned training samples (pixels) were manually and visually delineated every year [33,69].Next, a two-phase supervised strategy to balance errors of commission and omission [61,62] was employed on the RF probability images.In the first phase, seed pixels with strong burned signals were extracted using a threshold based on the average probabilities of the burned training samples.In the second phase, a region-growing process was applied up to a 50% threshold outwardly from the seed pixels to capture the entire burned patch using Rook's case contiguity, where only the burned patches containing at least one seed pixel were labeled as burned.
Finally, each annual BA map was exported as four 2 • × 2 • tiles at the 30 m resolution in the ESRI Shapefile format, and a thorough fire-by-fire interactive visual quality check was performed in a GIS environment to ensure the highest quality of the final burned patches.This allowed a few anomalies to be eliminated, such as very few continuous strip lines (one to two pixels wide), often starting from the edge of detected burned patches.These errors are caused mainly by the Landsat-7 ETM + SLC failure.Each burned polygon has a single attribute: indicating either the approximate burn detection date or null values for unobserved areas (i.e., clouds and cloud shadows).
At this step, it is important to note that the supervised BA contouring process under the BA Cartography tool may not be accomplished because of GEE's scaling errors [43].This issue was specifically encountered with years holding large amounts of BA.To avoid this shortcoming, we opted for a balanced and conservative approach when training the burned/unburned samples by avoiding very large polygons and complex geometries.The procedural steps for the generation and validation of the BA product (Phase I and II) and the intercomparison analysis (Phase III) are given in Figure 2.

Spatial Validation
Validation plays a pivotal role in ensuring the integrity and reliability of satellitederived BA products.The Committee on Earth Observation Satellites (CEOS) Working Group on Calibration and Validation (WGCV) defines validation as "the process of assessing, by independent means, the quality of the data products derived from those system outputs" [73].The spatial validation of our BA product was performed using two validation years: 2017 and 2021.These years were selected because they hold the largest annual BA during the Sentinel-2 era (2015-present).Sentinel-2 data were used to create the independent reference data.To assess the temporal fire reporting accuracy of the BA product, the burn detection dates from the BA Cartography tool were compared to the active fire dates from MODIS and VIIRS.
Five validation areas were sampled for 2017 and 2021 by the VA tool following a stratified random sampling methodology [74,75].At each sampling unit, this tool incorporates two key stratification criteria: the predominant Olson biomes [76] and fire activity according to FireCCI51 or MCD64A1.To ensure comprehensive and robust validation, long sampling units (LUs) lasting 100 days, and consisting of 110 × 110 km 2 tiles with a minimum temporal frequency of 30 days between consecutive Sentinel-2 images and cloud cover below 10%, were employed.This approach ensures a long data series with frequently available cloud-free images to create reliable reference data from continuous and frequent ground observations throughout the validation years.As the generation of high-quality reference perimeters in validation areas of 110 × 110 km 2 (~ 12,000 km 2 ) is time-intensive, the VA tool was designed to create smaller sampling units at the center of each Sentinel-2 tile [24,64,77].However, this subsampling resulted in sampling units not being representative across our study area.To solve this problem, we manually generated a tessellated grid with 168 cells of 20 × 20 km 2 (400 km 2 ) over the originally sampled validation areas [45].Thereafter, we selected the cells with the greatest fire activity based on

Spatial Validation
Validation plays a pivotal role in ensuring the integrity and reliability of satellitederived BA products.The Committee on Earth Observation Satellites (CEOS) Working Group on Calibration and Validation (WGCV) defines validation as "the process of assessing, by independent means, the quality of the data products derived from those system outputs" [73].The spatial validation of our BA product was performed using two validation years: 2017 and 2021.These years were selected because they hold the largest annual BA during the Sentinel-2 era (2015-present).Sentinel-2 data were used to create the independent reference data.To assess the temporal fire reporting accuracy of the BA product, the burn detection dates from the BA Cartography tool were compared to the active fire dates from MODIS and VIIRS.
Five validation areas were sampled for 2017 and 2021 by the VA tool following a stratified random sampling methodology [74,75].At each sampling unit, this tool incorporates two key stratification criteria: the predominant Olson biomes [76] and fire activity according to FireCCI51 or MCD64A1.To ensure comprehensive and robust validation, long sampling units (LUs) lasting 100 days, and consisting of 110 × 110 km 2 tiles with a minimum temporal frequency of 30 days between consecutive Sentinel-2 images and cloud cover below 10%, were employed.This approach ensures a long data series with frequently available cloud-free images to create reliable reference data from continuous and frequent ground observations throughout the validation years.As the generation of high-quality reference perimeters in validation areas of 110 × 110 km 2 (~12,000 km 2 ) is time-intensive, the VA tool was designed to create smaller sampling units at the center of each Sentinel-2 tile [24,64,77].However, this subsampling resulted in sampling units not being representative across our study area.To solve this problem, we manually generated a tessellated grid with 168 cells of 20 × 20 km 2 (400 km 2 ) over the originally sampled validation areas [45].Thereafter, we selected the cells with the greatest fire activity based on the BA estimates from FireCCI51 in the corresponding validation years.In all cases, ten cells were selected as new validation sites, with two cells assigned to each sampled Sentinel-2 tile, all of them lying in the Mediterranean forest biome (Figure 3).
the BA estimates from FireCCI51 in the corresponding validation years.In all cases, ten cells were selected as new validation sites, with two cells assigned to each sampled Sentinel-2 tile, all of them lying in the Mediterranean forest biome (Figure 3).To properly date the validation period, the CEOS Land Product Validation Subgroup (CEOS-LPVS) recommends generating reference data from multitemporal pairs of images rather than single pairs of pre-and post-fire images [77].Therefore, 227 Sentinel-2 consecutive images with cloud cover below 10% were identified for each validation site using the VA Dates tool.Validation periods were about 5-7 months for 2017 and 4-8 months for 2021 (Table A1).
Satellite-derived BA product validation requires independent reference data at a higher spatial resolution than the product to be validated [78].Accordingly, the RP tool was used to extract the 10 m resolution reference perimeters (RPs) from the selected pairs of consecutive Sentinel-2 images at each validation site.The RP tool performs the same BA mapping exercise as the BA Cartography tool, previously used in Section 3.1, albeit with some modifications.Instead of temporal composites, single-image pairs from the pre-and post-fire periods were employed.Additionally, the SLC image or QA60 and B1 bands were used to mask clouds and cloud shadows.Likewise, the RF classifier was iteratively trained by defining burned and unburned training pixels until the optimal delineation of burned patches was achieved.
It is important to note that the "cloud cover" information used by the VA Date tool to select the cloud-free Sentinel-2 images at each validation site is that of the entire Sentinel-2 scene (110 × 110 km 2 ) stored in the image metadata rather than being specific to the manually sampled validation sites (20 × 20 km 2 ).Consequently, very few consecutive dates were excluded in the reference data creation step owing to high cloud cover that obstructed ground observations.The resulting Sentinel-2 reference data (hereinafter referred to as "S2RD") were exported in the ESRI Shapefile format.A single attribute is given to the S2RD vector layers to indicate burned and unobserved polygons.Finally, a visual quality control of the short units of the RPs was performed before merging them into long temporal units at each validation site [69,77,79] (Figures A1 and A2).
At each validation site, the 2017 and 2021 BA maps were compared to their counterparts from S2RD.The commission error (CE), omission error (OE), overall accuracy (OA), To properly date the validation period, the CEOS Land Product Validation Subgroup (CEOS-LPVS) recommends generating reference data from multitemporal pairs of images rather than single pairs of pre-and post-fire images [77].Therefore, 227 Sentinel-2 consecutive images with cloud cover below 10% were identified for each validation site using the VA Dates tool.Validation periods were about 5-7 months for 2017 and 4-8 months for 2021 (Table A1).
Satellite-derived BA product validation requires independent reference data at a higher spatial resolution than the product to be validated [78].Accordingly, the RP tool was used to extract the 10 m resolution reference perimeters (RPs) from the selected pairs of consecutive Sentinel-2 images at each validation site.The RP tool performs the same BA mapping exercise as the BA Cartography tool, previously used in Section 3.1, albeit with some modifications.Instead of temporal composites, single-image pairs from the pre-and post-fire periods were employed.Additionally, the SLC image or QA60 and B1 bands were used to mask clouds and cloud shadows.Likewise, the RF classifier was iteratively trained by defining burned and unburned training pixels until the optimal delineation of burned patches was achieved.
It is important to note that the "cloud cover" information used by the VA Date tool to select the cloud-free Sentinel-2 images at each validation site is that of the entire Sentinel-2 scene (110 × 110 km 2 ) stored in the image metadata rather than being specific to the manually sampled validation sites (20 × 20 km 2 ).Consequently, very few consecutive dates were excluded in the reference data creation step owing to high cloud cover that obstructed ground observations.The resulting Sentinel-2 reference data (hereinafter referred to as "S2RD") were exported in the ESRI Shapefile format.A single attribute is given to the S2RD vector layers to indicate burned and unobserved polygons.Finally, a visual quality control of the short units of the RPs was performed before merging them into long temporal units at each validation site [69,77,79] (Figures A1 and A2).
At each validation site, the 2017 and 2021 BA maps were compared to their counterparts from S2RD.The commission error (CE), omission error (OE), overall accuracy (OA), Dice coefficient (DC) [80], and relative bias (RelB) were computed employing the assessment tool based on the confusion matrix [81].These accuracy metrics are widely used to assess the accuracy of satellite-derived BA products [44,79,82].Additional metrics were also considered to measure the total area correctly detected as burned (SurfBA) or unburned (SurfUB).Similarly, SurfCE and SurfOE represent the total committed and omitted burned surfaces, respectively.

Temporal Validation
The burn detection dates in the BAMTs-derived product were determined based on the most frequently occurring date (the mode) for all the pixels in each detected burned patch in the Landsat post-fire composite [43].Consequently, these dates do not imperatively correspond to the effective burn dates.To assess our BA product's temporal reporting accuracy, the temporal delay in days was computed between the burn detection dates and the active fire dates from the MODIS and VIIRS hotspots [27,28,69].Nevertheless, this analysis was performed by bearing in mind all the burned patches in the entire study area for a temporal validation period spanning from 2001 to 2023 rather than at the validation sites, because taking specific years as the number of burned patches to be compared was insufficient.To accomplish this, 1 × 1 km 2 windows were created around the yearly active fire hotspots from MCD14ML from 2001 to 2011, which represents the approximate MODIS pixels.Similarly, windows of 375 × 375 m 2 were defined around the VNP14IMGML hotspots from 2012 to 2023.Next, the most frequent dates from the active fire windows that intersected the individual burned patches in the corresponding years were extracted and then compared to the burn detection dates (Figure A3).This analysis was performed on 33.68% of the burned patches (N = 22,236), which represented the intersection rate with the active fire windows from MODIS and VIIRS.Burned patches were found in 49.14% of the active fire windows (N = 9030 from MCD14ML and N = 34,745 from VNP14IMGML) between 2001 and 2023.The remaining burned patches (66.32%) were not found in the active fire windows, which can be explained by active fire sensors' native coarse spatial resolution or satellite overpass.In addition, environmental conditions like thick smoke plumes and dense vegetation cover significantly reduce fire detection.Ideally, the MODIS sensor routinely detects active fires of ~0.1 ha and, under ideal observational conditions, fires up to ~50 m 2 in size can be captured, but those of lesser extents or with low emitted fire radiative power (FRP) may not be detected [83].Conversely, the VIIRS sensor has a minimum theoretical detectable fire size of ~20 m 2 at 1200 K fires burning in the daytime and one half at 900 K instantaneous fires at nighttime [58].

Intercomparison Analysis
In this section, an intercomparison was made between our validated BA product (hereinafter referred to as NEALGEBA; NE Algeria BA) and five existing BA products with different spatial resolutions and input data in terms of omissions and commissions.Furthermore, cross-correlation and temporal trend analyses were performed using the ground-based BA estimates from the DGF.

Spatial Accuracy
To highlight the differences in terms of burned patch delineation between NEAL-GEBA and existing BA products, the accuracy metrics (CE, OE, DC, RelB) obtained for all the BA products were compared.To this end, the same Sentinel-2 reference data at the same previously selected validation sites of 2017 were used to assess GABAM, C3SBA11, FireCCI51, MCD64A1, and EFFIS since this is the common year between all these BA products [44,69].The 2017 BA maps from FireCCI51 and MCD64A1 were already available as image collections in the GEE's data catalog and were directly evaluated for spatial accuracy against the S2RD in the assessment tool, whereas the GABAM, C3SBA11, and EFFIS maps were first adjusted to its required format.It is important to note that GABAM does not provide the burn date, as required by the assessment tool, but indicates only the year of burn.To overcome this situation, we assigned approximate burn detection dates to the burned patches on the 2017 GABAM map by overlaying the BA map from NEALGEBA.The burned patches with the largest overlap were assigned the same burn detection date, while those with no overlap were not dated.

Cross-Correlation with Ground-Based Fire Dataset
To assess NEALGEBA product's performance in BA estimation terms, a cross-correlation analysis was performed among the total annual BA estimates from our product, existing BA products, and the ground-based fire dataset from the DGF over the 13 available years during the 2004-2019 period using the "PerformanceAnalytics" R-package [84].However, the C3SBA11 product was not considered since only three years were available (2017, 2018, and 2019) among the overlapping years between the BA products.Here and in the subsequent step (Section 3.3.3),exclusively, land cover class 10 (tree cover), class 20 (shrubland), and class 30 (grassland) from the ESA WorldCover map 2021 [52] were assumed to represent the spatial extent of natural vegetation landscapes, which was used to mask out detected burned patches outside these classes from NEALGEBA and existing BA products.This was carried out following the definition of wildfires in Algeria [50], based on which the fraction of cropland fires is not taken into consideration when reporting BA estimates by the local forest services of the DGF.

Temporal Burned Area Trends
Here, the temporal trends of the annual BA estimates from the filtered NEALGEBA and DGF datasets for the 1985-2023 period (the DGF dataset lacks data for 1984) were analyzed using the Mann-Kendall trend test [85,86] and Sen's slope estimator [87].These two non-parametric methods are widely used to test the significance of trends in fire regimes [88].

Spatio-Temporal Patterns
The BA mapping phase resulted in the generation of spatially explicit annual BA maps covering all fires occurring in the natural landscapes and rural areas of the six wilayas on the NE coast of Algeria at the 30 m resolution for the 1984-2023 period (Figure 4a).In the last four decades, the total mapped BA was 1.16 M ha, with an average of 29,130.5 ha year −1 affecting forests (65.2%), shrublands (9.30%), grassland (22.2%), croplands (2.79%) and other types (<0.6%).This BA is double the official BA estimate from the DGF (0.65 M ha).The total annual BA showed wide temporal variability in the last 40 years, with 68.8% of the total BA occurring in 13 peak years above the mean.The most extreme years were 1994, 2000, 2017, and 2021, when more than 15.0%, 6.6%, 5.3%, and 8.4% of the total BA were observed, respectively.These years accounted for up to 35.3% of the total BA in the entire time series.Conversely, the lowest annual BAs were found in 1984, 1996, and 2018, which represent less than 1% of the total BA over the entire study period (Figure 4b).The study area showed a high fire recurrence rate with more than 52.7% of the total BA affected by two fire events or more in the last four decades (Figure 4c).Indeed, the total mapped BA (1.16 M ha) was twice the extent of the land affected by at least one fire (0.58 M ha).Only 47.3% of this land was burned once, while 25.9% was burned twice and 25.8% was affected by three fires or more (Figure 4c).

Fire-Size Distribution
The minimum detectable fire size equaled ~ 0.1 ha, which corresponds to the areal unit covered by a Landsat pixel on the ground.The maximum fire size was 26,300 ha, which occurred in the Tizi Ouzou wilaya in 2021.Fire size distribution using all the fire patches (N = 98,259) followed the self-organized criticality hypothesis [89] up to a fire size threshold of 1.5 ha, below which some small fires could be missing (Figure 6).

Fire-Size Distribution
The minimum detectable fire size equaled ~0.1 ha, which corresponds to the areal unit covered by a Landsat pixel on the ground.The maximum fire size was 26,300 ha, which occurred in the Tizi Ouzou wilaya in 2021.Fire size distribution using all the fire patches (N = 98,259) followed the self-organized criticality hypothesis [89] up to a fire size threshold of 1.5 ha, below which some small fires could be missing (Figure 6).

Fire-Size Distribution
The minimum detectable fire size equaled ~ 0.1 ha, which corresponds to the areal unit covered by a Landsat pixel on the ground.The maximum fire size was 26,300 ha, which occurred in the Tizi Ouzou wilaya in 2021.Fire size distribution using all the fire patches (N = 98,259) followed the self-organized criticality hypothesis [89] up to a fire size threshold of 1.5 ha, below which some small fires could be missing (Figure 6).

Fire Seasonality
Fire seasonality, based on fire number and BA, could be generated on 22,115 fire patches, which could be dated out of the 66,018 fire patches identified over the 2001-2023 period when the VIIRS/MODIS hotspots were available.NE Algeria displays typical Mediterranean fire activity, and July and August are the most burnable months, which prolong to September and October (Figure 7).
Remote Sens. 2024, 16, 2500 14 of 34 4.1.3.Fire Seasonality Fire seasonality, based on fire number and BA, could be generated on 22,115 fire patches, which could be dated out of the 66,018 fire patches identified over the 2001-2023 period when the VIIRS/MODIS hotspots were available.NE Algeria displays typical Mediterranean fire activity, and July and August are the most burnable months, which prolong to September and October (Figure 7).

Spatial Validation
The results of the accuracy assessment of the NEALGEBA maps of 2017 and 2021 at the 10 validation sites using S2RD are presented in Table 2.The corresponding validation images are found in Figures A4 and A5.For the validation year 2017, the OA was 98.22%, with a CE of 7.96%, an OE of 8.19%, and a DC of 91.92%.Notably, validation site 31SEA-E3 had the highest OE with 30.16%, which resulted in an underestimation of the total BA of 23.19% and the lowest DC (79%) of all the validation sites.The highest commissions (CE = 9.65%) occurred at 32SMF-Y3.The overall RelB indicated a slight underestimation of the total BA of 0.24%.The total area correctly detected as burned was 19,551 ha, while 170,000 ha of unburned surface was correctly detected at all the validation sites (covering 200,000 ha).Additionally, the product misclassified 1692 ha of the unburned surface as burned and omitted 1743 ha of the burned surfaces.For 2021, the product achieved comparable results to those for 2017 in terms of OA (98.15%) and commissions (7.92%), albeit with a twofold decrease in omissions (4.76%) that resulted in a higher DC (93.63%).Of all the validation sites, the highest commissions and omissions were for 32SKF-P3 with values of 25.38% and 14.88%, respectively, and a RelB value of 14.08%.The DC was the lowest at all validation sites with 79.52%.The total area correctly detected as burned was 27,070 ha, and 168,984 ha was correctly detected as unburned.There were instances where the product incorrectly classified 2553 ha of unburned surface as burned and omitted 1352 ha of the burned surface.Overall, cloud coverage was the primary source of BA omissions, while most commissions were caused by the BA algorithm over land covers spectrally

Spatio-Temporal Validation 4.2.1. Spatial Validation
The results of the accuracy assessment of the NEALGEBA maps of 2017 and 2021 at the 10 validation sites using S2RD are presented in Table 2.The corresponding validation images are found in Figures A4 and A5.For the validation year 2017, the OA was 98.22%, with a CE of 7.96%, an OE of 8.19%, and a DC of 91.92%.Notably, validation site 31SEA-E3 had the highest OE with 30.16%, which resulted in an underestimation of the total BA of 23.19% and the lowest DC (79%) of all the validation sites.The highest commissions (CE = 9.65%) occurred at 32SMF-Y3.The overall RelB indicated a slight underestimation of the total BA of 0.24%.The total area correctly detected as burned was 19,551 ha, while 170,000 ha of unburned surface was correctly detected at all the validation sites (covering 200,000 ha).Additionally, the product misclassified 1692 ha of the unburned surface as burned and omitted 1743 ha of the burned surfaces.For 2021, the product achieved comparable results to those for 2017 in terms of OA (98.15%) and commissions (7.92%), albeit with a twofold decrease in omissions (4.76%) that resulted in a higher DC (93.63%).Of all the validation sites, the highest commissions and omissions were for 32SKF-P3 with values of 25.38% and 14.88%, respectively, and a RelB value of 14.08%.The DC was the lowest at all validation sites with 79.52%.The total area correctly detected as burned was 27,070 ha, and 168,984 ha was correctly detected as unburned.There were instances where the product incorrectly classified 2553 ha of unburned surface as burned and omitted 1352 ha of the burned surface.Overall, cloud coverage was the primary source of BA omissions, while most commissions were caused by the BA algorithm over land covers spectrally akin to fire affected surfaces, such as harvested croplands, labored lands, water bodies, and shadowed areas.

Temporal Validation
The NEALGEBA product showed a significant fire detection delay compared to the active fire hotspots from MODIS and VIIRS.On average, burned pixels were detected 42.55 days after the hotspot detection: 21.07% on the first 15 days, 46.47% in the first month, and 75.67% in a 2-month window (Figure 8).This disparity is attributed to Landsat's low temporal revisit frequency (8-16 days), which could be worsened by the unavailability of cloud-free images for BA extraction purposes.
akin to fire affected surfaces, such as harvested croplands, labored lands, water bodies, and shadowed areas.

Temporal Validation
The NEALGEBA product showed a significant fire detection delay compared to the active fire hotspots from MODIS and VIIRS.On average, burned pixels were detected 42.55 days after the hotspot detection: 21.07% on the first 15 days, 46.47% in the first month, and 75.67% in a 2-month window (Figure 8).This disparity is attributed to Landsat's low temporal revisit frequency (8-16 days), which could be worsened by the unavailability of cloud-free images for BA extraction purposes.

Spatial Accuracy
The evaluation of the NEALGEBA product versus existing BA products showed its overall outperformance (Figure 9).Indeed, NEALGEBA had the lowest CE (7.96%), the

Intercomparison Analysis 4.3.1. Spatial Accuracy
The evaluation of the NEALGEBA product versus existing BA products showed its overall outperformance (Figure 9).Indeed, NEALGEBA had the lowest CE (7.96%), the highest DC (91.92%), and the best total BA estimation compared to S2RD (RelB = −0.24%).Despite being generated using the same input data (30 m Landsat images) and having a lower OE (3.89%), GABAM obtained a higher CE (13.33%) than NEALGEBA (7.96%), which resulted in an overestimated BA (RelB = 10.89%) compared to our product.As expected, coarser resolution products based on MODIS and Sentinel-3 OLCI: FireCCI51, MCD64A1, EFFIS, and C3SBA11 produced the highest CEs and Oes and the lowest DCs, with a significant lack of estimation of the total BA compared to S2RD.Specifically, FireCCI51, C3SBA11, and MCD64A1 had CEs ≥ 37.43%, which resulted in a marked overestimation of the total BA that reached 58.39%, 43.10%, and 36.03%,respectively.However, these three products had lower omissions ≤ 14.89% compared to the EFFIS product (≥33.26%).The last one had the highest OE of all the assessed BA products, which led to the highest BA underestimation rate (RelB = −12.69%).Detailed results of the accuracy assessment of the GABAM, FireCCI51, C3SBA11, MCD64A1, and EFFIS BA maps of 2017, along with their corresponding validation images, appear in Tables A2-A6, and in Figures A6-A10, respectively.
highest DC (91.92%), and the best total BA estimation compared to S2RD (RelB = −0.24%).Despite being generated using the same input data (30 m Landsat images) and having a lower OE (3.89%), GABAM obtained a higher CE (13.33%) than NEALGEBA (7.96%), which resulted in an overestimated BA (RelB = 10.89%) compared to our product.As expected, coarser resolution products based on MODIS and Sentinel-3 OLCI: FireCCI51, MCD64A1, EFFIS, and C3SBA11 produced the highest CEs and Oes and the lowest DCs, with a significant lack of estimation of the total BA compared to S2RD.Specifically, FireCCI51, C3SBA11, and MCD64A1 had CEs ≥ 37.43%, which resulted in a marked overestimation of the total BA that reached 58.39%, 43.10%, and 36.03%,respectively.However, these three products had lower omissions ≤ 14.89% compared to the EFFIS product (≥33.26%).The last one had the highest OE of all the assessed BA products, which led to the highest BA underestimation rate (RelB = −12.69%).Detailed results of the accuracy assessment of the GABAM, FireCCI51, C3SBA11, MCD64A1, and EFFIS BA maps of 2017, along with their corresponding validation images, appear in Tables A2-A6, and in Figures A6-A10, respectively.The visual analysis clearly shows that all the detected large, burned patches were consistent in all the BA products, although most relatively small fires were captured only by medium spatial resolution products (NEALGEBA and GABAM) (Figure 10a).GABAM demonstrated good spatial matching with S2RD but obtained a high CE (13.33%) at all the validation sites (Table A2).This resulted in relatively extended burned patch limits compared to NEALGEBA.In other words, GABAM tended to overestimate the BA and included other areas that may not have been burned.NEALGEBA showed higher OEs at all the validation sites, which resulted in slightly reduced burned patches compared to S2RD (Table 2).This indicates a tendency to exclude some burned patches, which leads to reduced BA compared to GABAM.These products slightly underestimate or overestimate burned areas, respectively, compared to the higher-resolution reference data from Sentinel-2.Furthermore, both omitted a few small, burned patches and misclassified a few unburned pixels as burned.It is important to consider the trade-off between the CEs and OEs when evaluating the performance of BA maps.The visual analysis clearly shows that all the detected large, burned patches were consistent in all the BA products, although most relatively small fires were captured only by medium spatial resolution products (NEALGEBA and GABAM) (Figure 10a).GABAM demonstrated good spatial matching with S2RD but obtained a high CE (13.33%) at all the validation sites (Table A2).This resulted in relatively extended burned patch limits compared to NEALGEBA.In other words, GABAM tended to overestimate the BA and included other areas that may not have been burned.NEALGEBA showed higher OEs at all the validation sites, which resulted in slightly reduced burned patches compared to S2RD (Table 2).This indicates a tendency to exclude some burned patches, which leads to reduced BA compared to GABAM.These products slightly underestimate or overestimate burned areas, respectively, compared to the higher-resolution reference data from Sentinel-2.Furthermore, both omitted a few small, burned patches and misclassified a few unburned pixels as burned.It is important to consider the trade-off between the CEs and OEs when evaluating the performance of BA maps.Conversely, the larger pixel size in FireCCI51 (250 m), C3SBA11 (300 m), and MCD64A1 (500 m) led to the size of burned patches being significantly overestimated.These products overlooked small fires and were unable to accurately delineate larger ones.They also displayed difficulty in discerning neighboring fire patches by perceiving them as single fire events, and they omitted the unburned islands in the main fire patches, which meant that the total BA was perpetually overestimated.EFFIS displayed a smoothing effect that oversimplified fire perimeters, with a loss of details inside burned patches (unburned islands), which led to significant BA commissions.Particularly in this product, no fires were captured at validation site 31SEA-E3; an important part of the BA extent was omitted at validation site 32SKF-V2; and relatively small and spatially fragmented fires, at a distance from large, burned patches, were completely omitted at the remaining validation sites.This resulted in significant omission errors.As Figure 10b depicts, the total detected amount of BA at the five validation sites in S2RD of 2017 was 21,295 ha.This was slightly reduced in NEALGEBA with 20,400 ha and slightly higher in GABAM with 23,140 Conversely, the larger pixel size in FireCCI51 (250 m), C3SBA11 (300 m), and MCD64A1 (500 m) led to the size of burned patches being significantly overestimated.These products overlooked small fires and were unable to accurately delineate larger ones.They also displayed difficulty in discerning neighboring fire patches by perceiving them as single fire events, and they omitted the unburned islands in the main fire patches, which meant that the total BA was perpetually overestimated.EFFIS displayed a smoothing effect that oversimplified fire perimeters, with a loss of details inside burned patches (unburned islands), which led to significant BA commissions.Particularly in this product, no fires were captured at validation site 31SEA-E3; an important part of the BA extent was omitted at validation site 32SKF-V2; and relatively small and spatially fragmented fires, at a distance from large, burned patches, were completely omitted at the remaining validation sites.This resulted in significant omission errors.As Figure 10b depicts, the total detected amount of BA at the five validation sites in S2RD of 2017 was 21,295 ha.This was slightly reduced in NEALGEBA with 20,400 ha and slightly higher in GABAM with 23,140 ha.Due to high commissions in FireCCI51, C3SBA11, and MCD64A1, the amount of BA was overestimated with 34,610 ha, 31,610 ha, and 29,690 ha, respectively, whereas EFFIS significantly detected less BA with only 18,880 ha due to high omission values.

Cross-Correlation with Ground-Based Fire Dataset
The evaluation of the different BA products compared to the DGF dataset was performed using the log10-transformed total annual BA estimates over the 13 overlapping years that comprised the 2004-2019 period between the different products.In general, the distinct BA products showed a strong correlation with the DGF data over the entire study area (Figure 11a).NEALGEBA, FireCCI51, and MCD64A1 obtained the highest correlations (r = 0.96), while GABAM yielded a slightly lower value (r = 0.94) even though it had the same spatial resolution as NEALGEBA (30 m).EFFIS had the lowest correlation with the DGF dataset.On the scale of the different wilayas, NEALGEBA demonstrated a better correspondence with the DGF estimates than existing BA products (Figure 11b).
Remote Sens. 2024, 16, 2500 18 of 34 ha.Due to high commissions in FireCCI51, C3SBA11, and MCD64A1, the amount of BA was overestimated with 34,610 ha, 31,610 ha, and 29,690 ha, respectively, whereas EFFIS significantly detected less BA with only 18,880 ha due to high omission values.

Cross-Correlation with Ground-Based Fire Dataset
The evaluation of the different BA products compared to the DGF dataset was performed using the log10-transformed total annual BA estimates over the 13 overlapping years that comprised the 2004-2019 period between the different products.In general, the distinct BA products showed a strong correlation with the DGF data over the entire study area (Figure 11a).NEALGEBA, FireCCI51, and MCD64A1 obtained the highest correlations (r = 0.96), while GABAM yielded a slightly lower value (r = 0.94) even though it had the same spatial resolution as NEALGEBA (30 m).EFFIS had the lowest correlation with the DGF dataset.On the scale of the different wilayas, NEALGEBA demonstrated a better correspondence with the DGF estimates than existing BA products (Figure 11b).

Temporal Trends of the Burned Area
The temporal evolution of total annual BA estimates from NEALGEBA and DGF displayed a consistent increasing trend over the 1985-2023 period (1984 is unavailable in the DGF dataset) yet differ in the annual rate, with NEALGEBA showing significant trends.

Temporal Trends of the Burned Area
The temporal evolution of total annual BA estimates from NEALGEBA and DGF displayed a consistent increasing trend over the 1985-2023 period (1984 is unavailable in the DGF dataset) yet differ in the annual rate, with NEALGEBA showing significant trends.For all the wilayas combined, NEALGEBA showed a significant increase rate of 532.4 ha year −1 at 10% (p-value = 0.053) in the total annual BA, while DGF showed a lower increase rate (197.1 ha year −1 ; Figure 12a).At the wilaya scale, the annual BA from the two products showed significant positive trends at 5% (p-value ≤ 0.05) for the two products in Béjaïa and was only significant for NEALGEBA in Jijel (p-value = 0.05) and El-Tarf at 10% (p-value = 0.056).In general, NEALGEBA showed higher increased BA rates (182.2,105.6, and 67.9 ha year −1 in Béjaïa, Jijel, Skikda, and El-Tarf, respectively; Figure 12b).For all the wilayas combined, NEALGEBA showed a significant increase rate of 532.4 ha year −1 at 10% (p-value = 0.053) in the total annual BA, while DGF showed a lower increase rate (197.1 ha year −1 ; Figure 12a).At the wilaya scale, the annual BA from the two products showed significant positive trends at 5% (p-value ≤ 0.05) for the two products in Béjaïa and was only significant for NEALGEBA in Jijel (p-value = 0.05) and El-Tarf at 10% (pvalue = 0.056).In general, NEALGEBA showed higher increased BA rates (182.2,105.6, and 67.9 ha year −1 in Béjaïa, Jijel, Skikda, and El-Tarf, respectively; Figure 12b).

Discussion
In this analysis, we reconstructed and validated 40 years (1984-2023) of historical fire events at fine spatial resolution in typical Mediterranean ecosystems of NE Algeria.The newly generated NEALGEBA product represents the first and most continuous time series of BA at fine spatial resolution in this part of the Mediterranean Basin, which faces a substantial fire occurrence threat.The BA product generation (phase I) proved the high

Discussion
In this analysis, we reconstructed and validated 40 years (1984-2023) of historical fire events at fine spatial resolution in typical Mediterranean ecosystems of NE Algeria.The newly generated NEALGEBA product represents the first and most continuous time series of BA at fine spatial resolution in this part of the Mediterranean Basin, which faces a substantial fire occurrence threat.The BA product generation (phase I) proved the high potential and reliability of the BA Cartography tool in generating spatially consistent annual BA maps based on a Random Forest supervised classification and a two-phased strategy [43,62].
Despite being labor-intensive and heavily relying on the visual interpretation of pre-and post-fire temporal composites, this semi-automatic procedure enabled high control over CEs and OEs and thus improved the BA product's quality.The analysts' expertise is more involved in selecting representative burned (seeds) and unburned samples with an iterative analysis of BA delineation [45]-a considerable asset that is not provided with fully automated methods [33].Additionally, the visual quality control and manual refinement of the generated fire perimeters allowed us to reduce potential anomalies such as those caused by the sensor.The RP, VA, VA dates, and assessment tools greatly facilitated the spatial validation exercise of satellite-derived BA products compared to previous studies [64,90], all in accordance with the BA assessment standardized protocol endorsed by the CEOS.These tools ensured the creation of high-quality reference data (RP tool) from consecutive 10 m cloud-free Sentinel-2 images (VA Date) located at the validation sites preselected by stratified random sampling (VA tool).The assessment tool allowed for a full-automatic comparison of the BA maps to the Sentinel-2 reference data and reported accuracy metrics (CE, OE, DC, and RelB) at each validation site.
The accuracy assessment of the 2017 and 2021 NEALGEBA maps showed remarkable results, with CEs of 7.96% and 7.92%, OEs of 8.19% and 4.76%, and DCs of 98.22% and 98.15%, respectively.These metrics are consistent with better performance than those obtained in the original case study in south-eastern Australia, in which a BA product for the 2019/2020 fire season was generated and validated using the same input data, with a CE of 11.80%, an OE of 8.90% and a DC of 89.60% [43].However, the larger pixel size (30 m) in NEALGEBA led to a subtle alteration in the extent of the burned patches, which meant that their boundaries slightly extended outwardly compared to the independent reference perimeters from the 10 m Sentinel-2 data.We also observed that almost all the spatially isolated small, burned patches were misclassified as burned, mainly in 2021.This was due to the algorithm's limitation of discriminating small, spectrally confusing surfaces with a similar spectral response to the burned surfaces.For the 2021 fire season, most of the BAs were in mountainous areas, which made it quite challenging to select representative and sufficient burned seed pixels to capture the entire burned patches and thus to reduce omission errors.We attempted to avoid burned pixels in shadowed areas to reduce commissions on the classification map.Moreover, the algorithm failed to ensure the continuity of some large, burned patches, and omissions occurred mostly on the edges of the main burned patches and on unburned islands, with very few small isolated burned patches that were completely omitted.Overall, the obtained accuracy metrics indicated the NEALGEBA product's spatial consistency.However, the temporal validation using the active fire hotspots from MODIS and VIIRS underlined its limitation for accurately reporting fire events over time.This is explained by the long revisit time of the Landsat satellites (8-16 days), atmospheric conditions (i.e., clouds), and temporal compositing, which could significantly delay BA detection.This is not uncommon in medium-resolution products [43,69] compared to MODIS-derived products that incorporate active fire information [24,29], and it underscores the need for further development to improve temporal uncertainty.Employing data from satellite sensors with a higher observation frequency, such as Sentinel-2 and the active fire hotspots from VIIRS, would reduce temporal reporting delays [27].
Coarse-resolution BA products were found to significantly overestimate the total BA on a finer spatial scale due to a larger pixel size (≥250 m), unlike the continental scale, on which the total BA was overly underestimated when compared to more accurate data from the Sentinel-2 MSI sensor [27,28,37,38].In addition, their limited temporal coverage (2001-present) prevents long-term fire studies compared to the BA products generated from the Landsat data archive dating back to 1984.The validation of the 2017 EFFIS BA map showed that the latter presented the highest omissions of all the assessed BA products, which resulted in a considerable underestimation of the total BA, plus a smoothing effect on fire perimeters that roughly delineated the burned patches.These inconsistencies have been previously reported in [41,91,92] and are attributed to the 250 m coarse-resolution input data from MODIS used to generate the EFFIS product.On the other hand, GABAM is, to date, the only available global high-resolution BA product to provide BA mapping at finer spatial resolution to reliably detect smaller burned patches [33].However, this product was generated in yearly composites by providing only the year of burn rather than the approximate burn detection date, which prevents a fire seasonal analysis.Moreover, while GABAM has the longest time span amongst global BA products , some years are still unavailable.In addition to the reported commissions over agricultural lands, significant systematic errors were observed when we examined the GABAM annual maps versus the corresponding Landsat post-fire image composites in a Long SWIR/NIR/Red color composition and NEALGEBA.The former represents BA commissions over water bodies, unburned forest areas, clouds, flares in oil/gas facilities, and Landsat strip errors (Figure A11a-c).Additionally, significant errors occurred in 2002 (Figure A11d) and were perhaps caused by the significant alteration of the primary functioning mode of Landsat-5 TM's scan mirror, known as the scan angle monitor (SAM), which led to internal synchronization problems.This failure caused diagonal patches of anomalous observations with remarkably high reflectance values in the long shortwave infrared (SWIR2) towards the scene footprint edges and led to false fire detections.The SAM system was then switched to an alternative one called the bumper mode [93], which overcame this anomaly.The semi-automatic BA extraction procedure, which uses the BA Cartography tool along with a thorough visual inspection of the mapped burned patches, allowed these anomalies to be mitigated and consistent results to be obtained on the NEALGEBA annual maps.Overall, these limitations highlight the challenges and complexities involved in using existing BA products to accurately characterize local and regional fire regimes, especially over lengthy periods.The newly generated BA product herein presented serves as a surrogate to existing BA products by offering improved spatio-temporal resolutions, allowing for a thorough assessment of fire impacts on forest ecosystems and, in turn, assisting in designing strategies and adapted action plans to mitigate their severity in NE Algeria.Additionally, BAMTs can be easily used by managers of forests and protected areas in Algeria to operationally extract burned perimeters and to integrate complementary field data, especially the day and time of ignition, which reduces temporal uncertainties.
Our first evaluation of the newly generated NEALGEBA could properly address the fire seasonal distribution that spans from July to October and peaks in August.This result is in accordance with the seasonal drought and fire weather index seasonality characterizing our region [5].A similar fire season length is reported in neighboring Tunisia [13] and in Portugal [94], and it is slightly longer than the usual July-September fire season reported in Morocco [95], Italy [96,97], Greece [98], Bulgaria [99], and the Iberian Peninsula [100].The length of the fire season in our region may be shaped by the early and late fire-prone weather conditions favored by climate change effects and socio-economic factors [5,9].Regarding affected vegetation, we obtained a fraction of burnable areas affected by fires reaching 2.95% year −1 , in fair agreement with Portugal (3.31%), making our BA estimates on the highest range of variability observed in the Mediterranean Basin.Only 0.19% was observed from Tunisia [13], Lebanon (0.58%), France (0.53%), Greece (0.57%), Spain (0.84%), and Italy (1.14%), as reviewed in [10].We also detected an increasing trend in BA (532.4 ha year −1 ) over the region.This trend is different from that in the northern part of the Mediterranean, where a general decreasing fire trend has been observed [101].More precisely, an abrupt decrease was observed in France in 1990 with increasing firefighting expenditures [60]; there was an increase in the 1980s and 1990s then a decrease after 2000 in Spain, owing to land abandonment [59,102]; and a decreasing trend was observed in Greece [14].In the Middle East and North Africa, no trends have been particularly observed since the 1980s [10,13], consistent with our study.However, the recent collapse of political regimes led to an abrupt increase in the BA in Tunisia [103], which we did not detect in Algeria, which was not affected by this political collapse.Throughout history, Algeria has been marked by significant political events that have led to heavy burning periods, such as the Algerian Civil War (the Black Decade) in the 1990s [6], which was reflected by the highest peak in BA observed in 1994.Algeria has remained quite stable since 2011, when the Arab Spring started in the southern Mediterranean Basin, and thus did not affect the BA during this recent period.Here we note the exceptional heat wave that hit the region in 2023 with record-breaking temperatures in April [104] and in July 2023, with some casualties during fire events, which were widely reported in the media but did not lead to the most extreme fire year compared to 1994 and 2021.Hence, NEALGEBA appears as a keystone database to provide accurate information on burned areas and temporal trends, and as a reference database to allow for the objective contextualization of fire years, to be further used in fire weather analyses, fire impacts assessments, fire model benchmarking or Euro-Mediterranean initiatives of fire-related issues.
In its current version, NEALGEBA covers all types of fires that have occurred across all landscapes in NE Algeria from 1984 through 2023.However, this product, as is the case with all satellite-derived BA products, exhibits specific limitations that necessitate reporting for future improvements.First, the burn dates indicated in our product do not match the effective fire dates, when the fires were actively burning.The BA Cartography tool computes the modal date from all pixels within each detected burned patch in the yearly Landsat post-fire composite.This results in significant fire detection delay, which impacts the product's temporal fire reporting accuracy, especially in large fires lasting several days.We highlight here that using data at a higher temporal resolution from Sentinel-2 MSI could significantly reduce this disparity [27,28,39].Second, the product's commission errors were primarily observed over agricultural lands (harvested or plowed croplands), which exhibit similar spectral features to burned areas, characterized by abrupt changes in the reflectance data, particularly in the near and shortwave infrared bands [28,33,36,61,62,105]. Here, we should also emphasize the high uncertainties associated with the detection of cropland fires [106].Third, very large fire events enduring several days may not be effectively captured as single burned patches due to the low revisit frequency of Landsat satellites and the availability of cloud-free images.Some large fire patches may not be spatially contiguous, owing to low burn signals over shadowed areas, sparse vegetation, or discontinuity in vegetation cover.Additionally, spotting fires can result in spatially isolated burned islands from the main fire patches.Fourth, Landsat sensor anomalies were one of the main challenges, especially the Landsat-7 ETM + SLC failure, which affected most of its time coverage.Post-processing procedures have been applied to mitigate commissions caused by these anomalies.However, some omissions or late fire detections (strips within an actual fire patch) should be acknowledged.Regarding future work, we envisage complementing the NEALGEBA product with detailed information contained in the ground-based fire inventories from local forest services of the DGF, especially for extreme fire events.For instance, burn detection dates can be corrected, thus reducing the product's temporal uncertainties.Possible attributes include forest name or locality, date and exact time of ignition/intervention/extinction, burned vegetation type and species, land ownership, cause of ignition, perpetrator of the fire, fire reporter, weather conditions, participating bodies in fire suppression, damage assessment, and investigation.Ongoing efforts involve expanding NEALGEBA to a country-level BA product with the continuous mapping of fire-affected areas for the upcoming years using higher-resolution imagery data from Sentinel-2 MSI.This aims to provide an accurate characterization of the spatiotemporal patterns of fires across a larger geographical scale.

Conclusions
This paper presents a new spatially explicit fire database (NEALGEBA) in the Mediterr anean-type ecosystems of NE Algeria, performed following a thorough quality checking of international standards.NEALGEBA provides a reconstruction of 40 years (1984-2023) of historical fire events from Landsat data to be used as a reference for the ongoing global processing of burned area datasets.We could demonstrate the higher performance of our approach to both the coarse resolution currently available and the global yet discontinuous Landsat-based GABAM, leading to a final dataset increasing the official annual BA provided by the national forest service DGF by 46.28%.
NEALGEBA also offers Algerian forest managers and policymakers valuable insights into the spatio-temporal distribution and magnitude of fire events at a fine resolution in the most fire-prone region of the country.Such databases represent the keystone prerequisite to further characterizing fire regime and identifying the key driving factors of fire occurrence.Furthermore, this dataset significantly contributes to further national and international fire hazard and impact assessments.It will act as a reference for contextualizing future fire extremes in Mediterranean scale synthesis, such as the 2023 exceptional heat wave, which we show not to have led to the most extreme fire year over the last 40 years.
Remote Sens. 2024, 16, 2500 8 of 34 procedural steps for the generation and validation of the BA product (Phase I and II) and the intercomparison analysis (Phase III) are given in Figure 2.

Figure 2 .
Figure 2. Flow diagram of the procedural steps of the methodology.(a) Generation of the BA product; (b) spatio-temporal validation of the generated BA product; (c) and intercomparison analysis.

Figure 2 .
Figure 2. Flow diagram of the procedural steps of the methodology.(a) Generation of the BA product; (b) spatio-temporal validation of the generated BA product; (c) and intercomparison analysis.

Figure 4 .
Figure 4. Spatio-temporal patterns of fires according to NEALGEBA in NE Algeria from 1984 to 2023.(a) Spatial extent of fires; (b) temporal distribution of the total annual BA; and (c) fire recurrence.

Figure 4 .
Figure 4. Spatio-temporal patterns of fires according to NEALGEBA in NE Algeria from 1984 to 2023.(a) Spatial extent of fires; (b) temporal distribution of the total annual BA; and (c) fire recurrence.

Figure 5 .
Figure 5. Kernel density distribution according to NEALGEBA in NE Algeria from 1984 to 2023.(a) Kernel density estimation (KDE); (b) BA-weighted kernel density estimation (WKDE); and (c) spatial distribution of fire size classes.

Figure 6 .
Figure 6.Fire size distribution: fire size (ha) and fire frequency on the log10 scale.

Figure 5 .
Figure 5. Kernel density distribution according to NEALGEBA in NE Algeria from 1984 to 2023.(a) Kernel density estimation (KDE); (b) BA-weighted kernel density estimation (WKDE); and (c) spatial distribution of fire size classes.

34 Figure 5 .
Figure 5. Kernel density distribution according to NEALGEBA in NE Algeria from 1984 to 2023.(a) Kernel density estimation (KDE); (b) BA-weighted kernel density estimation (WKDE); and (c) spatial distribution of fire size classes.

Figure 6 .
Figure 6.Fire size distribution: fire size (ha) and fire frequency on the log10 scale.Figure 6. Fire size distribution: fire size (ha) and fire frequency on the log10 scale.

Figure 6 .
Figure 6.Fire size distribution: fire size (ha) and fire frequency on the log10 scale.Figure 6. Fire size distribution: fire size (ha) and fire frequency on the log10 scale.

Figure 7 .
Figure 7. Fire seasonality based on monthly burned area and fire frequency during the 2001-2023 period.

Figure 7 .
Figure 7. Fire seasonality based on monthly burned area and fire frequency during the 2001-2023 period.

Figure 8 .
Figure 8. Distribution of the temporal difference in days between the burn detection dates from NEALGEBA and the active fire dates from MODIS and VIIRS across the study area from 2001 to 2023.M: mean and Mdn: median.

Figure 8 .
Figure 8. Distribution of the temporal difference in days between the burn detection dates from NEALGEBA and the active fire dates from MODIS and VIIRS across the study area from 2001 to 2023.M: mean and Mdn: median.

Figure 11 .
Figure 11.Pearson's correlation analysis of the total annual BA estimates from the BA products and the DGF dataset.(a) all wilayas combined and (b) for each wilaya.*** indicates significance at 0.1%.

Figure 11 .
Figure 11.Pearson's correlation analysis of the total annual BA estimates from the BA products and the DGF dataset.(a) all wilayas combined and (b) for each wilaya.*** indicates significance at 0.1%.

Figure 12 .
Figure 12.Temporal trends in the total annual BA estimates from NEALGEBA and DGF for the 1985-2023 period.(a) All wilayas combined and (b) for each wilaya.P: p-value, Tau: Kendall's Tau, Sen: Sen's slope.The red asterisk indicates statistically significant trends.

Figure 12 .
Figure 12.Temporal trends in the total annual BA estimates from NEALGEBA and DGF for the 1985-2023 period.(a) All wilayas combined and (b) for each wilaya.P: p-value, Tau: Kendall's Tau, Sen: Sen's slope.The red asterisk indicates statistically significant trends.

Figure A1 .
Figure A1.(a-g) Short units and (h) long units of the reference data at validation site 32SLF-V2 for the 2017 validation year.

Figure A2 .
Figure A2.(a-g) Short units and (h) long units of the reference data at validation site 31SFA-G3 for the 2021 validation year.

Figure A1 . 34 Figure
Figure A1.(a-g) Short units and (h) long units of the reference data at validation site 32SLF-V2 for the 2017 validation year.

Figure A2 .
Figure A2.(a-g) Short units and (h) long units of the reference data at validation site 31SFA-G3 for the 2021 validation year.

Figure A2 .
Figure A2.(a-g) Short units and (h) long units of the reference data at validation site 31SFA-G3 for the 2021 validation year.

Figure A3 .
Figure A3.Comparison of the burn dates from BAMTs-derived product and hotspot windows from VIIRS (a) and MODIS (b).

Figure A3 .
Figure A3.Comparison of the burn dates from BAMTs-derived product and hotspot windows from VIIRS (a) and MODIS (b).

Figure A4 .
Figure A4.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for NEALGEBA.Figure A4.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for NEALGEBA.

Figure A4 .
Figure A4.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for NEALGEBA.Figure A4.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for NEALGEBA.

Figure A5 .
Figure A5.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2021 for NEALGEBA.

Figure A6 .
Figure A6.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for GABAM.

Figure A5 . 34 Figure A5 .
Figure A5.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2021 for NEALGEBA.

Figure A6 .
Figure A6.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for GABAM.Figure A6.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for GABAM.

Figure A6 .
Figure A6.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for GABAM.Figure A6.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for GABAM.

Figure A7 .
Figure A7.Correctly classified BA, commission and omission errors, and unobserved areas in the validation sites of 2017 for FireCCI51.

Figure A8 .
Figure A8.Correctly classified BA, commission and omission errors, and unobserved areas in the validation sites of 2017 for C3SBA11.

Figure A7 . 34 Figure A7 .
Figure A7.Correctly classified BA, commission and omission errors, and unobserved areas in the validation sites of 2017 for FireCCI51.

Figure A8 .
Figure A8.Correctly classified BA, commission and omission errors, and unobserved areas in the validation sites of 2017 for C3SBA11.Figure A8.Correctly classified BA, commission and omission errors, and unobserved areas in the validation sites of 2017 for C3SBA11.

Figure A8 .
Figure A8.Correctly classified BA, commission and omission errors, and unobserved areas in the validation sites of 2017 for C3SBA11.Figure A8.Correctly classified BA, commission and omission errors, and unobserved areas in the validation sites of 2017 for C3SBA11.

Figure A9 .
Figure A9.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for MCD64A1.

Figure A10 .
Figure A10.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for EFFIS.

Figure A9 . 34 Figure A9 .
Figure A9.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for MCD64A1.

Figure A10 .
Figure A10.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for EFFIS.Figure A10.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for EFFIS.

Figure A10 .
Figure A10.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for EFFIS.Figure A10.Correctly classified BA, commission and omission errors, and unobserved areas at the validation sites of 2017 for EFFIS.

Figure A11 .
Figure A11.Examples of the BA detection anomalies observed on the GABAM annual BA maps (blue) vs. the corresponding post-fire image composites (displayed in a Long SWIR/NIR/Red color composition) and NEALGEBA (red).

Figure A11 .
Figure A11.Examples of the BA detection anomalies observed on the GABAM annual BA maps (blue) vs. the corresponding post-fire image composites (displayed in a Long SWIR/NIR/Red color composition) and NEALGEBA (red).

Table 1 .
Total land area and natural vegetation type areas and mean annual rainfall (P) in the six studied wilayas.

Table 1 .
Total land area and natural vegetation type areas and mean annual rainfall (P) in the six studied wilayas.

Table 2 .
Spatial validation results of the NEALGEBA maps of 2017 and 2021 at all the validation sites.CE: commission error, OE: omission error, OA: overall accuracy, DC: Dice coefficient, and RelB: relative bias, all expressed as percentages.SurfBA: surface correctly detected as burned, SurfUB: surface correctly detected as unburned, SurfCE: committed burned surface, SurfOE: omitted burned surface, all expressed as hectares.

Table 2 .
Spatial validation results of the NEALGEBA maps of 2017 and 2021 at all the validation sites.CE: commission error, OE: omission error, OA: overall accuracy, DC: Dice coefficient, and RelB: relative bias, all expressed as percentages.SurfBA: surface correctly detected as burned, SurfUB: surface correctly detected as unburned, SurfCE: committed burned surface, SurfOE: omitted burned surface, all expressed as hectares.

Table A3 .
Spatial validation results of the 2017 FireCCI51 map.

Table A4 .
Spatial validation results of the 2017 C3SBA11 map.

Table A2 .
Spatial validation results of the 2017 GABAM map.

Table A3 .
Spatial validation results of the 2017 FireCCI51 map.

Table A6 .
Spatial validation results of the 2017 EFFIS map.

Table A5 .
Spatial validation results of the 2017 MCD64A1 map.

Table A6 .
Spatial validation results of the 2017 EFFIS map.