Burned-Area Detection in Amazonian Environments Using Standardized Time Series Per Pixel in MODIS Data

Abstract

Fires associated with the expansion of cattle ranching and agriculture have become a problem in the Amazon biome, causing severe environmental damages. Remote sensing techniques have been widely used in fire monitoring on the extensive Amazon forest, but accurate automated fire detection needs improvements. The popular Moderate Resolution Imaging Spectroradiometer (MODIS) MCD64 product still has high omission errors in the region. This research aimed to evaluate MODIS time series spectral indices for mapping burned areas in the municipality of Novo Progresso (State of Pará) and to determine their accuracy in the different types of land use/land cover during the period 2000–2014. The burned area mapping from 8-day composite products, compared the following data: near-infrared (NIR) band; spectral indices (Burnt Area Index (BAIM), Global Environmental Monitoring Index (GEMI), Mid Infrared Burn Index (MIRBI), Normalized Burn Ratio (NBR), variation of Normalized Burn Ratio (NBR2), and Normalized Difference Vegetation Index (NDVI)); and the seasonal difference of spectral indices. Moreover, we compared the time series normalization methods per pixel (zero-mean normalization and Z-score) and the seasonal difference between consecutive years. Threshold-value determination for the fire occurrences was obtained from the comparison of MODIS series with visual image classification of Landsat Thematic Mapper (TM), Enhanced Thematic Mapper Plus (ETM+), and Operational Land Imager (OLI) data using the overall accuracy. The best result considered the following factors: NIR band and zero-mean normalization, obtaining the overall accuracy of 98.99%, commission errors of 32.41%, and omission errors of 31.64%. The proposed method presented better results in burned area detection in the natural fields (Campinarana) with an overall accuracy value of 99.25%, commission errors of 9.71%, and omission errors of 27.60%, as well as pasture, with overall accuracy value of 99.19%, commission errors of 27.60%, and omission errors of 34.76%. Forest areas had a lower accuracy, with an overall accuracy of 98.62%, commission errors of 23.40%, and omission errors of 49.62%. The best performance of the burned area detection in the pastures is relevant because the deforested areas are responsible for more than 70% of fire events. The results of the proposed method were better than the burned area products (MCD45, MCD64, and FIRE-CCI), but still presented limitations in the identification of burn events in the savanna formations and secondary vegetation.

1. Introduction

Human disturbance is the primary cause of burned events in the Amazon Forest, mainly in areas of deforestation, agriculture, and pastures [1,2,3,4]. Anthropogenic fire from farming and livestock usually escapes into neighboring forest areas, resulting in large-scale fire events [5]. Therefore, land-use determines the burning patterns of the Amazon forest, where areas with intense deforestation, fragmentation, and the presence of highways have a higher number of fires, whilst Conservation Units and Indigenous Reserves are essential barriers [6,7,8,9,10]. The Amazon forest has a high vulnerability to fire incidences due to the following factors: most species do not tolerate recurrent fire events, floristic and structural changes, and organic matter incineration that is essential for the maintenance of ecological processes [11,12,13].

The concentration of Amazon forest fires is along an area called the arc of deforestation, which presents an intensification of anthropogenic actions along the eastern and southern forest edges. Compared with the central zones of the forests, the arc of deforestation has lower biomass and drier climate resulting in larger burnings [14,15]. During severe drought periods, there is an increase in tree mortality on forest edges and the accumulation of dry leaves [16,17], which favor fire events mainly in degraded forests [18,19,20]. Models of climate change in the Amazon rainforest provide for the expansion of fires due to more frequent droughts and land use intensification [4,21].

In the monitoring of the Amazon Forest (5.5 million km²), the use of remote sensing data is necessary due to its synoptic view, rapidity, and high cost-effectiveness. The major impediment to optical remote sensing in the Amazon region is the high atmospheric interference from cloud cover and aerosols [22]. Therefore, high temporal resolution optical sensors are essential for the monitoring of Amazon forest fires, increasing the possibility of acquiring good quality images without atmospheric interference [23,24,25]. Among the high temporal resolution sensors, the Moderate Resolution Imaging Spectroradiometer (MODIS) sensor has been prominent in the mapping of fires in the different terrestrial biomes [26,27,28], offering global scale products of active fires [29,30] and burned areas [26,29]. Active fire detection captures the energy emission using thermal infrared imagery (3.6–12 µm range), and provides information only of burning pixels during the passage of the satellite, not allowing the burned area quantification [26,29,31]. In contrast, the post-fire mapping evidenced by ash, coal, and vegetation changes enable obtaining of the perimeter, area, and an estimate of damage caused by forest fires. Thus burned area products (MODIS-MCD45/MCD64) were extensively used in diversity studies about ecosystem structure change and biomass burning emissions [24,32,33,34]. Padilla et al. [35] performed a validation of the MCD45 and MCD64 products for several terrestrial biomes, obtaining estimated commission errors rates of 46% and 42%, and omission error rates of 72% and 68% for MCD45 and MCD64, respectively. Therefore, the burned area product tends to underestimate the extent of the burned area, where the omission errors were higher than commission errors [35,36,37,38]. Cardozo et al. [39] showed that the MODIS burned area product (MCD45) in the Amazon region at different dates, presents values higher than 90% of omission errors, which corresponds to low accuracy and a need for the development of new algorithms adjusted for this region.

Different remote sensing methods have been proposed for the burned area mapping. The most widely used processing using bi-temporal images combines two algebraic operations: (a) Spectral index calculation, where the most commonly used are the Normalized Burn Ratio (NBR) [40] and Normalized Difference Vegetation Index (NDVI) [41]; and (b) Seasonal differences between the pre- and post-burned indices, such as dNDVI [42] and dNBR [40]. Other indices have been proposed such as the Burned Area Index Modified (BAIM) [43], Global Environmental Monitoring Index (GEMI) [44], Mid-Infrared Burned Index (MIRBI) [45], and the variation of Normalized Burn Ratio (NBR2) [46].

The adaptation of the index differentiation method to long-term MODIS time series [47,48] considers the entire time series of one year with its antecedent, generating an image set rather than just a single difference image from a custom selection. However, the automatic difference over the time series shows a significant increase in signal-to-noise ratio, with respect to the bi-temporal images that were pre-defined by the analyst. Therefore, the application of the differencing method in continuous time series requires methodological improvements [47,48]. Lhermitte et al. [49] tried to solve this problem by selecting control pixels surrounding the unburned area based on time series similarity. This approach considers the spatial context (neighboring pixels) to minimize the external influences and phenological variations [50]. However, this method has severe limitations for complete automation, requiring supervision for not including neighboring pixels with different characteristics as compared to the analyzed pixel.

A proposed alternative is the use of normalized time series using the Z-score (mean zero and standard deviation equal to 1) or means (mean zero), which allow highlighting of the burned areas without changing the signal-to-noise ratio [47]. The normalization procedure is very promising because it does not require neighboring control pixels, considering only the time series of the pixel. This method is different from that proposed by the research of Diáz-Delgado et al. [50] and Lhermitte et al. [49], because it is fast and simple data processing restricted to the pixel, which generates an image that emphasizes the burning points in different environments.

The present paper aims to evaluate and compare different methods of burned area mapping from MODIS time series (BAIM, GEMI, MIRBI, NBR, NBR2, NDVI, and the near-infrared band) and normalized time series per pixel (Mean and Z-Score) in different land covers in the municipality of Novo Progresso, belonging to the Amazon region. This research also evaluates the spatial and temporal distribution of fire in different land use types.

2. Study Area

The study area is the municipality of Novo Progresso and surroundings (74,552.1518 km²), State of Pará, Southeast Brazilian Amazon (Figure 1). The climate is Humid Equatorial with annual precipitation ranging from 2050 mm to 2650 mm and average temperature above 25 °C [51,52]. The rainy season is from November to May and the dry season, from June to October, when most fire events occur [53].

The predominant vegetation is classified as Open Ombrophilous Sub-Montane Forest, but it also contains areas of Forest Savanna and Submontane Semideciduous Seasonal Forest [55]. This region has been inhabited by indigenous populations, currently restricted to the Baú and Menkragnoti Indigenous Reserves [60]. The construction of the federal highway BR-163 in the 1980s brought an intense migratory flow that increased the soybean cultivation, livestock, and infrastructure projects [61]. The Conservation Units curb deforestation growth within Integral Protection Units (Rio Novo National Park, Cristalino State Park, Nascentes da Serra do Cachimbo Biological Reserve, and Crystal Private Reserve Natural Heritage) and Sustainable Use Units (Jamanxim and Altamira National Forests, Iriri State Forest and Tapajós Environmental Protection Area) [56]. The Agrarian Reform Settlements and Military Areas complement the federal land types of the region [58,59]. The leading causes of fire are deforestation and pasture maintenance [2,9].

3. Materials and Methods

3.1. Methodology Flowchart

The methodology adopted the following steps, including MODIS 8-day composites data acquisition (MOD09A1 and MOD09Q1), calculation of spectral indices, filtering, time series normalization, burned area mapping, and validation (Figure 2). The research calculated seven spectral indices and applied in each index two types of time series normalization per pixel (Z-score and zero-mean) and the seasonal difference method, computing a total of 21 different procedures. Finally, we compared the burned areas by the different MODIS image processing with the following data: Visual interpretation of Landsat images (Thematic Mapper (TM)/ Enhanced Thematic Mapper Plus (ETM+)/Operational Land Imager (OLI)) and burned area products (MCD45 [31], MCD64 [26], and Fire CCI) [62,63].

3.2. MODIS Data

The MODIS sensor aboard the Terra and Aqua satellites was developed by the National Aeronautics and Space Administration (NASA) to provide continuous observations of the Earth’s surface [64]. This sensor has 36 spectral bands ranging from 0.4 μm to 14.4 μm at three spatial resolutions (250, 500 m, and 1 km) according to the band. We used the following MODIS products: (a) 8-day composite surface reflectance at 250 and 500-m resolution (MOD09A1 and MOD09Q1); (b) daily burned area at 500-m resolution (MCD45A1) [31]; (c) daily burned area at 500-m resolution (MCD64A1) [26], and (d) monthly burned area product at 300-m resolution derived from both MODIS and Medium Resolution Imaging Spectrometer (MERIS) sensors (FIRE CCI) [62,63]. The data download sites were the NASA Land Data Products and Services website (https://search.earthdata.nasa.gov/) and from the European Space Agency (ESA) Fire CCI website (http://www.esa-fire-cci.org/).

MODIS surface reflectance product is composed by bands 1–2 (0.6 µm–0.9 µm) at 250-m spatial resolution and bands 3–7 (0.4 µm–2.1 µm) at 500-m spatial resolution [65]. These images are available as daily data and 8-day composite product that select the best observation within the time interval and eliminate the presence of artifacts such as cloud or cloud shadow. The traditional temporal composite based on the maximization of values can cause loss of information of burned areas, which are characterized by low albedo [66]. In contrast, daily images have more noise and atmospheric effects that can also lead to the insertion of invalid information.

The acquisition of MODIS data covered the period 2000–2014, resulting in 683 8-day composite images. We discarded the poor-quality pixels with (a) cloud coverage, (b) related to the high and medium aerosol classes, and (c) relative to frontal scattering and with zenith angles greater than 60°, because of the influence of the spectral mixing of the targets with their shadows, especially in the near-infrared (NIR) band [66,67]. We converted the MODIS images from native sinusoidal projection to the Universal Transverse Mercator (UTM, WGS 84 ellipsoid) projection. Spectral bands with the spatial resolution of 500 m were resampled at 250 m using the nearest neighbor interpolation method to make the data compatible.

3.3. Spectral Indices

The data compared to the fire detection were: (a) Only the NIR band (band 2); (b) NDVI, NBR, NBR2, BAIM, GEMI, and MIRBI index; and (c) temporal differencing (Δ) (

Table 1. Spectral indices to highlight the burned areas adopted in this research.

Index	Equation	Bibliographic Reference
NBR/ΔNBR	$NBR=\frac{NIR-SWIR2}{NIR+SWIR2}$	[40,79]
NDVI/ΔNDVI	$NDVI=\frac{NIR-RED}{NIR+RED}$	[41]
NBR2/ΔNBR2	$NBR2=\frac{SWIR1-SWIR2}{SWIR1+SWIR2}$	[46]
BAIM/ΔBAIM	$BAIM=\frac{1}{{\left(0.05-NIR\right)}^2+\left(0.2-SWIR2\right)2}$	[43]
GEMI/ΔGEMI	$GEMI=n\left(1-0.25 n\right)-\frac{RED-0.125}{1-RED}$	[44]
MIRBI/ΔMIRBI	$MIRBI=\left(10\ast SWIR2\right)-\left(9.8\ast SWIR1\right)+2$	[45]
Δ	$prefire-postfire$

Notes: In GEMI $n=\frac{2(NI{R}^2-RE{D}^2)+{1.5}_{NIR}+{0.5}_{RED}}{NIR+RED+0.5}$.

). The NIR band allows a high separation of the burned areas compared to other spectral channels [68,69]. The NDVI uses the RED and NIR bands, and is widely used to highlight vegetation changes [70]. The NBR index developed to highlight burned areas adopts the normalized difference between the NIR and Short Wave Infrared (SWIR) [40]. The NBR2 and MIRBI indices are variations of the NBR and NDVI indices, considering the SWIR1 and SWIR channels to highlight the burn areas [45,46]. The BAIM index maximizes the spectral deviation between the burned area and other land covers from the NIR and SWIR channels [43]. The GEMI index was determined to evaluate the vegetation conditions, reducing the content of both the NIR and RED channels [45].

Seasonal differencing is a change detection method with extensive use in the analysis of the severity of a fire [47,48,71,72,73]. This method consists of subtracting the values of spectral indices or bands acquired before and after a fire event (pre-fire minus post-fire) to emphasize the changes related to the fire. The advantage of the method is the possibility of defining the fire severity, apart from indicating the areas with vegetation regeneration [46]. Despite the proven utility of the technique in many different environments, there are few details about its application in the Amazon region or other areas of tropical forests. Seasonal differencing techniques have been used with different spectral indices including: ΔNBR [48,74], ΔNDVI [75], ΔNBR2 [76], ΔBAIM [74], ΔGEMI [77], and ΔMIRBI [78].

3.4. Time Series Standardization per Pixel

The burned areas show differences in spectral behavior according to the vegetation type and the period after the fire [80]. In the Amazon, burning occurs mostly in pasture and agriculture, but there are significant impacts on forest and savanna areas [81,82]. Therefore, the detection of areas burned by remote sensing should be adaptable for environmental differences. The normalized time series allows to equate the behavior of distinct vegetation and highlight the burned areas [83]. The normalization procedures consider the temporal data per pixel. The two evaluated methods of time series normalization were: Zero-mean (mean equal to 0) and Z-score (mean equal to 0 and standard deviation equal to 1) expressed by the following formulations:

$$Zero-mean=x_t-\mu$$

(1)

$$Z-score=\frac{x_t-\mu }{\sigma }$$

(2)

where $x$ is the pixel value, $t$ is the time, μ is the mean of the remote sensing data over time ($t$) for a pixel, and σ is its standard deviation. Therefore, the two normalization methods present negative values when the remote sensing values are smaller than the temporal mean and positive values when they are above. The Z-score approach was successful for burned area detection in regions of savanna vegetation [83]. Normalization makes the burned-area features of different environments compatible and eliminates confusion with low-albedo targets (e.g., water). As an example, we compared the MODIS NIR temporal series of water bodies with burned areas in natural fields (Campinarana), grasslands, and wooded savannas, showing an overlap of values and the inability to distinguish different targets by a single threshold value (Figure 3a). The time series normalized by the mean allows separating the burned area and water bodies by a single threshold (Figure 3b). Moreover, the Ombrophilous forests have reflectance values in the NIR range that are higher than the seasonal forests and forested savannas, which remain higher even after a fire event and prevent the use of a single threshold for the detection of burning among these types of vegetation (Figure 4a). The method of time series normalization equalizes the fire points. Thus, forest time series normalized by the mean decrease false fire events significantly, but the omission of some points occurs mainly in the seasonal forest and forested savanna (Figure 4b).

In addition to the relevant results in the time dimension, normalization per pixel also causes a significant digital enhancement of the burned areas in the spatial dimension. In images without normalization, the savanna and watercourses have similar values to the burned areas, being confused in the application of a threshold value (Figure 5a). The normalized image allows an apparent distinction between the burned areas compared to water bodies, since they have significantly lower values (Figure 5b).

3.5. Seasonal Differencing

The seasonal difference is one of the most commonly used methods for the detection of burned areas, but studies focus on bi-temporal analysis based on preselected images. Therefore, most studies adopt a specific fire event, where selection by the analyst of two images near the birthday dates considers the low presence of clouds or noises. Few studies detect areas burned by seasonal differences in long-term time series of satellite images, using tens or hundreds of images. The seasonal difference algorithms in long time series assume fixed time intervals and present a certain amount of random noise, since it is impossible to get multiple cloudless and noise-free images on birthday dates in different years [47,48,49,83]. Thus, some proposals have been made to mitigate interference, considering spatial and temporal dimensions. A spatial approach involves selecting control pixels for pre-fire data by examining the contextual neighborhood around the focal pixel and averaging two or more similar pixels within spatial windows [84,85]. A temporal approach involves the use of the average of the birthday data of the complete time series rather than the simple data from the previous year [47,83].

In this research, the seasonal differencing considered two strategies: (a) Automatic method with a fixed annual interval, where a pre-fire image was subtracted from the post-fire with a range of one year exactly; and (b) Method to select a good quality image closer to the anniversary date when the previous year’s image had poor quality or cloud coverage. The latter option may have non-fixed subtraction intervals, performing the temporal search for a better-quality pixel within a time window.

3.6. Landsat Reference Data and Burned Area Mapping

The threshold definition for separation between burned and unburned areas used the method developed by Carvalho Júnior et al. [83]. This procedure compares a reference map of the burned areas, previously obtained from an image with better spatial resolution, and the burned areas acquired by a sequence of thresholds from a spectral band or spectral index. Thus, the method generates an overall accuracy curve from the comparison of the reference map and the classifications using threshold values (Figure 6). The fire detection in other time series images used the threshold value with the best overall accuracy. In this research, the definition of the best threshold used seven reference maps of the burned area produced by the visual interpretation of the Landsat TM and ETM+ images on the following dates: 15/08/2001, 19/09/2002, 28/07/2003, 30/07/2004, 03/09/2005, 10/08/2008, and 01/09/2010 (

Table 2. Low-quality pixel percentage of Landsat scenes used to define the threshold value.

Image Date	Sensor	Percentage of Low-Quality Pixels
		Normalization Method	Seasonal Difference (Fixed Interval)	Seasonal Difference (Closest High-Quality Image)
15/08/2001	ETM+	34%	-	-
19/09/2002	ETM+	25%	77%	19%
28/07/2003	TM	13%	82%	24%
30/07/2004	TM	44%	51%	13%
03/09/2005	TM	26%	51%	12%
10/08/2008	TM	47%	57%	45%
01/09/2010	TM	7%	49%	22%

, Figure 7). The acquisition of Landsat data came from the United States Geological Survey (USGS) website (https://earthexplorer.usgs.gov/). The visual image interpretation methodology [35,86] used at least two different dates to identify burned areas and exclude false positives. The fires in the region occur mainly during the dry season, where there is more availability of cloudless images and Landsat/MODIS pairs with matching dates, avoiding the possible omission of fires due to different dates.

The tests covered an area corresponding to 20% of the study area (14,689.638 km²), containing the highest rate of deforestation and burning in the analyzed image set. The images included three classes: burned area, unburned area, and low-quality pixels (clouds and noise) from the quality mask available with Landsat images on the EarthExplorer site (Table 2). The seasonal difference with the fixed interval presupposes the image selection with exactly one previous year, being able to have a high percentage of pixels with atmospheric interference. The seasonal difference with the closest high-quality image selects data with the lowest cloud coverage within a user-stipulated period.

We applied the average of the seven best thresholds in the complete time series of 14 years, generating a time series of burned area mask. The best threshold definition uses the equation $x\le y$, where x is the value of the burned pixel and y is the threshold value for band 2, GEMI, NBR, NBR2, NDVI, ΔBAIM, and ΔMIRBI, while the equation $x\ge y$ for BAIM, MIRBI, ΔBand2, ΔGEMI, ΔNBR, ΔNBR2, and ΔNDVI.

The overlap of the same-year polygons resulted in burned areas during the year. We used an algorithm that eliminated the burned areas that did not occur on at least two consecutive dates (8-day), which were considered noises because the ash of burns in savanna and forest requires a longer time to be extinguished [69,87,88]. The errors of areas burned in water bodies were corrected from the land use/land cover mapping produced by the Amazon Surveillance System (SIVAM) [55].

3.7. Dataset Comparison and Accuracy Analysis

Accuracy analysis used visual interpretations of independent areas of those adopted in the determination of optimal threshold values (Figure 8). The methodology assumed was proposed by Padilla et al. [89], which considers the confusion matrices between the classification of MODIS and Landsat images (

Table 3. Confusion matrix.

		True Condition (Landsat Classification)
		Burned Area	Unburned Area
Predicted Condition (MODIS Classification)	Burned area	True Positive (TP)	False Positive (FP)
	Unburned area	False Negative (FN)	True Negative (TN)

) and its indices: Overall Accuracy (OA) (Equation (3)), Commission Errors (CE) (Equation (4)), and Omission Errors (OE) (Equation (5)) [90].

$$OA=\frac{\sum TP+\sum TN}{Total Population}\ast 100$$

(3)

$$CE=\frac{\sum FP}{\sum FP+\sum TP}\ast 100$$

(4)

$$OE=\frac{\sum FN}{\sum FN+\sum TP}\ast 100$$

(5)

We also used the Dice coefficient (DC), described in Reference [91], which combines CE and OE into a single category according to the below Equation, as in Reference [38]:

$$DC=\frac{2\ast a}{2\ast a+b+c}$$

(6)

where “a” refers to the burned area mapped correctly in both methods, “b” the commission errors, and “c” the omission errors. Given DC is a measure of aggregate precision, it obscures differences in the individual measurements of CE and OE. Therefore, CB has limitations if CE and OE are not equally significant, in which case the accuracy analysis should be directed to the largest specific error [38].

The reference image corresponded to about 25% of the study area (19,486.3664 km²) obtained from the visual interpretation of TM, ETM+, and OLI images on the following dates: 10/10/2001, 10/05/2002, 09/22/2003, 10/02/2004, 09/19/2005, 09/22/2006, 09/09/2007, 09/27/2008, 09/30/2009, 09/09/2010, 09/12/2011, 09/30/2012, 10/03/2013, and 09/28/2014. In 2012, we used Landsat 7 ETM+ images (even with the presence of errors) because it was the only sensor in operation of the Landsat program [92], the scan-line issues were masked and disregarded in that year’s accuracy analysis.

The non-parametric McNemar’s test [93] assessed the independence of the results obtained by pair-wise comparisons of the classifications [94]. The McNemar test calculates the chi-square distribution (x²) considering the matrix elements of

Table 4. McNemar test between two classifications.

	Classification 2
	Correct	Incorrect	Total
Classification 1
Correct	f11	f12	f11 + f12
Incorrect	f21	f22	f21 + f22
Total	f11 + f21	f12 + f22	f11 + f12 + f21 + f22

in the following Equation:

$$x^2=\frac{{\left(f_{12}-f_{21}\right)}^2}{f_{12}+f_{21}}$$

(7)

The McNemar test analysis considered 4000 random samples for each year between 2001 and 2014. The significance assessment compared the value of x² obtained with the tabulated value, indicating the null hypothesis when the x² value was less than 5% of significance (3.841), or rejection of the null hypothesis when the x² value was higher.

The validation of the MCD45 [36], MCD64 [26], and Fire CCI [62,63] products used the same data set to facilitate comparison with the proposed methodologies. The MCD45 uses an algorithm that identifies burned areas daily, eliminating false events such as cloud shadows. The MCD64 identifies persistent changes in the time series of vegetation indices and includes the localization of active fire to define boundary values for the burned area and unburned area. The Fire CCI product combines active fire pixel information and spectral changes of the MERIS sensor time series. These products are available in monthly data with 500 m of resolution, where each pixel can indicate the exact day of the burn on Julian days (1–365), unburned area, cloud area, or water area [31,62].

3.8. Analysis of the Spatial Relationships between Land Use/Land Cover Classes and Burned Area

Land use dynamics is a crucial factor to understand the spatial pattern of fires. In this spatial analysis, we used the land use/land cover mapping produced by the Amazon Surveillance System (SIVAM) [55] and deforestation provided by the National Institute for Space Research by Project for Monitoring Deforestation in the Legal Amazon (PRODES) [54]. In the SIVAM mapping, land use/land cover classes considered those established by the Brazilian Institute of Geography and Statistics [55]: Ombrophylous Forests, Seasonal Forests, Campinarana, Savanna, Secondary Vegetation, Annual Deforestation, and Water. As the PRODES methodology considers one year beginning in August, we used Landsat images to update deforestation by the end of each year from 2001 to 2014, resulting in annual land cover maps. The intersection of information allowed us to evaluate the interactions between deforestation, land use, and fire. Moreover, we assessed the influence of the BR 163 highway on the occurrence of fires, calculating the percentage of burned area in a buffer zone with distance intervals of 20 km from the highway.

4. Results

4.1. Determination of the Best Threshold Value

The use of both the seasonal differencing and the zero-mean normalization improves the overall accuracy of the original images (

Table 5. Average of best thresholds and overall accuracy (OA) for different tested procedures in the years of 2001, 2002, 2004, 2005, 2008, and 2010.

Original Data
	BAIM	Band 2	GEMI	MIRBI	NBR	NBR2	NDVI
Threshold	17.48	0.206	0.560	1.312	0.148	0.230	0.480
OA	97.1%	97.0%	97.0%	97.0%	96.3%	96.7%	96.0%
Seasonal Difference (Δ = Pre-fire − Post-fire)
	ΔBAIM	ΔBand 2	ΔGEMI	ΔMIRBI	ΔNBR	ΔNBR2	ΔNDVI
Threshold	−31.2	0.04	0.11	−0.25	0.22	0.12	0.22
OA	98.7%	98.2%	98.3%	98%	98.2%	98.4%	98.3%
Seasonal Difference (Δ = Pre-fire − Post-fire) with selected images by quality next to birthday
	ΔBAIM	ΔBand 2	ΔGEMI	ΔMIRBI	ΔNBR	ΔNBR2	ΔNDVI
Threshold	−28.32	0.041	0.156	−0.406	0.278	0.138	0.172
OA	97.2%	97.1%	97.0%	97.0%	96.7%	96.5%	96.5%
Z-score normalization
	BAIM	Band 2	GEMI	MIRBI	NBR	NBR2	NDVI
Threshold	1.235	−1.744	−1.825	1.430	−1.925	−1.428	−1.918
OA	97.1%	97.4%	97.2%	97.0%	96.5%	96.5%	96.0%
Zero-mean normalization
	BAIM	Band 2	GEMI	MIRBI	NBR	NBR2	NDVI
Threshold	11	−0.109	−0.200	0.250	−0.284	−0.160	−0.203
OA	97.8%	97.9%	97.7%	97.3%	96.9%	96.9%	96.8%

). The zero-mean normalization obtained better precision results (overall accuracy) than the Z-score normalization (increasing false positives) and seasonal differencing (lower overall accuracy).

The BAIM presented the best accuracy index from the original data and the difference index (Table 5). The normalized NIR band by Z-score or by the mean presented the best fit for the burn area detection, with higher overall accuracy values than all the spectral indices evaluated. The BAIM, NIR band, GEMI, and MIRBI normalized by the mean had similar overall accuracy, while the NBR and its variation NBR2 presented the fifth-best performance, despite being a widely-used index for the mapping of burned areas. The NDVI showed the worst value for overall accuracy.

The annual difference already implies the loss of one year of the time series [83]. Moreover, the use of fixed intervals without the image selection by quality causes the reduction of the validation area due to the presence of clouds in the differenced images (Table 2). Although the seasonal differences highlight the burned areas, there is confusion between newly cleared areas with or without burning and omission of pixels burned in consecutive years.

The lower accuracy of the seasonal difference method compared to the normalization method was due to the following factors: (a) Pixel quality, recurrent fires, and continuous deforestation. Pixel quality is one of the main limitations of seasonal difference in the tropical region due to intense cloud cover. The selection of high-quality pixels in two consecutive years is unlikely in most of the Amazon region, requiring algorithms for the selection of images close to the anniversary or temporal composition methods. In contrast, the normalization method does not depend on pixel quality in consecutive years, since it uses the entire time series available, facilitating the delimitation of fires in regions with high cloud presence.

The seasonal difference in areas with high fire recurrence tends to omit burned areas in consecutive years or regenerating regions due to the decrease in the difference values (Figure 9). The continuous deforestation of the tropical regions generates areas burned in new deforestation with different values of the fires in old deforestation, preventing the use of a single threshold value. Differently, the normalization method allows the grouping of burned areas in new and old deforestation from only one threshold value.

Figure 10 compares the application zero-mean normalization and the seasonal difference in the BAIM image. The RGB composition of the 2005 and 2006 dates (R-2005/G-2006/B-2005) highlights the changes between the two methods. The RGB combination of zero-mean images showed deforested areas of 2006 in magenta color, burned deforestation areas in white, with burned areas in old deforestation in green (Figure 10a). In this case, a single threshold identifies most of the fires in the white and green areas. The seasonal difference method has inverse values in relation to the zero-mean normalization method, where deforested areas in 2005 are in greens and “burned areas” in 2006 are in magenta (Figure 10b). However, burned deforestation areas also have green tones, mixing with deforested areas and preventing their separation by a single threshold.

4.2. Validation and Data Comparison

The MODIS mapping validation indicated that the NIR band normalized by the mean achieved the best result (Dice coefficient of 0.647) (

Table 6. Accuracy indices of the burned areas using MODIS times series in the study area. “DC” is the Dice coefficient and “OA” is the overall accuracy.

Original Data
	BAIM	Band 2	GEMI	MIRBI	NBR	NBR2	NDVI
DC	0.588	0.527	0.337	0.326	0.422	0.430	0.399
OA	98.56%	98.71%	98.36%	98.31%	98.30%	97.71%	96.81%
Commission Errors	48.02%	34.21%	20.75%	27.44%	34.40%	61.99%	46.06%
Omission Errors	38.20%	53.75%	79.96%	76.95%	62.60%	45.49%	38.02%
Seasonal Difference (Pre-fire − Post-fire) with fixed interval
	ΔBAIM	ΔBand 2	ΔGEMI	ΔMIRBI	ΔNBR	ΔNBR2	ΔNDVI
DC	0.506	0.485	0.384	0.334	0.428	0.412	0.400
OA	98.15%	97.87%	97.11%	96.59%	97.51%	97.50%	97.15%
Commission Errors	29.82%	71.67%	40.71%	42. 37%	51.22%	54.21%	69.47%
Omission Errors	62.60%	42.85%	73. 57%	74.02%	61.23%	62.01%	59.41%
Seasonal Difference (Pre-fire − Post-fire) with selected images
	ΔBAIM	ΔBand 2	ΔGEMI	ΔMIRBI	ΔNBR	ΔNBR2	ΔNDVI
DC	0.592	0.606	0.481	0.357	0.537	0.487	0.425
OA	98.69%	98.74%	98.57%	98.49%	98.62%	98.52%	97.23%
Commission Errors	14.86%	32.30%	17.98%	25.56%	23.14%	35.46%	40.22%
Omission Errors	48.43%	38.71%	64.98%	76.45%	57.50%	58.50%	46.72%
Z-score normalization
	BAIM	Band 2	GEMI	MIRBI	NBR	NBR2	NDVI
DC	0.591	0.618	0.476	0.501	0.429	0.452	0.397
OA	98.38%	98.80%	98.54%	98.44%	98.69%	97.73%	96.38%
Commission Errors	45.42%	38.97%	34.16%	35.80%	45.77%	47.23%	56.81%
Omission Errors	39.36%	38.70%	62.58%	55.90%	45.45%	45.74%	45.59%
Zero-mean normalization
	BAIM	Band 2	GEMI	MIRBI	NBR	NBR2	NDVI
DC	0.625	0.647	0.566	0.539	0.508	0.503	0.479
OA	98.36%	98.99%	98.44%	98.59%	98.66%	98.31%	97.54%
Commission Errors	38.01%	32.41%	28.29%	33.21%	35.55%	46.61%	48.60%
Omission Errors	31.64%	34.75%	55.20%	49.60%	42.69%	42.48%	30.96%

). The zero-mean normalization obtained higher accuracy than the other methods in most of the tests performed. Z-Score normalization presented a significant improvement in the detection of burned areas in Band 2, GEMI, and MIRBI. The seasonal differencing with image selection by quality next to birthday improved the classification of all indexes, exhibiting better performance than the Z-Score normalization in the BAIM, GEMI, NBR, NBR2, and NDVI indices, and higher than the zero-mean normalization for the NBR index. The seasonal differencing using image selection showed the lowest commission error rate, and with the zero-mean normalization or Z-score had the smallest error of omission.

The McNemar test confirmed the statistical independence between most of the classification methods, presenting values of x² > 3.841 (Figure 11). The most accurate method (Band 2 with zero-mean normalization) was statistically different from the other classifications. The procedures with the highest degree of agreement between them were those derived from the GEMI, MIRBI, and NBR2 indices, while the methods based on Band 2, BAIM, and Seasonal Difference were independent.

Considering the same reference data (Landsat-TM, ETM +, and OLI Images), the burned area products (MCD45, MCD65, and FIRE CCI) presented less than 25% of commission errors, like the seasonal difference method, but due to the high error rate of omission, the Dice coefficient was lower than all the methods used. The MCD45 and FIRE CCI products had a Dice coefficient lower than 0.11 and omission errors higher than 90%. The MCD64 product obtained the Dice coefficient of 0.24 in detriment of the reduction of the omission errors, surpassing the products MCD45 and FIRE_CCI (

Table 7. Dice coefficient and overall accuracy of the burned areas for MCD45/MCD64/FIRE CCI data in the study area during 2001–2014 period.

Burned Area Product	Overall Accuracy	Dice Coefficient	Commission Errors	Omission Errors
MCD45	98.58%	0.11	24.34%	93.55%
MCD64	98.61%	0.24	23.83%	84.51%
FIRE_CCI	98.48%	0.09	24.89%	94.67%

). The MODIS products still need improvements [24], validations in Africa using Landsat images obtained 34% commission errors and 59% omission errors [86], and globally containing a commission error of 46% and omission error of 72%, where specific biomass showed more satisfactory results than others such as Borealis Forestry and Tropical and subtropical Savanna [38]. The imprecision of the proposed method occurred in small burn areas due to the spatial resolution difference between the sensors used for validation, showing commission errors in forest and savanna areas. The proposed method had a lower commission error rate than MODIS / MERIS products, which focused on the resolution difference with the validation image (Figure 12).

4.3. Mapping Validation by Land Cover Type

The accuracy analysis of the burned area mapping by different land cover used the result of higher accuracy (NIR band and zero-mean normalization). The overall accuracy was higher than 98% in all classes, due to the high proportion of unburned area correctly identified. The Campinarana and Deforestation classes presented the highest Dice coefficients (0.778 and 0.755, respectively), and lower commission and omission errors. Burnings in the Ombrophilous Forest obtained low commission errors, but significant errors of omission. Burnings in savanna vegetation and seasonal forests presented the highest rates of commission and omission errors, respectively (

Table 8. Estimated accuracy of burned area for each land use/cover classes. The water and secondary vegetation were not evaluated by the low level of representation of fires in the selected area.

Land Cover	Overall Accuracy	Dice Coefficient	Commission Errors	Omission Errors
Campinarana	99.25%	0.778	9.71%	27.60%
Deforestation/Pasture	99.19%	0.755	18.57%	34.76%
Ombrophilous Forest	98.68%	0.506	23.40%	49.62%
Savanna	98.62%	0.466	59.25%	37.42%
Seasonal Forest	98.57%	0.371	23.23%	76.85%

).

4.4. Land Use/Land Cover Classes and Burned Area Patterns

Fire events and deforestation showed a high correlation (

Table 9. Annual burned area (km² and percentage within each class) from MODIS daily images for the following classes of land use and land cover: Campinarana (CP), Deforestation/Pasture (D/P), Ombrophilous Forest (OF), Savanna (S), Seasonal Forest (SF), and Secondary Vegetation (SV).

Area (km2)
Year	CP	D/P	OF	S	SF	SV
2001	14.16 (0.5%)	280.66 (6.9%)	49.56 (0.1%)	2.55 (0.03%)	12.18 (0.2%)	0.28 (0.1%)
2002	12.84 (0.4%)	773.68 (14.4%)	84.08 (0.2%)	14.61 (0.18%)	12.77 (0.2%)	0.78 (0.2%)
2003	26.82 (0.9%)	497.05 (8.2%)	92.58 (0.2%)	9.06 (0.11%)	18.24 (0.3%)	1.05 (0.3%)
2004	51.60 (1.7%)	1203.22 (16.8%)	235.74 (0.5%)	15.58 (0.19%)	20.14 (0.3%)	0.12 (0.0%)
2005	33.10 (1.1%)	771.25 (9.9%)	111.60 (0.2%)	15.38 (0.19%)	23.97 (0.4%)	0.86 (0.3%)
2006	26.85 (0.9%)	781.16 (9.4%)	244.92 (0.5%)	22.15 (0.28%)	70.80 (1.1%)	1.52 (0.5%)
2007	18.24 (0.6%)	677.03 (7.6%)	143.28 (0.3%)	5.21 (0.06%)	22.24 (0.3%)	0.20 (0.1%)
2008	41.48 (1.4%)	625.33 (6.6%)	130.02 (0.3%)	28.77 (0.36%)	44.52 (0.7%)	1.37 (0.4%)
2009	15.19 (0.5%)	233.87 (2.3%)	121.21 (0.3%)	4.53 (0.06%)	52.62 (0.8%)	0.43 (0.1%)
2010	7.21 (0.2%)	749.36 (7.2%)	196.83 (0.4%)	11.10 (0.14%)	23.60 (0.4%)	1.62 (0.5%)
2011	23.66 (0.8%)	214.61 (2.0%)	90.71 (0.2%)	5.95 (0.07%)	15.54 (0.2%)	0.66 (0.2%)
2012	4.49 (0.1%)	395.61 (3.5%)	87.07 (0.2%)	2.78 (0.03%)	13.08 (0.2%)	0.40 (0.1%)
2013	10.77 (0.4%)	166.01 (1.4%)	47.75 (0.1%)	3.30 (0.04%)	6.06 (0.1%)	0.25 (0.1%)
2014	5.65 (0.2%)	748.24 (6.2%)	102.40 (0.2%)	3.18 (0.04%)	16.60 (0.3%)	2.06 (0.6%)

). Deforestation had the highest fire representativeness in the years studied, representing more than 70% of the total area burned annually. This attested to the anthropogenic influence on current burning patterns, both in deforestation and in agricultural management. In the study area, the Ombrophylous Forest had intense deforestation, decreasing from 70% in 2000 to 59% in 2014, which resulted in a high number of burnings in the years 2004, 2006, and 2010. More than 16% of deforested areas burned in 2004, while less than 2% burned in 2013. The other types of vegetation (Savanna, Seasonal Forests, and Secondary Vegetation) presented a low proportion of annual burned area, with average values smaller than 1% of the affected areas.

The highest occurrence of burned areas was concentrated near the BR-163 highway, especially in the buffer zone up to 40 km (

Table 10. Burned area percentage in buffer zones with 20 km distance intervals from the BR-163 highway.

Distance	20 km	40 km	60 km	80 km	100 km	120 km	140 km	160 km
Burned Area (%)	12.5%	10.1%	3.7%	1.6%	0.9%	2.7%	2.5%	0.7%

). More than 85% of the forest area showed no fires. In the 14 years examined, the burned area was slightly more than 11% of the study area (Figure 13). Most of the area burned only once or twice during the study period (7% and 3%). Approximately 1% of the pixels burned in more than three years in the time series, concerning the agricultural areas. The indigenous territories and the military base showed a low incidence of fires in the temporal series, in contrast, the Jamanxim National Forest (highly degraded conservation unit) exhibited a high frequency of fires. The private or unlicensed properties (without the direct intervention of public authorities) totaled 24,358.8 km² of the study area, and corresponded to the highest burning rates, followed by the settlement projects (Figure 13).

5. Discussion

Consistent with other research, seasonal differences and standardization techniques improved the burned area mapping [47,83,95]. However, seasonal differencing considering continuous data in time has the following limitations [83]: (a) The subtraction operation highlights the noise present in the two images (pre- and post-fire); (b) Fire detection is susceptible to interference of phenological variations, oscillation of cropping cycles, and land use changes. Therefore, the main difficulty of the seasonal difference is to ensure that other changes independent of the fire event are not confused with burned areas, such as deforestation or atmospheric variations. The seasonal difference in continuous data over time intensifies this limitation, as opposed to the use of discrete data chosen by the analyst.

Another alternative is the use of time series standardized by pixel, which has been used to identify burned areas in the savanna region of Central Brazil [83]. This procedure converts the data set to Z-Score values (with the mean of zero and standard deviation of one), facilitating the determination of a threshold value of the burned areas from a reference image. This technique does not alter the relative values of the temporal signature and avoids the errors described by the seasonal difference. However, the application of this method to forest areas generated false positives mainly related to the presence of cloud shadows not detected by the MODIS sensor quality mask due to the resolution difference between the products. The MODIS cloud mask has a 1-km resolution, while the images of this research had 250- and 500-m resolutions. The time series of spectral indexes in forest areas shows low standard deviation in comparison to the other targets, always having a photosynthetically active behavior with little change over time. Therefore, the Z-score normalization that has the standard deviation as the denominator causes an overestimation of the normalized forest values, differently from what occurs in the vegetation with high seasonal variation such as the savannah. An alternative proposed in this research was the use of zero-mean normalization, with only the subtraction by the mean, without the division by the standard deviation. This normalization approach allowed an adaptation of the method [83] for these vegetation types, and more accurate results in the detection of burned areas. The main advantage of the zero-mean normalization is that it equalizes the values of the burned area in the different scenarios of the Amazon region (forest, pasture, and deforestation), facilitating their detection from the single threshold.

Except for the NBR index, all the data tested (band of Near Infrared, BAIM, GEMI, MIRBI, NBR2, and NDVI) presented better accuracy using zero-mean normalization than the other procedures (original image, seasonal differencing, seasonal differencing with the selected images, and Z-score), attesting to the superiority of the zero-mean in detriment of the z-score normalization for the mapping of burning areas in forest regions. The NBR index presents a higher accuracy using the seasonal difference with the selected images.

Within the normalization procedures, the input data that had the best Dice coefficient were two: (a) BAIM in the original image and seasonal difference with fixed interval; and (b) Infrared band in the seasonal differencing with selected images, Z-score, and zero-mean normalization methods. Therefore, the best result used the infrared band and zero-mean normalization. Pleniou et al. [96] and Alonso-Canas [63] emphasized the ability to separate the burned areas from the near infrared. Although the NBR and dNBR indices are widely used in the literature for the detection of burned areas and fire severity [27,97], there are limitations in the study area due to the intense dynamics of deforestation. Roy et al. [97] warned about the ineffectiveness of the NBR index to measure the severity of burning in tropical forests.

Some limitations persist in the mapping of small areas burned due to MODIS spatial resolution, as verified by other authors [81,98]. The burned areas on small farms and agrarian reform settlements result in an underestimation of the total area burned. Moreover, the daily images may present better results than the MODIS products from the multitemporal compositing techniques [66,99]. Future research efforts could use high spatial and temporal resolution satellite data, which should reduce the omission of undetected small fires.

The comparison between the MCD45, MCD64, and FIRE CCI products and the reference burned area mapping demonstrated proved unfeasible for use in the Amazon region due to the low detection of burned areas. Libonati et al. [100] also highlighted a high percentage of omission errors in the MCD45 product. Even with MCD64 product enhancements, the omission errors may exceed 60% [35]. Despite the limitations of the burned area products in the Amazon, these were developed for the world surpassing in many cases other techniques of burned-area detection. The spatial resolution of 500 m from MODIS and 300 m from FIRE CCI products influenced the difference between the two methods.

The methodology used presents better fire detection capacity in pasture/deforestation areas and Campinarana vegetation. Deforested regions were the most affected by fires, where it is used at regular intervals for pasture maintenance, combating invasive plants, pests, and improving low soil productivity [13]. The lack of access to alternative agricultural management is also an essential factor to explain the high rates of fires in the region [101]. Locations with a single occurrence of burning in the time series are often fires in the native forest that originated externally in agricultural areas and became out of control. In this context, fragmented forests are more vulnerable to fires because of the increased border with agrarian regions [2,8].

The Novo Progresso region is located between a conservation unit and Indigenous Parks that are essential barriers to deforestation and fires [7,102]. However, the Jamanxim National Forest is among the most deforested and burnt protected areas of the Brazilian Amazon due to the proximity of the Cuiabá-Santarém highway. Currently, the Brazilian Congress is negotiating a 26% reduction proposal in the Jamanxim National Forest area to allow the construction of a railroad for the transport of grain from Mato Grosso by the port of Miritituba, Pará. This railway will undoubtedly cause impacts increasing the occurrence of deforestation and fires.

Although the most widespread methods for detecting burned areas are based on specific thresholds of spectral indices, the temporal normalization per pixel should be evaluated in the future with other methods. Instead of considering spectral indices, several studies use spectral mixing techniques to highlight burned areas. [85,103,104]. Thus, the temporal normalization techniques for the pixel can be further tested in the time series of fraction images. Another approach used is the supervised classification, either by pixel or object-based, which depends on the collection of samples of the different targets present in the scene [105,106]. Thus, the supervised classification on normalized data per pixel should be tested.

6. Conclusions

This study evaluated different procedures for the mapping of burned areas in the Amazon Forest, considering different image processing in the MODIS time series, such as seasonal difference, normalization (Z-score, zero-mean normalization, and without normalization), spectral indices and bands (BAIM, GEMI, MIRBI, NBR, NBR2, NDVI, and NIR band). The zero-mean normalization and seasonal difference showed improvements in burned area detection in comparison with the original images. Although the seasonal difference method was widely used in other biomes, the zero-mean normalization achieved better accuracy, confirmed by the McNemar statistical test. The use of the NIR band presented the best results in burned area detection. The zero-mean normalization allowed comparison of the burning behaviors of different environments in an integrated way, defining a single threshold value for the different vegetal formations. This approach to time series normalization per pixel is an innovation in the processing of remotely sensed data for detection of burned areas that must be improved for other environments and monitoring studies. In addition to the methodology presented, greater accuracy in the detection of burned areas in Campinarana and pastures, burned areas occur predominantly in deforested areas, which account for more than 70% of annually burned areas. The adopted methodology presented an overall accuracy (98.99%) greater than the burned area products (MCD45, MCD64, and FIRE CCI), and low errors of omission, despite having higher commission errors. The improvement of the monitoring of burned areas in the Amazon enables better environmental management and the estimations of their impacts.