Remote Sensing the Implications of Fire Management in the Andean Paramo: a Preliminary Assessment Using Satellite Remote Sensing

The upper ranges of the northern Andes are characterized by unique Neotropical, high altitude ecosystems known as paramos. These tundra-like grasslands are widely recognized by the scientific community for their biodiversity and their important ecosystem services for the local human population. Despite their remoteness, limited accessibility for humans and waterlogged soils, paramos are highly flammable ecosystems. They are constantly under the influence of seasonal biomass burning mostly caused by humans. Nevertheless, little is known about the spatial extent of these fires, their regime and the resulting ecological impacts. This paper presents a thorough mapping and analysis of the fires in one of the world's largest paramo, namely the " Complejo de Páramos " of Cruz Verde-Sumapaz in the Eastern mountain range of the Andes (Colombia). Landsat TM/ETM+ and MODIS imagery from 2001 to 2013 was used to map and analyze the spatial distribution of fires and their intra-and inter-annual variability. Moreover, a OPEN ACCESS 11062 logistic regression model analysis was undertaken to test the hypothesis that the dynamics of the paramo fires can be related to human pressures. The resulting map shows that the burned paramo areas account for 57,179.8 hectares, of which 50% (28,604.3 hectares) are located within the Sumapaz National Park. The findings show that the fire season mainly occurs from January to March. The accuracy assessment carried out using a confusion matrix based on 20 reference burned areas shows values of 90.1% (producer accuracy) for the mapped burned areas with a Kappa Index of Agreement (KIA) of 0.746. The results of the logistic regression model suggest a significant predictive relevance of the variables road distance (0.55 ROC (receiver operating characteristic)) and slope gradient (0.53 ROC), indicating that the higher the probability of fire occurrence, the smaller the distance to the road and the higher the probability of more gentle slopes. The paper sheds light on fires in the Colombian paramos and provides a solid basis for further investigation of the impacts on the natural ecosystem functions and biodiversity.


Introduction
Fire management practices are part of the traditional land use activities in the high mountains of the Neotropical region [1].However, scientifically little attention has been paid to the dynamics and potential impact of these fires in the paramo, a zonal grassland ecosystem [2] located primarily in the high mountains of the Neotropics [3].Paramo is a broad term employed in Latin America to describe the high altitude grasslands that cover the areas between the timberline (above 3200-3600 m a.s.l.) and the snow line (below ca.4500-5000 m a.s.l.) of the northern and the Equatorial Andes [4].
With an area of approximately 35,000 km 2 [2], these Neotropical ecosystems of Pliocene origin [5,6] are mostly located in Colombia, Ecuador, Peru, Venezuela and Costa Rica [7].Their unique landscapes and endemic flora provide shelter and habitat for a variety of mammals, birds, insects, amphibians and reptiles [4].Paramos are well-known for their extensive water storage and sophisticated water regulation capacity [8], which is a result of the combination of plants with low evaporation characteristics (tussock grasses and xerophytic herbs [9]) and the physical properties of the non-allophanic Andisols that dominate the region.The high water surplus feeds the rivers providing fresh water downstream [10].Paramos are also extraordinary reservoirs of carbon and contain high levels of organic carbon related also to their mineralogical background, climate and elevation and water saturation characteristics.These ecosystems play an important role in the global carbon balance with up to 520.9 t C ha −1 in the soil according to reports for belowground organic C found in the undisturbed paramos of the Chingaza National Protected Area in Colombia [11].For the paramos, these factors are of extremely high ecological, genetic and scientific importance [4].The maintenance of the hydrological and carbon storage functions is essential to ensuring a good water quality to the numerous human communities that depend on the paramo water supply, including the Bogota Metropolitan Area [10,11].
For a long time, the paramos remained a natural ecosystem largely unaffected by human activities [12].Recently, however, the growing population and the progressive land expansion for agricultural purposes into the neighbouring lowlands have led to increasing pressure on these fragile ecosystems, although the population density is still relatively low in these high mountain ecosystems [13].The fresh grasses of the paramo provide ideal places for grazing [14].To foster the growth of fresh grasses, tussock grasses and xerophytic shrubs are usually burned before being allocated for grazing.As a result, the combination of grazing and regular vegetation burning activities has become a common management practice to support the growth of palatable young grasses that are used to feed the livestock [15].At the same time, the most favourable areas for cropping are tilled, predominantly to grow potatoes and beans [1].
Burning, intensive grazing, tilling and replacement of the natural grassland with more nutritive grass species significantly affect the water balance of the paramos areas.Phenomena typically accompanying pasture farming and tillage, such as soil compaction and soil crusting, additionally alter the infiltration rates, water storage and regulation capacity of the paramos.This seriously compromises its water supply function.Some scientists also state the effects of human activities in terms of accelerated soil erosion given the properties in the paramos [8,16].Very few studies have been conducted to quantify the impact of these land management changes on the hydrodynamic properties of paramos [10], their floristic composition, the vegetation structure [4], the morphological evolution of the soils [8] or on carbon storage in the soils [11].
The aim of this study is to map the most recent burning activities, covering the time span from 2000 to 2013 through remote sensing data.In doing so, the spatial impact of the fire clearing processes will be assessed diachronically.The "Complejo de Páramos" of Cruz Verde-Sumapaz (Colombian Andes) was selected as a study area.It constitutes one of the largest paramo areas in the world.In addition, a logistic regression model analysis was undertaken to test the hypothesis that the dynamics of the paramo fires can be related to human pressures.Research results from this remote sensing approach provide a solid knowledge base for further field activities with the aim to observe the on-site effects of the burning and grazing practices on the soil quality, biodiversity and hydrological dynamics of paramos.

Study Area
The majority of the paramo ecosystems are located in the Colombian Andes.With an area of 266,750 ha and an altitudinal range from 3250-4230 m a.s.l., the "Complejo de Páramos" of Cruz Verde-Sumapaz in the Cordillera Oriental is one of the largest paramos known [17,18] (Figure 1).Much of its surface is located within the protected National Park area of Sumapaz covering ca.123,794 ha.In 1977, this region was declared national park [19,20] because (i) it is a high mountain area of Colombia extremely rich in flora; (ii) it is a biodiversity hotspot and (iii) it is a major source of freshwater supply for the Bogotá Savannah, the most densely-populated area of the country [4,21,22].)).The rainfall distribution during the year is mainly bimodal with all-season humid conditions.Higher precipitation occurs from March-May and from October to November, while the dry season occurs from December-February and from June to September (Figure 2).January is usually the driest month of the year.Both the temperature and the temperature range highly depend on the elevation.The average annual temperature is 7.3 °C with a daily temperature range between 3.3 and 11.4 °C (World Clim-Global Climate Database, 1950-2000) [24].

Paramo Delimitation
The "Complejo de Páramos" of Cruz Verde-Sumapaz was outlined using the classification proposed by Morales et al. (2007) [4].Accordingly, the paramo regions were subdivided into 22 ecosystems.For this study, special attention was paid to the four major units that characterize the Cruz Verde-Sumapaz Paramos: The very humid paramo is the most extensive paramo ecosystem (40%), followed by the humid paramo (30%), the dry paramo (7%) and the very humid sub-paramo (3%).Other ecosystems cover about 20%.

Landsat Imagery Acquisition and Pre-Processing
The multispectral characteristics, spatial resolution and temporal range and frequency of the Landsat sensors make the imagery highly suitable for mapping fire scars at a regional scale [25,26].The Global Visualization Viewer (http://glovis.usgs.gov)was used to preview the Landsat imagery and to select the best cloud-free scenes.Three Landsat Thematic Mapper (TM) and twenty-seven Enhanced Thematic Mapper (ETM+) satellite images were downloaded from the Earth Resources Observation and Science Center (EROS) of the United States Geological Survey (USGS) (Table 1).These images ensured the best coverage of the study area over a time period of thirteen years (from 2001 to 2013, excluding the years 2008 and 2009 where no cloud-free Landsat images were available).Images acquired during the main dry season (December-March) were preferred because of (i) lower cloudiness, (ii) more frequent occurrence of fires [27] and (iii) the lack of issues related to seasonal variations in vegetation conditions.Moreover, the images acquired during the dry season were free from atmospheric haze and, therefore, had less inter-image differences between the Sun's angle and azimuth, atmospheric transmissions and soil moisture [28].
The EROS Landsat images for the study area are provided as standard level-one terrain-corrected (L1T) [23].This image type has not been atmospherically corrected, but is rather geometrically adjusted to remove systematic geometric errors related to the positions of the sensor, the satellite and the Earth [23].Therefore, within the image pre-processing stage, the Landsat imagery was corrected for atmospheric disturbances using the ENVI 5 software.This operation included surface reflectance computation [28] and normalization based on dark object subtraction (DOS, [29]) to improve the accuracy of the vegetation change detection procedure [30].

Burned Area Mapping
An image differencing technique [31] based on image subtraction [32] was applied to assess the burned areas for the "Complejo de Páramos" of Cruz Verde-Sumapaz.This method subtracts the spatially-registered dates of two Landsat images following a pixel-by-pixel procedure.The burned area was calculated based on the delta (Δ) between the vegetation indices (VIs) calculated from a pre-fire image (VIspre-fire) and a post-fire image (VIspost-fire) (Equation ( 1)).Negative values indicated areas that had been subjected to vegetation cover changes caused by fires.
The burned areas were detected for all years by subtracting the vegetation index values of a Landsat image from those of the preceding year.However, if the pre-fire images presented spatial information gaps due to cloud, shadow or noise interference, repeated change detection operations were performed using other images to reduce these gaps.
The potential of single-band information (i.e., red, near-infrared and short-wave infrared), as well as for multiband vegetation indices (i.e., Normalized Difference Vegetation Index (NDVI, [33]), Normalized Difference Water Index (NDWI, [34]), Normalized Burn Ratio (NBR, [35])) for distinguishing burned areas in the paramo regions was tested using the fire that burned ca.2200 hectares of paramo near Granada (Meta Department) in 2003.The NBR yielded the best results in discriminating the pixel value differences between burned and unburned paramo.However, the NDVI showed marked falls in its values recorded on burned surfaces, resulting in a slightly better discrimination of some areas not severely affected by the fires.Accordingly, the image differencing procedure was carried out using both indices.The NBR combines information on the near-infrared and the mid-infrared bands (Equation ( 2)), while the NDVI makes use of the vegetation's radiometric information in the red and near-infrared bands (Equation ( 3)): Annual difference layers for each index were obtained (NBR, NDVI), exported from ENVI into a Geotiff format and, subsequently, imported into ArcMap 10.2.The Raster Calculator tool of the spatial analyst extension was used to isolate the paramos areas that were subject to change from the areas that did not experience change.A threshold value for each change map was defined after the analysis of mean and standard deviation values of the burned areas.The annual VIs difference data, i.e., the rough burned areas, were converted from a grid format into a polygon shape file.Vegetation changes that were incorrectly shown as positive, i.e., clouds, shadows and noise that survived the masking operations, were identified and filtered out through an on-screen visual interpretation of the Landsat RGB colour composition, the annual VIs change detection rasters, Google Earth and rough burned areas.
Considering the small difference observed in this study between burned areas mapped by the two vegetation indices (±4.5%), the test area outcomes (Granada) and the results of the accuracy assessment (Section 4.2), the burned area detected by NBR and NDVI was merged into a single layer (Figure 3).This operation enabled coverage of a wider region of the spectrum and detection of the burned areas better than if each index had been used separately.

Accuracy Assessment of the Burned Areas
The accuracy of the burned area detected from Landsat images rested on a confusion matrix by a per-pixel analysis [36] (geometric accuracy) and a linear correlation (thematic accuracy).A set of burned areas was mapped via onscreen visual interpretation on Google and Bing high resolution imagery (acquired via Web Map Service (WMS)).These changes were used as reference data, while the burned areas resulting from the Landsat image analysis represented the classified data.

Statistical Analysis
A logistic regression model analysis was used to test the hypothesis that the dynamics of the paramo fires can be explained by human pressures measured as the accessibility of the area [37].Researchers have applied this logistic regression model to forest fire contexts [38].It estimates the relationship between the dichotomous dependent variable (presence/absence of a fire event) and a set of independent variables.The statistic fit between the dependent variable and the independent variables is expressed by the ROC (receiver operating characteristic) indicator.The ROC value is defined as an interval from 0 to 1.The threshold value is 0.5 (i.e., slope = 1) and indicates no predictive power of the model [39].Higher predictive values indicate a predictive relevance of the model and a good model fit.
In this study, the dependent (predicted) variable is represented by the centroid point of each burned area polygon and by an equal number of random points, which represent the absence of a burned area.The distance from the road [40], the distance from urban areas (MODIS land cover product MCD12Q1 with a spatial resolution of 500m [41]), the SRTM data V4 Digital Elevation Model [42] and the derived slope gradient were used as independent (predictive) variables.The data were selected based on the following criteria: (i) the DEM's power to represent the microclimatic condition, as well as the human pressure at different altitudes; (ii) the ability to identify the land suitability for agricultural uses based on the degree of slope; and (iii) the road distance and urban distance strictly considered in relation to the fragmentation and accessibility of the territory.

Data Acquisition and Processing of the MODIS Imagery
The Moderate-Resolution Imaging Spectroradiometer (MODIS) was employed to enhance the fire information obtained by Landsat imagery.Thanks to a quasi-daily revisiting cycle, the MODIS imagery has a high frequency of cloud-free days; however, the spatial resolution is relatively coarse (250 m).The following MODIS time-series products (2001-2013) were acquired through NASA's Earth Observing System (http://reverb.echo.nasa.gov/reverb/)and processed in ENVI 5.0 and ArcGIS 10.2: (i) 16-day MOD13Q1 NDVI composite imagery; (ii) monthly MODIS burned area imagery (MCD45A1) [43]; and (iii) daily active fire point information (Fire Information of Resource Management System Active Fire Data, Collection 4, FIRMS [44]).
The 16-day composite MOD13Q1 NDVI [41] imagery was employed to run a new series of pixel-oriented image differencing operations in order to cover the spatiotemporal gaps left by the Landsat imagery.The FIRMS [44] active fire point information was employed as supporting knowledge for the image differencing operations and further analyses.

Burned Area
Figure 4a shows the accumulated burned paramo land mapped based on the Landsat TM/ETM+ imagery from January 2001 overlaid on the shaded relief map. Figure 4b presents the fires mapped from the MODIS imagery.During the observed period (2001-2013), the detected fires burned a total area of 57,179.8hectares (ha), which is equal to 21.2% of the total area of the high altitude grassland in the 'Complejo de Páramos' of Cruz Verde-Sumapaz.Approximately 77.4% (44,277.8ha) of the burned areas was recorded by Landsat imagery, while 21.5% (12,273.4ha) and 1.1% (628.6) were measured by MODIS MOD13Q1 NDVI data and MODIS MCD45A1 (Figure 5).
The annual burned areas range from 413 (in 2012) to 14,861.3 ha (in 2001) (average 4399 ha y −1 ) (Table 2).The most extreme vegetation changes occurred during 2001 (equal to 5.5% of the paramo area).From January 2001-mid-2007, a large number of fire scars were observed (equal to 89.5% of the total burned area).After this period (mid 2007-2013), the burned paramo areas rarely exceeded 1000 ha y −1 with a peak in 2010 (4026.2ha), which remains substantially smaller than the ones observed in the period 2001-2007.
The estimation of the burn severity [45] was based on the results of the ∆NBR index (10 3 ), which can rank between −2000 and 2000.Four different levels of severity were defined following the general severity classification provided by the United States Geological Survey (USGS) [46,47] ((unburned (∆NBR −100 to 99), low severity (∆NBR −500 to 101, 100 to 269), moderate severity (∆NBR 270 to 439) and high severity (∆NBR 440 to 1300)).Most cases fell into the moderate severity burn equivalent to 40.6% of the burned area, followed by the low severity (31.4% of the burned area) and the high severity burn (21.8% of the burned area).
The very humid paramo ecosystem shows the largest vegetation disturbance caused by fires.With 30,710.1 ha, this ecosystem type accounts for 56.6% of the total burned area mapped (29% of the of the very humid paramo ecosystem).However, 14,299.5 ha (18%), 4536.9 ha (24%) and 640.3 ha (9%) were mapped for the humid paramo, dry paramo and very humid sub-paramo ecosystem types, respectively.Statistics reveal that ca.28,604.3ha, equal to 50% of the total detected burned area, occurred in the Sumapaz National Park.This means that ca.22.6% of the alpine grassland protected by the Sumapaz National Park was affected by fires between 2001 and 2013.The burned land in the Cruz Verde paramo, in contrast, was much smaller, totalling 45.6 ha (0.6% of the Cruz Verde paramo area).From an administrative perspective, Meta was the administrative department most affected by fires with 45.6% of the total burned area (22% of the administrative department), followed by Distrito Capital de Bogotá , Cundinamarca and Huila with 27.2% (19% of the administrative department), 25.1% (20% of the administrative department) and 7.5% (38% of the administrative department), respectively.Most of the fires were located in the flatter paramo areas.Approximately64% of the fires occurred on hillslopes with slope gradients ranging from gentle to moderate slope (<15°).Fires on moderately-steep slopes and steep slopes equalled 26% and 10%, respectively.
Information on the temporal distribution of the fires during the year was obtained by analyzing the FIRMS [44] active fire data.The months with the highest number of records throughout the period 2001-2013 are February (n = 382), January (n = 135) and, to a lesser extent, also December (n = 96).February 2004 is the month that showed the largest number of fire events during the thirteen years observed by MODIS, with 96 active fires detected in total (FIRMS [44]).

Accuracy Assessment
The accuracy assessment analysis was performed using a dataset of 21 burned areas mapped from high resolution imagery.The thematic accuracy analysis shows that 90.5% (n = 19) of the visually identified burned areas were also highlighted in image differencing burned data.Only two burned areas detected by visual interpretation were not detected at all by the image differencing.The linear correlation results in a coefficient of determination (r2) of 0.92 (a < 0.01).It was observed that the burned areas derived by image differencing (n = 19) tend to underestimate the fire area (−9.2%) compared to the visual interpretation validation dataset (Figure 6).The geometric accuracy was carried out using a per-pixel analysis (confusion matrix).The producer accuracy of the detected burned areas was 90.1% with a Kappa Index of Agreement (KIA) of 0.746.If only those burned areas present in both detected and validation datasets are considered (n = 19), the producer accuracy rises to 91.5% with a KIA of 0.771.These later values refer to the real ability of the proposed approach to detect burned areas in the paramo, excluding from the analysis those burned areas (n = 2) that may remain unmapped because of temporal (different acquisition date) or spatial (clouds, gaps) inconsistencies between the datasets.

Statistical Analysis
The results of the logistic regression model illustrated in Table 3 and Figure 7 suggest a significant predictive relevance of the predictive variables road distance (0.55 ROC), DEM (0.56 ROC) and slope gradient (0.53 ROC).The urban distance, in contrast, appeared not to be a significant predictor for the presence of burned surfaces (ROC 0.51).The road and the slope gradient showed indirect behaviours.This indicates that the probability of fire occurrences is higher the shorter the distance to the road is.It also increases the more gentle the slope gradient becomes.For the DEM, the analysis suggests a direct probability behaviour.The ROC curve for the distance from the urban areas did not show a significant slope.

Discussion
By combining vegetation change detection techniques using both Landsat, as well as MODIS imagery, a substantial vegetation disturbance caused by fires could be revealed for the study area, the "Complejo de Paramos" Andino of Cruz Verde-Sumapaz.The findings show that ca.57,180 hectares of the paramo alpine grassland were burned in the course of the study period, between 2001 and 2013.This is tantamount to 21.2% of the total area, of which 28,604 hectares were located in Sumapaz National Park.Interestingly, as the fires are predominantly caused by human [1], these findings show that the actual situation in the paramo does not meet the requirements laid out in law No. 153-1977 and the latest conservation plans according to which the flora, fauna, landscapes and water ecosystems of the Columbian National Park areas must be well preserved for scientific, educational and recreational purposes [4].Based on the maps and information provided by this study, policy makers may need to take a closer look at the situation.
With respect to the applied methodology, a broad series of surface conditions may have caused an unexpected drop in the vegetation index values (VIs) while performing the image differencing operations.This may have influenced the quality of the results.This is often the case for crops, deciduous forests, as well as for grasslands with seasonal dominance of dried vegetation [48,49].However, the rather homogenous vegetation cover of the paramo regions significantly reduced most of the sources of falsely valuated VI changes, thereby facilitating the mapping operations.False-positive burned areas appeared only in some agricultural areas located in the lower paramo.The on-screen visual interpretation and rectification procedures helped to identify and correct these false-positives to map only the effective burned areas (Figure 8).Nevertheless, the extremely high cloudiness characteristic for the area challenged the comprehensive fire monitoring.High clouds and cloud shadow covers come along with a restricted number of Landsat scenes in the United States Geological Survey database [23].Two-thirds of the images were acquired from the Landsat ETM+ after May 2003 due to the gaps of the Landsat 7 Scan Line Corrector-Off (SLC-Off) malfunction [50].
Clouds, shadows and gaps in the cover of the Landsat scenes ranged from 6.1%-95% (mapped based on supervised classification operations).Much of these image disturbances were corrected by partially closing the information gaps in the Landsat scenes (SLC-Off data, as well as cloud, shadow or noise interference) by means of a repeated change detection procedure that used images of preceding years (e.g., the NDVI values of 2003 were subtracted from those of 2001 instead of 2002).This resulted in a considerable reduction of the spatial gaps.As frequently stressed in the literature [27], temporal gaps mainly between March and April remained a challenge and called for further analysis.In line with the analysis of the FIRMS [44] active fire data, the majority of fire records, ca.96%, occurred between January and March [27].This corresponds also to the period with the best cloud-free scenes (January-February) and, therefore, provided confidence that the vast majority of the burned areas were detected and adequately mapped.A temporal analysis of the MODIS active fire data [44] indicated that a 100% complete mapping of fire scars in the area can only be achieved based on cloud-free Landsat scenes acquired during March and April.This, however, appears to be a rather unrealistic expectation in a study area that is known for its extensive cloud cover in March and April.The analysis of the FIRMS [44] active fire data supported the hypothesis of an underestimation of the burned areas.Figure 9 shows a sizeable fire in 2007 where about 2528 ha of alpine grassland were burned.The image was acquired by the Landsat ETM+ sensor on 23 February 2007 (3 p.m.).FIRMS [44] reports twenty-eight active fires across the three day period that also spatially coincided with this fire.Five of these fires were detected some hundreds of meters away from the burned area derived from the Landsat scene base.This example, like many others, indicates an even further expansion of the northern fire front after the Landsat image was taken.The attempt to map fires of a given year (e.g., 2002) based on detected fire scars in the following year (e.g., 2003) did not prove to be successful.Although these operations can provide useful information for some ecosystems [32], this is not the case for paramo, where the soil conditions and the high rainfall levels favour a rapid succession of parts of the vegetation cover [1] (Figure 10).However, it implies that a full restoration of the pre-fire floral and faunal conditions does not take place within a year [51,52].Currently, there is a knowledge gap concerning the post-fire vegetation dynamics in these Neotropical regions, and studies that observe the vegetation development in areas affected by fire events at paramo sites are extremely rare [1,53].This lack of knowledge limits the understanding of how these fires threaten the biodiversity of the paramo ecosystem [54,55] and leaves room for contrasting hypotheses [1].It also prevents modelling of fires based on detected fire scars with a regrowing vegetation cover.
To create a thorough map of the burned areas in the "Complejo de Páramos" Andino of the Cruz Verde-Sumapaz area, further image differencing operations were undertaken employing MODIS NDVI data (MOD13Q1).The quasi-daily data acquisition of the MODIS sensors allowed us to temporally extend the burned area detection processes.This increased the captured information about the burned area by 22% (12,273 ha) and became possible because the fires in the area often reach considerable dimensions (about 32% of the fires mapped by Landsat were larger than 6.25 ha) (Figure 11).
The analysis of FIRMS [44] active fires indicated a clear concentration of fires during the primary dry season (January-March).These insights into the temporal occurrence of fires are mainly consistent with the results reported by other researchers [55,56,57,58].Nevertheless, Amaya and Armenteras (2012) [27] observe the distribution of active fires across Cundinamarca and noted a bimodal distribution with relevant peaks during the secondary dry season (from July-September).During the same period, the active fires observed in the paramos of Cruz Verde-Sumapaz were only marginal (0.4%).A further analysis of the active fires across Cundinamarca revealed different temporal fire patterns in the region.There was a considerable occurrence of fire events in the lowlands during the second dry season (mostly in August), which was not observed in the paramos.In the years 2001, 2004 and 2007, the paramo area showed both the largest burned area size and the highest number of FIRMS [44] active fires.These results for the paramo are consistent with regional and national observations [27,58].The high congruence between the thermal anomalies observed by Pinzón and Tobiá n (2009) [59] and the spatial occurrence of fires highlights the close relationship between these factors at a seasonal and inter-annual level.The years most affected by fires were remarkably arid and coincided with the presence of El Niño (i.e., 2002, 2003, 2004, 2006 and 2007 [60,61]), which seems to influence the occurrence, as well as the extent of the fires.Interestingly, the detected annual burned areas show a sharp decrease for the second half of the study period (2008)(2009)(2010)(2011)(2012)(2013).Despite a relevant lack of information for 2008 and 2009 due to the unavailability of cloudfree Landsat scenes, the low counts of FIRMS [44] active fires in 2011, 2012 and 2013 suggest a structural break compared to previous developments.Statistical analyses between the detected fires and the distribution of the natural ecosystems, climate data (i.e., long-term monthly temperature, monthly precipitation and other bioclimatic variables; World Clim data base) excluded a greater occurrence of fires in the dryer sections of the paramo.The burned area ratio in the very humid paramo (0.29) was noted to be greater than in the dry paramo (0.24).Future research is invited to look further into the reasons for the occurrence of the structural break.The fire history of the paramo and the current flora and fauna indicate that it belongs to the fire-dependent ecosystems [1].The results of the logistic regression model do not suggest the grassland fire dynamics that are typical for areas with large Wildland Urban Interfaces zones [62].However, the insights obtained provided further evidence for the human-induced nature of the fires.
The insights of this study, as well as local literature provide a reason to assume fires as the result of clearing processes carried out by the local communities primarily to create grazing areas [1].Fires are a practical way of removing tussock grass and dried vegetation, thereby facilitating the growth of fresh and nutritious grass.The low occurrence of burned areas left fallow and the low number of FIRMS [44] active fires in the Cruz Verde compared to the Sumapaz paramo support this hypothesis.Paramo of Cruz Verde is located in the proximity of urban areas.Therefore, it appears easier to control and regulate the fires compared to the more remote Sumapaz paramo.The overall climate conditions and the appearance of El Niño events act as catalysts [63].The varying seasonal fire distribution noted for the paramo area and the adjacent lowlands leave room to speculate about differences in the land management practices.This would also be in line with the diminished number of FIRMS [44] active fires observed from 2008-2013.Policy makers are therefore encouraged to examine the situation and to discuss the introduction of regulations and controls on the rural agricultural practices to adequately counteract the outlined threatening developments.

Conclusions
Fires are a serious threat to the integrity of the high-Andean ecosystems.At the same time, livestock grazing, agriculture and other modern practices are essential for some of the rural societies inhabiting the complex paramo regions.Their current and future living conditions, development and economic well-being strongly depend on the exploitation of the resources that the mountain area offers.Simultaneously, these mountain areas fulfil important ecological and economic functions for the adjacent lowlands [9].The paramos area ensures water supply for the population, agriculture, food production and industry at lower altitudes [2].A population of ca.7.5 million people in the metropolitan area of Bogotá depends on this water supply, as more than 80% of the city's water comes from the Chingaza catchment (Chingaza and Sumapaz National Parks), while only 7% comes from the Tunjuelo River Basin [64].Accordingly, Bogotá almost completely depends on the water provided by the mountain wetlands, especially since the forested areas in the region have widely been subject to deforestation.The degradation of the paramo vegetation is likely to trigger processes of accelerated soil erosion [8] from which serious negative impacts on the system of water storage and water quality for the Bogotá metropolitan area can be expected.
The complexity of the human structures and the environmental abundance make the sustainable resource management in these mountainous ecosystems particularly challenging.The needs of the mountain population, the interests of the lowland population and the demands of the biodiverse environment often seem to have conflicting goals.With regard to paramo, a well thought-out management concept may ensure mutual benefits.A possible approach could be to target the preservation of the high biological diversity and endemicity of the area, which would ensure the provision of an adequate water quality and quantity.With a careful management change, the rural practices of the mountain population could be adapted to target the protection of the biological heritage.This may provide opportunities for the development of ecotourism, which may act as an income substitute for the rural communities that reduces their reliance on unsustainable fire management practices.The preservation of the natural vegetation cover in the paramos ensures the vital water supply from this "Water Tower of the Andes" [65].This would also include the installation of groundwater filter systems and of soil conservation measures to keep soil erosion rates low to prevent the siltation of channels and lowland flooding.
The results presented shed light on the fire dynamics in the Colombian paramos and provide a solid basis for further investigation of the impacts on the natural ecosystem functions and biodiversity.Based on the maps and information provided by this study, policy makers may need to take a closer look at the situation.Further research should also investigate the natural and human impact on the spatio-temporal distribution of fires in the paramos, including paleoenvironmental studies.

Figure 2 .
Figure 2. Annual distribution of temperature and precipitation at the study site (based on World Clim-Global Climate Database, 1950-2000).

Figure 3 .
Figure 3. Graphic scheme of the GIS analysis to detect the burned areas.The polygons deriving from the ΔNDVI and ΔNBR analysis were merged.As illustrated in the right-hand image, burn severity was assessed only for the ΔNBR data.NBR, Normalized Burn Ratio.

Figure 4 .
Figure 4. Location of the fire scars detected in the 'Complejo de Páramos' of Cruz Verde-Sumapaz from January 2001-March 2003 through Landsat (a) and MODIS (b) imagery (background data: SRTM data V4).

Figure 5 .
Figure 5. Descriptive statistics of the detected burned areas and FIRMS [44] active fires.

Figure 6 .
Figure 6.Examples of the dataset used during the accuracy assessment analysis.

Figure 7 .
Figure 7. Probability curve modelled as the paramo burned data as a function of the road distance (a), altitude (b), slope gradient (c) and urban distance (d).

Figure 8 .
Figure 8. Example of the detected burned area (a) Polygons obtained as the final outcome of the overall procedure where the false-positives were filtered out by on-screen visual interpretation.(b) The burn severity obtained from the NBR.

Figure 11 .
Figure 11.Landsat image showing multiple scars and active fires in the paramo ecosystems (data base: Landsat TM/ETM; the black line indicates the study area limits).

Table 1 .
The Path 008 Row 057 Landsat time series.

Table 3 .
Logistic regression model of the centroid points of paramo burned areas and the predictive variables road distance, altitude, slope gradient and urban distance.