Fire Data as Proxy for Anthropogenic Landscape Change in the Yucat á n

Fire is one of the earliest and most common tools used by humans to modify the earth surface. Landscapes in the Yucatán Peninsula are composed of a mosaic of old growth subtropical forest, secondary vegetation, grasslands, and agricultural land that represent a well-documented example of anthropogenic intervention, much of which involves the use of fire. This research characterizes land use systems and land cover changes in the Yucatán during the 2000–2010 time period. We used an active fire remotely sensed data time series from the Moderate Resolution Imaging Spectroradiometer (MODIS), in combination with forest loss, and anthrome map sources to (1) establish the association between fire and land use change in the region; and (2) explore links between the spatial and temporal patterns of fire and specific types of land use practices, including withinand between-anthromes variability. A spatial multinomial logit model was constructed using fire, landscape configuration, and a set of commonly used control variables to estimate forest persistence, non-forest persistence, and change. Cross-tabulations and descriptive statistics were used to explore the relationships between fire occurrence, location, and timing with respect to the geography of land use. We also compared fire frequencies within and between anthrome groups using a negative binomial model and Tukey pairwise comparisons. Results show that fire data broadly reproduce the geography and timing of anthropogenic land change. Findings indicate that fire and landscape configuration is useful in explaining forest change and non-forest persistence, especially in fragmented (mosaicked) landscapes. Absence of fire occurrence is related usefully to the persistence of spatially continuous core areas of older growth forest. Fire has a positive relationship with forest to non-forest change and a negative relationship with forest persistence. Fire is also a good indicator to distinguish between anthrome groups (e.g., croplands and villages). Our study suggests that active fire data series are a reasonable proxy for anthropogenic land persistence/change in the context of the Yucatán and are useful to differentiate quantitatively and qualitatively between and within anthromes.


Introduction
This study examines anthropogenic land cover change and persistence across the Yucatán Peninsula of Mexico by incorporating fire occurrence as a geographic context-specific indicator of human land practices.Land use and land cover change (i.e., land change) is an International Council for Science (ICSU) millennial research challenge [1] and has been addressed by various research communities [2][3][4][5].A consensus has emerged that human-induced land changes have increased at unprecedented rates throughout the last century [5][6][7] to such an extent that the long-term well-being of human-environment relationships now requires more sustainable management of land uses [8,9].Without disagreeing with the previous assertion, other researchers have argued that humans have altered ecosystems to such a degree that present landscapes would be best characterized in terms of anthropogenic biomes or anthromes, as they occupy the largest portion of the terrestrial land surface [10].These concerns have guided numerous efforts to detect, map, and monitor land change, especially with the aid of remote sensing technology (e.g., [11,12]).Complementary mapping efforts have been undertaken to understand the level of direct human influence, as well as efforts to model the causes, constraints, mechanisms, and trajectories of land change and persistence (e.g., [2,10,[13][14][15][16]).
While much has been learned from these research efforts, certain research challenges remain unaddressed [4,16,17].Two of these challenges are pertinent for highly dynamic and spatially heterogeneous mosaicked humanized landscapes: the need for additional explicit studies seeking to untangle region-specific human-environment interactions and the need to incorporate context and temporal information from satellite data sources into the analysis [4,[16][17][18].
The accuracy and spatial explicitness of land use and land cover maps has improved in the past 20 years with the refinement of the spatial and spectral resolution of remotely sensed data [12,16,19].For local-to-regional scale mapping, land change research has relied heavily on cross-sectional assessments (e.g., two or three map dates; [20]).Although valuable, cross-sectional approaches often overlook high frequency land dynamics occurring between land cover map dates [21] and can potentially misrepresent recurrent changes as unidirectional trends, leading to over-or under-estimation of change [22].In the last decade, high temporal frequency remotely sensed datasets having1 km spatial resolution and improved calibration such as the Advanced Very High Resolution Radiometer (AVHRR), the Global Inventory Modeling and Mapping Studies (GIMMS) NDVI data and the MODIS data product suite [23] have been made available.While AVHRR and GIMMS data are ideal for monitoring dynamic, long-term environmental processes at continental and global scales [11,18,24,25], their application is limited at sub-national scales [20].The MODIS product suite in contrast, although limited by its shorter time-span, offers not only better or equivalent spatial, spectral and temporal capacity [26], but also the potential use of the MODIS archive of derived products (e.g., sea ice, active fire) as ancillary data to yield improvements in land cover mapping and change detection (e.g., [27,28]).
In addition to improvements in spatial and temporal details of land cover mapping, researchers have made significant progress to specify and even re-conceptualize traditional land cover/land use categories.Specifically, the re-framing of categories in terms anthropogenic land use intensity (e.g., anthromes) places the focus on human activities.Many of these efforts have relied on population density to characterize differences between areas that are highly or sparsely intervened by humans.Heterogeneity is another criterion commonly used to characterize anthropogenic landscapes.In those instances, landscapes are often measured on a homogeneity continuum from highly heterogeneous areas of "built" landscapes combining many different land uses in relatively small areas, to others presenting large swaths of homogenous features [10].
Fire is one of the earliest and most significant tools used by humans to modify landscapes and to extract benefits from ecosystem goods and services (e.g., land clearing for agriculture or urban development, forest thinning).At the same time, human land use and management (e.g., fuel load quantity and distribution, regulation and ignition sources) have altered naturally occurring fire regimes [29][30][31].The "pyric transition" model has been proposed to describe the relationship between human intervention and the intensity of fire use across time and space (e.g., none in dense urban contexts, high in agricultural extensive systems, low or background in relatively "untouched" ecosystems) [29].This relationship, however, is highly context-specific, in terms of the biophysical and the cultural characteristics of a region.
The present study takes advantage of that potential by using MODIS active fire data time-series in combination with categorical land cover/use maps to examine the link between land change, human activities and fire in the Yucatán Peninsula.

Study Area
The Yucatán Peninsula (henceforth, the Yucatan), defined in this study by the Mexican states of Campeche, Quintana Roo and Yucatán (Figure 1), is a region characterized by a mostly flat, karstic, quaternary plain.The region is a highly dynamic and complex mosaic of subtropical dry forest in different states of ecological succession driven by human and environmental disturbance (e.g., fire, tropical storms, agriculture) [32][33][34][35].Unless interrupted by specific disturbances such as fire-friendly plant invasion or management practices, post-fire vegetation grows back at a rapid rate after a burn.In much of the Yucatán Peninsula, secondary vegetation can reach a height of 10 m within-five years of clearing [36][37][38].As a consequence, there is typically only a short window of time to detect land change. .This relationship, however, is highly context-specific, in terms of the biophysical and the cultural characteristics of a region.The present study takes advantage of that potential by using MODIS active fire data time-series in combination with categorical land cover/use maps to examine the link between land change, human activities and fire in the Yucatán Peninsula.

Study Area
The Yucatán Peninsula (henceforth, the Yucatan), defined in this study by the Mexican states of Campeche, Quintana Roo and Yucatán (Figure 1), is a region characterized by a mostly flat, karstic, quaternary plain.The region is a highly dynamic and complex mosaic of subtropical dry forest in different states of ecological succession driven by human and environmental disturbance (e.g., fire, tropical storms, agriculture) [32][33][34][35].Unless interrupted by specific disturbances such as fire-friendly plant invasion or management practices, post-fire vegetation grows back at a rapid rate after a burn.In much of the Yucatán Peninsula, secondary vegetation can reach a height of 10 m within-five years of clearing [36][37][38].As a consequence, there is typically only a short window of time to detect land change.Fire in the Yucatán is not only an environmental disturbance, but also the main land-clearing tool in the region.It precedes most land conversion (e.g., urbanization, clear-cut deforestation) and land modification events such as shifting cultivation, pasture maintenance, and secondary vegetation regrowth [32,[37][38][39][40]. Persistence and maintenance of non-forest land cover types (e.g., pasture, various crops, early secondary vegetation) are also controlled by fire.In the Yucatán, two types of fire are commonly distinguished by local population and regulation: quemas and incendios.Quemas (~agricultural burns) are, in theory, legally registered and strictly regulated and scheduled for the purpose of swidden (slash-and-burn) agriculture.This includes subsistence maize-based cropping (milpa), pasture maintenance and mechanized agriculture (e.g., sugar cane, commercial corn).Incendios forestales (wildfires) are, in theory, accidental fires that are ignited by cigarette butts, broken glass, or lightning.They typically occur during the dry season (March-May), can affect large areas and last for several days.These fires, if a consequence of human action, are illegal and usually generate government agency response, especially if they grow out of control and affect farming, forest extraction and natural protected areas (NPA) [40,41].Although the distinction can be easily Fire in the Yucatán is not only an environmental disturbance, but also the main land-clearing tool in the region.It precedes most land conversion (e.g., urbanization, clear-cut deforestation) and land modification events such as shifting cultivation, pasture maintenance, and secondary vegetation regrowth [32,[37][38][39][40]. Persistence and maintenance of non-forest land cover types (e.g., pasture, various crops, early secondary vegetation) are also controlled by fire.In the Yucatán, two types of fire are commonly distinguished by local population and regulation: quemas and incendios.Quemas (~agricultural burns) are, in theory, legally registered and strictly regulated and scheduled for the purpose of swidden (slash-and-burn) agriculture.This includes subsistence maize-based cropping (milpa), pasture maintenance and mechanized agriculture (e.g., sugar cane, commercial corn).Incendios forestales (wildfires) are, in theory, accidental fires that are ignited by cigarette butts, broken glass, or lightning.They typically occur during the dry season (March-May), can affect large areas and  [40,41].Although the distinction can be easily blurred (e.g., escaped burns), quemas tend to be more ubiquitous in the region, are smaller in area, and shorter in duration than incendios.In contrast, incendios are less frequent, cover large areas, and temporally discrete events that are often amplified in contact with hurricane blow downs [42].
Given the previously described links between various human activities and fire practices, we expect that the spatio-temporal pattern of fire can serve as a proxy for certain anthropogenic land use practices in the Yucatán (e.g., milpa, mechanized agriculture, sugar cane cultivation and new urban development).Naturally occurring wild fires are less frequent in this region than fire of anthropogenic origin.This characteristic represents an opportunity to evaluate fire as a proxy for various land uses and land cover changes in the region.
The primary goal of this study is to assess land cover change and persistence across the Yucatán Peninsula by incorporating MODIS active fire occurrence data as a site-specific indicator of human-induced land change.The specific objectives are to: (a) describe the spatio-temporal dynamics of fire and its association with land use practices and land change; (b) determine the importance of fire presence to explain types of land persistence and change using a spatial multinomial logit model; and, (c) explore the ability of fire frequency to discriminate between anthropogenic landscape types (anthropogenic biomes or anthromes).
Beside the specific goal stated, lessons learned in this study may be of wider interest in other regions of the world.Spatial driver variables used for land-change modeling and prediction are typically expensive to produce or not frequently updated as most countries and regions do not possess such information at the level of detail required.Similarly, land cover maps are rarely produced annually.Most land cover maps are commonly produced using a census like approach with 5-10 years gaps at the country or regional levels.For an increasing number of applications, mapping rapid changes within-these temporal gaps across countries boundaries may be crucial.It is in this context that active fire data becomes useful as a predictor.It has both high temporal frequency and is readily available worldwide.To evaluate the usefulness of fire as predictor in both data poor and rich areas, we modeled land change with and without the presence of additional spatial drivers.Results reported in this specific study should therefore be of potential interest to a wider audience of researchers and land management practitioners.

Data
A comprehensive active fire dataset was constructed using NASA's 1 km MOD14A2 (from TERRA) and MYD14A2 (from AQUA) 8-day composite product, using various time ranges within-the time period 2000-2010 [43] in order to match fire data and land cover data availability.Both data products provide detection confidence levels [44] of which only high and medium were used in order to minimize potential false positive detections [27,45] (Figure 2).The 1 km fire grid cell locations were converted to centroid points to match and extract the 30 m grid cell land cover attributes [42].
Land change/persistence information was obtained from Southern Mexico Forest Cover and Clearance from c. 1990 to c. 2000 and to c. 2007 digital maps developed by the Center for Applied Biodiversity Science (CABS) at Conservation International [46].The original maps are based on the interpretation of (~30 m resolution) Landsat-5 Thematic Mapper satellite imagery, and are composed of five categories: forest, non-forest, mangrove, water, and cloud/shadow.Forest cover is defined as mature forest >3 m closed canopy [44].Thematic map accuracy was reported to have a Kappa index of 86% with an error of commission of 11.27% and an error of omission of 4.86% for the dominant tropical dry forest class prevalent in the study area [44].Forest and mangrove classes were merged, whereas water, clouds and shadow pixels were masked and removed from further analysis to produce a map of land persistence and change for 2000-2007, with three resulting categories: forest persistence (FP), non-forest persistence (NFP), and forest-to-non-forest (CH) (Figure 3).It must be noted that in contrast with FP ("old growth" forest persistence), NFP by definition, allows for land changes between sub-classes within-a broad non-forest category.For example, a transition from agriculture to pasture or from pasture to secondary forests would be included in NFP.Although lacking detailed land categories, the resulting three-category map has the advantage of being more consistent than most multi-temporal assessments because it was produced using a consistent methodology by the same research team explicitly for the purpose of forest change mapping [46,47].
Land 2017, 6, 61 5 of 19 lacking detailed land categories, the resulting three-category map has the advantage of being more consistent than most multi-temporal assessments because it was produced using a consistent methodology by the same research team explicitly for the purpose of forest change mapping [46,47].Additional variables commonly used in land-change modeling were included in the analysis to control for environmental and socio-economic factors.Elevation and slope were derived from a 90 m SRTM digital elevation model.Precipitation data were derived from Tropical Rainfall Measuring Mission (TRMM).Protected Areas from CONANP (Mexican Natural Protected Area Commission) were processed to match the study area.Communal property lands (ejidos), population density change, as well as the distance to road variable were derived from the Mexican Census Bureau [49].Cattle density by municipality was obtained from the 2007 Censo Agropecuario [49].Even though MODIS fire data are constantly updated, in this study, we only include data to the year 2010.We consider the 2000-2010 ten-year period sufficient to describe the salient points about the spatio-temporal patterns of a fire in the region.More importantly, the land cover data and anthrome data that we use do not go beyond 2007 and 2008, respectively.For this same reason, the multinomial modeling component of this study does not use any fire data from after 2007 (see methods section).

Methods
The methodology has three components: (a) a characterization of the spatio-temporal patterns and co-location of land change/persistence and fire; (b) an assessment of the effect of fire presence, landscape patch size and other covariates land-change outcomes (FP, NFP and CH), using a spatial multinomial logit model for the period 2000-2007; and, (c) an exploration of the ability of fire frequency to distinguish human landscape categories defined in anthrome groups and anthrome classes.

Characterization of Land Change and Fire Patterns
A cross-tabulation-based change analysis was performed on the simplified CABS land cover map pairs 1990-2000 and 2000-2007, resulting in one map with two persistence classes (FP and NFP) and one change class (CH) (Figure 3).Spatio-temporal dynamics of fire occurrence were derived by manipulating the MODIS fire archive via map algebra operations (e.g., overlay, reclassification) into map and graph summaries of annual and multi-year fire frequency.After accounting for the different time-spans of the TERRA and AQUA series and adjusting for double counting, three separate fire occurrence sets were created: a TERRA ten-year summary map (2000-2010) (Figure 2) and a TERRA-AQUA combined annual time-series graph (2000-2009) (Figure 4).Global (Moran'I) and local (Getis-Ord Gi* and LISA) measures of spatial autocorrelation were calculated for all sets of maps [50][51][52].Finally, a second cross-tabulation analysis was performed for the land change/persistence map (FP, NFP and CH) (Figure 3) and the binary fire (reclassification of Figure 2) maps in order to describe and explore the spatial correspondence between them and the relationship between fire occurrence and the different land cover categories.The result of this cross-tabulation map is shown in Figure 5.

Effect of Fire and Patch Size on Land Change/Persistence: A Spatial Multinomial Logit Model
A spatial multinomial logit model using 2003-2007 fire "density" or "intensity" (fire count/five years) as independent variable (FIRE), and 7 control variables was performed.This model specification is used in the literature (e.g., [53,54]) and is preferred because the dependent variable consists of three possible outcomes that are categorical/nominal in nature [55][56][57].To account for the positive spatial autocorrelation present in remotely sensed data [58,59] spatial lags variables based on rook-case adjacency were calculated.This was done by generating, for each observation, the

Effect of Fire and Patch Size on Land Change/Persistence: A Spatial Multinomial Logit Model
A spatial multinomial logit model using 2003-2007 fire "density" or "intensity" (fire count/five years) as independent variable (FIRE), and 7 control variables was performed.This model specification is used in the literature (e.g., [53,54]) and is preferred because the dependent variable consists of three possible outcomes that are categorical/nominal in nature [55][56][57].To account for the positive spatial autocorrelation present in remotely sensed data [58,59] spatial lags variables based on rook-case adjacency were calculated.This was done by generating, for each observation, the

Effect of Fire and Patch Size on Land Change/Persistence: A Spatial Multinomial Logit Model
A spatial multinomial logit model using 2003-2007 fire "density" or "intensity" (fire count/five years) as independent variable (FIRE), and 7 control variables was performed.This model specification is used in the literature (e.g., [53,54]) and is preferred because the dependent variable consists of three possible outcomes that are categorical/nominal in nature [55][56][57].To account for the positive spatial autocorrelation present in remotely sensed data [58,59] spatial lags variables based on rook-case adjacency were calculated.This was done by generating, for each observation, the average sum of its neighbors of the same class (CFP, CNFP, and CCH), and including the lag variable in the model (see Equation ( 1)) [59,60].Lag values range from 0 (no like neighbors) to 1 (all four like neighbors), while average Moran's I values range from virtually random for CCH (0.0035) to weak-moderate positive for CNFP (0.392) and CFP (0.391).In this particular case, spatial lags also serve as an indicator of landscape configuration and heterogeneity (i.e., similar to patch size and cohesion).
Pseudo-likelihood estimates were calculated for the main effects and spatial lag parameters of the multinomial model as follows (Equation ( 1)): where π hij is the probability that a pixel h in the map of any class i for which fire has been detected will become one of land change/persistence category j (FP, NFP or CH), j = r, with land change/persistence category r being the reference category.There are separate sets of intercept parameters αj and regression parameters β j for each of the logit models.The matrix x hi is the set of explanatory variables: spatial lag (CY), fire (FIRE), and the commonly used control variables described.Two logits equations are modeled for each FIRE and CY population: the logit comparing non-forest persistence (NFP) to forest persistence (FP) and the logit comparing change (CH) to forest persistence (FP) [56].

Fire Frequency and Anthromes Characterization
Anthromes map datasets were processed to match them to the study area, and temporal fire count data were aggregated to generate a fire frequency variable (Figure 6).In order to explore the association between fire counts and anthrome groups graphically, we used boxplots of mean fire frequency by each category present in the study area: dense settlements, villages (in effect, village and surrounding areas), croplands, rangelands, forested, wildlands.We measured the differences in mean fire frequency across anthromes using a negative binomial model and a Tukey pairwise comparison.A binomial model was preferred over Poisson model because the overall (as well as per category) mean and variance are noticeably different (by a factor 100), which suggests over-dispersion in the Poisson parameter.We also conducted Tukey comparisons for fire frequency within the cropland anthromes group to assess the potential of using fire count for finer definition of the human imprint of landscape.1)) [59,60].Lag values range from 0 (no like neighbors) to 1 (all four like neighbors), while average Moran's I values range from virtually random for CCH (0.0035) to weakmoderate positive for CNFP (0.392) and CFP (0.391).In this particular case, spatial lags also serve as an indicator of landscape configuration and heterogeneity (i.e., similar to patch size and cohesion).Pseudo-likelihood estimates were calculated for the main effects and spatial lag parameters of the multinomial model as follows (Equation ( 1)): where πhij is the probability that a pixel h in the map of any class i for which fire has been detected will become one of land change/persistence category j (FP, NFP or CH), j ≠ r, with land change/persistence category r being the reference category.There are separate sets of intercept parameters αj and regression parameters βj for each of the logit models.The matrix x′hi is the set of explanatory variables: spatial lag (CY), fire (FIRE), and the commonly used control variables described.Two logits equations are modeled for each FIRE and CY population: the logit comparing non-forest persistence (NFP) to forest persistence (FP) and the logit comparing change (CH) to forest persistence (FP) [56].

Fire Frequency and Anthromes Characterization
Anthromes map datasets were processed to match them to the study area, and temporal fire count data were aggregated to generate a fire frequency variable (Figure 6).In order to explore the association between fire counts and anthrome groups graphically, we used boxplots of mean fire frequency by each category present in the study area: dense settlements, villages (in effect, village and surrounding areas), croplands, rangelands, forested, wildlands.We measured the differences in mean fire frequency across anthromes using a negative binomial model and a Tukey pairwise comparison.A binomial model was preferred over Poisson model because the overall (as well as per category) mean and variance are noticeably different (by a factor 100), which suggests overdispersion in the Poisson parameter.We also conducted Tukey comparisons for fire frequency within the cropland anthromes group to assess the potential of using fire count for finer definition of the human imprint of landscape.
For the data processing, we used the raster (version 2.5-8), sp (version 1.2-4) and rgdal (version 1.2-7) packages to perform some of the preprocessing for the spatial inputs.The raster package was also used for aggregation of the fire occurrences, cropping and reclassification of anthromes data.We used the implementation of multinomial model in the R "nnet" package (version 7.3-12) for analyses.We also ran the negative binomial model using the "MASS" package (version 7.3-12) in R version 3.4 and used the glht method from the multcomp R package (1.4-6) to perform the Tukey comparison.For the data processing, we used the raster (version 2.5-8), sp (version 1.2-4) and rgdal (version 1.2-7) packages to perform some of the preprocessing for the spatial inputs.The raster package was also used for aggregation of the fire occurrences, cropping and reclassification of anthromes data.We used the implementation of multinomial model in the R "nnet" package (version 7. [3][4][5][6][7][8][9][10][11][12] for analyses.We also ran the negative binomial model using the "MASS" package (version 7. [3][4][5][6][7][8][9][10][11][12] in R version 3.4 and used the glht method from the multcomp R package (1.4-6) to perform the Tukey comparison.

Spatio-Temporal Patterns of Fire and Land Change
As a proxy indicator of human land use, fire frequency can be interpreted as a quantitative measure of "density" (in space and time) of land use in this case.It corresponds to how actively "used" or "under-used" a 1 km parcel of land is in relation to the types of land use that require burning (e.g., high for pasture and cultivation land).
Fire occurrence for the period 2000-2010 presents a distinct non-random spatial pattern across the study area (Moran's I = 0.62) (Figure 2).The individual observations presenting the highest detection frequencies (>3 per year) appear in spatial clusters and are concentrated in the northern and western section of the Peninsula in the states of Campeche and the Yucatán.Quintana Roo has comparatively fewer and less clustering fire detections.Areas with less fire frequency are observed in Campeche and Quintana Roo, especially in the presence of federal and state run natural protected areas (e.g., Sian Ka'an and Calakmul bioreserves) (Figure 1), and in locations where road and settlement density is low.This road-settlement-fire co-location is more explicit in the Yucatán and northern Campeche, and less so in southern Campeche, where road density decreases.
Locations affected by fire (Figure 5 and Table 1) cover roughly 45% of the study area.More than half of the fire occurrences concentrate in NFP areas, followed by 40% of FP, and 2% in CH areas (1% of the total).Locations with no fire occurrence cover more than 50% of the study area.No-fire occurrences are two times more common in FP than in NFP, and they are minimally present in CH areas.Land change/persistence figures (total rows percentage in Table 1) indicate that roughly 2/3 of the FP shows no fire detections, whereas 60% of NFP presents fire occurrence.CH shows 86% of its area as having burned at some time (although it is hardly visible in Figure 5, because CH covers only 1% of the total area).
Table 1.Cross-tabulation of fire occurrence and land change/persistence: row, column and total percentages per class.Preliminary interpretation matrix (extended/detailed legend for Figure 5).

Fire
No Appendix A).FP is also 1.07 times more likely than NFP.In general, CH has a negative relationship with CY.This relationship can be attributed to CH occurring more often in isolated pixels corresponding to fragmented patches.Conversely, FP and NFP occur in larger continuous areas (Figure 3).

Comparison of Fire by Anthromes
Boxplots of mean fire frequency by anthromes groups suggest that fire can be used to distinguish among anthromes groups thereby serving as an additional potential variable to examine human imprints in the context of the Yucatán region (Figure 6).The croplands anthromes group has a higher fire frequency on average compared to most other anthromes categories.We note that rangelands display a very large variance and thus may be difficult to distinguish from croplands and villages.Finer anthropogenic categories within the cropland groups also appear to be contrasted by the fire frequency variable (see boxplots Figure 7).

Comparison of Fire by Anthromes
Boxplots of mean fire frequency by anthromes groups suggest that fire can be used to distinguish among anthromes groups thereby serving as an additional potential variable to examine human imprints in the context of the Yucatán region (Figure 6).The croplands anthromes group has a higher fire frequency on average compared to most other anthromes categories.We note that rangelands display a very large variance and thus may be difficult to distinguish from croplands and villages.Finer anthropogenic categories within the cropland groups also appear to be contrasted by the fire frequency variable (see boxplots Figure 7).
Tukey pairwise comparison tests between anthromes groups (applied to the negative binomial model Table 3) highlight that six and seven pairwise comparisons are significantly different at 5% and 10% level respectively.Coefficients indicate relative increase or decrease of the expected log of fire frequency for categorical change in relation to the reference variable.For instance, there is a −1.37 decrease in the log of fire frequency when an observation is found in the dense settlements rather than in croplands.We find that fire frequency variable is useful to differentiate croplands with almost all categories with the exception of Villages and Rangelands.Both overlap substantially on the boxplots (Figure 7).Tukey pairwise comparison tests between anthromes groups (applied to the negative binomial model Table 3) highlight that six and seven pairwise comparisons are significantly different at 5% and 10% level respectively.Coefficients indicate relative increase or decrease of the expected log of fire frequency for categorical change in relation to the reference variable.For instance, there is a −1.37 decrease in the log of fire frequency when an observation is found in the dense settlements rather than in croplands.We find that fire frequency variable is useful to differentiate croplands with almost all categories with the exception of Villages and Rangelands.Both overlap substantially on the boxplots (Figure 7).
Tukey comparisons were also performed for fire frequency within the cropland anthromes group to assess the potential of using fire frequency for finer definition of the human imprint of landscape.We found that within croplands two out 10 and four out of 10 comparison results in significant differences in the log of fire frequency (Table 4) at 5 and 10 percent significance level, respectively.This suggests that fire frequency is less useful at distinguishing within the cropland anthromes.The strongest contrasts are found between "residential rainfed mosaics-remote croplands mosaic" and between "residential irrigated cropland-remote croplands".Both indicate decreases in the log of fire frequency for locations that are found in residential rainfed and irrigated cropland compared to remote croplands.

Discussion and Conclusions
Analysis of the spatial patterns of fire in the Yucatán using fire data visually reproduces the geography of human land-based activities during 2000-2010, including the influence of roads, settlements, and natural protected areas.The intra-annual temporal dimensions of fire data approximate key moments of seasonal agriculture and other land-based activities (Figure 4).These findings are consistent with studies finding higher fire frequencies in the proximity to roads and settlements region [40], and with research about fire and deforestation in other tropical settings (e.g., [28,[62][63][64]) and more generally, about the relationship between land change, road development and market accessibility [53,54,64,65].
The multinomial model confirms the significance of fire intensity (defined as count of fire flagged pixels over time) and landscape configuration (represented by the spatial lag variable CY) for the peninsula for all three types of land change/persistence outcomes.Notably all change is strongly affected by the presence of fire.The odds of a forested area becoming non-forested multiply by six in the presence of fire.This is not surprising given that the cross-tabulation analysis showed that nearly 87% of pixels that changed from 2000 to 2007 had been flagged as fire by MODIS at some point.The dominance of no-fire observations over fire occurrence, the (pseudo) absence of fire, makes for a better predictor of forest persistence than fire occurrence per se, especially for large patches (high CY values).Predicting non-forest persistence proved more difficult, signaling the diversity of ecological and land use practices and transitions between them that define the NFP category as discussed in the data section.What exactly accounts for the remaining 13% is difficult to determine.We can speculate however that possible factors include data error and accuracy mismatches.As some have noted, MODIS products do omit fires for various reasons [40].We can also speculate on some detected fires being non-agricultural.Finally, we can also imagine some fires starting as agricultural burns but becoming large multiday forest fires as they grow out of control [40].
Fire also proves to be a reasonable indicator to distinguish between anthromes groups (Figure 7) but less among finer anthromes categories.Our analysis shows that, for example, with the aid of a fire frequency croplands are separable from the rest of categories and that dense settlements are clearly distinct from villages.Other categories like rangeland prove more difficult to separate since they present larger variance.Some within-variation within anthropogenic categories is also distinguishable.For example, populated irrigated is clearly separable from populated rainfed and from remote croplands.These two former categories, however are not as separable from each other with the aid of fire frequency alone.
Perhaps some insight about differences distinguishing between-and within-anthropogenic landscapes can be gained by looking into the distinction between quemas (legally regulated agricultural fires) from incendios (large escape non-regulated forest fires).Both can be addressed through the use of MODIS active fire data, but cannot be thoroughly distinguished with these data alone.Within-variation exists in both categories.On the one hand, quemas can be due to subsistence and commercial agriculture, such as sugar cane cultivation along the border with Belize (i.e., either rain fed or irrigated anthrome categories), or Mennonite corn fields in Bacalar,( i.e., the populated irrigated anthrome category), as their spatial and temporal signatures can differ quite noticeably (Figures 2 and 4).On the other hand, incendios can be truly accidental or be the result of purposeful burning with nonagricultural goals (e.g., forcing land zoning code changes for urban development, hunting practices).Future research should exploit these differences in order to account for a more nuanced description of the relationships between fire and human-induced land change in the region.
We conclude that in subtropical forest settings, in which burning is a fundamental landscape tool, fire is a reasonable proxy for land change.Use of fire data improves the estimation of land change and expands our understanding of human-environment relationships that produce complex mosaics of highly dynamic landscapes, such as those in the Mexican Yucatán.

Figure 1 .
Figure 1.Study Area: Mexico's Yucatán Peninsula.Major federal roads, main cities and towns, state boundaries, and protected areas.

Figure 1 .
Figure 1.Study Area: Mexico's Yucatán Peninsula.Major federal roads, main cities and towns, state boundaries, and protected areas.

Figure 2 .
Figure 2. Fire frequency map, 2000-2010 MODIS TERRA-only series.Significant Getis-Ord Gi statistic (0.05 level) hot spots and LISA-based clusters are shown in inset (upper left). 1 = 2006 post-hurricane Wilma fires, and urban development related fires; 2 and 5 = Mennonite land clearing for commercial maize cultivation; 3 = PEMEX oil processing plant; and, 4 and 6 = high yield commercial corn and bean cultivation.

Figure 2 .
Figure 2. Fire frequency map, 2000-2010 MODIS TERRA-only series.Significant Getis-Ord Gi statistic (0.05 level) hot spots and LISA-based clusters are shown in inset (upper left). 1 = 2006 post-hurricane Wilma fires, and urban development related fires; 2 and 5 = Mennonite land clearing for commercial maize cultivation; 3 = PEMEX oil processing plant; and, 4 and 6 = high yield commercial corn and bean cultivation.

Figure 2 .
Figure 2. Fire frequency map, 2000-2010 MODIS TERRA-only series.Significant Getis-Ord Gi statistic (0.05 level) hot spots and LISA-based clusters are shown in inset (upper left). 1 = 2006 post-hurricane Wilma fires, and urban development related fires; 2 and 5 = Mennonite land clearing for commercial maize cultivation; 3 = PEMEX oil processing plant; and, 4 and 6 = high yield commercial corn and bean cultivation.

Figure 3 .
Figure 3. Land change/persistence classes based on the CABS (2009) land cover maps.Figure 3. Land change/persistence classes based on the CABS (2009) land cover maps.

Figure 3 .
Figure 3. Land change/persistence classes based on the CABS (2009) land cover maps.Figure 3. Land change/persistence classes based on the CABS (2009) land cover maps.

Figure 4 .
Figure 4. Cumulative annual and 8-day frequency (counts) of fire occurrence (FIRE), 2000-2009 in the three-state Yucatán Peninsula: blue = morning purple = afternoon.Inset shows detail of intra-annual fire season for 2005.

Figure 4 .
Figure 4. Cumulative annual and 8-day frequency (counts) of fire occurrence (FIRE), 2000-2009 in the three-state Yucatán Peninsula: blue = morning purple = afternoon.Inset shows detail of intra-annual fire season for 2005.

Figure 4 .
Figure 4. Cumulative annual and 8-day frequency (counts) of fire occurrence (FIRE), 2000-2009 in the three-state Yucatán Peninsula: blue = morning purple = afternoon.Inset shows detail of intra-annual fire season for 2005.
These fires, if a consequence of human action, are illegal and usually generate government agency response, especially if they grow out of control and affect farming, forest extraction and natural protected areas (NPA)

Table 3 .
Between-anthromes Tukey pairwise comparison for the negative binomial models with fire count and anthromes groups.Note that the reference variable is in italic.

Table 3 .
Between-anthromes Tukey pairwise comparison for the negative binomial models with fire count and anthromes groups.Note that the reference variable is in italic.

Table 4 .
Within-anthrome Tukey pairwise comparisons using the negative binomial models with fire count and anthromes croplands categories.Note that the reference variable is in italic.