Publisher Correction: Identifying priority conservation areas in a Saharan environment by highlighting the endangered Cuvier’s Gazelle as a flagship species

An amendment to this paper has been published and can be accessed via a link at the top of the paper.

www.nature.com/scientificreports www.nature.com/scientificreports/ decisions establishing national and global priorities in biodiversity conservation are lacking observational data. Researchers, funders and journals have been urged to conduct, fund and disseminate relevant field research 6 .
The conservation status and distribution of the Sahara's megafauna is a paradigmatic case to illustrate this lack of fundamental and key data for well-designed conservation actions in remote areas. Field studies are becoming increasingly necessary because of the collapse in Saharan wild ungulate and mammalian carnivore populations over the last century 7 . In this remote region of northern Africa, any field work faces great logistical and safety challenges, due to remoteness, lack of infrastructure, extreme environmental conditions and, in many cases, ongoing armed conflicts 8 . One of the species that inhabits this desert is the endangered Cuvier's Gazelle (Gazella cuvieri Ogilby 1840), a medium-sized ungulate endemic to North Africa. Over the last century, the Cuvier's Gazelle population has undergone major fragmentation and its numbers have declined dramatically to 2,360-4,560 individuals due to overhunting and habitat loss [9][10][11][12][13] . Currently, the distribution of Cuvier's Gazelle is limited to the Atlas Mountains and the neighbouring mountain ranges in Morocco, Algeria and Tunisia 11,[13][14][15] . It is included within the subset of Sahelo-Saharan antelopes along with dama gazelle (Nanger dama), dorcas gazelle (Gazella dorcas), slender-horned gazelle (G. leptoceros), scimitar-horned oryx (Oryx dammah) and addax (Addax nasomaculatus). These species are well-adapted to life in the extreme conditions of this region, such as high temperatures, drought and a seasonal lack of food 7,10 . Like other Sahelo-Saharan antelopes, Cuvier's Gazelle is poorly studied due mainly to its elusiveness and remote habitats 7,16 . Therefore, the limited knowledge about its biology and actual conservation status may lead to a lack of protection and even local extinctions. Cuvier's Gazelle is defined as "Endangered" in Morocco, Algeria and Tunisia 10 but globally classified as "Vulnerable" to extinction 3 . Given its importance in North Africa, the IUCN has elaborated an international strategy that includes several actions focussed on protection from illegal hunting, management of habitat, monitoring and environmental sensibility to achieve better preservation and recovery of populations 13 . Only through increasing knowledge regarding its natural history and habitat requirements will it become possible to identify the areas that need protection for more effective conservation of the species.
Changes are currently being proposed that would reclassify Cuvier's Gazelle as a mountain ecotype of slender-horned gazelle 17 , a species strictly associated with the great sand desert and ergs of the Sahara 18 . Covering a wide variety of habitats from sea level to 2600 metres, Cuvier's Gazelle mainly lives in mountain ranges and associated plateaus, such as semi-arid open Mediterranean forests of cork oak (Quercus suber), holm oak (Quercus ilex), Pinus spp., sandarac (Tetraclinis articulata), Atlas cedar (Cedrus atlantica), argan (Argania spinosa) and Phoenicean juniper (Juniperus phoenicea), but also in steppes, maquis or scrubland areas, and even in cereal fields in Algeria and Morocco 9,11,14,15,19 . In contrast, Cuvier's Gazelle appears to avoid areas covered in heavy snow at high altitudes in winter and ranges to the Sahara Desert on rocky mountains and desert plateaux, limited to areas with argan and thorn trees Acacia spp. forests 9 . In this environment, recent surveys in the Bas Drâa-Aydar region (Morocco) show promising distribution and numbers of Cuvier's Gazelle 20 , since this population has been described as one of the largest populations of the species, with 935 individuals (95% CI 597-1607) 20 . This population of Cuvier's Gazelle inhabits the extreme north-western portion of the Sahara Desert and is probably the most important in terms of numbers and genetic diversity, making it essential to the species' longevity and conservation 20 . In this region, with apparently suboptimal conditions, the topography, food availability and extreme climatic conditions may be crucial factors determining the species' presence, but human settlements may also play a role and may restrict Cuvier's Gazelle presence and usable habitat 9 .
In this study, we applied species distribution models, also known as habitat suitability models, a technique that combines information on species occurrence or abundance with environmental estimates and/or spatial characteristics [21][22][23] . Species distribution models are mainly used to address the potential effects of climate change on species distribution [24][25][26] but are also used to improve the understanding of ecological factors for conservation planning and to detect unknown potential distribution areas for rare species [27][28][29] . Our goal is to identify the main ecological features of the southernmost subset of the Cuvier's Gazelle population in the Sahara; for this purpose, we integrated data from a broad-scale field survey carried out from 2011 to 2014 and species distribution models to model habitat selection at the landscape scale, using a predictive distribution map. The final purpose of this study is to contribute to the biological conservation management of Cuvier's Gazelle by identifying the extent of its potential habitat and high priority areas to ensure the species' survival in the southernmost area of its distribution.

Methods
Study area. The survey was conducted in the extreme north-western part of the Sahara Desert between 28°30′-26°50′N and 11°40′-9°25′W in the Guelmim-Es Semara region (Fig. 1). The area is delimitated by two important geographical features: the lower Drâa River to the north and the upper basin of the Sequiat Al Hamra to the south. It is a typical Saharan landscape with a subtropical desert and a low-latitude hot, arid climate (Köppen-Geiger classification 30  . The area contains diverse terrain with rough and hilly areas (jbels), flat areas with saline depressions (sebjas), plateaux (hammadas), clay plains (dayas), stony plains (regs) and some small dune areas (ergs). There are large mountainous reliefs: Aydar Mountain (western zone); three main jbels: Zini, Rich and Ouarkziz, the last with an approximate total length of 400 km (250 km in Morocco), and at the eastern edge, the Hammada of Tindouf. The whole region is a contact point between two ecoregions: the Sahelo-Arabian and the Macaronesian regions. Thus, is it possible to find features typical of Mediterranean, tropical vegetation and Macaronesia [31][32][33] . Woody vegetation is scarce and is mainly located in ravines and oueds (ephemeral sandy rivers) with open savannah-like forests of horn trees such as acacia (Acacia raddiana), sometimes along with balanites (Balanites aegyptiaca) and calotrope www.nature.com/scientificreports www.nature.com/scientificreports/ or sodom apple (Callotropis procera), and abundant African tamarisk (Tamarix africana) bushes along the gueltas. Other typical species are argan (an endemic tree of Morocco that reaches its southernmost limit here), Periploca laevigata, Launaea arborescens, sumac (Rhus tripartitum), Maerua crassifolia and Euphorbia officinarum, this latter being a key species from a more continental Macaronesian region. There is also important extremophilic vegetation that survives in the hammadas, regs and ergs, such as rose of Jerico (Anastatica hierochuntic), Panicum turgidum, Nucularia perrinii, colocynth (Citrullus colocynthis) and Mesembryanthemum cryptanthum. field surveys and environmental predictor variables. Presence-absence data were collected from April 2011-April 2014 by completing walking surveys and sampling direct sightings and indirect signs, such as tracks, isolated dung piles and dung middens (accumulations of dung piles in latrines) [34][35][36][37] . After 194 transects on foot counting a total of 2169 km, 615 dung points were georeferenced and 61 gazelles sighted at 21 different sites. The mean distance covered by walking surveys was 12.08 ± SE 0.72 km (range 3.8-22.5 km). The correct identification of indirect signs was tested by genetic analyses at the Research Centre for Biodiversity and Genetic Resources at Porto University, Portugal, using methods described in Silva et al. 17 . For more detailed descriptions on the field data, see Gil-Sánchez et al. 20 .
To develop models that predict habitat suitability and determine the influence of variables on Cuvier's Gazelle presence, we selected a set of previously published metric variables ranked in five categories (Table 1). Mean values for topographic morphology metrics (altitude, slope and topographic ruggedness index) and two temperature-moisture metrics (heat load index and compound topographic index) were obtained from Shuttle Radar Topography Mission (SRTM) with 1 arc-second (~30 m 2 ) spatial resolution 38 . For topographic distance and human factor metrics, we used the Euclidean distance computed in ArcGIS 10.4 39 , measuring the distance to the coast and a two-tiered settlement size classification (DISTSETT1 and DISTSETT2). To describe vegetation, we used a primary productivity parameter: the average annual maximum green vegetation fraction or AMGVF 40 . AMGVF data are based on 12 years (2001-2012) of Collection 5 MOD13A standardised difference vegetation index (NDVI) data, with a spatial resolution of 30 arc-seconds (~1 km 2 ). Finally, from 19 climate variables in the WorldClim database 41 , we selected variables that captured the main environmental gradients in the study area: mean annual temperature (BIO1) and annual precipitation (BIO12). WorldClim data were downloaded at a spatial resolution of 30 arc-seconds (~1 km 2 ). The original spatial resolution of predictor variables was maintained to calculate the average of pixel values in each grid-cell by focal analysis. We then resampled the pixel size to 1 × 1 km to equal data scale for data processing. All spatial analysis were conducted using ArcGIS 10.4 42,43 . We analysed the pairwise correlation between explanatory variables: first removing those strongly correlated with coefficients of correlation above 0.80 44 and then testing for multicollinearity in the data with the Variance Inflation Factor (VIF), by using stepwise elimination of highly inflating variables with a threshold of 10 23 . After these analyses, we removed slope, heat load index and distance to coast (Supporting Information, Appendix S1), resulting in eight variables selected for modelling: topographic ruggedness index, altitude, compound topographic index, annual www.nature.com/scientificreports www.nature.com/scientificreports/ maximum green vegetation fraction, annual mean temperature, annual precipitation, distance to cities and urban villages, and distance to rural villages.
Species modelling strategy. The study area was enclosed in a 29 UTM (Universal Transverse Mercator) grid-cell measuring 5 × 5 km (Fig. 1). Due to the high mobility of gazelles, all records were gathered in each 5 × 5 km plot 45 with a measure of occurrence (presence/absence). Presence was determined by dung middens and direct observations. To avoid imperfect detection in the absence of data (i.e. false negative surveys), we only included grids with transects of more than 5 km repeated at least twice 20 . Of the 84 surveyed grids, we removed those sites within a radius of 5 km to avoid spatial correlation, thus obtaining 51 grid-cells in all (41 positive cells and 10 negative cells). We then computed the nearest neighbouring index in ArcGIS 10.4 to prevent using clustered data 42 (Supporting Information, Appendix S2).
To build a species distribution model, we used the R package 'biomod2' 46 for R 3.6.2 47 . We used two algorithms common in distribution modelling: (1) Generalised Linear Model 48 (GLM) and (2) Maximum entropy modelling of species distributions implemented in MAXENT 3.0.4 beta software 49 , to ultimately obtain an assembled model based on the variables we had selected. The use of an ensemble model or consensus algorithms in species distribution model is a powerful tool that prevents the selection of one single best model and thus eliminates or limits model selection bias while improving predictions of the current range of a species 23,50,51 . GLM is a flexible method and an extension of the linear regression models, which allows the response variable to follow a non-linear distribution and non-constant variance function. We fit GLMs for binary responses (the presence-absence data) using a logic link function, quadratic terms for each predictor and an automatic stepwise procedure with AIC. MAXENT is a reliable, effective technique based on a machine learning algorithm known as maximum entropy, which can use only presence, but also performs well when compared with presence-absence procedures that utilise both real and pseudo-absence data 52 68,69 Flow accumulation by catchment size (wetness index) described as a function of both the slope and the upstream contributing area per unit width orthogonal to the flow direction and a quantification of catenary topographic convergence. www.nature.com/scientificreports www.nature.com/scientificreports/ quadratic, product, threshold and hinge features. We ran 10 sets of pseudo-absences (PA) and equalled the number of pseudo-absences as available presences to prevent a sampling bias 53 . To build a consensus ensemble forecasting model, we chose a mix of both algorithms that had shown good performance and therefore used them to project the potential spatial distribution in the region (Fig. 2).

Topographic distance and human factor
Model evaluation and spatial forecasting. The predictive power of occurrence in the models was tested using a split-sample cross-validation approach. To do this, the dataset was randomly split into two subsets: 75% for training data and 25% for testing data. This process was repeated 20 times for each model (× 10 pseudo-absence samplings × 2 algorithms). Therefore, 400 models were built in total. Models obtained through cross-validation were compared using a predictive accuracy metric: the area under the receiver operating characteristic curve (AUC). Only predictive models with a value of over 0.8 AUC were selected to create a final ensemble model. Differences in the predictive accuracies of different algorithms were tested using a Wilcoxon signed rank test. The models with a threshold that optimises AUC (above 0.8) were transformed into one single model by binary transformation, using the committee average approach in the training area 23,51,54 . In this process, the selected models are given equal probability levels by averaging the results of binary transformation and allowing direct comparisons across models 46,50 . To determine which environmental variables were most important, we selected 10 permutations to estimate variable importance in two ensemble methods: using algorithms and with all models. Then, the three ensemble model approaches (GLM, MAXENT and the consensus model) were analysed to assess the contribution and variable response of the predictive factors. To integrate results of different methods and avoid model selection bias 50,51 , we created a committed averaging map using the ensemble model of all algorithms to forecast the potential spatial occurrence in the southernmost known distribution area of Cuvier's Gazelle. The threshold value below which the species was considered absent was <0.5 55 with a range for probability of occurrence: low (0.5) to high (1). Then, the potentially suitable area was calculated by adding together all of the 5 × 5 km presence grid-cells, considering the species' geographic Atlantic Saharan population border to be the Drâa River.

evaluation of models and habitat inferences.
Only 97 GLM and 43 MAXENT models with good accuracy (above 0.8 AUC) were selected to build three ensemble models. Significant differences were found between the two algorithms (W = 2682, p-value = 0.006), though GLM models performed better (AUC mean value: GLM, AUC = 0.86 SD = 0.048 and MAXENT, AUC = 0.83 SD = 0.039). However, the committee averaging ensemble model for MAXENT had the highest AUC value ( Table 2). The eight variables selected for modelling showed different response curves, but the shape was similar for each variable in all ensemble models, except for distance to rural villages in MAXENT (Fig. 3a). The topographic ruggedness index and distance to cities and urban villages www.nature.com/scientificreports www.nature.com/scientificreports/ were key variables determining gazelle presence for the three ensemble models, though its respective contribution showed some differences (Fig. 3b): topographic ruggedness index was the most important variable for the consensus model and MAXENT, whereas distance to cities and urban villages was most important in GLM. In both cases, the probability of the Cuvier's Gazelle occurrence increased with higher rugged terrain index and longer distances to cities and urban villages. Variables such as annual mean temperature in GLM and in the combined model, and distance to rural villages in MAXENT played a secondary role. In this instance, higher temperatures and a decrease in distance to rural villages affected negatively to gazelle presence. The rest of the variables in the different approaches showed low contribution determining gazelle presence (Fig. 3b).
predicted distribution range and priority conservation areas. Considering the northern limit of our study area to be the Drâa River, the predicted suitable Cuvier's Gazelle habitat resulting from the consensus ensemble forecasting model yielded an area of 15820 km 2 (Fig. 2). The forecasted potential spatial occurrence at its southernmost limit of distribution revealed two well-demarcated areas: "A" and "B, " with estimated areas of 4385 km 2 and 11315 km 2 , respectively (Fig. 2). These areas are regions lacking formal protection, with only a relatively low percentage (5%) included in Khnifiss National Park. Presence of Cuvier's gazelles was higher in the eastern area of Khnifiss National Park, on the steep slopes of the Hammada of Tindouf, in the mountainous terrain of Aydar and across the large reliefs in the area (Jbels Zini, Rich and Ouarkziz). The consensus ensemble forecasting model showed a continuous suitable habitat strip to the north of the Drâa River, in the Jbel Bani and beyond the limits of the ecoregion.

Discussion
Our study illustrates an example of how to successfully deal with a lack of key data in developing well-designed wildlife conservation actions in remote areas, such as the hard environment of the Sahara Desert, through combining large scale field surveys and distribution modelling. This approach has clear applications for any type of study on the population distribution and dynamics of threatened species that live in these types of regions, not only wild ungulates but also carnivores or other mammals of medium to large size. The present study provides robust ensemble models that achieved a good estimate of the importance of the variables and their response curves 23,46,56 . Our analyses confirm that the species presence and distribution of Cuvier's Gazelle in the study area are primarily influenced by the variability or complexity of the terrain and the increasing distance to large human settlements. As a result, we found a higher projected probability of Cuvier's Gazelle presence in zones with a complex network of hills and ravines and in remote areas with a low human population density and less accessibility (e.g. paved roads). These results are consistent since Cuvier's Gazelle is well-adapted to rugged mountain areas 57 , where they can find shelter from predation and poaching, but also where environment conditions are more favourable against the extreme Saharan climate. In such areas, Cuvier's Gazelle can feed in ravines and dry river basins, where vegetation remains productive throughout the year. However, we were able to confirm that they avoid the more productive oueds, which are located in flat areas, probably because of greater pressure from livestock and human presence, factors that may hold great importance in habitat selection in the area 20 . This is consistent with the annual maximum green vegetation fraction and the compound topographic response curves, both positively correlated to the size of river basins, and therefore more productive oueds. The human factor represented by the distance to the nearest large human settlement is one of the many factors constraining the potential distribution of the species. A similar phenomenon was described in Chammem et al. 35 , in which the species presence and distribution for dorcas gazelle in Tunisia seemed to be affected negatively by human presence and land use, rather than habitat characteristics. The projected species distribution model area probably reflects a shelter range that allows survival in apparently suboptimal conditions, due to both poaching and the extreme Saharan environment, where the species reaches the southernmost limit of its range. The impact of illegal wildlife hunting is one of the causes or the main cause of decline and extinction of ungulate species throughout the Sahara Desert 7,8,10,16,58 , as it is in Morocco 59,60 . Therefore, the better-preserved Cuvier's Gazelle populations in Morocco may currently be related to remote territories that have a low human population density. Furthermore, the species shows high environmental heterogeneity throughout its range in North Africa 18 . Cuvier's Gazelle and slender-horned gazelle were suggested to inhabit two different ecotypes of mountain and lowland sand dunes, respectively, with a parapatric population structure 17 , which may explain the ecological flexibility of Cuvier's Gazelle and thus its survival under extreme environmental conditions. In the study area, other threats such as overgrazing and habitat loss must be studied, but at least do not seem to have had major effects thus far, in principle. Nevertheless, this situation appears to be changing. Recently, we noticed an increase in the number of domestic grazers and local pastoralists. This appears to be related to the opening of new roadways (some paved) and the use of large flexible plastics for water storage, allowing for more stable livestock exploitation in areas far from water sources.
implications for conservation policies. At Table 2. Predictive accuracy of ensemble models by binary transformation using the committee average approach in the training area for: a) GLM (EMbyGLM), b) MAXENT (EMbyMAXENT) and c) both algorithms (EMbyAll).
www.nature.com/scientificreports www.nature.com/scientificreports/ www.nature.com/scientificreports www.nature.com/scientificreports/ populations are located in the Western High Atlas (north of Agadir) and on the south side of the central High Atlas and East Atlas 60 . We extended the species distribution model forecast area to the south of its known distribution area, identifying two well-delimited key areas for Cuvier's Gazelle (A-B; Fig. 2), a result consistent with the previous scattered and opportunistic records inside and outside of the study area, which resemble our species distribution model 9,57,60 . Within both areas, our species distribution model identified sectors with high Cuvier's Gazelle probability of presence (close to 1). These locations hold a suitable habitat where action aimed at conservation is urgent. The whole region is considered public land under the Moroccan government's management, with nomads moving about in temporary camps with herds of goats, sheep and dromedaries. During the field surveys, we observed that these people opportunistically try to hunt gazelles, whereas groups of poachers coming from the nearest large cities operate at will. Moroccan authorities planned to declare one national park (Bas Drâa NP), one Biological and Ecological Interest Site (Oued Tirhzer) as well as two hunting reserves (Messeied-Abeteih and Oued Chbeyka) 61,62 . However, the surface area covered by these planned protected areas is insufficient to protect a key space for Cuvier's Gazelle (i.e. a viable population), since it would cover just 25% of the range we have estimated and, moreover, these projects have been stalled for more than two decades. Moreover, local Cuvier's Gazelle densities are quite low (0.08 individuals/km 2 ) and the total population has been estimated at 935 (95% CI 597-1607) individuals 20 . Our results show the need to extend the protected areas designed in the Moroccan strategy for the preservation of endangered species 60 , as well as the urgency in the implementation of this national strategy to protect the current natural populations of Cuvier's Gazelle. New protection areas forming a backbone for the region's sustainable development with proper protection (in terms of poaching) and a low human impact will be crucial to preserving the species. The Cuvier's Gazelle population inhabiting this region will have a better chance against any harsh conditions (water shortages and extremely high temperatures) in a changing climate scenario 7 . Furthermore, the region has been proposed as a dispersal corridor and fauna shelter for the Sahara-Sahel region 16 . These features make this part of Morocco a crucial area not only for the Cuvier's Gazelle preservation, but also for other endangered ungulate species still present in the region such as dorcas gazelle and Saharan barbary sheep (Ammontragus lervia sahariensis), and carnivores such as striped hyaena (Hyaena hyaena), sand cat (Felis margarita), caracal (Caracal caracal) and, in all likelihood, the critically endangered Saharan cheetah (Acinonyx jubatus hecky). The latter could find a strategic area for recovery here, the area could also aid in the recovery of other extinct ungulate species like mhorr dama gazelle (Nanger dama mhorr), scimitar-horned oryx and addax, through reintroduction projects. Some of these species have begun to be reintroduced in the Souss-Massa region 63,64 . The current work can provide a basis for the definitive design of protected areas and reintroduction projects.

Data availability
The datasets generated and/or analysed during the current study are available from the corresponding authors upon reasonable request.