Threatened Habitats of Carnivores: Identifying Conservation Areas in Michoac á n, M é xico

: The present study contributes to bridging the gap in research related to the presence and distribution patterns of carnivore mammals in western M é xico and identiﬁes priority areas for biodiversity conservation in western Michoac á n, M é xico. The distribution of 11 carnivore species ( Canis latrans ; Urocyon cinereoargenteus ; Herpailurus yagouaroundi ; Leopardus pardalis ; Leopardus wiedii ; Puma concolor ; Panthera onca ; Conepatus leuconotus ; Bassariscus astutus ; Nasua narica ; Procyon lotor ) in western M é xico was modeled through the application of a two-scale approach, including a large modeled region that corresponded to the western part of the country, for which consensus models were obtained that represent the species’ bioclimatic envelopes (historic occurrence records); and the second modeled study area that includes only the western portion of the state of Michoac á n in which compounded models of the species’ habitat suitability (ﬁeld occurrence records) for this region were proposed. Using species’ habitat suitability models as biodiversity units, prioritization exercises were carried out on important areas for the conservation of these species, as well as the comparison and analysis of the existing natural protected areas (NPA) and existing proposed conservation areas in the study area. The different exercises for prioritizing areas for conservation yielded similar results and show the potential percentages of the landscape that can be subjected to conservation programs. The highest conservation priority values were mainly located in the Costas del Sur and Cordillera del Sur provinces. This study signiﬁes a ﬂexible basis from which future studies on planning and designing a network of natural protected areas can be carried out in this region.


Michoacán's Mammals
The state of Michoacán is located in the central-western part of México and includes a variety of physiographic regions [1,2] (Figure 1, Table A1). Michoacán's heterogeneity in climate, physiography and lithology promotes ecological conditions that allow for the presence of a wide variety of vegetation communities [3] associated with a high fauna diversity. A high diversity of vertebrate species is found in the transition zones between physiographic regions in the state of Michoacán [4].

Carnivores and Biodiversity Indicators
Carnivores play an important ecological role; these can belong to many trophic levels due to their variance in size, their separation between temporal and spatial niches, and their food habits [6,7]. Their presence and population conditions represent valuable ecological metrics because this is a group indicator of the conservation status in different ecosystems. Species of carnivores are considered to be umbrella species, since the conservation of their distribution aids in protecting many other species inhabiting the same habitats [6]. As large predators, carnivores control their prey populations, affecting the trophic web dynamics, and therefore the ecosystem's energy flows [8].
Issues related to the global conservation of carnivores are widely documented [9][10][11][12][13][14]. Their specific food habits, their need for large areas and their low tolerance to human presence make carnivores highly vulnerable to extinction processes [15]. The strongest threats are the destruction and fragmentation of habitats and the intense hunting of them and their prey [16].

Carnivores and Biodiversity Indicators
Carnivores play an important ecological role; these can belong to many trophic levels due to their variance in size, their separation between temporal and spatial niches, and their food habits [6,7]. Their presence and population conditions represent valuable ecological metrics because this is a group indicator of the conservation status in different ecosystems. Species of carnivores are considered to be umbrella species, since the conservation of their distribution aids in protecting many other species inhabiting the same habitats [6]. As large predators, carnivores control their prey populations, affecting the trophic web dynamics, and therefore the ecosystem's energy flows [8].
Issues related to the global conservation of carnivores are widely documented [9][10][11][12][13][14]. Their specific food habits, their need for large areas and their low tolerance to human presence make carnivores highly vulnerable to extinction processes [15]. The strongest threats are the destruction and fragmentation of habitats and the intense hunting of them and their prey [16].

Species Distribution Models
Analyzing patterns of species distributions during the last few decades has been a research tool used in important fields such as conservation planning. Therefore, many species distribution algorithms and procedures have been developed, which aim to improve the prediction of such models [17]; such is the case for the application of geographic information systems (GIS) and the development of statistic and machine learning algorithms oriented to generate robust species distribution models [18][19][20][21][22][23][24][25][26][27][28][29][30]. Considering the wide assortment of algorithms [31], an ensemble of consensus of models has been proposed as an approach to generate more robust species distribution predictions [32].
For some authors, ecological niche modeling (ENM) is not synonymous with species distribution modeling (SDM) [33]. However, the niche concept is key to differentiating between potential (e.g., environmental envelopes) and actual distribution models (i.e., habitat suitability). Peterson and Soberón [34] recommend being careful when distinguishing between the environmental envelope models generated under the ENM approach and their existing or realized manifestations; the reconstructed distributions should necessarily consider the historical landscape transformations to obtain models of the species' actual distributions. In this study, we followed Grinell's [35] concept of habitat suitability to refer to species distribution models. Therefore, we differentiate between those models determined by only bioclimatic variables and topography (bioclimatic envelopes) and those determined by both land-use/land-cover suitability and bioclimatic envelopes (habitat suitability).

Conservation Areas
Historically, the establishment of natural protected areas has been based on opportunity ad-hoc criteria, which is different from systematic biodiversity knowledge [36,37]. However, currently, species distribution models are used as key input for identifying and selecting natural areas with maximum representation of species richness, meaning that computing methods have been developed for the selection of natural areas [38][39][40][41]. Such methods have focused on the following three quantitative approaches [42,43]: (a) identification of hot spots in species richness, (b) selection of areas with rare species, and (c) complementarity of natural areas.
Systematic conservation planning requires setting clear objectives from which it is possible to propose measurable and explicit conservation goals [44]. Conservation goals consist of explicitly quantifying a minimum amount of biodiversity elements, such as species distributions and vegetation types, which may be protected at local, regional, or continental scales [45].
Several natural protected areas exist in western Michoacán, México, covering approximately 328,000 ha, comprised mostly of the Zincuirán-Infiernillo biosphere reserve [4,46]. Several proposals to protect this region s natural areas also exist, including Terrestrial Priority Regions (TPR) [47], Important Areas for Bird Conservation (IABC) [48], and Priority Terrestrial Sites for Biodiversity Conservation (PTSC) [49] at the country scale ( Figure 2), and Priority Areas for Conservation in Michoacán (PACM) [50] and the Conservation Areas System in Michoacán (CASM) [51] at the state scale.

Statement of the Problem
Elaborating national schemes for prioritizing mammal conservation in México has made evident the need for generating reliable species distribution models at a higher spatial resolution [52]. Most of the studies on Michoacán's mammals have focused on estimating species richness [5,[53][54][55][56][57][58]. These studies suggest that the region is crucial for the connectivity of populations, but there is no spatial analysis of how biodiversity patterns are distributed across the region.
Western Michoacán is representative of western México's biodiversity and endemism. The Faja Transvolcánica and Sierra Madre del Sur are identified as relevant regions because of their high concentrations of endemic plants and animals [59]. There are various aspects that explain this region's research and biodiversity conservation challenges. For example, the region is representative of the entire biodiversity in western México's slope, including various species of carnivores. There is a lack of information about the regional patterns of biodiversity and the biosphere reserve Zicuirán-Infiernillo (ZIBR) is the only large natural protected area, as the remaining areas are not significant in terms of size and care. There are several proposed areas defined as priorities for conservation, with some of them comprising very large areas. Finally, the region has experienced significant changes in land-use/landcover, reducing and altering natural habitats for wildlife, including carnivore species.

Statement of the Problem
Elaborating national schemes for prioritizing mammal conservation in México has made evident the need for generating reliable species distribution models at a higher spatial resolution [52]. Most of the studies on Michoacán's mammals have focused on estimating species richness [5,[53][54][55][56][57][58]. These studies suggest that the region is crucial for the connectivity of populations, but there is no spatial analysis of how biodiversity patterns are distributed across the region.
Western Michoacán is representative of western México's biodiversity and endemism. The Faja Transvolcánica and Sierra Madre del Sur are identified as relevant regions because of their high concentrations of endemic plants and animals [59]. There are various aspects that explain this region's research and biodiversity conservation challenges. For example, the region is representative of the entire biodiversity in western México's slope, including various species of carnivores. There is a lack of information about the regional patterns of biodiversity and the biosphere reserve Zicuirán-Infiernillo (ZIBR) is the only This study aims to determine the level of spatial correspondence between natural protected areas (NPA) and the existing proposed areas for conservation, with actual hot spot areas identified based on the distribution of carnivore species. One research question relates to the ZIBR's role in protecting areas with high carnivore richness and the specific priority sites within the large regions proposed for conservation.
This study considers a group of carnivores as a biodiversity indicator, whose spatial patterns of richness and rareness would suggest the location of areas that should be considered as priority for conservation. This analysis or selection of areas also includes vulnerability conditions resulting from deforestation and fragmentation processes. We propose models of suitable habitats of carnivore species in western Michoacán, México, generated from bioclimatic envelopes and habitat/species associations, corresponding to western México. We also propose to obtain estimates of vulnerability in areas where carnivore species are distributed, due to deforestation and fragmentation processes (i.e., land-use/land-cover changes), and identify important areas for conservation, by applying systematic conservation planning tools (e.g., prioritization and complementarity), based on the distribution of suitable habitats of carnivore species and the potential transformation in the land cover.

Carnivore Species Occurrence Data
Occurrence data of carnivore species in western Michoacán were obtained from field work carried out during the years 2010-2014, which was part of a project named "Jaguar Conservation" (unpublished study). These data were collected by a photo-trap system. We installed between 15 and 30 photo-trap stations per sampling site. The stations were uniformly distributed according to a grid and the landscape characteristics, separated 1.5-2.5 km from each other. The digital cameras used were the Cuddeback brand (PO Box 10447, Green Bay, WI 54307) and the capture and ambush model. The cameras were active for 24 hrs for 1.5-2.5 months. The stations were situated on trails used mainly by wild felines and waterholes. Sampling was mainly carried out during the dry season, between February and July, without using lures. The location of stations is shown in Figure 1.
On the other hand, historical occurrence data of carnivore species in western México were acquired from accessing the following different databases: Global Biodiversity Information Facility (GBIF), Mammal Networked Information System (MaNIS), Unidad de Informática para la Biodiversidad UNAM (UNIBIO), and Comisión Nacional para el Conocimiento y Uso de la Biodiversidad (CONABIO). These databases are public and they were cleaned and integrated into a single database, avoiding duplicates and records without complete taxonomic and georeferenced data.
A preliminary species list of carnivores in Michoacán [5] and the recorded species by Nuñez (unpublished data) were used to verify their taxonomic status following the work of Ramirez-Pulido et al. [60,61] and supported by the work of Ceballos et al. [62,63], and CONABIO. Eleven species were finally included in this study, as shown in Table 1. Table A2 provides general information about the species' biological and ecological characteristics. This study used a total of 263 occurrence records that belonged to 11 carnivore species found in western Michoacán during the years 2010-2014. Nasua narica was the species with the highest number of records (81), while Bassariscus astutus, Conepatus leuconotus and Herpailurus yagouaroundi had the lowest with six records. On the other hand, a total of 2146 historic occurrence records, for the same 11 carnivore species but across western México, were integrated from querying public databases. The highest number of records belonged to Nasua narica (376), while the lowest was Leopardus wiedii with 77 records (Table 1).

Species Distribution Modeling
Species distribution models were generated by applying a consensus approach [32]. The software MODECO [64] was used to run the following different prediction algorithms: support vector machines (SVM), generalized linear models (GLM) and artificial neural networks (ANN). MaxEnt (maximum entropy) was also applied, but it was run independently because its software provides a more complete analysis.

Prediction Variables
The variables for predicting the species bioclimatic envelopes in western México included topography (elevation and aspect) and 19 bioclimatic variables (see Table 2 for a description) from the WorldClim project [65] (https://www.worldclim.org/ (accessed on 19 July 2020). This data set (21 variables) was used to generate preliminary bioclimatic envelope models, allowing the selection of more significant prediction variables. The selected variables differed in number and composition for each species, and they accounted for ≥95% of the total percentage. All of these environmental variables have a spatial resolution of about 1 × 1 km 2 in geographic coordinates. Table 2. Variables used in the elaboration of bioclimatic envelope models for carnivore species in western México. Selected variables per species represented ≥95% in prediction importance. A complete set of variables was used to select the important prediction variables for each species.

Bioclimatic Envelope Models
Maxent [66] was the algorithm applied to generate bioclimatic envelope models for selecting important prediction variables, using the complete set of climate (n = 19) and topography (n = 2) variables. We used MaxEnt's default settings with linear and quadratic feature types and 30% of the samples for accuracy assessment. The selected prediction variables for each species are shown in Table 2 and the response curves of the most important prediction variables appear in Figure A1. These selected variables were then used to run four different prediction algorithms per species, whose general settings are as follows: MaxEnt as previously described; support vector machines (C-SVC as SVM types and radial basis function kernel); generalized linear models (logit link function); artificial neural networks (prob-ANN trained based on backpropagation, 0.3 momentum and 0.1 learning rate). The algorithms applied were chosen because these were considered to be promising prediction tools [67,68] and are included in a group of algorithms in ModEco that run with presence data vs. backgrounds.
Four models for each species, for a total of forty-four bioclimatic envelope models, were constructed by applying the algorithms MaxEnt, support vector machine, generalized linear model and artificial neural networks. Afterwards, a consensus approach [32] was applied, so that the median was the representative statistic for combining four prediction models for each species.

Habitat Suitability Models
Considering that occurrence data for the carnivore species distributed in western Michoacán were recently collected (2010-2014), generalized linear models (GLM) were run to estimate habitat associations, a type of analysis that required pseudo absence data. These data for each carnivore species were generated from those bioclimatic envelope models obtained for western México. By applying this approach, it was expected to reduce the risk of generating pseudo absences in places that were potentially suitable for the species, in contrast with the random approach that is often suggested [69,70]. The use of preliminary habitat suitability maps as a spatial guide for generating pseudo absences is an efficient means to reduce the possibility for obtaining false negatives [70,71].
GLMs are an extension of linear models and allow us to manage data sets that lack normal error distributions and constant variances. Binary response variables, such as presence/absence data, show these characteristics, meaning that the association between carnivore species occurrence and habitat types in western Michoacán was described by the GLM run in the software R (R Core Team, 2013). The independent variables consisted of land-use/land-cover data obtained for the years 2013-2014 (INEGI-Serie V) [72,73] and the topographic position index [74]. INEGI's data consisted of vector coverage and a 1:250,000 scale, which was rasterized to 250 × 250 m 2 resolution. Therefore, the independent variables consisted of land-cover categorical variables; temperate forest, tropical semi-evergreen forest, tropical dry forest, introduced grassland, cropland, bare soils, water bodies and urban settlements. The only independent numeric variable was the topographic position index (TPI), which consists of differentiating landforms, based on elevation gradients; flat or near continuous slopes are shown by values near to zero, ridge or hill areas by large positive values, and the bottom of valleys or gullies by large negative values [74]. The significance of variables in the GLM was obtained by applying the Akaike information criterion.
Habitat suitability models for western Michoacán were obtained by adjusting the bioclimatic envelope models for western México by using the species/habitat association (probability) models obtained from the GLM analysis. A weighted sum of both types of models-the bioclimatic envelope and habitat associations-was conducted, assigning a weight of 2 to the former and 1 to the latter, and then the models were re-scaled to a 0-1 range of values. The results yielded continuous models that represented the estimates of habitat suitability for 11 carnivore species distributed in western Michoacán, determined not only by the climate and topography, but also by the land-cover conditions at the time when the species occurrence data were gathered in the field. We resampled the higher resolution models (i.e., species/habitat-type associations, 250 × 250 m) to the lower resolution models (i.e., bioclimatic envelope models, 1 × 1 km).

Model Validation
The evaluation of the models consisted of applying the partial ROC test (receiver operating characteristic) [75], which provides the partial and ratio values of the area under the curve (AUC), among other statistics. This approach has been proposed because it removes the emphasis on absence data, emphasizes the role of omission error when evaluating niche model predictivity and analyzes limited sectors of the ROC space that are not directly relevant [75]. The application of partial ROC was carried out by employing the NicheToolBox [76] (http://shiny.conabio.gob.mx:3838/nichetoolb2/ (accessed on 3 March 2021), using 30% of the training records with 500 iterations and calculating the AUC average [76].

Land-Use/Land-Cover Change
A land-use/land-cover modeling exercise was conducted to obtain the deforestation vulnerability estimates (i.e., deforestation likelihood) to be incorporated into the identification of priority areas for the conservation of carnivore species in western Michoacán. Vector coverages on a 1:250,000 scale, generated by the Mexican government (Instituto Nacional de Estadística y Geografía, INEGI) and known as Serie II [77] and Serie V [78], were used as two different dates to determine the changes in the region's land use/land cover. Indeed, Serie II was created using Landsat TM images gathered in 1993, while Serie V was elaborated based on SPOT images gathered in 2012-13.
The Land Change Modeller module in the IDRISI Selva software was used to generate probability transition models related to the following three main submodels: deforestation (change from conserved temperate and tropical forests to cropland and introduced grassland); regeneration (change from cropland and introduced grassland to conserved and secondary temperate and tropical forests); disturbance (change from conserved temperate and tropical forests to secondary forests). This study focused on the transition probability estimates that corresponded to the deforestation submodel, because these estimates of deforestation vulnerability were included as key variables in the identification of priority areas for the conservation of 11 carnivore species in western Michoacán. Given the nature of the transition probability output, it was necessary to calculate the inverse and to rescale values from a 0-<1.0 range to 0-1.0, so that the model represented the probability of persistence of the different land-cover/land-use types. For instance, high values represented areas with low probability of change (high persistence) and low values represented high deforestation or low persistence probabilities.

Identification of Priority Areas for Conservation
Individual habitat suitability models by species were summated (1 model per species, 11 models in total) to identify areas with a high concentration of suitable species habitats or richness estimates. Continuous models were converted to binary models (i.e., presence/absence) based on the application of a prediction threshold. In this case, such a threshold consisted of the sensibility/specificity maximization criteria [79] (see Table A3).
The software Zonation v1.0 [80] was used to prioritize the analysis of conservation areas. Zonation is based on the complementarity of areas and serves as a balance to rank conservation priorities across the landscape through the iterative elimination of cells or planning units to direct at a minimum aggregated loss of conservation values. The cell elimination order depends on an elimination rule, which can be selected according to conservation objectives [80][81][82].
The suitable habitat models (continuous) for each species and the inverse deforestation model were used as inputs, while the additive benefit function (ABF) and the core area zonification (CAZ) were applied as removal rules. In general terms, the ABF rule exclusively prioritizes exclusively high richness areas, whereas the CAZ rule determines important high richness areas, while at the same time identifies areas with unique biodiversity elements (e.g., rareness areas).
Different weights were assigned to each species, depending on their conservation status (see Table 3) and different modeling parameters were used. Smoothing and boundary length penalty (BLP) are two aggregation methods that generate solutions more compactly. The former conserves well-connected areas and the latter is the most frequently used method. Additionally, BLP does not consider situations such as species fragmentation, but rather impacts the boundary length to produce more compact conservation areas. Table 3. Weights assigned to each species in the prioritization of conservation areas. A higher weight represents the highest level of concern considered by the prioritization algorithm (zonation).

SPECIES Conservation Status (NOM059) Weight
Canis Different conservation priorities are classified in the resulting conservation models, allowing for the specification of different thresholds in the percentage of landscape. For instance, 15% of the landscape is sometimes chosen as a threshold to describe the conservation areas identified.

Comparing Identified and Already Proposed Priority Areas
The conservation areas identified in this study were compared (overlap) with other existing areas that have been proposed for conservation. Particular attention is paid to a recent version of the priority terrestrial conservation areas (RTP) proposed by CONABIO [49] (Priority Terrestrial Sites for Biodiversity Conservation (PTSC); see Figure 2), which consist of 256 km 2 hexagons distributed across the country and with a ranked priority level (medium, high and very high).
A flow diagram, as shown in Figure 3, helps us to explain and clarify the different steps taken in the applied methods.

Variable Importance
The most important prediction variables varied among species (Table 2); 13 variables predicted 5-8 species and 8 variables predicted <5 species (see Table 2). The b10 variable (mean temperature of warmest quarter) was not important for any species, while b2 (mean diurnal range (Mean of monthly values (max temp-min temp)) was important for eight species. Table 4 shows a summary of the estimated association probability between habitat types and each species' presence (GLM analysis) in western Michoacán. All models had a partial ratio AUC value above 1.2, which indicates good-fitting models. The results indicate that the lowest probability areas for all the species correspond to bare soil and cropland areas. Some species also showed low occurrence probabilities in other land-cover types; for instance, Urocyon cinereoargenteus, Herpailurus yagouaroundi and Leopardus pardalis in temperate forests, and Conepatus leuconotus and Bassariscus astutus in tropical semi-evergreen forests ( Figure A2).

Associations Species/Habitat Types
The high occurrence probability areas for Urocyon cinereoargenteus, Leopardus pardalis and Nasua narica consisted of almost exclusively tropical semi-evergreen forests for the species Canis latrans, Leopardus wiedii and Procyon lotor. The highest presence probability areas were areas with tropical semi-evergreen forests and temperate forests, whereas Panthera onca was highly associated with tropical semi-evergreen forests and tropical dry forests. Temperate forests are areas where Conepatus leuconotus and Bassariscus astutus had the highest occurrence probability, while tropical dry forests were the land-cover types with the highest probability for Herpailurus yagouaroundi and Puma concolor (see Figure 4).

Habitat Suitability Models
The habitat suitability models for each species were obtained by the weighted sum of the bioclimatic envelope models for western México and habitat/species associations models for western Michoacán. Seven species (Urocyon cinereoargenteus, Herpailurus yagouaroundi, Leopardus pardalis, Panthera onca, Conepatus leuconotus, Nasua narica and Procyon

Habitat Suitability Models
The habitat suitability models for each species were obtained by the weighted sum of the bioclimatic envelope models for western México and habitat/species associations models for western Michoacán. Seven species (Urocyon cinereoargenteus, Herpailurus yagouaroundi, Leopardus pardalis, Panthera onca, Conepatus leuconotus, Nasua narica and Procyon lotor) showed suitable habitats almost exclusively in the Michoacán coast region, corresponding to the Costas del Sur physiographic province, including counties ("municipios") such as Aquila, Coahuayana, Chinicuila, Lázaro Cárdenas, Coalcomán and Arteaga, zones covered by tropical dry forest and tropical semi-evergreen forest ( Figure A3).
The suitable habitats of these seven species are partially included in small natural protected areas (NPA), such as "Lagunas Costeras y Serranias Aledañas a la Costa Norte de Michoacán", "Playa Maruata y Colola" and "Playa Mexiquillo" sanctuaries; the Terrestrial Priority Region (TPR) "Sierra de Coalcomán"; the Important Area for Bird Conservation (IABC) "Coalcomán-Pomaro"; polygon 6 in the Priority Areas for Conservation in Michoacán (PACM); and areas 01 and 02 in the Conservation Areas System in Michoacán (CASM). Nevertheless, recent records exist for most species within the "Zicuirán Infiernillo" biosphere reserve (ZIBR); only Panthera onca y Puma concolor showed high occurrence probability within this area ( Figure A3).
The other four species (Canis latrans, Leopardus wiedii, Puma concolor and Bassariscus astutus) showed suitable habitats, expanding from the coastal region towards inland Michoacán, including the Costas del Sur and Cordillera Costera del Sur physiographic provinces ( Figure A3). Temperate forest is the main vegetation type in this latter province.
Bassariscus astutus is the species that showed the most restricted suitable habitats in western Michoacán; small areas with high occurrence probability are found in the Cordillera Costera del Sur province ( Figure A3), within the Coalcomán and Aguililla counties. On the other hand, Puma concolor was the species with the most widespread suitable habitats. It is important to mention that the Depresión del Tepalcatepec province did not include species with high occurrence probability because it is largely covered by cropland.

Land-Use/Land-Cover Change
The highest values in the inverted deforestation model correspond to temperate forest areas, while the lower values are associated with tropical forest areas ( Figure 4A). Apparently, temperate forests in western Michoacán have been transformed at lower rates than tropical formations within the time-period 1993-2014. Inference can be made that tropical forests have higher vulnerability to the deforestation processes that occur in western Michoacán.

Identification of Conservation Priority Areas
A species richness model was created by adding together the individual binary habitat suitability models for western Michoacán. This richness model shows that more species (e.g., 10-11 species) co-occur in the Costas del Sur physiographic province, involving the counties of Coahuayana, Aquila and Lázaro Cárdenas. However, there are small high richness areas scattered across the region. The Cordillera Costera del Sur province included 7-9 richness values, while the Depresión de Tepalcatepec province included mainly low richness values (1-3 species; Figure 4B).
One small NPA ("Lagunas Costeras y Serranías Aledañas de la Costa Norte de Michoacán") and two sanctuaries demonstrated spatial correspondence with the highest richness areas within the Costas del Sur province. The ZIBR, located in the Cordillera Costera del Sur province, also included high species richness areas (7-9 species).
Currently, the proposed priority areas and regions for conservation have not been actualized due to issues of size (exceptionally large or exceedingly small). Thus, important portions of the TPR, including "Sierra de Coalcomán" and the IABC "Coalcomán-Pomaro", as well as polygon 6 in PACM, areas 01 y 02 in CASM, the sanctuaries "Playa de Maruata y Colola" y "Playa Mexiquillo", and a small area of NPA "Lagunas Costeras y Serranías Aledañas de la Costa Norte de Michoacán", all show spatial correspondence with the high-priority conservation areas identified by the 18 priority models. The southern portion of ZIBR also included high-priority areas. To a lesser extent, other areas also contain identified priority areas, such as the NPA "Pico de Tancíntaro", "El Jabalí", "Volcán Nevado de Colima", "Las Huertas", "Bosque Mesófilo Nevado de Colima", and the TPR and IABC "Tancíntaro".
In general, each of the generated 18 models are located in the most important conservation areas within the western Michoacán coastal region. However, the models generated by applying the CAZ (core areas zoning) identified small highest priority areas towards the inland and northern sections of the study area ( Figure A4).
The models generated in this study allow for flexibility in selecting different landscape area threshold percentages for conservation. For example, as shown in Figure A4, a 15% threshold percentage corresponds to the area in yellow, pink, cherry and red; a 5% threshold percentage corresponds to the areas in pink, cherry, and red, and finally, a 2% threshold percentage is represented by the areas in red ( Figure A3).
The two aggregation methods in the study resulted in small homogeneity differences in the models. In general, the BLP method tends to generalize and slightly increase the size of priority areas, reducing fragmentation, independently of the elimination rule. On the other hand, the smoothing method does not seem to significantly affect the priority areas' spatial configuration ( Figure A5).
The weights of species, depending on their conservation status, seem to have an apparent effect on the priority areas' connection. Assigning a species' weight generates a continuum in the coastal region, independently of the aggregation method and reduces the size of priority areas (the CAZ rule) located north of the study area ( Figure A5).
The inclusion of the deforestation model does not significantly affect the identification of priority areas. The only apparent effect is the presence of a set of priority areas located in the northeastern study area (Figure A4a-c).
Despite these variations among the models, we found that all the models show a high concentration of high priority areas in the coastal region, which coincides with the species richness model, meaning, in other words, similar spatial patterns. These areas correspond to tropical forest formations, which show high vulnerability to degradation and fragmentation processes that occur in the region.
[49] identified 54 conservation priority hexagons (256 km 2 each) within the study area, with 10 being considered as high priority and 6 as extreme priority. Eight hexagons of high and extreme priority are distributed along the Michoacán coast, in correspondence with the high conservation priority areas identified by this study (2-15% percentage of landscape) ( Figure 5).

Modeling Scales
This study uses a two-scale approach for modeling the habitat suitability of carnivore species in western Michoacán. The large scale corresponded to western México, whose modeling was based on using climate prediction variables following the hypothesis that

Modeling Scales
This study uses a two-scale approach for modeling the habitat suitability of carnivore species in western Michoacán. The large scale corresponded to western México, whose modeling was based on using climate prediction variables following the hypothesis that species distributions are primarily determined by this type of environmental factors at regional scales [84][85][86]. These models represent bioclimatic envelope models because climate and topography variables are used to build fundamental species ecological niches, which are projected to a geographic space [87]. The second scale consisted of western Michoacán; an area centered in the larger region, which is 1/5 of the size of the latter. This smaller area allowed for the collection of the current carnivore species occurrence data by means of applying a photo camera sampling method, and to perform land-use/land-cover change modeling, necessary data for generating actual species distribution models and applying analytical procedures for conservation planning and NPA network design [87,88].

Suitable Species Habitat Models
Previous studies have attempted to generate species distribution models from refining environmental envelope models in the field of conservation planning [87,[89][90][91][92][93]; these have consisted of using a digitally documented suitable habitat map to intersect an environmental envelope model. In this study, this masking procedure is replaced by using a habitat/species associations model to restrict a bioclimatic envelope model. This model combination makes it possible to include the landscape's habitat transformations in the climate/topography predicted model at local scales [94]. Indeed, the use of species/habitattype associations, particularly on mammal species, has been supported in research on the use and selection of habitats that are widely abundant in ecological studies. Such studies apply statistical tests (e.g., logistic regression) to determine the levels of associations between species occurrence and habitat types [95][96][97]. However, these studies are conducted at local scales, limiting their application to conservation goals [98].
The GLM results showed that seven species had a significant variable that predicted species occurrences. This species occurrence probability is an objective and quantitative means to determine the levels of habitat suitability for each of the carnivore species. Four species that did not have significant predicting variables had the smallest occurrence sample size, which would be the first factor to be considered in explaining the GLM model performance [99][100][101].
The identification of suitable habitats where species are likely to occur has demonstrated to be a powerful conservation tool [102][103][104][105][106]. Moreover, this information is considered as a measurement of biodiversity units that serve as inputs for applying algorithms and geographic information systems (GIS) to identify priority conservation areas [92,107]. The efforts made towards prioritizing conservation areas in México have revealed the need to generate species distribution models at more reliable and finer scales [52].
The geographic description of most of this study's carnivore species' actual distribution models show similar spatial patterns as those described in the previously cited literature. On the other hand, this study's habitat suitability models contrast with those distribution models proposed in studies with a national scope, such as the "Project DS006 Modeling mammal species distributions in México, a GAP approach" [108]; the former are significantly more restricted than the latter. For instance, even though Bassariscus astutus is reported to occur in most of Michoacán [109], and some historic records do exist [110,111], very few individuals have been recorded in the last two decades [58,112]. The Bassariscus astutus model showed the most restricted distribution among the 11 species (see Figure A3). The restricted nature of the previously mentioned habitat suitability models is closely related to the land-cover/land-use changes observed across the region [113]. A significant number of historic records related to Michoacán's carnivore species are currently located in sites with agricultural and human settlement land uses [106].

Selection of Conservation Priority Areas
Different characteristics of biodiversity elements and the application of modeling tools and parameters affect the results in prioritization exercises for the selection and complementarity of conservation areas [81,114]. This study did not generate drastic changes when varying the zonation modeling parameters to identify conservation priority areas. There was a consensus for locating priority areas mainly in the Costas del Sur province and Cordillera Costera del Sur to a lesser extent. One important difference consisted of locating priority areas in the study area's northern portion; when the prioritization rule was CAZ and with no species conservation status weighting, high priority areas appeared in "Pico de Tancíntaro" and "Volcán Nevado de Colima" NPA ( Figure A4a-c).
There are two sanctuaries, "Playa de Maruata y Colola" and "Playa Mexiquillo", which contain the highest priority importance areas (95-100%). These sanctuaries consist of very reduced areas (220 and 74 ha, respectively), with the original goal of protecting three species of marine turtles (Lepidochelys olivacea, Dermochelys coriacea, and Chelonia agassizii) [115,116]. Although these sanctuaries are marine NPAs, they can serve as a geographic reference for inland areas for biodiversity conservation.
The ZIBR attracts particular interest because of its location, biological diversity, endemism, size (265,000 ha) and legal status. The ZIBR is one of the largest biosphere reserves in México and was created in 2007. This NPA is occupied by high species richness (7-9 species), but has not been given the highest conservation importance priority status. However, this study's results were obtained from a reduced number of species, and it is likely that other wider taxonomic groups would produce higher prioritization estimates within the ZIBR, based on the maximization of its biodiversity. This study used a group of umbrella species as indicators of biodiversity and could be considered as a preliminary approach for evaluating conservation priority areas [6,7].
The deforestation model in this study was inverted so that it represented the probability of natural areas to persist with their native land-cover type ( Figure 4). Therefore, the prioritization procedure selected, with high prioritization importance, conserved areas of tropical dry forests, tropical semi-evergreen forests, and temperate forests. However, the prioritization procedure without the deforestation model resulted in the same selection of areas. This could be since both high species richness areas coincide with the areas of land-cover types with a high probability of persistence, given the land-cover/land-use changes that occurred between 1993 and 2014.
A comprehensive Marxan proposed prioritization study [49] was compared to this study's results. Such a comprehensive study included the optimization analysis of 1450 biodiversity elements across México, including critical vegetation types, plant species at risk, the richness of plants and vertebrate species, resident birds, reptiles, amphibia, and mammals, along with 19 layers of threatening factors [49]. Despite the different scopes applied in each study, there were coincidences detected in the location of priority sites and areas. Our study included half of CONABIO's high-priority conservation sites, which were included within the study area ( Figure 5). In other words, both studies identify the Costas del Sur province in western Michoacán as an important region for the conservation of biodiversity. This kind of concurrence supports the idea that a small group of umbrella species may be used as biodiversity indicators, as biodiversity units are useful for prioritizing biodiversity conservation areas.

Study Limitations and Future Research
The main limitation of this study is the small sample sizes, particularly for the four carnivore species that occur in western Michoacán, despite the fact that photo-trap sampling was carried out for several months during the years 2010-2014. Such a small sample size resulted in the lack of robustness of some of the GLM models.
This study proposes the identification of conservation areas for 11 carnivore species distributed in western Michoacán, México, which represent the basis for further research on the following: (1) suitable sites to conduct species monitoring across the region; (2) conservation of Michoacán's coastal region, mainly corresponding to the Aquila and Lázaro Cárdenas counties; (3) issues about the extent of the current proposed regions for biodiversity conservation, given their larger sizes and location; (4) connectivity of areas identified as a priority for biodiversity conservation; and (5) comparison of prioritization and complementarity exercises based on analyzing multiple taxonomic groups, including animal and plant taxa of different taxonomic hierarchy and functional groups, such as key and umbrella species.

Conclusions
Bioclimatic envelope models in western México and habitat/species associations in the western Michoacán state were combined to obtain habitat suitability models for 11 carnivore species distributed in western Michoacán. The habitat suitability models were considered as biodiversity units for identifying conservation areas for 11 species of carnivores. Natural protected areas (NPA) and the existing proposed conservation areas were spatially related to our identified conservation areas, highlighting the prioritization of Michoacán's coastal region for biodiversity conservation. This study contributes a new method for refining environmental envelope models into actual habitat suitability models by including the actual landscape configuration, which has historically been severely modified. By using 11 carnivore species as biodiversity surrogates, we identified the priority conservation areas that have not been included in the current NPA system and have not been specifically selected by the existing conservation proposals. The high-priority conservation areas identified in this study spatially corresponded to the high-priority areas proposed by a national study that used multiple taxa and prediction variables. This study signifies a flexible basis from which future studies on planning and designing a network of natural protected areas can be carried out in the region.  Medium-size canid with a 500-600 mm total length and 3-5 kg weight. Opportunistic feeding habits, consuming rodents, lagomorphs, fruit, and insects. Home range varies from 1 to 8 km 2 , depending on the season, population density and habitat quality. They live mainly in forests and shrubland, but can also live in every vegetation type, from sea level to 3500 masl. Jaguarundi As a felid, it is a small animal (888-1372 mm total length and 3.5-9 kg weight). It can be found in tropical evergreen and dry forests, but can also live in mangroves, cloud forests, desert scrubs, and occasionally in temperate forests, from sea level to 2000 masl. It is well adapted to arboreal environments and its food habits consist of a wide range of prey, including invertebrates, reptiles, birds, and small mammals. Its activity area is 88.8-99.9 km 2 for males and 20.1 km 2 for females.

Leopardus pardalis (Linnaeus, 1758)
Ocelot This is a medium-size felid (920-1367 mm total length and 6-15 kg weight). Lives mainly in tropical evergreen forests, tropical forests and mangroves, but also in cloud forests and desert scrubs occasionally. It can be found from sea level to 2000 masl. Prefers dense land-cover habitats and feeds on small and medium-size rodents, but also on invertebrates, reptiles, birds and other mammals, including the temazate deer. Its activity area varies from 3.5 to 17.7 km 2 for males and 0.7-14.6 km 2 for females.

Ringtail
This is a medium-size carnivore with a 616-811 mm total length and 870-1100 g weight. Lives in desert scrubs, temperate forests, tropical semi-desert, scrubland, and even urban parks, from sea level to 2880 masl. It is an omnivore and feeds mainly on small mammals, insects, fruits, birds, reptiles and occasionally nectar. They are solitary animals with nocturnal habits. Its activity area is variable and depends on habitat type, season, and sex.

Whitenosed coati
Procyonids with a size similar to a medium size dog (850-1.340 mm total length and 4-6 kg weight). It is a highly social species, gathering in groups of up to 20 individuals, made up of almost exclusively adult females and juveniles. It is an omnivore collector animal that feeds mainly on fruits and invertebrates such as insects, myriapods, arachnids, crustaceans, and annelids. To a lesser extent, they hunt terrestrial vertebrates, such as rodents, amphibians and reptiles. They can be found mainly in tropical evergreen and tropical dry forests, but they are also common in temperate forests and desert scrubs, from sea level to 2900 masl.

Raccoon
Medium-size procyonid with a 603-950 mm total length and 3-9 kg weight. It has twilight and nocturnal habits. They have adapted to a wide variety of habitat types with permanent water bodies. They are more abundant in tropical evergreen and tropical dry forests, mangroves, and areas with aquatic and semi-aquatic vegetation associated with marshes, swamps, and wetlands, from the sea level to almost 3000 masl. It is an omnivore and feeds on a wide variety of plant and animal foods. It is considered as a solitary species and its activity area widely varies from 0.2 to 4946 ha, although a typical annual average is 65 ha for males and 39 ha for females.  Figure A1. Response curves of the most important prediction variables, when using only the corresponding variable for creating a maxent model.  A land-cover map is included to visualize its spatial configuration, along with the different habitat association models for each species, for the whole study area.