Fire effects in the northern Chihuahuan Desert derived from Landsat-5 Thematic Mapper spectral indices

Abstract Fire effects on desert ecosystems may be long-lasting based on ecological impact of fire in these environments which potentially is detected from multispectral sensors. To assess this, we analyzed changes in spectral characteristics from 1986 to 2010 of pixels associated with the location of fires that occurred between 1986 and 1999 in Big Bend National Park, USA, located in the northern Chihuahuan Desert. Using Landsat-5 Thematic Mapper (TM) data, we derived spectral indices including the simple ratio (SR), normalized difference vegetation index (NDVI), soil-adjusted vegetation index (SAVI), and normalized burn ratio (NBR) from 1989, 1999, and 2010 from the TM data and compared changes in spectral index values for sites with and without observed fire. We found that the NDVI and SAVI had significantly different values over the time for burned sites of different fire sizes. When differences of the spectral indices were calculated from each time period, time since fire was correlated with the SR and NBR indices. These results showed that large fires potentially had a persistent and long-term change in vegetation cover and soil characteristics which were detected by the extraordinary long-data collection period of the Landsat-5 TM sensor.


Introduction
The historical importance of fire in deserts has been questioned because the return intervals and plant community effects of past fire are relatively unknown. 1 Fire return intervals for grass dominated desert ecosystems have been estimated at 7 to 10 years to sustain dominance. 2 Where shrubs dominate, fire return intervals of 20 years may be required to allow time for shrubs to reach reproductive maturity. 3 The ecological role of fire as a process in deserts is to promote species diversity by reducing dominance of woody shrubs through changes in site water and nutrient availability. 4 Without fire, woody plants establish to the exclusion of other species. In sagebrush-steppe ecosystems found in the western US, pine and juniper species expand into unburned grasslands replacing fine herbaceous fuels with woody fuels that foster intense wildfires and recruitment of invasive grass species. 5 Loss of fire in prairie ecosystems with resultant woody plant encroachment may increase herbivory of herbaceous vegetation favoring woody plant dominance. 6 In desert communities, fire effects may last decades with few studies having collected data for the necessary timeframe. 7 The need to better understand fire effects on desert ecosystems is important for management decisions, especially regarding the use of prescribed burning and wildfire firefighting investment.
Remote sensing in semiarid regions may be utilized to monitor broad-scale plant community change over multiple years where differences in vegetation seasonality and soil condition can be minimized. [8][9][10][11] Multispectral remote sensing detection of fire impact of vegetation is based on the reduction of foliage display 12 and soil charring and ash persistence surrounding sparsely vegetated areas. 13 Over longer periods, differences in vegetation response to fire and regrowth potentially may be detected as latent fire effects as a function of fire size and intensity. [14][15][16] Spectral vegetation indices such as the simple ratio (SR), normalized difference vegetation index (NDVI), and the soil-adjusted vegetation index (SAVI) correlate well with vegetation changes such as leaf area index and biomass following fire. 11,[17][18][19][20] The normalized burn ratio (NBR) was also derived to detect fire effects in vegetation based on relative changes in near-infrared reflectance (NIR). 21 In addition, the NBR also includes the short-wave infrared (SWIR) region (2.08 to 2.35 μm) to account for soil reflectance changes associated with charring and mineral oxidation. 22 The calculated change in these index values over the time is also related to fire effects, such as severity, 23 which can be an ecological indicator of postfire vegetation recovery. However, because vegetation composition and cover may respond to fires related to site productivity potential and fire severity, these indices have varying response to fires in landscapes depending on time between observations. Satellite-derived spectral indices may be used to monitor fire-disturbed sites if the postdisturbance recovery rate of the ecosystem is greater than the time of repeat observation by the satellite sensor. In addition, sensor properties such as spatial and spectral resolutions influence detection of specific vegetation type and canopy changes. 24 Because desert vegetation has generally slower growth rates and with high mortality from fire, the influence of past fires in these ecosystems may be detected from changes in surface reflectance that are apparent over long time scales. Given the recent loss of data from the Landsat-5 TM sensor, the purpose of this study is to show that historical satellite-derived spectral indices are important for detecting wildfire effects in desert vegetation due to the slow recovery rate of vegetation in this environment. We also assess different spectral indices to show that some are more sensitive to the size of past fires and time-since-fire in a desert ecosystem potentially as a function of how vegetation and soil were affected by fire.

Methods
For this study, we assessed satellite-spectral responses associated with past fire events observed for Big Bend National Park (BBNP), USA (29.3°N, 103.2°W) area. This park is composed of vegetation types representative of northern Chihuahuan Desert with a total area of approximately 2724.5 km 2 . The park's vegetation includes lowland desert shrub (e.g., creosote, Larrea tridentata), grassland (e.g., grama-grass, Bouteloua sp.), yucca-ocotillo (e.g., yucca, Yucca sp., ocotillo, Fouquieria splendens), riparian woodland (e.g., salt cedar, Tamarix ramosissima), and montane-forested areas (e.g., pine, Pinus sp.). Riparian vegetation is found along the Rio Grande River, which forms the southern border of the park and the international boundary between the U.S. and Mexico. Montane forests grow in the Chisos Mountains located in the center of the park which forms a unique "sky island" ecosystem. Elevation within the park ranges from 550 m on the Rio Grande to 2388 m in the Chisos Mountains. Average annual precipitation ranges from 24 to 35 cm, depending on elevation, with most precipitation occurring as rainfall between May and October also referred to as the summer monsoon.
A database of fires in BBNP was obtained from the BBNP Resource Management Division that contained information on the size, location, and date of 138 fires known to have occurred from 1986 to 1999 (Fig. 1). This figure shows that most of the fires burned in the northern portion of the park and vicinity. However, this geographic distribution may be misleading as detection of fires in the southern portion of the park is limited by visibility and accessibility to the south. For our analysis, these locations, referred to as fire sites, were grouped into the following three size classes: ≤0.1 ha (Class 1; n ¼ 67), >0.1 ha to ≤5 ha (Class 2; n ¼ 40), and >5 ha (Class 3; n ¼ 31).
Based on a paired-plot experimental design, we also randomly located 138 addition sites, referred to as no-fire sites by first creating a 2500 ha circular polygon centered around the geographic location of each known fire location. The area of the polygon was based on the maximum area of any individual fire that occurred during the analysis period. Geographic locations of the points were then randomly generated outside of these polygons and within the BBNP boundary to potentially exclude pixels influenced by known fire locations. To meet the assumptions of our paired-plot design, we also assigned a random fire size class value to each of the nofire sites so that the number per fire size classes designations were also equivalent between fire and no-fire sites.
We attempted to reduce confounding effects of different vegetation types on our analysis of spectral reflectance changes over the time associated with the fire and no-fire sites based on the following method. First, we determined the potential vegetation type of the fire sites by extracting the thematic value from a supervised classification produced from 1986 TM data for a related study. 25 From this, we found that 83%, 9%, 1%, 3%, and 4% of the fire sites were located in desert shrub, grassland, yucca-ocotillo, riparian, and montane forest vegetation types, respectively. We then attempted to match the percentage of sites with the same potential vegetation type for the no-fire sites. This was accomplished by also extracting the thematic value of 1986 classified data for the initial location of the no-fire site. Once vegetation type was established for all 138 no-fire sites, we compared the percentage of pixels representing vegetation types of the no-fire and fire sites. If the numbers were not equal, such that one vegetation type category was over-or under-represented, we adjusted geographic location of the site from an over-represented type to the nearest pixel of an under-represented vegetation type outside of the 2500 ha polygons. This process was repeated until the percentage of 1986-based vegetation types associated with the fire and no-fire sites were fairly equivalent. Most of the pixels (>70%) were located in sites classified as desert shrub vegetation (Table 1). However, the distance to known fire locations and scarcity of some vegetation types, e.g., ocotillo, resulted in minor differences in the distribution of vegetation by fire type and size.

Satellite Data Processing
We obtained Landsat-5 Thematic Mapper (TM) spectral data for 1986, 1999, and 2010 encompassing BBNP with dates chosen based on the availability of fire records for our study site between 1986 and 1999. The 2010 data provided additional observation of fire-affected sites to assess Class 1 fires were defined as fires with an estimated size ≤0.1 ha, Class 2 fires defined with size >0.1 to ≤5 ha, and Class 3 fires defined with size >5 ha. Also shown are locations without known fire occurrence between 1986 and 1999 used as controls for our analysis. longevity of fire effects and potential postfire recovery of the ecosystem. These TM data acquired for late March were constrained by data availability with minimal cloud coverage (<10%).
Data were horizontally georegistered to the Universal Transverse Mercator (UTM) coordinate system in the northern section of zone 13 by first identifying road intersections within BBNP in the satellite with reference coordinates derived from U.S. Geological Survey (USGS) 1:100,000 DLG transportation data. The transportation data are in the North American Datum for 1983 datum which uses the Geodetic Reference System 1980 spheroid model. Due to significant topographic relief within the park due to their Chisos Mountains, digital elevation model data at a 1:24,000 scale from the USGS were used to provide elevations of selected reference locations. Using 30 points to derive a projective transformation matrix, we found a cumulative root-mean-square error value of <10 m for all data. From this rectification process, we assumed that pixel geographic locations analyzed over the time were accurate for all years of satellite observation used in this study.
All satellite data for each reflective band were converted from DN values to radiance values using gain and offset values for bands 1-5 and 7 based on the equation where L i is the at-sensor radiance (mW · cm −2 · μm −1 · sr −1 , G is the gain, DN is the digital number, and O is the offset. To account for path radiance from atmospheric sources that likely differed between data collection years, we used a dark target subtraction method. With the absence of dark, clear lakes in this desert region, we instead found mean L i values for the south wall of the Santa Helena canyon near the western edge of the park. This canyon is oriented west to east and has near vertical walls which extend 300 m from the river's surface. During March, the southern walls are highly shaded with minimal direct or indirect reflectance. Mean L i values were derived from a 20 × 20 pixels area visually identified within the canyon's shadow form the data for each year. The resultant L i values were then subtracted from all image pixels for each band for each year to approximate atmospheric scattering effects. These modified L i values were then converted to apparent reflectance using the equation 26 where ρ λ is the reflectance by wavelength region, L λ is the measured spectral radiance, β s is the solar elevation angle, and E sλ is the solar top of the atmosphere irradiance by wavelength. Gain and offset values were obtained from individual Landsat TM header files, while E sλ values were obtained from literature sources. 26 The solar elevation angles for each set of data were provided as part of the satellite metadata.

Spectral Index Analysis
We calculated four vegetative indices (SR, NDVI, SAVI, and NBR) from the derived ρ λ data for all pixels. The SR was calculated using the ratio of ρ NIR and ρ red (Ref. 27), where the spectral band extent of ρ NIR and ρ red is 0.76 to 0.90 μm and 0.63 to 0.69 μm for TM, respectively. The NDVI (a logarithmic stretch of SR) was calculated as ðρ NIR − ρ red Þ∕ðρ NIR þ ρ red Þ. 28 The SAVI was calculated using the equation ½ðρ NIR − ρ red Þ∕ðρ NIR þ ρ red þ LÞð1 þ LÞ, where L ¼ 0.05 for arid environments. 29 The NBR was calculated from ðρ NIR − ρ SWIR Þ∕ðρ NIR þ ρ SWIR Þ where ρ SWIR is the SWIR channel (TM band 7) with a spectral sampling extent of 2.08 to 2.35 μm. 21 To meet the normality assumptions of the parametric statistical tests, the index values were transformed with a log e ðx þ 1Þ function. Because spectral indices such as the NBR are sensitive to prefire vegetation conditions, 30 we assessed potential differences in the 1986 fire and no-fire pixels transformed index values with a one-way analysis of variance (ANOVA). Next, potential changes in the transformed spectral index values over the time were assessed with a repeated measure ANOVA in which the presence of fire and size were assessed independently and as cofactors over the time.
To assess potential correlation between time since fire and spectral index values, we first calculated the difference between nontransformed spectral index values for each fire pixel for the each pair of satellite data acquisition dates (e.g., 1986 to 1999 dSR ¼ SR 86 − SR 99 ). The derived differenced values were grouped by size class and found to be normally distributed based on the Kolmorgorov-Smirnoff test. Each set of differenced spectral index values were assessed by fire size class using a one-way ANOVA with significant differences among fire size classes identified posthoc using a Tukey's test. We then calculated the parametric Pearson's correlation coefficient (r xy ) for each differenced spectral index value and the number of days following fire from March 22, 2010, for each fire size class. This date was chosen because it was the day following satellite acquisition of the 2010 TM data.

Results
Analysis of the spectral index values showed that the mean, untransformed SR, NDVI, and SAVI values were lowest for 1999 and highest in 2010 for all fire sites (Fig. 2). For no-fire sites, the highest mean SR, NDVI, and SAVI values were found for 1999 with lower values for 1986 and 2010. For the NBR index, both mean fire and no-fire site values decreased in value by year of analysis. The ANOVA results of the transformed spectral index values for the 1986 fire and no-fire were not significantly different indicating that all sites had relatively similar reflectance prior to the fires the occurred from 1986 to 1999.
From the repeated measures ANOVA, we found both NDVI and SAVI had within the index value differences over the time with both fire occurrence and size as significant cofactors ( Table 2). Both NDVI and SAVI were significant when both the fire occurrence and size were included regardless of changes over the time. The NDVI values were also significantly different for fire occurrence and size as individual variables. The SAVI values were significantly different over the time and when compared to size. Finally, the NBR values were significantly different over the time.
For the analysis of the differenced spectral indices, we found that the dSAVI and dNBR had significantly different values for the 1999 to 2010 (Table 3). Based on Tukey's posthoc test, Classes 1 and 2, as a group, were different from Class 3 site index values. Only the dNBR was significant for the 1986 to 2010 period which we found differed among Classes 1 and 2, as a group, and Class 3 site values. For our analysis of time since fire, we found significant negative correlation for Class 2 fires for dSR and dNBR 1986 to 1999 and 1986 to 2010 periods (Table 4). No significant correlation was found for the 1999 to 2010 period for any of the differenced spectral index values; however, the r xy values were highest during this period for Class 3 fires.

Discussion
Our analysis showed that NDVI and SAVI were related most to changes in reflectance associated with vegetation cover that may be ecologically linked to fire occurrence and size. Because of this finding, we expect other vegetation-sensitive indices such as the enhanced vegetation index 31 to be similarly responsive. The lower SR, NDVI, and SAVI values we found for the 1999 fire site compared confirms the assumed effect that fire would result in vegetation removal at each point (Fig. 2). Because the no-fire site SR, NDVI, and SAVI values were higher in 1999 than 1986, this supports that fire was the main cause for reflectance changes during this time instead of other environmental factors. These overall results are somewhat surprising given that the average time between the fire occurrences and the date of the Landsat TM for 1999 was 7.3 years. The rate of vegetation recovery detected by vegetation spectral indices is related to inherent ecosystem growth rates and fire severity. 22 For 2010, we attributed higher values of SR, NDVI, and SAVI for the fire sites to fire effects on soil resources and plant community composition. Mineral nitrogen, for example, may be released from burned vegetation where fire temperatures <160°C, which if present, would potentially increase vegetation production. Fire also affects certain plant species differently based on survivorship of reproductive tissue. In our study area, black grama (Bouteloua eriopoda) is locally extirpated by burning as most of its vegetative reproductive structures are aboveground. 32 For our study area and period of analysis, any grasses negatively impacted by fire were likely replaced by woody shrub species. 25 This change in plant type may have resulted in a higher ρ NIR as green foliage is mostly present on shrubs during spring months. In contrast, the lower SR, NDVI, and SAVI values we found for no-fire sites in 2010 compared to 1999 are potentially due to vegetation community maturity with concomitant resource limitation and productivity decline.
Significant differences found for NDVI and SAVI values indicated that the combined changes in index values over the time, the fire occurrence, and size of fire all influenced spectral reflectance ( not significantly different between treatments, this shows that fire was a dominant mechanism that influenced the spectral response over a 25-year period. The sensitivity of NDVI to variation in vegetation density is well documented. [33][34][35] The lower values for 1999 for both NDVI and SAVI showed that the main effect of fire was to reduce the vegetation cover for an extended fire period of time. Because both NDVI and SAVI showed similar results, this also indicated that loss of vegetation cover was an important indicator of fire occurrence. Changes in NDVI and SAVI associated with time and fire were also due to differences in fire size based on the MANOVA analysis. Fire size may be related to difference in NDVI and SAVI as compared to no-fire sites as a function of varying fire effects on vegetation and soil. The smaller fires were likely to be less severe as they were characterized in the database as spot fires with only single plants or small patches of vegetation burned. These types of fires may have lower broadscale plant mortality and potentially less effect on ρ red and ρ NIR . 22 For Class 3 fires, size of fire generally is associated with prefire moist climate conditions that favor higher fuel loadings to sustain a fire spread over large areas. 36 In addition, these larger fires tend to have broad plant mortality with lower productivity than areas without fires or with smaller fires. When the spectral index values were expressed as differences between time periods, a slightly different result was found. Despite the occurrence of fires occurring between 1986 and 1999, none of the difference index values for fire pixels showed significant differences among fire size classes for this period (Table 3). Differences were found for the 1999 to 2010 period, however only for dSAVI and dNBR which differentiated the small (1 and 2) and large size (3) class fires. Because both SAVI and NBR were significant, we propose that soil reflectance was most influenced by fire size during 1999 to 2010. For SAVI, the dSAVI values for 1986 to 1999 showed small positive changes in the mean SAVI values for Classes 1 and 2 pixels compared to Class 3 pixels. This is contrasted by larger negative dSAVI values for Classes 1 and 2 pixels compared to the Class 3 pixels for 1999 to 2010. We interpret this as a direct effect of fire on vegetation recovery, with smaller fires having less plant mortality, and higher postfire regrowth. For the dNBR, the positive and significant Table 1 Percentage of sites based on the 1986 vegetation classification for fire and no-fire sites and by fire size indicated by classes. Class 1 fires were defined as fires with an estimate size ≤0.1 ha, Class 2 fires defined with size >0.1 to ≤5 ha, and Class 3 fires defined with size >5 ha.  result for the Class 3 fires from 1999 to 2010 analysis indicates lower ρ NIR and higher ρ SWIR than the previous period. This is interpreted as continued low rate of vegetation recovery coupled with potentially long-lasting change in soil surface reflectance. The slight negative dNBR for the Classes 1 and 2 fire pixels for the 1999 to 2010 period was likely related to vegetation recovery  as any soil changes were likely minimal and easily removed by erosion for these smaller fires. The significant dNBR results for the 1986 to 2010 by fire size class supports other evidence that fire can have lasting ecological impacts in this desert region. [37][38][39] Our results indicated that smaller fires had less different NBR values between 1986 and 2010 compared to larger fires. In other ecosystems, changes in spectral index values may be directly related to burned area and show a significant time component ranging from months to years associated with growth rates of postfire vegetation. 17,30,40 Our study also demonstrated significant pixel-level detection of fire effects based on the dNBR, however, mapping actual fire perimeters may require biannual satellite data observations to discriminate burned areas from naturally exposed soil within this type of landscape. 15 Time since fire was found to correlate with dSR and dNBR values, however, only for fire size Class 2 ( Table 4). The uniqueness of the Class 2 fires as a significant result can be explained by considering the potential variance in fire effect and recovery for this intermediate fire size class. For Class 1 fires, vegetation mortality due to fire is likely to be minimal with a higher recovery rate. For the Class 3 fires, fire effects were likely uniform with persistent low vegetation recovery 36 and soil impacts. For the Class 2 fires, vegetation mortality was likely spatially heterogeneous resulting in time since fire as an important variable associated with dispersion and regeneration of different desert plants. Most of the fires observed in the park for this study were from the late 1980s with second peak in 1999 (Fig. 3). These pulses of fire occurrence influenced the variety of pixels with differing vegetation cover and soil conditions indicated by the SR and NBR, respectively. However, the lack of significance found for 1999 to 2010 indicates that the postfire recovery period was not influenced by the timing of fires during the 1986 to 1999 period. This suggests that, although fires having long-lasting effect for large-intense fires, intermediate size fires may reach an earlier recovery "climax" a decade following a smaller fire. This study highlights the importance of long-term remote sensing observation of land surfaces to better monitor potential ecological changes associated with disturbance. The Landsat-5 TM sensor is one of the longest continuous space platforms from a single instrument in the history of remote sensing. Comparable long-term observations from coarse resolution sensors such as the Advanced Very High Resolution Radiometer require extensive intersensor calibration. 41,42 Though some changes in sensor calibration have occurred, 43 the 28 years of historical Landsat-5 TM data will continue to be important for research in natural resource science particularly where ecological process rates are slow or disturbance effects are long-lasting.

Conclusion
Changes in spectral reflectance detected by satellite sensors associated with disturbance in ecosystems are related to intrinsic rates of ecological change constrained by climate, plant species diversity, and disturbance intensity. We have shown that typical spectral indices such as NDVI, SAVI, and NBR to detect the occurrence of fires in deserts, even after nearly of decade following a burn event. Because fire effects may vary with fire size, differences in spectral indices that include ρ NIR and higher ρ SWIR may be more sensitive to differences in vegetation loss and recovery or soil color alteration, respectively. Satellite information may be especially useful for studying fires in deserts because of the large areas observed as density of burned areas in deserts may be sparse. In addition, some satellite sensor missions that extend over multiple decades are important for adding new insight into desert ecosystem processes because their operation is at an ecologically relevant time scale.