Testudo graeca ibera, commonly referred to as the “spur-thighed tortoise” or “İber Kaplumbağası” in Turkish, is a terrestrial chelonian belonging to the Testudinidae family. [1] This subspecies is naturally distributed across the eastern Mediterranean region, including Anatolia, the Balkans, and parts of the Caucasus, where it occupies a wide range of Mediterranean and semi-arid habitats. [2] Its presence in these regions reflects long-term ecological adaptation rather than recent human-mediated introduction. Notably widespread in temperate and Mediterranean climate zones, Testudo graeca ibera exhibits considerable ecological plasticity, allowing it to persist under varying temperature and precipitation regimes. [3] The species commonly inhabits scrublands, open woodlands, grasslands, and agricultural landscapes, demonstrating tolerance to habitat fragmentation and seasonal climatic fluctuations. Due to its sensitivity to environmental change and its broad yet structured distribution, this subspecies represents an informative model organism for ecological niche modeling studies aimed at understanding climate–habitat relationships. # Methods
To acquire the occurrence data for Testudo graeca
ibera, we used the Global Biodiversity Information
Facility (GBIF) database, which provides globally curated species
occurrence records. Occurrence data were obtained by querying the GBIF
database for records corresponding to Testudo graeca ibera and filtering
the results to include only observations with available geographic
coordinates. The downloaded dataset was imported into R as a
comma-separated values (CSV) file. To prepare the data for ecological
niche modeling, only the geographic coordinate information—specifically
the decimalLongitude and decimalLatitude
columns—was retained. Records containing missing or incomplete
coordinate values were removed using the na.omit function
to ensure data quality and reliability. The resulting cleaned dataset
consisted of georeferenced occurrence points representing the observed
distribution of Testudo graeca ibera. These occurrence records formed
the basis for subsequent ecological niche modeling analyses conducted in
this study.
To format the occurrence data for ecological niche modeling using the
biomod2 package, the dataset was prepared by combining
species presence records of Testudo graeca ibera with
environmental variables. The following steps outline the data
preparation process:
Preparing the Species Occurrence Data: The
geographic coordinates (decimalLongitude and
decimalLatitude) of Testudo graeca ibera presence
records, obtained from GBIF, were formatted as a DataFrame. This
DataFrame was then converted into a spatial object using the
vect function from the terra package, with the
Coordinate Reference System (CRS) defined as
EPSG:4326.
Preparing the Environmental Variables:
Bioclimatic variables were downloaded from the WorldClim database at a
spatial resolution of 10 minutes. The variables were stored as
.tif files and loaded into R as a SpatRaster object using
the rast function from the terra package.
These variables were used as environmental predictors for modeling the
distribution of Testudo graeca ibera.
Selected bioclimatic variables include:
Integrating Species Data with Environmental
Variables: Environmental variable values associated with the
species occurrence locations were obtained using the
extract function. This procedure produced a dataset in
which geographic coordinates were linked with the corresponding
bioclimatic values for each occurrence point.
Formatting Data for Modeling: Occurrence records
were prepared for analysis with the biomod2 package through
the BIOMOD_FormatingData function. To obtain a balanced
dataset, pseudo-absence locations were randomly sampled within the
defined study area. Species presence records were coded as
1, whereas pseudo-absence records were assigned a value of
0. The main parameters applied during the data formatting
process are outlined below:
resp.var: Species presence data.
expl.var: Environmental predictor variables
(bioclimatic variables).
resp.xy: Coordinates of species presence
points.
PA.nb.absences: Number of pseudo-absence points
generated (set to 10,000).
This procedure prepared the dataset for species distribution modeling
with the biomod2 package. By linking species presence
records with environmental predictors and incorporating pseudo-absence
points, a comprehensive dataset suitable for ecological niche modeling
was obtained.
Historical and future bioclimatic datasets were obtained to characterize the ecological niche of Testudo graeca ibera and to assess its potential distribution under present and projected climate conditions. The data acquisition process followed the steps outlined below:
Historical Bioclimatic Variables:
Past bioclimatic variables were obtained from the WorldClim database at a spatial resolution of 10 minutes.
The variables, provided in .tif format, were
imported into R as a SpatRaster object using the rast
function from the terra package.
The following bioclimatic variables were selected for this study:
BIO1 Annual Mean Temperature
BIO5 Max Temperature of Warmest Month
BIO6 Min Temperature of Coldest Month
BIO12 Annual Precipitation
BIO15 Precipitation Seasonality
BIO17 Precipitation of Driest Quarter
These variables were selected based on their established influence on habitat suitability and environmental tolerance in terrestrial species.
Future Bioclimatic Variables:
Future bioclimatic variables were obtained from the WorldClim database for the 2081–2100 time period. The data were derived from the UKESM1-0-LL climate model under the SSP585 scenario, which represents a high greenhouse gas emissions pathway.
Future bioclimatic variables were acquired at the same spatial resolution as the historical dataset (10 minutes) to ensure consistency across analyses.
The future climate .tif files were imported into R
and aligned with the historical variable names to ensure compatibility
during modeling and projection steps.
Variable Preparation:
Both historical and future bioclimatic datasets were organized as
SpatRaster objects using the rast function from the
terra package, enabling their efficient use as
environmental predictors in the modeling process.
The Coordinate Reference System (CRS) of the future climate variables was harmonized with that of the historical datasets to maintain spatial consistency across analyses.
Selection of Bioclimatic Variables:
By incorporating both historical and future bioclimatic variables, the model was able to characterize the current ecological niche of Testudo graeca ibera and assess potential shifts in habitat suitability under projected climate conditions. This approach is essential for evaluating the potential effects of climate change on the species’ spatial distribution.
## Modeling
Species distribution models were fitted using the
BIOMOD_Modeling function from the biomod2
package, applying Generalized Linear Models (GLM). Model calibration
employed a random cross-validation approach with three replicates,
allocating 80% of the data to model training. Model accuracy was
assessed using the True Skill Statistic (TSS) and the Area Under the
Receiver Operating Characteristic Curve (ROC). A total of 10,000
pseudo-absence points were generated using a random sampling
strategy.
Current habitat suitability was projected using historical
bioclimatic variables (BIO1, BIO5, BIO6, BIO12, BIO15, and BIO17)
through the BIOMOD_Projection function.
Future climate projections for the 2081–2100 period under the SSP585 scenario were prepared for analysis. Following the harmonization of layer names with the historical dataset, habitat suitability projections were generated using the same GLM model.
Model performance was evaluated using the
get_evaluations() function, with mean TSS and ROC values
calculated across all cross-validation replicates.
The projected current habitat suitability for Testudo graeca ibera reveals a spatial pattern largely concentrated in the eastern Mediterranean basin and adjacent regions. The GLM-based model identifies high suitability across Anatolia, the southern Balkans, and parts of the eastern Mediterranean coastline, with additional suitable areas extending into the Caucasus region and surrounding lowland habitats. Moderate suitability is also observed in parts of the Middle East and western Asia, reflecting the species’ tolerance to semi-arid and Mediterranean climatic conditions. Overall, the model highlights a relatively constrained but ecologically coherent distribution, closely aligned with the known native range of Testudo graeca ibera. Areas characterized by extreme climatic conditions, such as northern Europe and highly arid desert regions, are consistently predicted as unsuitable. These results suggest that temperature extremes and precipitation patterns play a key role in limiting the current distribution of the species, with suitable habitats primarily restricted to regions exhibiting mild winters, warm summers, and moderate seasonal precipitation variability.
Under the SSP585 scenario for the 2081–2100 period, future habitat suitability projections for Testudo graeca ibera indicate a pronounced shift in suitable areas toward higher latitudes and more northern regions. Compared to current conditions, suitability across the eastern Mediterranean basin and parts of Anatolia is markedly reduced, with many coastal and southern areas becoming increasingly unsuitable. In contrast, regions in the northern Balkans, parts of Eastern Europe, and higher-altitude areas show relatively higher suitability under future climate conditions. Overall suitability across North Africa and the Middle East declines substantially, reflecting the species’ sensitivity to increased temperatures and altered precipitation regimes. The projected distribution suggests a contraction of suitable habitats from traditionally occupied southern and lowland regions, accompanied by a relative expansion toward cooler and more temperate zones. These patterns indicate that future climate change may substantially constrain the geographic range of Testudo graeca ibera, potentially forcing the species toward climatically marginal but cooler habitats.
The following maps illustrate the spatial distribution of the bioclimatic variables used in the species distribution model for Testudo graeca ibera. Each map represents the geographic variation of a specific climatic factor, such as temperature and precipitation, across the study area. Color gradients indicate differences in variable intensity, providing an environmental context for understanding the species’ predicted habitat suitability.
The map below represents the predicted current distribution of Testudo graeca ibera*, derived from species distribution modelling using present-day bioclimatic variables.
Based on future climate projections for the 2081–2100 period under the SSP585 scenario, Testudo graeca ibera is expected to undergo a pronounced contraction and spatial shift in habitat suitability. Compared to current conditions, suitability markedly decreases across much of the eastern Mediterranean region, including large parts of Anatolia and surrounding lowland areas. While some relatively suitable habitats persist in higher-latitude regions of Eastern Europe and the northern Balkans, overall suitability within the species’ traditional Mediterranean range declines substantially. Regions characterized by increasingly warm and dry conditions, particularly North Africa, the Middle East, and southern coastal zones, show a strong reduction in habitat suitability. In contrast to expectations of broad range expansion, the projections indicate only limited potential persistence in cooler and more continental environments rather than widespread northward expansion. Overall, the results suggest that future climate change may restrict Testudo graeca ibera to a narrower range of climatically suitable areas, leading to a loss of habitat in regions that currently support stable populations, especially within Mediterranean and semi-arid ecosystems.
In light of the updated future projections, the bioclimatic preferences of Testudo graeca ibera—which favor temperate and Mediterranean conditions—are only partially maintained, while notable departures from earlier expectations become evident. Although previous assumptions often implied that warming climates might allow persistence or gradual shifts toward cooler regions, the present model suggests that altered temperature extremes and changing precipitation patterns, including increased drought intensity and seasonal variability, may substantially limit the species’ future distribution. Rather than retaining suitability across southern and coastal Mediterranean regions, the projections indicate a pronounced decline in habitat suitability throughout much of Anatolia, the Middle East, and North Africa. Limited areas of relative suitability persist in more northern and continental regions, such as parts of Eastern Europe and the northern Balkans, where climatic conditions remain comparatively moderate. This pattern contrasts with earlier expectations of stable suitability within traditional Mediterranean habitats and highlights a potential contraction toward climatically marginal but cooler environments. These discrepancies may be attributed to differences among climate models, the choice of emission scenarios, and the exclusion of additional ecological factors such as habitat fragmentation, land use change, and biotic interactions, which can strongly influence real-world species persistence. Consequently, while Testudo graeca ibera continues to exhibit an affinity for temperate climatic conditions, its projected future distribution appears more restricted than previously anticipated, reflecting the complex and region-specific impacts of climate change on Mediterranean and semi-arid ecosystems.
## full.name PA run algo metric.eval
## 1 Testudo.graeca.ibera_PA1_RUN1_GLM PA1 RUN1 GLM TSS
## 2 Testudo.graeca.ibera_PA1_RUN1_GLM PA1 RUN1 GLM AUCroc
## 3 Testudo.graeca.ibera_PA1_RUN2_GLM PA1 RUN2 GLM TSS
## 4 Testudo.graeca.ibera_PA1_RUN2_GLM PA1 RUN2 GLM AUCroc
## 5 Testudo.graeca.ibera_PA1_RUN3_GLM PA1 RUN3 GLM TSS
## 6 Testudo.graeca.ibera_PA1_RUN3_GLM PA1 RUN3 GLM AUCroc
## 7 Testudo.graeca.ibera_PA1_allRun_GLM PA1 allRun GLM TSS
## 8 Testudo.graeca.ibera_PA1_allRun_GLM PA1 allRun GLM AUCroc
## 9 Testudo.graeca.ibera_allData_allRun_GLM allData allRun GLM TSS
## 10 Testudo.graeca.ibera_allData_allRun_GLM allData allRun GLM AUCroc
## cutoff sensitivity specificity calibration validation evaluation
## 1 654.0 98.993 95.588 0.946 0.938 NA
## 2 781.0 98.322 96.025 0.985 0.988 NA
## 3 693.0 98.322 95.975 0.943 0.923 NA
## 4 780.5 97.651 96.400 0.986 0.982 NA
## 5 782.0 97.651 96.150 0.941 0.921 NA
## 6 780.5 97.987 96.138 0.986 0.984 NA
## 7 625.0 98.660 95.670 0.943 NA NA
## 8 709.0 98.391 95.900 0.985 NA NA
## 9 14.0 97.319 95.330 0.925 NA NA
## 10 18.5 97.051 95.840 0.984 NA NA
- [1] Fritz, U., & Havaš, P. (2007). Checklist of chelonians of the world. Vertebrate Zoology, 57(2), 149–368.
- [2] Vamberger, M., Stuckas, H., Sacco, F., D’Angelo, S., Arculeo, M., Cheylan, M., & Fritz, U. (2017). Differences in gene flow and genetic structure in the spur-thighed tortoise (Testudo graeca) across the Mediterranean region. Zoology in the Middle East, 63(2), 113–128.
- [3] IUCN. (2023). Testudo graeca. The IUCN Red List of Threatened Species. International Union for Conservation of Nature.
- [4] Hijmans, R.J., Cameron, S.E., Parra, J.L., Jones, P.G., & Jarvis, A. (2005). Very high resolution interpolated climate surfaces for global land areas. International Journal of Climatology, 25(15), 1965–1978.
- [5] Fick, S.E., & Hijmans, R.J. (2017). WorldClim 2: New 1-km spatial resolution climate surfaces for global land areas. International Journal of Climatology, 37(12), 4302–4315.
- [6] Thuiller, W., Georges, D., Engler, R., & Breiner, F. (2023). biomod2: Ensemble platform for species distribution modeling. R package version 4.2-5.
- [7] Araújo, M.B., & Peterson, A.T. (2012). Uses and misuses of bioclimatic envelope modeling. Ecology, 93(7), 1527–1539.
- [8] Urban, M.C. (2015). Accelerating extinction risk from climate change. Science, 348(6234), 571–573.