Prescribed Burning Reduces Large, High-Intensity Wildﬁres and Emissions in the Brazilian Savanna

: Brazil has recently (2014) changed from a zero-ﬁre policy to an Integrated Fire Management (IFM) program with the active use of prescribed burning (PB) in federal Protected Areas (PA) and Indigenous Territories (IT) of the Brazilian savanna (Cerrado). PB is commonly applied in the management of ﬁre-prone ecosystems to mitigate large, high-intensity wildﬁres, the associated emissions, and high ﬁre suppression costs. However, the effectiveness of such ﬁre management in reducing large wildﬁres and emissions over Brazil remains mostly unevaluated. Here, we aim to ﬁll the gap in the scientiﬁc evidence of the PB beneﬁts by relying on the most up-to-date, satellite-derived ﬁre datasets of burned area (BA), ﬁre size, duration, emissions, and intensity from 2003 to 2018. We focused on two Cerrado ITs with different sizes and hydrological regimes, Xerente and Araguaia, where IFM has been in place since 2015. To understand ﬁre regime dynamics, we divided the study period into three phases according to the prevalent ﬁre policy and the individual ﬁre scars into four size classes. We considered two ﬁre seasons: management ﬁre season (MFS, which goes from rainy to mid-dry season, when PBs are undertaken) and wildﬁres season (WFS, when PBs are not performed and ﬁres tend to grow out of control). Our results show that the implementation of the IFM program was responsible for a decrease of the areas affected by high ﬁre recurrence in Xerente and Araguaia, when compared with the Zero Fire Phase (2008–2013). In both regions, PB effectively reduced the large wildﬁres occurrence, the number of medium and large scars, ﬁre intensity, and emissions, changing the prevalent ﬁre season from the WFS to the MFS. Such reductions are signiﬁcant since WFS causes higher negative impacts on biodiversity conservation and higher greenhouse gas emissions. We conclude that the effect on wildﬁres can still be reduced if effective ﬁre management policies, including PB, continue to be implemented during the coming decades.


Introduction
Savanna ecosystems encompass around 20% of the world's land surface, one of the most fire-prone landscapes [1,2], contributing to approximately 60% of gross global mean fire emissions [1]. Over the past two decades, savanna biomes have been subjected to remote sensing techniques have helped identify fire patterns, characteristics, and changes worldwide, including tracing burning patterns resulting from fire management efforts [37]. However, few studies have focused on such satellite-derived databases to analyze PB patterns in Brazil. Emergent literature seeks to show how PB has been contributing to decreasing wildfires in Brazilian savannas, mainly based on field campaigns or local observation [25,26,[38][39][40][41][42][43][44]. The majority of those studies rely on statistical summaries of active fire counts and total values of BA [6,39,43,45,46]. Therefore, a systematic analysis considering additional fire characteristics, such as fire intensity, date of burn occurrence, emissions, and patchiness, is urgently needed to enable further characterization of PB benefits.
Here, we aim to analyze how PB implementation affects the occurrence of large wildfires, relying on the most up-to-date satellite-derived fire datasets of BA, fire size, fire danger, fire intensity, and associated emissions from 2003 to 2018. We focus on two ITs over the Cerrado, namely Xerente and Araguaia, where PB has been conducted regularly within the IFM program since 2015. We examine variations in fire regime defining three distinct fire policy and management phases: No Fire-policy Phase (NFP)-from 2003 to 2007, a period characterized by an absence of fire policies; Zero Fire Phase (ZFP)-from 2008 to 2013, when the zero-fire approach (fire use prohibition to manage the landscape) was prevalent; Integrated Fire Management Phase (IFMP)-from 2014 to 2018 when the IFM was implemented.

Study Region
The Cerrado is a mosaic of water dynamics and different ecosystems [47] that are subject to high intra-annual variability in precipitation [48] and a well-defined dry season in winter, generally from May to October [49]. According to Silva et al. [50], the fire intraannual trends Cerrado is spatially heterogeneous with marked north-south fire activity gradient, with increasingly large and less intense fires in the latest agricultural frontier in the north. More intense fires are found in the Amazon frontier biome due to uncontrolled deforestation and different histories of land conversion. Throughout the year, the infrequent large fires occur in late-dry-season and are responsible for the majority of burned area, whereas smaller fires, albeit more frequent in early-dry-season, have been decreasing in number [50].
The two study areas encompass the Indigenous Territories of Xerente and Parque Araguaia, located in the state of Tocantins within Cerrado (Figure 1a), hereafter Xerente and Araguaia, respectively. These two territories were homologated by the federal government in 1989 and 1998, respectively [51]. These areas have a tropical wet-dry climate, according to Köppen-Geiger's definition, with annual precipitation above 1600 mm, with 85% concentrated from November to April, i.e., mean rainfall in the dry season months is below 80 mm/month. The yearly average temperature is 24.9 C, with the hottest months being August, September, and October [52].
Xerente covers 162,542 ha, with the indigenous population [53] distributed in small settlements throughout the territory and an economy based on social programs, public jobs, itinerant agriculture, hunting, gathering, and handcrafts. They possess a vast traditional knowledge of how to manage the Cerrado with fire, which is used to plan PB activities. Araguaia covers 1,358,499 ha of the largest river island in the world (Bananal Island) [54]. The villages are distributed along the Araguaia and Javaé rivers, and the island's interior is sparsely populated. These areas are flat, subject to flooding, and dominated by native grasslands used for grazing in partnership with non-indigenous ranchers. Biomass production is very high in these areas, allowing extreme and frequent fire events, often associated with cattle activities [55].  to 2020 for Xerente (black) and Araguaia (grey). Dashed lines represent the trendline obtained using nonparametric Theil-Sen Robust Linear Regression (TS) and a two-tailed Mann-Kendall (MK) test of significance. NFP, ZFP, and IFMP are highlighted in blue, yellow and green, respectively. (c) FWI intra-annual variability for Xerente (black) and Araguaia (grey).
Among the ITs served by the IFM program, Xerente and Araguaia have the longest history of official PB implementation [46]. Xerente received the first fire control activities in 2008 based on zero-fire policies through the Tocantínia municipal brigade, composed mainly of indigenous people. From 2014, fire began to be increasingly used in land man-agement. Xerente became one of the model brigades, and its brigades have been frequently requested to help in management and firefighting operations around the country [56]. Araguaia also received the first fire protection programs in 2008, generally associated with brigades located in nearby municipalities, with little indigenous participation. Fire use was prohibited, and major combat operations were registered in the savannas. Since 2014, with the indigenous brigade's implementation and the IFM program beginning, firefighting operations in savanna areas have reduced, but they still occur in forest areas. In these two regions, the IFM program's implementation included the professionalized indigenous brigade's training, vehicles, infrastructure, and safety equipment. The aim is to use local traditional knowledge in defining the fire regimes to be implemented [57,58].

Fire, Climate, and Emissions Datasets
We used 5 datasets, which were as follows: (i) PB locations and dates provided by the IBAMA/PREVFOGO (the Brazilian Environmental Agency responsible for the IFM program in TIs); (ii) individual fire scars from the Global Fire Atlas (GFA) database [31]; (iii) Fire Radiative Power (FRP) information derived from MCD14ML active fire product [59]; (iv) Fire Weather Index (FWI), derived from the ERA5 climate reanalysis [60]; and (v) fine particulate matter (PM 2.5 ) concentrations derived from the Copernicus Atmosphere Monitoring Service (CAMS) [61]. These are described in detail below and summarized in Table 1. PB information, derived by IBAMA/PREVFOGO, encompasses the geographical coordinates (latitude/longitude in degrees) and dates (in Julian days) carried out from 2015 to 2018. This period was when the IFM began its implementation over the two study regions; however, the extent of PB is not available.
GFA is a global fire database derived from the MCD64A1 BA collection 6 product [31] that provides individual fire characteristics based on the information generated by 8-day images from the MODIS satellite [62]. Here, we relied on individual fire scars over the study region containing scar size (ha) and date of occurrence (Julian day) from 2003 to 2018. Since GFA is based on the MCD64A1 500-m product, the smallest BA mapped is around 25 ha, precluding the identification of very small scars [31]. As PB may be small, there is no guarantee that this satellite-derived database can detect all burned areas from the IFM program. Thus, for the sake of simplicity, we assumed that GFA scars account for all fires, indistinguishably, both planned and unplanned.
FRP represents a measure of the instantaneous release of combustion energy and is associated with the fire intensity of the burning process [50]. We used FRP (MW) data derived from 1 km MODIS/Aqua Thermal Anomalies/Fire locations from 2003 to 2018, and only pixels with presumed vegetation fires (type 0) were used to minimize false alarms. The MCD14ML FRP product has limited detection above the threshold of 9-11 MW [63]. However, very low FRP mainly occurs away from the diurnal peak of fire activity in Cerrado (between 15-18h local time) [64], thus not impairing our assessments.
FWI is based on the ERA5 reanalysis variables from the European Centre for Medium-Range Weather Forecasts (ECMWF) [65], namely, the daily values of air temperature and relative humidity at 2 m-height, wind speed at 10 m-height, and total precipitation, from 1979 until the present. FWI provides a complete historical reconstruction of meteorological conditions favorable to the start, spread, and sustainability of fires and is based on two indices-namely, Initial Spread Index and Build-Up Index-which are the fire spread rate and fuel available for combustion, respectively [60]. This is one of the most reliable and globally applied fire weather indices and has shown to be particularly suitable to rate meteorological fire danger for the Cerrado ecosystems [49,[66][67][68]. Here, we used information from 1980 to 2020 at 0.25 spatial resolution.
We used the monthly PM 2.5 values from the fourth-generation global CAMS reanalysis dataset (EAC4) [61] with 0.75 spatial resolution. PM 2.5 refers to particles with a diameter of less than 2.5 µm and remain suspended longer in the air than heavier particles, where fires are one of their main sources [69]. It is worth mentioning that the EAC4 dataset can ensure good temporal consistency [61].

Pre-Processing and Statistical Analysis
We analyzed temporal patterns and changes of PB and climate in each region. To this end, the number and location of PBs conducted since 2015 were analyzed. Moreover, climate trends in each region were analyzed from 1980 until 2020 through FWI. We searched for significant trends using the nonparametric Theil-Sen Robust Linear Regression (TS) and a two-tailed Mann-Kendall (MK) test of significance [74,75]. The TS estimator is robust against outliers and suitable for application in time series with significant interannual variations, such as annual fire counts, precipitation, and deforestation [76,77]. The null distribution of the MK test statistics follows a normal distribution when n 8 [74,75]. The p-value of an MK statistic was then determined using the normal cumulative distribution function [78]. The trend analysis was conducted to analyze the FWI values dataset for each region through pyMannKendall, a package in the python programming language [79].
Then, we performed a comprehensive analysis of BA and the total number of scars during the last two decades using the GFA dataset, according to 4 scar size classes: very small (VS, 0-50 ha); small (S, 50-100 ha); medium (M, 100-1000 ha); large (L, >1000 ha), based on previous works [50,[80][81][82]. For each scar size class, we conducted separate analyses for the 2 fire seasons, determined by management practices: management fire season (MFS) encompassing the months when PBs are conducted, and wildfire season (WFS), corresponding to the months when PBs are not conducted, and fires are considered "bad-fires" by local inhabitants. The MFS was determined in each IT based on local ecological knowledge, as the periods when fire impacts on biodiversity reduce due to biophysical conditions, plant phenology (especially flowering and fruiting periods), and animal reproduction/nesting/nursing periods [16,46,83,84]. In Xerente, PBs are carried out from January to June (MFS), whereas any fires between July and December are considered wildfires (WFS). Due to soil flooding dynamics, the MFS in Araguaia extends from January to July, and the WFS from August to December [46]. As mentioned before, the IFM program implemented in 2015 over Cerrado has 3 main objectives: (I) to decrease late-dry-season wildfires; (II) to distribute fires over the year, decreasing extreme wildfire damages in the late-dry-season; (III) to protect fire-sensitive vegetation. Here, we aim to evaluate IFM objectives (I) and (II). Specifically, we addressed IFM objective (I) using 3 different approaches: We divided the study period (2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016)(2017)(2018), where satellite-derived information is available) into 3 distinct phases: NFP consists of the period from 2003 to 2007 when there was almost a lack of consistent fire policies over the region [44], hereafter the "No Fire-policy Phase". ZFP includes the period from 2008 until 2013, when a total fire suppression policy, hereafter the "Zero Fire Phase", was applied [16,17,85]. Finally, IFMP goes from 2014 to 2018, thus representing the IFM implementation including "controlled" and "prescribed" burns, hereafter the "Integrated Fire Management Phase". To analyze the differences between phases, for each fire season, we first considered the changes between (i) NFP and ZFP and (ii) ZFP and IFMP, considering annual BA and scar number (SN) frequency in each IT.
First, we examined variations in annual fire patterns along with the 3 distinct fire policy phases and 2 fire season periods employing boxplots. Second, for the WFS, we compared the Kernel probability density function (PDF) of (i) daily BA for large scars, (ii) daily fire intensity, and (iii) monthly fire emission during the 3 phases, highlighting extreme events using the 90th percentile (hereafter BAp90, FRPp90, and PMp90, respectively) of each distribution. For those extreme-related indexes, we computed the percentage changes between (i) NFP and ZFP and (ii) ZFP and IFMP. Changes between phases were confirmed based on 1000 samples generated by the bootstrapping approach, which can be used to assess variation in confidence between independent variable samples [86].
Third, a spatial analysis was conducted to evaluate changes in fire recurrence of large fire scars along with the 3 phases. For each pixel, fire recurrence was defined as the average number of years with fire occurrence during each phase. Fire recurrence was classified into 3 classes: zero (no fire occurrence), low (less or equal than 0.6 fire occurrence per year), and high (greater than 0.6 fire occurrence per year). We analyzed the proportion of pixels (i.e., area) in each fire recurrence class and evaluated changes in this amount between the 3 phases.
To evaluate IFM objective II, we analyzed seasonal patterns and changes of BA distribution according to each scar size class, fire period, and phase.

Overview of PB, Fire and Climate Trends
The number and geographical distribution of PB have increased in both study areas from 2015 until 2018 (Table 2 and Figures S1 and S2), according to the dataset derived by IBAMA/PREVFOGO. When the IFM program began, the PB number was low in both regions. However, PB is widespread around both regions after three years, increasing 8 and 70 times in Xerente and Araguaia, respectively. Annual FWI values from 1980 to 2020 revealed a systematic increase in both regions since 1980, in response to higher regional temperatures and lower relative humidity and precipitation. FWI has increased by 31% and 47% over Xerente and Araguaia during the last four decades, respectively. In addition, we observed a positive and significant (p < 0.001) trend for both study areas with higher FWI values in Araguaia compared with Xerente ( Figure 1b). The mean annual cycles of temperature, relative humidity, accumulated precipitation, and wind during the 40-year period are also represented through FWI (Figure 1c). Both regions present strong seasonality; the austral summer months are associated with very low values of FWI, reflecting high values of both relative humidity and precipitation (the rainy season). The higher values of FWI are observed between May and October, with maximum values in August, and consistent with the increasingly warmer and drier conditions over the region. Araguaia shows higher FWI values than Xerente, especially during the dry season.
During the MFS in Xerente, large scars were infrequent but account for around onequarter of the total BA (Figure 2a). In contrast, medium scars were more frequent and account for a large portion (65%) of the total BA. In the WFS (Figure 2c), large scars represent just one-quarter of the total number of scars but account for more than 75% of total BA in the region. Conversely, very small and small, together represent around half and one-quarter of the fire scars in the MFS and WFS, respectively; however, their contribution to the total BA is negligible in both fire seasons.
In Araguaia, very small and small fires together account for less than 25% of scars and represent less than 1% of the total BA in the MFS (Figure 2b) and WFS (Figure 2d). Medium scars prevail in both seasons, particularly during the MFS, though their contribution to the total BA is less than 25% (5%) in the MFS (WFS). In contrast, although less frequently, large scars correspond to more than 75% (95%) of the total BA in the MFS (WFS) (Figure 2b,d).

Observed Changes
We compared annual BA (Figure 3a (Table S1). During the same period of change, in Araguaia, the reduction reached 47.4% for SN, but only 0.5% for BA. In contrast, median values of BA and SN during IFMP were 270% (660%) higher than ZFP for Xerente (Araguaia) (Table S1). In the WFS, the medians for annual BA are all at the same level during the three phases for both regions. However, the boxplots show very different distributions of patterns. We observed that there is a greater variability for BA, as well as larger differences in extreme values in all phases. In those extreme cases, the reduction from ZFP to IFMP was 15.1% in Xerente and 23.5% in Araguaia (Table S1). Regarding SN, the median values decreased from NFP to ZFP and increased from ZFP to IFMP, particularly for Araguaia. Again, looking for differences between the spreads of the groups, we observed that some phases are more variable than others.
Overall, these results demonstrated that (i) for the MFS, median annual values for BA and SN increased from ZFP to IFMP; and (ii) for the WFS, extreme annual values for BA decreased. Figure 4 displays BA (for large scars), FRP, and PM 2.5 probability density functions, and identifies the 90th percentile BAp90, FRPp90, and PMp90 for the three phases during the WFS. Most fires display lower BA and intensity in all periods, and therefore, median values provide very limited discrimination among the three periods. However, those distributions have heavy right tails; therefore, high percentiles, such as the 90th percentile, improve the discrimination between each period. The exception being PM 2.5 distribution for Xerente, which did not show heavy tails. WFS values of BAp90 were much higher during ZFP compared with NFP and IFMP ( Figure 4). In Xerente (Araguaia), BAp90 increased 26.3% (43.5%) from NFP to ZFP and decreased 17.9% (20.2%) from ZFP to IFMP (Table S2). Regarding fire intensity, in Xerente and for the WFS, both NFP and IFMP displayed lower FRPp90 values than ZFP, when fires are much more intense (with a heavier right tail than that of IFMP). The increase in fire intensity during the zero-fire period is also notable in Araguaia during the WFS. In Xerente (Araguaia), FRPp90 increased from 96 MW (157 MW) in NFP to 106 MW (167 MW) in ZFP and showed a reduction of up to 105 MW (150 MW) in IFMP. This represents a slight reduction in extreme FRP values due to the IFM program compared with the zero-fire period in Xerente. However, in Araguaia, this reduction reached 10.3% (Table S2). During the WFS, we observed an increase of 7.7% and 22% in median values of PM 2.5 between ZFP and IFMP in Xerente and Araguaia. However, extreme PM 2.5 emissions values increased from NFP to ZFP (16.3% and 10.1%) and decreased to IFMP (10% and 18.2%) in Xerente and Araguaia.
In IFMP, not only did the WFS burn at a lower intensity than during ZFP, but extreme BA and emissions tended to be lower. Relative changes between phases are visually confirmed based on 1000 samples generated by the bootstrapping approach. Based on each sample (for instance, FRP values for IFMP), we create new 1000 samples using a bootstrapping approach, then for each one of the 1000 samples, we compute percentiles 10th to 90th, in steps of 10, in a way to calculate confidence intervals for bootstrapped samples [86]. Results of this analysis are provided in Supplementary Material. These boxplots show values of BA ( Figure S3), FRP ( Figure S4), and PM 2.5 ( Figure S5) percentiles 10th to 90th for Araguaia and Xerente and for the entire year, MFS and WFS for each group of 1000 samples generated by bootstrapping approach. Values for each phase are sequentially depicted for each percentile from 10th to 90th. For the WFS, these analyses reveal that there are no overlaps between the medians of each phase. In general, there is an increase in BA, FRP, and PM 2.5 from NFP to ZFP in the higher percentiles, followed by a conspicuous decrease from ZFP to IFMP. This behavior indicates that changes between these phases are statistically significant.
Fires during each phase have different seasonal cycles for all scar sizes ( Figure 6). Table S4 presents the modal months of high fire activity, highlighting the first and second highest peak in the case of bimodal patterns. Fire scar distribution varies within each phase ( Figure 6, Table S4), especially in Araguaia. In Xerente, with fewer scars, the different patterns between phases were not clear. However, we can see a flare on fire distribution on very small scars on IFMP compared with ZFP (Figure 6a-c, Table S4a). We observed a slight decrease on very small and small scars with a bimodal pattern, where the higher values are concentrated in July and September (Figure 6a). Medium and large scar numbers have similar values in the three phases, with higher values in mid-August and September (Figure 6e-g, Table S4a). In Araguaia, very small and small scar numbers have the lowest frequency in IFMP, showing a bimodal pattern, where the highest values are concentrated in July and September for very small scars (Figure 6b, Table S4b) and in July for small scars (Figure 6d, Table S4b), contrasting with the peak in September during the ZFP. The bimodal pattern is also observed for medium and large fires (Figure 6f,h, Table S4b). Medium scars have a displacement from mid-August (NFP) to mid-September (ZFP), shifting to July (highest first peak) during IFMP. Large scars showed an accentuated shift during IFMP, with the highest peaks in July and September, contrasting with the August peak during the ZFP.  Panels (a,b) show the spatial variation between the three phases. Panels (c,d) depict alluvial diagrams repreScheme 0 (0.6 fire occurrence per year, orange) and high frequency (>0.6 fire occurrence per year, red). Accordingly, fire seasonality changed between phases, specifically in Araguaia, where the scars are more distributed throughout the year, allowing better resources to prevent and combat wildfires.  (VS, a,b), small (S, c,d), medium (M, e,f) and large (L, g,h) for the three phases. The grey column represents WFS.

Evaluating PBs
The main focus of our research is to evaluate whether PBs contribute to reaching the main goals of the IFM program [6,44]. Our findings reveal that PB in both IT contributed to (i) the average rise of BA and SN from management season and the decrease from extremes wildfires season; (ii) lower-intensity and smaller-extension extreme wildfires; (iii) lower emissions during extreme events; (iv) the increase of low-fire-recurrence regions; and (v) changes in fire seasonality, decreasing extreme wildfires in the late-dry-season. By revealing changes in fire patterns across different fire policy phases, our results reinforce the PB concept, within the IFM program, as an essential strategy to properly reduce emissions without compromising local biodiversity and control fuel loads in the WFS to prevent large and severe wildfire occurrences [4,42,82]. The observed changes in the fire regime due to the implementation of PB, namely, reduced and a more fragmented BA with more scars in the WFS, are in agreement with other studies that analyzed PB programs worldwide [13,18,[87][88][89]. In both ITs, the total BA decreased (increased) in the WFS (MFS) with IFM, supporting PB as a potent land management tool to reduce the large extension and number of wildfires or unplanned fires. The reduction in fire intensity and PM 2.5 extreme values in the first four years of the IFM program and the decrease in high-recurrence areas are promising and crucial for fire-prone ecosystem maintenance. Low-intensity fires have low impacts on fire-resistant vegetation and help create mosaics in the landscape with different fire histories, fundamental to conserve the biodiversity in tropical savannas, and the ecosystem services linked to fire events [42]. Similarly, other studies have shown CO 2 emissions reductions in the western United States (18-25%) [90] and Australia (38%) [8] linked to PB actions.
The distribution of scar size class showed a distinct dynamic within both study regions and seasons. In the MFS, Xerente shows a predominance of medium scars that account for most of the total BA, whereas in Araguaia, few large fires correspond to most of the BA in this season. In the WFS, few large fires represent the greatest total BA for both ITs, as observed by other areas over the Cerrado biome [50,80,81,91] and are similar to those of tropical savannas in Africa and Australia [31,32,34,92,93]. Those patterns are in agreement with a recent study [50], which showed that few scars (7.1%) were responsible for more than 71.8% of the total BA in the last 20 years in Bananal ecoregion, where Araguaia and Xerente are located.
Previous studies have discussed the long-term fire regime patterns over Cerrado [50,82], however, few considered the role of the fire management policies in these patterns. Most current information was obtained from in situ studies that cover relatively small portions and use simple statistics or categorical measures of fire. Despite differences in study regions and datasets/treatment characteristics, our results about the constraints of key fire characteristics (i.e., intensity and area burned) under WFS and the shifts of fire season were consistent with those previous works, which found that PB successfully alters fire behavior relative to untreated areas in Cerrado. For instance, Mistry et al. [94] explored the traditional use of fire as a management tool in the Krahô IT and verified that fires for various reasons throughout the dry season produced a mosaic of burned and unburned patches in the landscape. Schmidt et al. [6,44] showed that the fire season changed over protected areas with a reduction of around 50% in late dry fires over the program's first three years. Falleiro et al. [46] verified a 17% of BA reduction in comparison over 15 ITs that have IFM plans. Alvarado et al. [95] compared protected areas over Brazil and South Africa and observed that active fire suppression did not alter the total BA but modified the seasonal fire distribution. Batista et al. [41] observed the fire regime in Canastra National Park, and found the wildfires in the late-dry-season more severe and extensive, supporting the necessity of PB to manage fire-resistant Cerrado ecosystems. Eloy et al. [43] analyzed the two first years of the IFM program in the Jalapão region, recognizing a BA reduction of (4-21%) in the late-dry-season, revealing that burning mosaics can be beneficial to wildfire prevention. Finally, Melo et al. [45] investigated nine ITs in northeastern Cerrado regarding fire recurrence patterns, fire-dependent characteristics of Cerrado, and the tendency to extend the BA due to climate impacts.
Increases in the effectiveness of wildfire suppression are expected with PB [96], and investments in IFM are effective for reducing BA [40]. However, Araguaia and Xerente territories are embedded in anthropogenic land covers (as shown in Figure 1), including mainly pastures and agriculture fields, where burning is a routine procedure [43]. This fact highlights the importance of surrounding conservation plans for those protected areas. An important aspect of fire management in the Brazilian savanna is that IFM is still only implemented in federal PA, including ITs, and wildfires during the WFS are still very common in most remnant vegetation in the Cerrado [50]. In private areas, fire policies are badly implemented, with bureaucratic burning permits that are usually not granted to (or even demanded by) farmers who commonly used illegal, and therefore untanned fires to manage their lands, resulting in frequent wildfires, which can even reach PAs and ITs [17]. This fact highlights the importance of a national fire management policy implementation that allows for the effective use of prescribed and controlled fires over different land ownership contexts. Such a national policy has been under discussion in the Brazilian parliament since 2018 (Law project 11.276/2018), but it has not advanced in the past years due to an environmentally unfriendly governmental agenda.
Finally, we observed through FWI that a decrease in precipitation combined with high temperatures and low humidity over the two regions contributes to a positive trend of fire risk, especially in the last years of the temporal series. These climate variations lead to a decrease in soil moisture, affecting vegetation flammability and creating ideal conditions for widespread biomass burning [18,92]. Worldwide, climate change has been increasing the likelihood of aggravated fire seasons and higher BA, and Cerrado is no exception [97]. Until the end of this century, an increasing potential fire weather risk to the Cerrado now seems almost inevitable due to global warming [49]. Based on our results, in the presence of extensive prescribed fire management, dry and hot trends do not inevitably increase wildfire activity. Accordingly, we argue that the impact of wildfires can still be potentially reduced if effective fire management policies, including PB, continue to be implemented over the coming decades.

Limitations of Current Datasets
Although the present study addressed only two indigenous territories in Cerrado, the approach could be adapted and applied in a variety of regions, as it is the first to consider other factors of the fire regime besides BA or active fire counts, such as frequency, fire intensity, seasonality, and emissions, to evaluate PB effectiveness.
It is noteworthy that we cannot distinguish the fire origin with limited PB information and BA product size mapping (>25 ha) due to sensor limitation [27]. PB fires are low-intensity and usually small to keep them under control, making satellite detection difficult [98], even by medium-spatial-resolution satellite imagery such as Landsat [99]. In addition, the frequent presence of clouds in the tropical region during the active burn season obscures the satellites' view [100], contributing to large uncertainty in estimating the BA of those small fires [101,102]. This is especially true for MCD64 BA detections, where accuracy is limited in regions with small and fragmented scars [103][104][105][106][107]. In the Brazilian savannas, omission errors were also reported, in particular, for very small scars (less than 100 ha) [82]. Indeed, scar size is the largest driver of accuracy across all land covers over Cerrado [91].
There are several challenges in the current PB database in Cerrado, which probably impact some of our analyses. First, PBs were recently implemented in Brazil; thus, a longterm time series is not yet available, in contrast to other countries [89]. Second, there is an absence of systematic information regarding the extent of PB and the perimeters in both ITs, enabling a proper distinction between management and wildfires within satellite-derived databases. Hence, an official, accurate, and comprehensive PB database, including fire extension, perimeter, and duration, could be important to help understand and better plan fire management regimes. For that to be possible, such an endeavor requires high and continuous investments in equipment, human resources, and technology. Accordingly, sufficient funding should be assigned to IFM, looking forward to generating more effective, systematic, and long-lasting results on controlling wildfire's destructive impact.
Finally, the analysis presented here is founded on the historical variability of fire patterns and the evolution of PB activities and climate trends, not taking into account other possible fire drivers. As stated early [7,108], it is unrealistic to expect that PB alone will reduce the occurrence of wildfires, and for this reason, further analysis should be undertaken in terms of other drivers of fire change. The effectiveness of PB in the long run is a challenging task, and the consideration of the dynamic influence of social and economic change is required [7,108]. For instance, fire dynamics in both regions are affected by several other factors rather than by management, such as land ownership conflicts, land-use changes, fires from surrounding unprotected areas, and the weather conditions at the time of the fire occurrence [44]. Moreover, variations in fire suppression efforts and in the availability and quality of equipment used may also introduce additional influence on the observed results obtained here. In addition, anthropogenic pressures from surrounding unprotected areas are known to be responsible for fire occurrence in protected areas in Cerrado, with human negligence, high road availability, and expansion of agricultural and pasture fields [109]. However, this confounding factor of fires from surrounding unprotected areas was not considered here. Therefore, although challenging, the interactive effects of such variables and prescribed burning should be taken into account in future studies. Nevertheless, we are confident that our results provide novel and useful insights to fill current gaps in the knowledge about PB effectiveness in Cerrado.

Final Remarks
Previous studies have shown that fire management practices, including PB, should be employed to further protect fire-prone ecosystems by creating burning mosaics-a strategy already used in West Africa, Australia, and by the US forest services-to limit the risk of large fire events [50,[110][111][112]. Brazil was one of the last countries with tropical savannas to adopt an effective fire management policy. This resistance to change from a zero-fire policy to the presently implemented IFM program can be explained by a historical colonial legacy [6,16,17,55], but also by considerable uncertainty related to the effects of different fire regimes on very diverse ecosystems [4,42]. Since local communities are deeply involved in fire management activities and in a landscape where fire ignitions are predominantly human, the exchange of knowledge among environmental managers, researchers, and these traditional communities is essential to develop strategies to better conserve tropical savannas [108,113]. Associated with other ecological variables, our results may also provide crucial insights about regional ecological impacts and the effects of differing fire-management practices within each region, thus supporting the development of locally adjusted fire-management practices.
Based on our results, we believe that the IFM program can potentially reduce wildfires and disseminate successful PB applications to other areas, creating a more reliable, unique, and continuous national planned fires and wildfires database to serve other fire management plans. Future work should devote more attention to the socioeconomic, biodiversity, and emissions implications of PB and, to this end, expand to encompass long-term experiments and monitored management programs [114]. Thus, we expect that our results could be useful to inform and support a better allocation of financial and human resources in managing fires as well as in decision-making strategies under current and future conditions. Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/fire4030056/s1, Figure S1: Prescribed burning (PB) over the Xerente IT between 2015 and 2018. Basemap: Google Satellite; Figure S2: Same as Figure S1, but for Araguaia; Figure S3: Boxplot showing values of burned area (BA) percentiles 10th to 90th for Araguaia (left columns) and Xerente (right columns) and for the entire years (upper row), MFS (middle row) and WFS (bottom row) for each group of 1000 samples generated by bootstrapping approach. Values for Phases 1, 2 and 3 are sequentially depicted for each percentile from 10th to 90th. The medians are represented by black circles; Figure S4: Same as in Figure S3, but for FRP; Figure S5: Same as in Figure S3, but for PM2.5. Table S1: Relative difference (%) of the percentile 50th and 90th (in parenthesis) from annual accumulated burned area (BA) and total scar number (SN)(according to Figure 3) between NFP and ZFP, and between ZFP and IFM in Xerente (a) and Araguaia (b) for each fire season: management fire season (MFS) and wildfires season (WFS); Table S2: Relative difference (%) of the percentile 50th and 90th (in parenthesis) from burned area (BA) of large scars, fire radiative power (FRP) and pyrogenic emission (PM2.5) distribution (Figure 4) between NFP and 2 and between ZFP and 3 in Xerente (a) and Araguaia (b) for wildfires season (WFS); Table S3: The proportion of pixels (i.e., area) in each fire recurrence class: zero, low (0.6 fire occurrence per year) and high (>0.6 fire occurrence per year) in Xerente (a) and Araguaia (b), according to Figure 5; Table S4: Modal month in each class for all three phases on Xerente (a) and Araguaia (b), according to Figure 6.