Ecological niche modeling of genetic lineages of the great gerbil, Rhombomys opimus (Rodentia: Gerbillinae)

Great gerbil (Rhombomys opimus Lichtenstein, 1823) is distributed in Central Asia and some parts of the Middle East. It is widely found in central and northeast parts of Iran with two distinct genetic lineages: R. o. sodalis in the northern slopes of the Elburz Mountains and R. o. sargadensis in the southern slopes. This large rodent acts as the main host of natural focal diseases. No study has surveyed the ecological niche of the lineages and how their distribution might be influenced by different climatic variables. To examine the distribution patterns of this murid rodent, we aimed to determine the habitat preferences and effects of environmental variables on the ecological niche. Using a species distribution approach for modeling of regional niche specialization, suitable habitats predicted for R. o. sodalis were mainly located in Golestan province in northern Iran, along the northern slope of Elburz, while R. o. sargadensis, showed great potential distribution along the southern slope of Elburz and around the Kavir Desert and the Lut Desert. Despite the widest potential distribution of R. o. sargadensis from northeast to northwest and through Central Iran, the geographic range of R. o. sodalis was smaller and mostly confined to Golestan province. The results support the presence of the two genetic lineages of Rhombomys in Iran and confirm that there is no significant niche overlap between the two subspecies. Furthermore, it provided several perspectives for future taxonomic studies and prevention hygiene programs for public health.


Introduction
Great gerbil (Rhombomys opimus Lichtenstein, 1823), which is generally considered to be a monotypic species based on coloration and size [1], is known to be distributed throughout Central Asia and some parts of the Middle East, including Kazakhstan, Uzbekistan, Kyrgyzstan, Tajikistan, Turkmenistan, North China, South Mongolia, Iran, northern Afghanistan, and southwest Pakistan [2,3]. Generally, this species is abundant in harsh climatic conditions with hot and dry summers and cold winters characterized by low average annual precipitation and relative humidity [4,5]. Rhombomys opimus has a social structure, and individuals often burrow close to one another; the family group (adult male, several adult females and young individuals of several generations) lives in one complex burrow system [3,6]. This large gerbil has a wide distribution range in central and northeastern Iran and occupies sandy or clay deserts, usually in foothills and mountainous areas with scattered shrubby vegetation, especially Haloxylon ammodendron (Amaranthaceae), succulent plants such as Salsola spp., Climacoptera spp., and Suaeda spp. (Amaranthaceae) [3,7]. The patchy distribution of saline ecological biotopes in the north, northwest and center of the country leads to isolated populations of the great gerbil, and high molecular variability within Iranian populations has been reported [1,8]. Ellerman [9] listed seven nominal subspecies for this species, whereas the real number of described forms is 12 [2]. An accurate number of subspecies critically need to be estimated in a future taxonomic revision, but provisionally can be accepted as six:  [2].
In Iran, Rhombomys opimus sodalis is only found at elevations of approximately 600-1000 m in Golestan province on the northern slope of the Elburz Mountains in addition to small patches in North Khorasan (Bojnourd) and Razavi Khorasan (Sarakhs and Dargaz) provinces in northeast of Iran, whereas R. o. sargadensis, which is considered the more widely distributed subspecies, could be found at higher elevations than R. o. sodalis on the southern slope of the Elburz Mountains and around the Kavir Desert in Central Iran [1,8,10,11] (Fig 1).
These two subspecies are also morphologically distinct from each other with known differences in coloration and size [1]. Due to their sympatric presence in Iran, based on the report of sympatric haplotypes of R. o. sodalis in a population of R. o. sargadensis in the Shahrood district (Semnan province), near Kavir Desert [8], it has been suggested that migration of R. o. sodalis from Golestan province to the territories of R. o. sargadensis in Semnan province is possible. However, Oshaghi et al. [8] and Bakhshi et al. [10] attributed these "sympatric haplotypes" to either common ancestry or migration. Nevertheless, mtDNA markers are more likely to reflect shared ancestral populations rather than recent migration. The reciprocal crosses between individuals of the two subspecies attempted by Oshaghi et al. [8] show that there are no pre or post-zygotic barriers between R. o. sodalis males and R. o. sargadensis females. However, the reciprocal cross of R. o. sodalis females with R. o. sargadensis males resulted in the death of the paired individuals, suggesting the possibility of some prezygotic isolation mechanism. Bakhshi et al. [10] obtained offspring from this reciprocal cross, although he provides no details of the number of crosses attempted or the number of successful crosses. Moreover, the Elburz Mountain Chain in northern Iran acts as a natural barrier between populations of R. o. sodalis and R. o. sargadensis [10], and consequently, future speciation is predicted due to the large intraspecific variation among populations distributed in different localities across Iran [8].
Great gerbils can damage crops and irrigation canals, destroy the vegetation and cause the die-off of plants over their colonies by biting the main roots. Hence, they are in direct competition with livestock and have been recognized as carriers of several zoonotic diseases (e.g., plague, leishmaniasis, leptospirosis, and chronic respiratory disease) [12][13][14]. The presence of flea (Siphonaptera) assemblages and sandflies (Diptera) (zoonotic cutaneous leishmaniasis transmitters) in great gerbil burrow systems increases their ability to sustain the pathogen and the probability of infectious emergence of rodent-borne diseases [15][16][17][18][19][20].
With respect to the lack of information on current habitat delimitation and the potential distribution pattern of great gerbil genetic lineages (or subspecies) in Iran, we seek to elucidate their environmental niches through species distribution models. Accordingly, this study aimed to examine (i) whether the genetic lineages of R. opimus have different habitat suitability and niche specialization throughout the Iranian Plateau and/or in the contact zone, where they are found in sympatry, and also (ii) to understand the role of environmental variables on their regional spatial distributions.

Environmental and species data
To use a machine learning model, 19 bioclimatic variables were selected from the WorldClim database (www.worldclim.org) at a spatial resolution of 1 km 2 and then processed to model the target genetic lineage distribution. Since multicollinearity among bioclimatic variables may result in misunderstanding the contribution ratio of the most important variables to the model [21], highly correlated variables (Pearson correlation coefficient: r� 0.75) were ignored for final analysis. The final set contained 13 climatic variables, as shown in Table 1.
Distribution data were gathered from the Global Biodiversity Information Facility (GBIF database), VertNet as a publicly accessible database of vertebrate biodiversity data from natural history collections around the world, and published papers and books [e.g., 3,22,23]. A total of 102 distribution records within the country were obtained (R. o. sodalis: n = 56 and R. o. sargadensis: n = 46) (Fig 2). The points were screened in ArcGIS 10 (ESRI, Redland, USA) with nearest neighbour analysis to assess spatial autocorrelation filtered using SDMTools [24,25]. This analysis discovered a low clustering among presence records. An Excel file including point localities for R. opimus subspecies is available at S1 Appendix. No ethical approval is required for this study because techniques performed here do not involve animals.

Modeling potential occurrence
To model the current geographic distribution range of great gerbil subspecies, the maximum entropy modeling algorithm (MAXENT v. 3.3.3. program, www.cs.princeton.edu) was used [26]. MAXENT has been found to perform better than many other modeling methods for occurrence data to predict a species' distribution [27][28][29]. MAXENT is proficient of calculating species distribution using presence-only records; it can consider both continuous and discrete variables in the model to identify the important environmental variables affecting species distribution [27,30].
MAXENT was applied, and 70% of the occurrence records were used as training data and the remaining 30% to test it (as test data). To determine the model performance, the calculated value of the area under the curve (AUC) of the receiver operating characteristic (ROC) curve on the training and testing data was considered. The AUC evaluator indicates the power of the model in distinguishing presence from absence records. A calculated value close to 1 indicates the high predictive ability of the model, while a value of 0.5 suggests that the model lacks sufficient power to predict the species distribution range [25,[31][32][33][34]. Spatial prediction maps of habitat suitability for any given location, as the model output, range from 0 (very low) to 1 (very high) relative habitat suitability of species presence [26,28]. Jackknife analyses used to estimate the importance of each of the variables that reduce the model reliability when omitted. Moreover, ENMTools [35] was used to test the percentage of niche overlap between the two predicted models using Schoener's D [36] and Hellinger's based I [37] indices.

Results
According to the obtained AUC values, the model predictive accuracy of both subspecies indicated high performance (AUC = 0.77 for R. o. sargadensis; AUC = 0.92 for R. o. sodalis), so the MAXENT approach seems to perform well for modeling ecological niche segregation. Based on our models, the most current suitable habitats of R. o. sodalis were mainly projected in Golestan province, along the northern slope of the Elburz Mountains, and penetrated further into North Khorasan provinces as well as small patches in the western parts of Iran (Fig 3). Under our projection, the other subspecies, R. o. sargadensis, showed great potential general distribution in most parts of central, eastern and northern Iran, along the southern slope of the Elburz Mountains and around Kavir and Lut Deserts in Central Iran, as well as northwesternmost parts of Iran (Fig 4). Species distribution models obtained by MAXENT showed that precipitation of the driest month (Bio14) and precipitation seasonality (Bio15) are the most important predictor variable determining the current distribution of R. o. sodalis, while the most important variables for R. o. sargadensis are precipitation of the coldest quarter (Bio19) and precipitation of the wettest month (Bio13) ( Table 1). The results also showed that the annual mean tempreture (Bio1) and mean temperature of the driest quarter (Bio9) are the most important variables in explaining the distribution of both currently assumed mitochondrial lineages. In addition, most suitable habitats of great gerbils of the genetic lineage R. o. sargadensis were situated in northeastern, northwestern and Central Iran, while the geographic range of

Discussion
As a desert-adapted rodent, the distribution of great gerbil is significantly associated with temperature, precipitation, terrain, vegetation and other ecological environmental factors [5]. The results showed that temperature shapes the ecological niche of both R. o. sargadensis and R. o. sodalis, followed by the annual precipitation amount. However, Gholamrezaei et al. [38] indicated slope as the main variable affecting the distribution pattern of R. opimus. Furthermore, Gao et al. [5] showed that R. opimus is distributed in the area of elevation between 200 and 600

PLOS ONE
m with a slope of 0-3 degrees, an average annual temperature from 6 to 10˚C and an annual precipitation of 120-200 mm in Xinjiang, northwest China.
According to the models obtained by MAXENT, the potential distribution of R. o. sargadensis was strongly constrained in an isolated habitat patch in the marginal part of the species distribution range; a patch of habitat in the extreme northwest of Iran might contain genetically isolated populations (Fig 4).
With regard to different karyological statuses that have been reported for R. o. sodalis from Gonbad and Bandar Torkaman, and Gorgan, considering that all studied regions are located in Golestan province and in the range of R. o. sodalis [1,39], we expected differentiation of the ecological niche of the species. However, previous chromosomal analysis has indicated that both genetic lineages of R. opimus in Iran have 2n = 40 chromosomes [40,41].
Haplotypes belonging to two different mtDNA lineages are found together in some populations of northern Iran, which makes it possible to hypothesize that there is hybridization of these two genetic lineages (but it can be proven only by analysis of nuclear genes(. While Oshaghi et al. [8] made a prediction regarding future speciation, it is dependent on retaining geographic or ecological barriers between these two subspecies that have not been identified. Some hybrid zones at a regional scale were observed in which the two genetic lineages may be in contact (Kopet-Dag Mountains in the northeast and Ghaflankooh Mountains in northwest Iran) (Fig 1). These indicated a potential niche overlap in the distribution range of both. Subsequently, the hypothesis of future speciation within the species based on geographic variation of haplotypes among localities [8] need to be assessed.
Habitat suitability for hosts and/or vectors can have effects on the spreading of pathogens and may influence both their abundance and movements [42,43]. Suitable areas that may lead to a higher contact rate of hosts and vectors can act as corridors for the transmission of pathogens through a larger landscape, and unsuitable habitats may act as barriers because they prevent the transmission of pathogens by hosts or vectors [33,44]. For each of the host or vector species, these corridors and barriers can show specification, depending on their movement abilities. In the present study, these are areas along the northern and southern slopes of the Elburz Mountains and the eastern parts of the Zagros Mountains around the Kavir Desert, which have been inferred as suitable areas for inhabitation of the great gerbil. The Sabalan, Sahand, and Ghaflankooh Mountains in northwestern Iran and the Saridash Mountains on the border of Iran-Turkey are other possible functional corridors for the distribution of this species. The southern shores of the Caspian Sea, throughout Golestan province and northwest of North Khorasan province, which are covered with dense forest and trees, with developing cultivation and agriculture activities that can provide food and shelter for rodent populations, may be considered as structural corridors (See Fig 1).
Harsh climatic conditions (arid or semi-arid climates) in some parts of the Kavir Desert act as another structural barrier for great gerbil to expand its distribution throughout Iran. For the great gerbil, the Elburz Mountains in the north and the Zagros Mountains in the west of the country play major barriers against dispersal.
Lizhi et al. [45] showed that the distribution of Chinese populations of great gerbil has been altered due to human activities. A study in the southern Kyzylkum Desert (western Uzbekistan) demonstrated that the population density and the size and structure of social groups varied yearly in response to changing conditions of precipitation and temperature; years with considerable rain and snowfall during winter seem to produce enough and available succulent vegetation cover to facilitate subsequent breeding and expanding family groups [46]. The great gerbil limits reproduction to periods of rainfall and the subsequent growth of green vegetation [47][48][49][50]. Thus, breeding and survival in this species seem to greatly depend on environmental conditions rather than any other factors, such as social organization [46]. Several studies reported fluctuations in rodents and their associated parasite populations with changes in climate, habitat structure and feeding resources. For example, Ari et al. [50] noted that fluctuation in the population of fleas harbored by great gerbils can occur due to rainfall, relative humidity, and temperature; in warm moist weather, rodent hosts are more available for the growth of bacteria Yersinia pestis, and hence, the transmission rates of plague infection may increase.
Great gerbil is considered as the most important reservoir host of Leishmania major in Iran, which is transmitted by sand flies of the genus Phlebotomus [16,17,19]. Zoonotic coutaneous leishmaniasis (ZCL) due to L. major is known as one of the zoonoses increasing in Iran [17]. As an example, Rassi et al. [51] stated that the rate of infection of great gerbils to this parasite is high and may reach to 92.5% at endemic areas of ZCL in Kalaleh, Golestan province, north of Iran. In another study [16], seasonal variations of natural infection with Leishmania in population of great gerbils in Badrood district, Esfahan province, central Iran were surveyed. The lowest and highest infection rates were observed in summer and fall, respectively. Gerbils were found to be infected with three species of L. major, L. turanica and L. gerbilli, which transmit in the population of R. opimus in central part of Iran. Leishmania major infection is generally accompanied by L. turanica in infected great gerbils, showing the highest rate in fall.
The distribution of L. turanica in rodents showed coincidence with the distribution patterns of sandflies [52]. Deep burrows of great gerbils, which may extend to three meters in depth depending on the stability of the soil [53], will have more stable temperatures inside and hence are likely to increase the prevalence of ectoparasites such as fleas, which are less prone to survive low humidity and extreme temperatures [54]. This rodent species plays a role as a reservoir host for the fleas of Xenopsylla and Nosopsyllus genera [20,55], tick Hyaloma [4,56], and the mite Ornithonyssus bacoti [55,[57][58][59]. Moreover, the oxyurid Dentostomella translucida is a nematode parasitizing the digestive system of great gerbils [60,61]. Kamranrashani et al. [62] reported that great gerbils are host for several species of cestodes and nematodes in Golestan province. Furthermore, dwarf tapeworm Hymenolepis nana, which is a cyclophyllidean zoonotic enteric parasite, known to be occurred in different rodents, including gerbils, of Golestan and Razavi Khorasan provinces in north of Iran [63].
Contact zones create candidate habitats where future researches could examine the extent of ecological divergence between the known lineages [64,65]. Therefore, increasing knowledge regarding the rodent distribution and effective environmental variables, the patterns of rodent populations, colonization and movements, and finally, the preparation of risk maps for making decisions against zoonotic diseases and development are of great concern.

Conclusion
The present study provided projections of the potential geographic distribution of the two subspecies of great gerbil, with a fundamental role in the epidemiology of zoonoses. For R. o. sargadensis distributed across most areas of Iran, southern Afghanistan and western Pakistan, the sampling area represents more than half of the geographic range of the subspecies, and the resulting model can potentially be of relatively good quality. However, this is not the case for R. o. sodalis. This subspecies is distributed across northeastern Iran, Turkmenistan and northern Afghanistan. Cumulatively, future studies should consider covering the whole geographical distribution range, especially identifying contact zones, population structure, and comprehensive distribution samplings for genetic studies are required to clarify the taxonomic status of this species. Using a more comprehensive dataset for ecological niche modeling and habitat evaluation will increase the precision of the models and estimate the probable future distribution of the species considering the role of climate changes affecting the environmental variables.