Estimation of Burned Area in the Northeastern Siberian Boreal Forest from a Long-Term Data Record ( LTDR ) 1982 – 2015 Time Series

A Bayesian classifier mapped the Burned Area (BA) in the Northeastern Siberian boreal forest (70◦N 120◦E–60◦N 170◦E) from 1982 to 2015. The algorithm selected the 0.05◦ (~5 km) Long-Term Data Record (LTDR) version 3 and 4 data sets to generate 10-day BA composites. Landsat-TM scenes of the entire study site in 2002, 2010, and 2011 assessed the spatial accuracy of this LTDR-BA product, in comparison to Moderate-Resolution Imaging Spectroradiometer (MODIS) MCD45A1 and MCD64A1 BA products. The LTDR-BA algorithm proves a reliable source to quantify BA in this part of Siberia, where comprehensive BA remote sensing products since the 1980s are lacking. Once grouped by year and decade, this study explored the trends in fire activity. The LTDR-BA estimates contained a high interannual variability with a maximum of 2.42 million ha in 2002, an average of 0.78 million ha/year, and a standard deviation of 0.61 million ha. Going from 6.36 in the 1980s to 10.21 million ha BA in the 2010s, there was a positive linear BA trend of approximately 1.28 million ha/decade during these last four decades in the Northeastern Siberian boreal forest.


Introduction
Boreal forests contain as much as five times the carbon in temperate forests and almost double that in tropical forests [1][2][3][4].Boreal forests act as giant storehouses of carbon emissions to play a vital role in the global carbon cycle and in regulating climate change [5,6].However, large and recurrent wildfires, one of the most important disturbances that affect boreal forests, have just the opposite effect, increasing gas emissions and changing the forests' biophysical characteristics [7][8][9].Therefore, the estimation of Burned Areas (BA) is fundamental to quantifying the role of wildfire in the climate system.In fact, research efforts are increasing in attempts to fully understand this role.At a regional scale, the Canadian-based Boreal Ecosystem-Atmosphere Study (BOREAS) seeks to improve predictions of the CO 2 stored in boreal forests and analyzes satellite data to adjust simulations and to anticipate the impact, interactions, and its effects on climate change [10][11][12].At a global scale, the programs that stand out in this regard are the Global Observation of Forest and Landcover Dynamics (GOFC/GOLD) under the Global Terrestrial Observing System (GTOS) and the International Geosphere-Biosphere Program (IGBP).
Currently, several satellite missions provide researchers with daily imagery and higher-level land and atmosphere products required to implement accurate algorithms to map BA globally [13][14][15][16].The Moderate-Resolution Imaging Spectroradiometer (MODIS) onboard the National Aeronautics and Space Administration's (NASA) Terra and Aqua satellites, the Advanced Very High Resolution Radiometer (AVHRR) onboard the series of National Oceanic and Atmospheric Administration (NOAA) satellites, and VEGETATION (VGT) onboard Systeme Pour I'Observation de la Terre (SPOT) 4 and 5 satellites are examples of such initiatives.
Multiple algorithms have mapped BA from remote sensing data since the mid-1980s, both locally and globally.The most widely used methods detect spectral changes in land cover from pairs of preand post-fire images.In a recent study, Huang et al. [17] found that near-infrared (NIR), red-edge, and shortwave infrared (SWIR) bands are the most sensitive to the change in reflectance caused by wildfires.These results are in line with the findings of other similar studies that apply vegetation indices specifically designed to map BA [18][19][20][21].Other algorithms develop strategies based on the temperature difference between active fires and their cooler surroundings [22][23][24].However, it is very difficult to accurately estimate the actual BA with this methodology, unless the radiometric, spatial, and temporal resolutions of the remote sensor are in accordance with the temperature, size, and duration of the fire [25][26][27][28].Finally, hybrid algorithms combine both of the previous methods and can successfully map BA [22,25,[29][30][31].
There are several global BA products based on satellite data [32].For their reliability and period covered, the best BA products come from MODIS Collections 5 and 5.1 (MCD45A1) and Collection 6 (MCD64A1) [22,[33][34][35] from the year 2000 onwards; the Global Fire Emissions Database version 4 (GFED-4) BA available from mid-1995 through the present [36] and considered the longest global BA product currently available; and the most recent MODIS product, Fire_cci v5.0, that will be available for years 2001 to 2016 [37].
The enhanced temporal, spatial, and spectral resolutions of the products previously defined and the availability of supporting geospatial datasets for validation have rendered many studies that mapped BA in boreal forest in North America [19,[38][39][40][41][42][43][44].However, BA satellite remote sensing in Siberia is scarce [45][46][47][48].This region lacks validation datasets [49] because wildfire records are poorly documented [13,50,51].In addition, BA estimations from satellite data show strong discrepancies from those reported in the official records [14,16,47,48,52,53].Goldammer et al. [54] indicate that these disagreements are due to the lack of financial resources that limit the authorities from monitoring wildfires, the large extent of Siberia, and the diverse vegetation types.
Wildfires are a natural phenomenon in the boreal forests of Siberia.Most of these fires occur during the fire season that runs from May to August.Wildfire frequency correlates with air temperature anomalies and the Standardized Precipitation Evapotranspiration Index (SPEI) [51].Historically, human activities have altered that natural fire regime [55].Official statistics in Russia indicate that humans cause more than 65 percent of fires and lightning strikes cause about 17 percent.The remaining fires reported as unknown cause were also mostly due to lightning.Nevertheless, lightning is not such a critical fire cause in other fire-prone regions like the Mediterranean Basin.For example, Vazquez and Moreno [56] reported lightning as the cause of only 3.3% of fires in Spain from the year 1974 to 1994.
Despite the absence of accurate data, it is clear that the boreal forest fires in Siberia cause smoke pollution in large areas.The consequent risk to human health and other negative impacts to this ecosystem lack a corresponding protection and prevention plan.Nevertheless, its worst contribution refers to the increase in global CO 2 emissions, which in turn can significantly increase the occurrence of fires and cause a positive feedback for CO 2 and global temperatures.The Intergovernmental Panel on Climate Change predicts that Siberia will be one of most-affected regions from climate change [57].Therefore, it is essential to have reliable information on the number of fires and the area burned to assess the conversion of Siberian boreal forests from a carbon sink to a carbon source.
The goal of the present study is to map BA in the Northeastern Siberian boreal forest.The proposed method applies a Bayesian network algorithm to the Long-Term Data Record (LTDR) from NASA.
The LTDR is the longest global daily dataset available for use in climate studies.The LTDR produces a consistent long-term dataset at a spatial resolution of 0.05 • (~5 km) based on daily data from the AVHRR on board NOAA satellites, and daily data acquired by MODIS onboard NASA's Terra and Aqua satellites [58].This unique remote sensing data time series can help us to monitor and understand wildfires [59], climate, and vegetation interaction.The BA product generated here, named LTDR-BA, is compared to other BA products such as MCD45A1, MCD64A1, and GFED-4.Furthermore, Landsat-TM scenes evaluated the spatial accuracy of these products in 2002, 2010, and 2011.Siberia lacks validation data for the pre-MODIS era (1982)(1983)(1984)(1985)(1986)(1987)(1988)(1989)(1990)(1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998)(1999).Hence, we also tested the algorithm's performance in neighboring Alaska to ensure consistency between the pre-and post-MODIS eras.

Study Area
The Northeastern Siberian boreal forest study region in Figure 1 contains boreal forest similar to its neighbor Alaska.It consists mostly of rather small Larix Mill.(larch) conifers due to the severe subarctic continental climate.The mean annual precipitation runs from 1800 mm in the driest interior valleys to 3000 mm in the mountains.The maximum precipitation occurs in mid to late summers.The low evapotranspiration rates support the forest vegetation throughout the region.Temperatures at valley locations can reach 20 • C in the summer and −30 • C in the winter, decreasing with elevation and increasing with coastal proximity [60].

Study Area
The Northeastern Siberian boreal forest study region in Figure 1 contains boreal forest similar to its neighbor Alaska.It consists mostly of rather small Larix Mill.(larch) conifers due to the severe subarctic continental climate.The mean annual precipitation runs from 1800 mm in the driest interior valleys to 3000 mm in the mountains.The maximum precipitation occurs in mid to late summers.The low evapotranspiration rates support the forest vegetation throughout the region.Temperatures at valley locations can reach 20 °C in the summer and −30 °C in the winter, decreasing with elevation and increasing with coastal proximity [60].

LTDR Data Preprocessing
A NASA server provides the daily global LTDR version 4 (1981-2016) and version 3 continuation dataset (2000-2008) (https://ltdr.nascom.nasa.gov/cgi-bin/ltdr/ltdrPage.cgi).The files in hierarchical data format (HDF) contain separate Scientific Data Sets (SDS).All SDS arrays cover the globe at a 0.05 • spatial resolution as a Climate Modeling Grid (CMG) with 7200 × 3600 pixels.Further processing included the conversion from HDF to a binary sequential (BSQ) format (Table 1).Transformation to physical values considered the quality assessment (QA) fields to filter clouds.The Clouds from AVHRR Phase I (CLAVR-1) algorithm generated a cloud mask from the cloud flag in QA fields of the Daily Surface Reflectance product (AVH09C1) [61].The CLAVR-1 is a pixel-level cloud mask that combines all of the spectral information from AVHRR.Missing and invalid values found in LTDR v4, from the NOAA-16 and 18 satellites from 2000 to 2008, made the MODIS-based LTDR v3 (continuation) the preferred choice for this period.The maximum brightness temperature T criterion constructed 10-day composites from the LTDR-BSQ files for the entire study period.The objective of this step was mainly to eliminate atmospheric effects, such as residual clouds and cloud shadows, present in the images that could interfere with the discrimination of the burned pixels [62,63].This composition criterion achieved the highest performance in discerning burned/unburned area in the study region [63].

BA Algorithm
The BA algorithm (Figure 2) identifies potential fire dates for a given year using the highest burned boreal forest index criterion (BBFI) from the 10-day composite LTDR-BSQ files [42].Following Moreno-Ruiz et al. [64], the algorithm calculates a set of twelve statistical variables based on the surface reflectance bands ρ 1 and ρ 2 , the brightness temperature T, BBFI, and Global Environmental Monitoring Index (GEMI), whose expressions are BBFI ), for the pre-fire and post-fire 10-day composite of potential fire dates for the fire event year, the year before, and the year after.The twelve parameters are grouped in three categories.The first one detects burned pixels with vegetation in the year of the fire.The algorithm selects the variables that measure the maxima of BBFI, GEMI, and T; the differences of the medians of GEMI before and after the fire; and the minimum of ρ 2 in the year of the fire.The second category detects nonburned pixels in the year before the fire.The variables in this category measure the differences in the GEMI, BBFI, and ρ 2 indices' medians between the post-fire period of the year the fire occurred in and the previous year, as well as the maximum value of the BBFI index for the year before the fire.The third category finds the remaining burned pixels in the year after the fire.The variables in this category measure the maximum of BBFI and T in the year after the fire, together with the difference of the medians of the GEMI index between the year of the fire and the year after.
A Bayesian network classifier calculated the normalized probability for unburned/burned classes using a training set based on unburned/burned reference data for the study region in the year 2010.The higher positive pixels indicated a higher classification probability as a burned pixel.This normalized probability is computed by multiplying the individual conditional probabilities for each selected statistical variable and burned/unburned class based on the training set.The Bayesian network classifier is available within the Waikato Environment for Knowledge Analysis (WEKA) machine learning package (http://www.cs.waikato.ac.nz/~ml/index.html, last accessed on 3 December 2017).The algorithm also applies neighborhood conditions for spatial coherence to minimize commission and omission errors.Isolated burned pixels are downgraded to an unburned class whereas low-probability burned pixels surrounded by burned neighbors are upgraded to a burned class.The resulting probability maps label pixels as burned when their probability value is greater than 0. Annual burned/nonburned maps combined the output from the 10-day composites, after transforming to an Albers conical equal area projection with 5 km pixels.
Remote Sens. 2018, 10, x FOR PEER REVIEW 5 of 15 network classifier is available within the Waikato Environment for Knowledge Analysis (WEKA) machine learning package (http://www.cs.waikato.ac.nz/~ml/index.html, last accessed on 3 December 2017).The algorithm also applies neighborhood conditions for spatial coherence to minimize commission and omission errors.Isolated burned pixels are downgraded to an unburned class whereas low-probability burned pixels surrounded by burned neighbors are upgraded to a burned class.The resulting probability maps label pixels as burned when their probability value is greater than 0. Annual burned/nonburned maps combined the output from the 10-day composites, after transforming to an Albers conical equal area projection with 5 km pixels.

BA Time Series
The LTDR-BA product requires variables from the year, and those before and after.Therefore, it generated BA for the Northeastern Siberian boreal forest from 1982 to 2015 with data from 1981 to 2016.A comparison to other products, such as MCD45A1, MCD64A1, and GFED-4, assured the LTDR-BA's performance.These products were not always available for the entire LTDR-BA period (Table 2).The University of Maryland FTP site stores these datasets (ftp://fuoco.geog.umd.edu/).These monthly BA images needed merging to generate annual BA maps and reprojection to Albers conical equal area, to minimize distortions and maintain equivalent square shapes [39].The MCD45A1 (Collection 5.1) contains the burn date and detection confidence levels 1 through 4. The MCD64A1 (Collection 6) has five layers: the burn date, the burn date uncertainty, the quality assurance, and the first and last burning days.The GFED-4 includes seven layers: the BA (ha/100), the BA uncertainty, the source, the tree cover distribution, the land cover distribution, the fire persistence, and the peat fraction.The LTDR-BA comparison to GFED-4 included only from 1995 to 2000, since GFED-4 uses MCD64A1 after year 2000.

BA Time Series
The LTDR-BA product requires variables from the year, and those before and after.Therefore, it generated BA for the Northeastern Siberian boreal forest from 1982 to 2015 with data from 1981 to 2016.A comparison to other products, such as MCD45A1, MCD64A1, and GFED-4, assured the LTDR-BA's performance.These products were not always available for the entire LTDR-BA period (Table 2).The University of Maryland FTP site stores these datasets (ftp://fuoco.geog.umd.edu/).These monthly BA images needed merging to generate annual BA maps and reprojection to Albers conical equal area, to minimize distortions and maintain equivalent square shapes [39].The MCD45A1 (Collection 5.1) contains the burn date and detection confidence levels 1 through 4. The MCD64A1 (Collection 6) has five layers: the burn date, the burn date uncertainty, the quality assurance, and the first and last burning days.The GFED-4 includes seven layers: the BA (ha/100), the BA uncertainty, the source, the tree cover distribution, the land cover distribution, the fire persistence, and the peat fraction.The LTDR-BA comparison to GFED-4 included only from 1995 to 2000, since GFED-4 uses MCD64A1 after year 2000.

Accuracy Assessment of LTDR-BA
The existence of BA products with higher spatial resolutions in the MODIS era split the analysis into before and after the year 2000 to assess the LTDR-BA accuracy.This assessment plotted annual LTDR-BA in comparison with MCD45A1 and MCD64A1 for the MODIS era.The test also included correlation analysis between the products.The lack of official BA data in the study region made it necessary to build a reference BA for years 2002, 2010, and 2011 with a large set of 152 Landsat-TM images obtained from the United States Geological Survey (USGS)-provided pre-and post-fire information at 30 m spatial resolution (Table 3).The year 2010 was selected because it trained the LTDR-BA algorithm and the other two years because they rendered the highest discrepancies between LTDR-BA and MCD64A1.Expert visual interpretation identified and delineated fires over 1000 ha.The interpreter visualized pre-and post-fire false color Landsat-TM composites with bands 7 (2.09-2.35µm), 4 (0.77-0.90 µm), and 3 (0.63-0.69 µm).The following step was to resample the results to 500 m and transform to Albers conical equal area projection.A linear regression between the annual BA products and the reference evaluated the BA proportions on 50 km × 50 km grids.The distribution of these grids was uniform over the entire region.The commission and omission errors gathered the data at this 50 km spatial resolution, to prevent errors due to co-georeferencing of the images [43].Landsat-TM images in the USGS database are scarce for the pre-MODIS era in the Northeastern Siberian region.Consequently, the LTDR-BA algorithm tested its consistency in Alaska since it is close to Siberia and similar in latitude, vegetation, and weather conditions.Cloud cover frequency could differ between Siberia and Alaska, but the selection of the best cloud-free radiometric values from 10-day composites minimizes the possible impact of this factor.Ground truth data from the Alaskan Fire Service (AFS) are available [39,42,[65][66][67].The 1982-1999 fire polygons registered in the AFS database served as the BA ground truth for the Alaskan region (70 • N, 168.5 • E, 60 • N, −141 • E).The maximum area within the pixel was assigned a burned/nonburned value after resampling the results to 500 m and transforming to Albers conical equal area projection.These reference annual BA maps used the percentage BA at the subpixel level to aggregate to 5 km.The accuracy assessment plotted annual LTDR-BA in comparison with the AFS dataset.The test also included correlation analysis between them.The comparison to GFED-4 was also possible since 1995, not only for Alaska but also for the Northeastern Siberian region.

Annual NE Siberian Boreal Forest LTDR-BA Maps for the Period 1982-2015
The LTDR-BA algorithm generated annual maps of the Northeastern Siberian boreal forest for the 1982-2015 period.Grouping BA by decade in Figure 3 demonstrates the largest BA in the 2000s.Furthermore, a linear regression model built with normalized BA suggests a positive trend in BA with time (R 2 = 0.8335, slope = 128,323.5,intersect = −247,719,953 and p-value < 0.1) (Figure 4).This normalization process assigns the average BA/year of its decade to the years with no data.Therefore, 1980 and 1981 receive the average BA/year from the 1980s decade and 2016, 2017, 2018, and 2019 take that from the 2010s.The LTDR-BA estimated approximately 27.49 million ha for the 1982-2015 period with an uneven annual distribution.The largest BA of 2.63 million ha occurred in 2002, whereas the lowest was in 2004, with only 0.14 million ha (Figure 5).MCD64A1 estimated the largest annual BA values, followed by LTDR-BA and MCD45A1.The highest discrepancy among the products happened in 2002, which was the year with the largest BA.The R 2 between MCD64A1 and MCD45A1, and LTDR-BA in the common period (2000-2015) were 0.929 and 0.890, respectively.GFED-4 estimated larger annual BA than LTDR-BA in the period of overlap (1995-2000), with higher discrepancies in 1996 and 1998 and a lower R 2 of 0.438.The LTDR-BA and GFED-4 estimated similar BA for the year 2000, larger than the other two BA products.

Annual NE Siberian Boreal Forest LTDR-BA Maps for the Period 1982-2015
The LTDR-BA algorithm generated annual maps of the Northeastern Siberian boreal forest for the 1982-2015 period.Grouping BA by decade in Figure 3 demonstrates the largest BA in the 2000s.Furthermore, a linear regression model built with normalized BA suggests a positive trend in BA with time (R 2 = 0.8335, slope = 128,323.5,intersect = −247,719,953 and p-value < 0.1) (Figure 4).This normalization process assigns the average BA/year of its decade to the years with no data.Therefore, 1980 and 1981 receive the average BA/year from the 1980s decade and 2016, 2017, 2018, and 2019 take that from the 2010s.The LTDR-BA estimated approximately 27.49 million ha for the 1982-2015 period with an uneven annual distribution.The largest BA of 2.63 million ha occurred in 2002, whereas the lowest was in 2004, with only 0.14 million ha (Figure 5).MCD64A1 estimated the largest annual BA values, followed by LTDR-BA and MCD45A1.The highest discrepancy among the products happened in 2002, which was the year with the largest BA.The R 2 between MCD64A1 and MCD45A1, and LTDR-BA in the common period (2000-2015) were 0.929 and 0.890, respectively.GFED-4 estimated larger annual BA than LTDR-BA in the period of overlap (1995-2000), with higher discrepancies in 1996 and 1998 and a lower R 2 of 0.438.The LTDR-BA and GFED-4 estimated similar BA for the year 2000, larger than the other two BA products.

Accuracy Assessment of LTDR-BA
The BA accuracy assessment with the Landsat-TM images for 2002, 2010, and 2011 in Table 4 revealed that MCD64A1 performed the best (88.80%), followed by LTDR-BA (64.68%) and MCD45A1 (43.51%).It is worth mentioning that the three products obtained the lowest total BA in 2011.This year had also a more fragmented BA, with smaller fire sizes.The R 2 between the MCD45A1, MCD64A1, and LTDR-BA, and the reference Landsat-TM in 50 km × 50 km grids also indicated higher R 2 and a slope closer to 1 for MCD64A1, followed by LTDR-BA and MCD45A1 (Figure 6).LTDR-BA incurred the highest commission errors, while MCD45A1 rendered the highest omission errors (Table 5).

Accuracy Assessment of LTDR-BA
The BA accuracy assessment with the Landsat-TM images for 2002, 2010, and 2011 in Table 4 revealed that MCD64A1 performed the best (88.80%), followed by LTDR-BA (64.68%) and MCD45A1 (43.51%).It is worth mentioning that the three products obtained the lowest total BA in 2011.This year had also a more fragmented BA, with smaller fire sizes.The R 2 between the MCD45A1, MCD64A1, and LTDR-BA, and the reference Landsat-TM in 50 km × 50 km grids also indicated higher R 2 and a slope closer to 1 for MCD64A1, followed by LTDR-BA and MCD45A1 (Figure 6).LTDR-BA incurred the highest commission errors, while MCD45A1 rendered the highest omission errors (Table 5).

Accuracy Assessment of LTDR-BA
The BA accuracy assessment with the Landsat-TM images for 2002, 2010, and 2011 in Table 4 revealed that MCD64A1 performed the best (88.80%), followed by LTDR-BA (64.68%) and MCD45A1 (43.51%).It is worth mentioning that the three products obtained the lowest total BA in 2011.This year had also a more fragmented BA, with smaller fire sizes.The R 2 between the MCD45A1, MCD64A1, and LTDR-BA, and the reference Landsat-TM in 50 km × 50 km grids also indicated higher R 2 and a slope closer to 1 for MCD64A1, followed by LTDR-BA and MCD45A1 (Figure 6).LTDR-BA incurred the highest commission errors, while MCD45A1 rendered the highest omission errors (Table 5).The AFS database registered approximately 5.56 million ha BA in the 1982-1999 period for the Alaskan boreal region, with an uneven annual distribution (Figure 7).The largest BA of 1.29 million ha occurred in 1990, whereas the lowest was in 1995, with only 0.02 million ha burned.The LTDR-BA product trained in the Northeastern Siberian boreal region followed the AFS BA pattern, with an overall underestimation.This was not the case for the GFED-4 product, since it reported a similar total BA from 1995 to 1999.The highest discrepancy between LTDR-BA and AFS occurred in 1997, when the overall R 2 was 0.898.
The AFS database registered approximately 5.56 million ha BA in the 1982-1999 period for the Alaskan boreal region, with an uneven annual distribution (Figure 7).The largest BA of 1.29 million ha occurred in 1990, whereas the lowest was in 1995, with only 0.02 million ha burned.The LTDR-BA product trained in the Northeastern Siberian boreal region followed the AFS BA pattern, with an overall underestimation.This was not the case for the GFED-4 product, since it reported a similar total BA from 1995 to 1999.The highest discrepancy between LTDR-BA and AFS occurred in 1997, when the overall R 2 was 0.898.

Discussion
Annual LTDR-BA estimations are comparable to other BA products with higher spatial resolutions, like MCD45A1 and MCD64A1 in the MODIS era (Tables 4 and 5 and Figures 5 and 6).The detailed analysis for the years 2002, 2010, and 2011 shows significant variability among the positional accuracy of the BA estimations, with the worst results in 2011, especially for the LTDR-BA product.The count, size, fragmentation, and severity of the fires registered that year could explain these higher errors.Moreno-Ruiz et al. [43] also found similar variability in their BA accuracy analysis for the North American boreal region.In addition, the analysis of the scatter plots of the different products versus the reference Landsat-TM data using 50 km × 50 km grids in Figure 6 reveals that the average determination coefficients of the linear regression model of the BA percentage versus the reference data are similar.
LTDR-BA rendered a 15% commission error (Table 5).This value is similar to that for the official MCD64A1 MODIS product.In contrast, the 43% omission error is relatively large in comparison to the 23% of the MCD64A1, but still lower than the 57% of the MCD45A1 product.Moreno-Ruiz et al. [43] analyzed the Pareto boundaries for the North American BA from Alaska and Canada reference data versus the LTDR-BA product.They demonstrated that the main contribution to the commission and omission errors came from the low data spatial resolution in relation to the BA size.This fact leaves little room for LTDR-BA algorithm improvement, whereas other algorithms applied to higher spatial resolutions like MODIS could benefit from further development.
LTDR-BA and MCD64A1 products have a coincident temporal pattern for the MODIS era.The 2000-2008 period shows the greatest discrepancy between them, corresponding to the use of LTDR

Discussion
Annual LTDR-BA estimations are comparable to other BA products with higher spatial resolutions, like MCD45A1 and MCD64A1 in the MODIS era (Tables 4 and 5 and Figures 5 and 6).The detailed analysis for the years 2002, 2010, and 2011 shows significant variability among the positional accuracy of the BA estimations, with the worst results in 2011, especially for the LTDR-BA product.The count, size, fragmentation, and severity of the fires registered that year could explain these higher errors.Moreno-Ruiz et al. [43] also found similar variability in their BA accuracy analysis for the North American boreal region.In addition, the analysis of the scatter plots of the different products versus the reference Landsat-TM data using 50 km × 50 km grids in Figure 6 reveals that the average determination coefficients of the linear regression model of the BA percentage versus the reference data are similar.
LTDR-BA rendered a 15% commission error (Table 5).This value is similar to that for the official MCD64A1 MODIS product.In contrast, the 43% omission error is relatively large in comparison to the 23% of the MCD64A1, but still lower than the 57% of the MCD45A1 product.Moreno-Ruiz et al. [43] analyzed the Pareto boundaries for the North American BA from Alaska and Canada reference data versus the LTDR-BA product.They demonstrated that the main contribution to the commission and omission errors came from the low data spatial resolution in relation to the BA size.This fact leaves little room for LTDR-BA algorithm improvement, whereas other algorithms applied to higher spatial resolutions like MODIS could benefit from further development.
LTDR-BA and MCD64A1 products have a coincident temporal pattern for the MODIS era.The 2000-2008 period shows the greatest discrepancy between them, corresponding to the use of LTDR version 3 (continuation) obtained from the MODIS.During these years, it was not possible to select LTDR version 4 obtained from AVHRR-NOAA N16 and N18, due to the lack of radiometric coherence and missing data.Given that version 3 (continuation) lacks the corresponding T3 channel of AVHRR, LTDR-BA had to process the MODIS channel equivalent to the T4 channel.This band is less sensitive to the active fire response, which resulted in a lower BA estimate [43,69].Additionally, it should also be noted that the Bayesian algorithm makes use of the year prior to the occurrence of fire to reduce false detections.The MODIS sensor was not operational in June of 2001 (https://modaps.modaps.eosdis.nasa.gov/services/production/outages_terra.html), which also has a negative influence on the results of the algorithm for the year 2002.
The Alaskan boreal region served to validate the LTDR-BA for the pre-MODIS era, since there are no reference values available in Northeastern Siberia for these years.There was a strong correlation between the LTDR-BA data and the reference AFS data (Figure 7).In contrast, GFED-4 data did not match the LTDR-BA data, nor the reference AFS.GFED-4 and LTDR-BA adjust particularly well for this period in Northeastern Siberia.Years 1995-1999 correspond to AVHRR-N14 in the LTDR data set.A poor adjustment of the calibration coefficients occurred due to the degradation of N14 channels 1 and 2 [70].
Decadal LTDR-BA projections for the Northeastern Siberian boreal region indicate a linear increase of 1.28 million ha/decade.Alaska rendered similar results from the AFS with a linear increase of 1.83 million ha/decade.The LTDR-BA time series delivered a high correlation (R 2 = 0.898) with the AFS database in Alaska in the period 1982-2015 (Figure 6).The estimated linear increase for Alaska would be 1.58 million ha/decade, similar to that obtained in Siberia.These results suggest that the effect of fire on boreal forests has a homogeneous incremental pattern of area burned over the last four decades.Numerous studies on global also change estimate growing fire activity in boreal regions as a result of global warming [71][72][73].Extrapolating this trend to the entire boreal region (North America and Eurasia) results in a total increase of about 50 million ha of BA per decade.

Conclusions
This study built a comprehensive BA time series  for the Northeastern Siberian boreal region using a Bayesian network algorithm from the LTDR v3 continuation and v4 datasets at 0.05 • spatial resolution.This time series is the longest developed for this region.It constitutes a unique data source for the pre-MODIS era (1980s and 1990s) for long-term studies of fire dynamics in this region.
The LTDR-BA product showed comparable results in terms of spatial and temporal accuracy with other BA products with higher spatial resolution such as MCD45A1 and MCD64A1, where the latter obtained the overall best results.The LTDR-BA, MCD64A1, and MCD45A1 products underestimated BA in comparison with the reference Landsat-TM scenes for 2002, 2010, and 2011 by 35%, 11%, and 56%, respectively.The annual BA distributions follow a similar temporal pattern, with over 0.95 R 2 .
The LTDR-BA algorithm trained in Northeastern Siberia underestimated BA in Alaska for the pre-MODIS era by approximately 45%, with respect to the ground truth data obtained from the AFS.The annual BA distributions followed a similar temporal pattern, with over 0.95 R 2 .Consequently, this study infers a BA underestimation for Northeastern Siberia of between 35% and 45%.LTDR-BA followed a similar temporal pattern to the GFED-4 for the available years 1995-1999.This product becomes a reliable source for the 1980s and 1990s, when other products with a higher spatial/spectral resolution like MCD64A1 and MCD45A1 are not available.It is also the longest continuous BA time series for this region.As a result, almost four decades of LTDR-BA data grouped by year and decade estimated a linear increase in BA of 1.28 million ha/decade.This increasing trend provides new evidence of how global warning affects fire activity in the Northeastern Siberian boreal forest.

Figure 2 .
Figure 2. Flowchart with the iterative design process of the Bayesian algorithm (a) and its application later to obtain annual maps of burned area (b).

Figure 2 .
Figure 2. Flowchart with the iterative design process of the Bayesian algorithm (a) and its application later to obtain annual maps of burned area (b).

Figure 3 .
Figure 3. BA (red) differentiated from nonburned (green) and water bodies (blue) for the Northeastern Siberian boreal region for the 1982-2015 period grouped by decade.Figure 3. BA (red) differentiated from nonburned (green) and water bodies (blue) for the Northeastern Siberian boreal region for the 1982-2015 period grouped by decade.

Figure 3 .
Figure 3. BA (red) differentiated from nonburned (green) and water bodies (blue) for the Northeastern Siberian boreal region for the 1982-2015 period grouped by decade.Figure 3. BA (red) differentiated from nonburned (green) and water bodies (blue) for the Northeastern Siberian boreal region for the 1982-2015 period grouped by decade.

Figure 4 .
Figure 4.Estimated and normalized decadal LTDR Burned Area for the Northeastern Siberia boreal region and the estimated linear trend.

Figure 4 .
Figure 4.Estimated and normalized decadal LTDR Burned Area for the Northeastern Siberia boreal region and the estimated linear trend.

Figure 4 .
Figure 4.Estimated and normalized decadal LTDR Burned Area for the Northeastern Siberia boreal region and the estimated linear trend.

Figure 6 .
Figure 6.Scattering plots of the BA proportions in 50 km × 50 km grids from the different BA products vs reference data for years 2002, 2010, and 2011 for the Northeastern Siberian boreal region.p-value < 0.001 for all cases.

Figure 7 .
Figure 7. 1982-2015 annual Alaskan BA from GFED-4 and LTDR BA products in comparison to the Alaskan Fire Service (AFS) database.
Furthermore, the LTDR-BA dataset provides the largest BA time series for the Northeastern Siberian boreal region from remote sensing data from 1982 to 2015.No other reliable BA sources for this study region are available in the pre-MODIS era.Over 20 million ha burned throughout the Russian Federation in the period 1982-1999, according to the official data provided by Avialesookhrana (https://www.fire.uni-freiburg.de/iffn/iffn_37/15-IFFN-37-Russia.pdf).Nevertheless, LTDR-BA detected 12 million ha in the Northeastern region of Siberia alone.Other studies also confirm high discrepancies with the official data.For example, Conard et al. [68] estimated 13.3 million ha burned in the year 1998 in Russia from AVHRR data, compared to the 3 million ha registered by Avialesookhrana.Additionally, Sukhinin et al. [48] also estimated 4.2 times the official area from AVHRR images in the period 1996-2002.

Figure 7 .
Figure 7. 1982-2015 annual Alaskan BA from GFED-4 and LTDR BA products in comparison to the Alaskan Fire Service (AFS) database.
Furthermore, the LTDR-BA dataset provides the largest BA time series for the Northeastern Siberian boreal region from remote sensing data from 1982 to 2015.No other reliable BA sources for this study region are available in the pre-MODIS era.Over 20 million ha burned throughout the Russian Federation in the period 1982-1999, according to the official data provided by Avialesookhrana (https://www.fire.uni-freiburg.de/iffn/iffn_37/15-IFFN-37-Russia.pdf).Nevertheless, LTDR-BA detected 12 million ha in the Northeastern region of Siberia alone.Other studies also confirm high discrepancies with the official data.For example, Conard et al. [68] estimated 13.3 million ha burned in the year 1998 in Russia from AVHRR data, compared to the 3 million ha registered by Avialesookhrana.Additionally, Sukhinin et al. [48] also estimated 4.2 times the official area from AVHRR images in the period 1996-2002.

Table 1 .
Long-Term Data Record (LTDR)-binary sequential (BSQ) and their respective bands from the original LTDR-hierarchical data format (HDF) files.

Table 3 .
Reference BA data from Landsat-TM images.

Table 4 .
Relative MCD45A1, MCD64A1, and LTDR BA percentages versus the reference Landsat-TM data for the Northeastern Siberian boreal region in 2002, 2010, and 2011.

Table 4 .
Relative MCD45A1, MCD64A1, and LTDR BA percentages versus the reference Landsat-TM data for the Northeastern Siberian boreal region in 2002, 2010, and 2011.

Table 4 .
Relative MCD45A1, MCD64A1, and LTDR BA percentages versus the reference Landsat-TM data for the Northeastern Siberian boreal region in 2002, 2010, and 2011.

Table 5 .
Commission Scattering plots of the BA proportions in 50 km × 50 km grids from the different BA products vs reference data for years 2002, 2010, and 2011 for the Northeastern Siberian boreal region.p-value < 0.001 for all cases.