Introduction

In Mexico, the bighorn sheep is an iconic species naturally distributed across the Baja California peninsula and the state of Sonora; populations maintained under semi-captivity exist in the states of Chihuahua and Coahuila (^{Sandoval and Espinoza 2001}; ^{Sánchez 2005}; ^{Uranga and Valdez 2011}). In these areas, this species contributes to the ecosystem production chain by consuming plant resources and as a source of energy for large predators and scavengers (^{Monson 1980}). In addition, it represents a key element for the conservation of its habitat, mainly due to its importance as a game species, since the prices of permits for legal hunting are high, ranging between $45,000 and $60,000 USD (^{Lee 2011}; ^{Ruiz 2014}).

The sustainable use of the species requires the development of a comprehensive and adaptive management plan, designed based on scientific information addressing the needs and interests of the managers of this resource (^{Walters 1986}; ^{Mandujano 1994}). However, the dimensions of the bighorn sheep habitat are far too extensive to propose feasible management strategies, thus making it necessary to identify the sites of higher probability of occurrence of individuals of this species to focus conservation efforts on them (^{Lara et al. 2001}; ^{Escobar et al. 2015}).

Studies on habitat preferences of bighorn sheep in the Baja California peninsula and Sonora, using use-availability methods, methodologies based on generalized linear models, and satellite remote sensing techniques, have shown that the species does not use habitat resources according to their abundance, as some are preferred while others are clearly overlooked. Elevation, slope, ruggedness and plant cover are the characteristics of the habitat that determine the sites preferred by the species (^{López et al. 1999}; ^{Lara et al. 2001}; ^{Guerrero et al. 2003}; ^{Escobar et al. 2015}). Hence, the potential distribution of the species in a particular locality can be determined using an ecological niche model, as it is a representation of the environmental suitability for the establishment of a particular species, calculated from field observations (presence/absence) and the availability of a series of environmental variables (climate, topography, vegetation, among others) that act as predictors (^{Soberón and Peterson 2005}; ^{Ferrier and Guisan 2006}; ^{Naoki et al. 2006}).

Although several methods to develop niche models are available, the one based on a statistical approximation named maximum entropy and run by the Maxent software (^{Phillips 2013}) offers some unique advantages relative to others. It works with presence data only, while other models require presence and absence data; it shows a consistently good performance versus other similar methods, in differentiating between sites that are either suitable or unsuitable for the establishment of a particular species, especially with a small number of samples (^{Elith et al. 2006}; ^{Phillips et al. 2006}; ^{Hernández et al. 2008}); also, it is freely available.

In the Baja California peninsula, the bighorn sheep is distributed from Sierra de Juarez on the border with the United States, to the Sierra de San Juan de la Costa adjacent to Bahía de La Paz (^{Sánchez 2005}; ^{Ruiz 2014}). Sierra de Juarez is listed as a high-priority terrestrial region for biodiversity conservation in Mexico due to its physical and biotic characteristics; for being a key area for the functional connectivity for wildlife populations across the landscape, given its location in an international border; and because this mountain range is associated with six hydrological basins that flow to the main cities in the state of Baja California (^{CONABIO 1998}; ^{Arriaga et al. 2000}; ^{Buchalski et al. 2015}).

The objective of this work was to construct a model of potential distribution for the bighorn sheep in Sierra de Juarez, to facilitate the identification of the areas of greatest importance for the species across its range. This information is key to design a management plan aimed at ensure the conservation of the species and its habitat.

Materials and Methods

Study Area. The habitat available for bighorn sheep in Sierra de Juarez is located in the area of ravines and canyons of the mountain range, which starts in the border with the United States of America in La Rumorosa and stretches southward to Valle de la Trinidad (31.300° N, -115.370° W and 32.630° N, -116.035° W), covering an area of approximately 145,105 ha. Land tenure is shared by the ejidos Emiliano Zapata, Cordillera Molina, Mission Santa Catarina, Jamau and 16 de Septiembre (^{RAN 2016}; Figure 1).

Figure 1 Study area (black dots) showing 2012 aerial census records (yellow) and sampling sites (white): 1) Cañón de Los Álamos; 2) Cañón de los Llanos; 3) Cañón del Ranchito; 4) Cerros Boludos; 5) Cerros de la Ponderosa; 6) La Rosa; 7) Cañón del Tajo; 8) Cerros de los Laureles del desierto; 9) Cañón de Guadalupe; 10) Cerros del Ejido Manatou; 11) Cañón del Alamar; 12) Cañón del Palomar. Ejidos located in the study area are: Emiliano Zapata (dark blue), Cordillera Molina (red), Misión de Santa Catarina (green), Jamau (light blue), 16 de septiembre (orange). 

Physiographically, the area is formed by an asymmetric mountain range that is steep-sloped to the east and with gentle slopes to the west. The predominant orientation in the area is to the east, with slopes between 0 % and 72 %, dominated by those lower than 20 %, and most ruggedness values range between 0 and 50 m, with a maximum of 228 m (^{INEGI 2013}).

Six vegetation types cover the study area: pine forest in the highest peaks, chaparral scrubs growing along the western slope, microphyll desert scrub to the east, rosette-shaped desert scrub patches, natural palm-tree patches in the bottom of canyons with water availability all year round, and gallery vegetation along riverbanks (^{INEGI 2013}).

Sierra de Juárez shows five types of climates in ravines: very arid, semi-warm [BWh(x’)]; very arid, temperate with rainfall throughout the year [BWk(x’)]; very arid, temperate with winter rainfall [BWks]; temperate [Cs]; and sub-humid semi-cold [CB’s] (^{García and CONABIO 1998}).

Construction of the database on bighorn sheep evidence. Twelve sites were selected for being recognized by Sierra de Juárez local inhabitants as bighorn sheep distribution areas between La Rumorosa and Cañón del Palomar, which were traveled from January to June 2016 (Figure 1). In each site, one transect was established, the length of which varied according to accessibility, ranging between 2 km and 24 km long. The presence of bighorn sheep was recorded and geo-referenced in each individual transect through indirect evidence (excreta and footprints), discarding those records separated from each other by less than 30 m, to avoid over-fitting.

The database was supplemented with records of the aerial census conducted by the San Diego Zoo in 2012, where four persons flew over Sierra de Juárez for 10 hours on board a Hughes 500™ helicopter with doors removed, from La Rumorosa to Valle de la Trinidad (^{García 2014}; Figure 1).

Environmental Characteristics of the Study Area. We used the free software QGIS 2.8.6 (^{QGIS Development Team 2016}) for geospatial data handling and processing. This software processed the information of a digital terrain model corresponding to the study area, with a spatial resolution of 30 x 30 m pixel (^{INEGI 2013}). We calculated the area covered by particular orientation, slope and ruggedness index data, the latter estimated with the method by ^{Riley et al. (1999)}. The vector information regarding the land-use and vegetation layer of the study area (^{INEGI 2013}) was rasterized to a 30 x 30 m pixel resolution, preserving the vegetation type field.

Plant cover was estimated from Landsat 8 OLI/TIRS satellite images (039 Path, 038 Row) for 21 April 2016, with a 30 x 30 m pixel resolution (^{USGS 2016}), using the calculated value of a Normalized Difference Vegetation Index (NDVI), where NDVI = [(Band 5 - Band 4) / (Band 5 + Band 4)] (^{Purevdorj et al 1998}; ^{Pettoreli 2013}). With this index, a value between -1 to 1 is assigned to each pixel in the matrix of the study area, where figures below 0.1 correspond to rocky outcrops with sparse vegetation, those between 0.2 and 0.5 to scrubs and grasslands, and those from 0.5 to 1 to forests and cropland (^{USGS 2015}).

The types of climates in the study area were updated using high-resolution climate areas for Mexico (^{Cuervo et al. 2014}). These were interpolated to a 30 x 30 m pixel resolution with the thin-plate interpolation technique using the software ANUSPLIN 4.3 (^{Hutchinson 2006}). The parameters of climate surface areas were combined to produced updated types of climates, by entering the respective data into the PATN.EXE 3.1.2.0 program (^{PATN Development Team 2013}), and conducting a cluster analysis with mixed variables using the Gower’s similarity coefficient. The resulting climate types were classified according to the nomenclature used by ^{García and CONABIO (1998)}.

Potential Distribution Modelling. An ecological niche model was constructed using the Maximum Entropy algorithm with the software Maxent version 3.3.3k (^{Phillips 2013}). The database of bighorn sheep evidence records were entered, along with the GIS information layers for orientation, slope, ruggedness index, climate types, vegetation types and NDVI.

The algorithm implementation was used to construct a model of average values, performing 20 replicates with 1,000 iterations each, where 80 % of occurrence records were used for constructing the model and 20 % for the test analysis, by activating the Random seed option in order to use different occurrence localities in each replicate for the construction of the models. A regularization factor β of 1 was used, disabling the Extrapolate and Do clamping options to avoid data overfitting. The selected output was Logistic, which classifies pixel values according to a zero-to-one scale, and interprets the scale as the probability of occurrence of the species (^{Phillips et al. 2006}).

The predictive accuracy of the model was determined by calculating the area under the curve (AUC) of the receiver operating characteristic (ROC). The fraction of sites misclassified as absences (omission errors) was determined by calculating the omission rate and the average predicted area (^{Phillips et al. 2006}). The model obtained was reclassified to produce a binary (presence-absence) map, where a cut-off threshold was applied for presence of training data points (^{Pearson et al. 2007}). Model variables were assessed through a Jacknife test, which yielded the information and percentage contributed by each variable. In addition, response curves were plotted to relate bighorn sheep presence data with each environmental variable used to elaborate the maximum-entropy model (^{Phillips et al. 2006}).

Results

A total of 196 bighorn sheep records of were obtained, with 116 corresponding to indirect evidence of its presence found during transect walks carried out from January to June 2016, and 80 to direct observations of individuals during the 2012 aerial census (Figure 1). The potential distribution range of the species covers 49,844 ha (34.3 %) of Sierra de Juárez, an area that is concentrated to the northern region of the cordillera and drops steadily toward the south of the mountain range, being finally scattered in patches in areas near Sierra de Las Tinajas (Figure 2).

Figure 2 Potential distribution of bighorn sheep in Sierra de Juárez (blue). 

The omission rate and predicted area showed that omission in test points did not fully match the predicted omission rate (Figure 3a). For its part, the area under the curve (AUC) of the receiver operating characteristic (ROC) yielded a model with a value of 0.877 and a standard deviation of 0.011 (Figure 3b). Taken together, the variables plant cover (26.7 %), vegetation type (24.1 %), and ruggedness (22.2 %) contributed with 73 % of the information necessary for the construct of the model (Table 1).

Figure 3 Plots of average omission rate and predicted area (a) and average receiver operating characteristic sensitivity (b), obtained from the modeling carried out for Ovis canadensis at Sierra de Juárez. In Figure a, the red line marks the average area; the blue strip, the standard deviation of the average area; the light blue line, the average omission rate; the yellow strip, the standard deviation of the omission rate; and the diagonal black line, the predicted omission rate. In Figure b, the red line marks the average; the blue strip, the standard deviation; and the black diagonal, the expected line if the model is not better than a random one. 

The response curves generated by the potential distribution model show that the highest probability of bighorn sheep encounter occurs in sites with the following characteristics: very dry semi-warm [BWh(x’)], very arid temperate [BWk(x’)] and very arid temperate [BWks] climates; gallery and natural palm-tree vegetation; NDVI from 0.05 to 0.07; orientation from 0 to 160°; slope from 0 to 65 %; and ruggedness from 35 to 160 m (Figure 4).

Figure 4 Response curves of bighorn sheep presence data vs each environmental variable used to elaborate the maximum-entropy model. The red line indicates the average, and the blue strip, the standard deviation. The vegetation types are: BP= pine forest; C= chaparral scrub; MDM= microphyll desert scrub; MDR= rosette-shaped desert scrub; PN= natural palm-tree forest; VG= gallery vegetation. 

Discussion

The potential distribution range of the bighorn sheep covers regions of northern and central Sierra de Juarez, as according to the niche model these have environmental conditions that provide shelter and resources to the species (^{Soberón and Peterson 2005}). However, the modeling reported here should be regarded as a partial instrument, as the closeness between the localities where the species was recorded led to an overfitting of the data caused by the dependency between data points (Figure 3a). The end result is that the distribution area identified for the bighorn sheep was smaller than the true area. However, the model functions properly to differentiate between sites with high vs. low probability of bighorn sheep occurrence (Figure 3b; ^{Phillips et al. 2006}).

In the model, vegetation was the primary component of the habitat that explained most of the bighorn sheep distribution in Sierra de Juarez (Table 1). It provides food resources for the species, is an indicator of water availability, and is closely related to its anti-predatory strategy (^{Monson 1980}; ^{Wilson et al. 1980}; ^{Holl 1982}; ^{Álvarez et al. 2009}; ^{Escobar et al. 2015}). The type of vegetation that predominates in regions identified as distribution areas of bighorn sheep is microphyll desert scrub (^{INEGI 2013}). This includes shrubs and herbs of importance as food sources for the bighorn sheep, such as Ephedra nevadensis, Eriodyction angustifolium, Senecio spartioides, Ceanothus gregii, Eriogonum fasciculatum, Simmondsia chinensis, Ditaxis lanceolata, Galium wigginsii, Cardiospermum corindum, Bromus ciliatus, Erioneuron pulchellum, Ferocactus sp. and Agave sp. (^{Reyes, 1976}; ^{Sánchez, 1978}). In addition, this plant community has a low vegetation index (NDVI = 0.05 to 0.07), providing more visibility or lesser visual obstruction for the detection of predators (^{Monson 1980}; ^{Wilson et al. 1980}; ^{Holl 1982}; ^{Álvarez et al. 2009}; ^{Escobar et al. 2015}).

The response of the distribution of organisms to the slope turned out to be indifferent (Figure 4), in contrast with the findings reported by ^{López et al. (1999)}, and ^{Guerrero et al. (2003)}, who concluded that the desert sheep evades sites with slopes below than 60 %. This situation is attributed to the fact that the suitability values of topography vary according to the scale used (^{Divine et al. 1996}). Differences may stem from the fact in the studies just mentioned, the slope was measured in the field, so that the true value was recorded on site. By contrast, in our study slope was calculated from a digital terrain model with a 30 x 30 m resolution; this approach involves intrinsic inaccuracies related to the original source of the data used to elaborate it, and to the interpolation of elevations that do not represent the true altitude of the terrain for a particular location. In addition, the accuracy of the elevation models is higher in flat versus rugged areas (^{Riley et al. 1999}). It has been recorded that slope is a structural element of minor importance in the elaboration of models of the potential distribution of wild ungulates that inhabit mountainous areas; in contrast, vegetation, climate and the anthropic component are the variables with the greatest influence on their distribution (^{Keya et al. 2016}; ^{Khan et al. 2016}).

Ruggedness was one of the variables with the highest contribution to the potential distribution model (Table 1). This variable is related to the proportion of land that facilitates the escape of the species, resulting from the calculation of relief heterogeneity (^{Riley et al. 1999}; ^{Escobar et al. 2015}). Relatively smooth land (35 to 160 m) was associated with the peninsular distribution of the bighorn sheep in Sierra de Juárez (Figure 4). These results are consistent with those recorded in the Santa Isabel and Mechudo mountain ranges by ^{Álvarez (2009)} and ^{Escobar et al. (2015)}. In this areas, the bighorn sheep is associated with ruggedness values between 21 m and 267 m, typical of mid- and high-altitude sierras with ravines that offer protection to sheep herds.

The distribution of the bighorn sheep was found to be related to natural palm-tree and gallery patches (Figure 4), both of which are associated with surface water bodies. Likewise, it was found that bighorn sheep individuals thrive along slopes facing east (Figure 4), consistent with the findings reported for Sierra del Mechudo by ^{Guerrero et al. (2003)}, as this is the prevailing slope orientation in both mountain ranges. As regard climate, the very arid semi-warm [BWh(x’)], very arid temperate [BWk(x’)], and very arid temperate [BWks] types were also associated with the distribution of the species (Figure 4).

The data obtained can be used as a baseline to elaborate an integrated management plan of the bighorn sheep habitat in Sierra de Juárez with the strategic aim to organize the resource owners in a monitoring group. Regular inspections to the bighorn sheep distribution area should be conducted in order to prevent poaching and habitat disruption. In addition, its populations should be monitored, and water bodies used by them should be kept in good condition. Another objective should be to promote the productive diversification of sites with tourism potential, equipping ranches adjacent to bighorn sheep areas to increase their capabilities to conduct nature-tourism activities in sites that are key for the bighorn sheep and that are currently used as rangeland. The above are measures that will contribute to the conservation of this iconic species and its environment.