Modelling the nesting-habitat of threatened vulture species in the caucasus: An ecosystem approach to formalising environmental factors in species distribution models

1.   Introduction

Located at the top of trophic pyramids, vultures play an important role in ecosystems as primary scavengers providing a vital ecological service in preventing pollution. Through their specialisation in the rapid carcass elimination, they reduce the potential risk of spreading disease to wildlife, livestock and humans (Aresu et al., 2021; Jha et al., 2022; Sharma et al., 2023). Sensitive to anthropogenic disturbance, vultures are considered indicators of undisturbed natural ecosystems as well as flagship or umbrella species (Guerrero-Casado et al., 2013; Mishra et al., 2018; Ortiz-Urbina et al., 2020). Significant causes of declining vulture population globally include decreased food availability, habitat loss and degradation, wind farms, poisoning by veterinary drugs (e.g., diclofenac), and electrocution (Guerrero-Casado et al., 2013; Kirazlı, 2016; Farashi and Alizadeh-Noughani, 2018; Xirouchakis et al., 2021; Sharma et al., 2023). The primary avian scavengers breeding in the Caucasus are the endangered Cinereous Vulture (Aegypius monachus), Egyptian Vulture (Neophron percnopterus), Griffon Vulture (Gyps fulvus), and Bearded Vulture (Gypaetus barbatus). Despite long ornithological observations of these birds in the region (Belik, 2004, 2014; Abuladze et al., 2011; Edisherashvili, 2011; Tilba, 2014; Ilyukh, 2017; Komarov, 2017; Karyakin et al., 2018; Karimov and Mammadov, 2019; Parfenov, 2019; Perevozov, 2020; Mnatsekanov, 2022; Shevtsov and Ilyukh, 2022; etc.), current knowledge about the localisation of their nesting sites remains extremely poor. Information on climatic, orographic and biotic characteristics of nesting habitats and vulture competition in the Caucasus is still lacking.

Nesting site availability largely determines the sustainability and growth of the vulture populations (Sandesh et al., 2022). Gaps in understanding the environmental predictors of bird nesting site distribution hinder the development of successful conservation and management strategies in areas of potential conflict between scavengers and humans (Guerrero-Casado et al., 2013; Farashi and Alizadeh-Noughani, 2018). In mountainous rugged landscapes, ground-based methods of detecting nesting sites are often difficult and ineffective. A proper solution could be the application of Species Distribution Modelling (SDM) or Ecological Niche Modelling (ENM) based on geographic information systems and statistical techniques. Information generated by SDM/ENM models from sparse occurrence points and environmental variables helps to investigate ecology and distribution of rare species and identify suitable nesting sites that can be used to create buffer zones (Guerrero-Casado et al., 2013; Farashi and Alizadeh-Noughani, 2018; Aresu et al., 2021; Jha et al., 2022). SDM/ENM models have been applied successfully in different parts of the world for vulture studies and conservation planning (Guerrero-Casado et al., 2013; Farashi and Alizadeh-Noughani, 2018; Brink et al., 2020; Ortiz-Urbina et al., 2020; Khwarahm et al., 2021; Vignali et al., 2021; Jha et al., 2022; Sharma et al., 2023). In the Caucasus, similar studies have previously only been conducted in Georgia (Gavashelishvili and McGrady, 2006; Gavashelishvili, 2012).

A current methodological challenge of species distribution modelling beyond species ‘bioclimate envelope’ modelling is to consider biotic interactions in SDM/ENM models. Positive and negative interspecific interactions significantly affect the species distribution and the sustainability of populations. Species distribution modelling and especially ecological niche modelling should therefore include analysis of biotic factors in order to be biologically meaningful (Soberón and Peterson, 2005; Peterson and Soberón, 2012; Simões and Peterson, 2018). A typical approach to considering biotic factors in SDM/ENM models includes the application of prey distribution/density shapefiles (food availability factor) and digital land cover maps (forest frequency, crown density, tree height, etc.) sourced from official government records or published literature (Farashi and Alizadeh-Noughani, 2018; Ortiz-Urbina et al., 2020; Khwarahm et al., 2021; Vignali et al., 2021; Jha et al., 2022). In our view, an effective way to formalise the biotic component of an ecological niche consists in including ENM/SDM models of some species (prey, predators, phorophytes, etc.) as biotic layers in models of the studied species (Pshegusov et al., 2022). This approach is particularly important in the absence of reliable information on associated species.

The widespread influence of human activity on natural processes makes it necessary for SDM/ENM models to include anthropogenic factors as critical predictors of species distribution. The common method of incorporating human activity into models involves estimating distances from the studied object (e.g., nesting sites) to infrastructure or protected areas (Ortiz-Urbina et al., 2020; Vignali et al., 2021; Sharma et al., 2023). The most popular tool in this process is the Euclidean distance, which does not consider one of the most important factors for mountainous landscapes, i.e., the altitude gradient. Therefore, in this study we used the Path Distance measure (path_landuse and path_setlment) calculated with surface distance, horizontal straight-line distance and vertical factor (McCoy et al., 2001).

The advantage of formalising abiotic, biotic and anthropogenic factors in SDM/ENM models lies in the ecosystem approach to ecological niche modelling, considering environmental, trophic, competitive and other relationships of the studied species. This approach makes it possible to analyse the ‘geographical extent of the realized ecological niche’, which corresponds more closely to the real spatial distribution of the species (Soberón and Peterson, 2005; Peterson and Soberón, 2012). The application of this method, however, requires careful biological sense control of the modelling results.

In this context, we analysed abiotic, biotic and anthropogenic factors that could affect the nesting site distribution of four threatened Caucasian vultures. In particular, the aims of this study were: (1) to investigate the ecological niches of vultures by formalising environmental factors and food resources (natural and anthropogenic) in SDM models; (2) to model the spatial distribution of the scavenger nesting sites in the Caucasus; and (3) to assess the similarity of vulture nesting-habitats and the overlap/differentiation of their ecological niches. We hypothesised that joint nesting and feeding caused high competition between Cinereous Vulture, Egyptian Vulture, Griffon Vulture, and Bearded Vulture in the region.

2.   Materials and methods

2.1   Study area and target species

The study covers an area of about 390,000 ​km² in the Caucasus region (between 38 and 47° N and 36–50° E), which includes seven large climate-orographic units (the North Caucasus and Transcaucasia (parts of the Greater Caucasus), the Ciscaucasia, the Colchis and Kura-Araks Lowlands, the Lesser Caucasus, and the Transcaucasian Highland) in the Russian Federation, Georgia, Azerbaijan, and Armenia (Fig. 1). Most of the area is mountainous with a warm summer continental climate (Dfb according to the Köppen-Heiger classification) from the foothills to the middle mountains, a cool summer continental climate (Dfc) and an alpine climate (ET) in the highlands (Pshegusov et al., 2022). A hot summer continental climate (Dfa) prevails in the Ciscaucasian plains, while the Western Caucasus, Western Transcaucasia and Colchis Lowland have a predominantly oceanic (Cfb) and humid subtropical (Cfa) climate. The prevailing climate in the southeast of the Kura-Araks Lowland and Transcaucasian Highland is cold semi-arid (BSk). In general, the dryness of the climate in the Caucasus increases from northwest to southeast (Pshegusov et al., 2022). Beech (Fagus orientalis), hornbeam (Carpinus betulus) and oak (Quercus spp.) forests are widespread in the foothills and middle mountains of the Caucasus. Spruce (Picea orientalis) and fir (Abies nordmanniana) forests distributed in the middle mountains and highlands of the western Caucasus and western Transcaucasia are replaced by Scots pine (Pinus sylvestris) and birch (Betula spp.) forests in the central and eastern Greater Caucasus. Oaks and junipers (Juniperus excelsa, J. oxycedrus, and J. foetidissima) compose the forests of the Transcaucasian Highland and south-eastern Lesser Caucasus. The plains, foothills and lowlands of the Caucasus are mainly developed for agriculture. The steppe and subalpine grasslands of the middle mountains and highlands traditionally serve as pastures.

Figure  1.  Geographical location and zoning of the study area. 1: Western Caucasus; 2: Central Caucasus; 3: Eastern Caucasus (parts of the North Caucasus); 4: Western Transcaucasia; 5: Central Transcaucasia; 6: Eastern Transcaucasia.

DownLoad: Full-Size Img PowerPoint

Cinereous Vulture is a large, well-recognised tree-nesting raptor, widely distributed in the Palearctic from Spain and Portugal in the west to China and Mongolia in the east (Costillo Borrego et al., 2011; Arkumarev et al., 2020; Lee et al., 2021). Despite its wide range, the IUCN Red List (BirdLife International, 2022) classified the species as ‘Near Threatened’. In European countries it is assessed as ‘Vulnerable’ (Costillo Borrego et al., 2011; Guerrero-Casado et al., 2013), in Turkey and South Korea as ‘Critically Endangered’ (Moreno-Opo et al., 2012; Kirazlı, 2016) and in Russia as ‘Rare’ (Red Data Book of Russian Federation, 2008). The surviving nucleus of Cinereous Vulture remains mainly in western Europe (Stoynov et al., 2019). Breeding sites also persisted in the Caucasus, where scavenger populations in Russia, Armenia, Georgia and Azerbaijan continue to decline due to habitat loss and food limitation (Gavashelishvili et al., 2006; Red Data Book of Russian Federation, 2008; Arslan and Kirazlı, 2022).

Egyptian Vulture is a medium-sized cliff-nesting scavenger widely spread in southern Europe, northern Africa, southern and western Asia (Farashi and Alizadeh-Noughani, 2018; Kabir et al., 2019). Species populations declined dramatically throughout the range (Arkumarev et al., 2020; BirdLife International, 2022) and it is categorized globally as ‘Endangered’ (BirdLife International, 2022). The North Caucasian population of Egyptian Vulture remains stable (Karyakin et al., 2018), whereas in western Transcaucasia, scavenger nesting sites actually disappeared (Belik, 2014; Tilba, 2014). It listed in the Red Data Book of Russian Federation (2008) as ‘Rare’.

Griffon Vulture is a colonial large cliff-nesting raptor distributed from the Iberian Peninsula to the Altai plateaus, with its core range in southern Europe (Xirouchakis et al., 2021). Although the scavenger population in Europe is now steadily increasing (BirdLife International, 2022), in the Caucasus the species still faces several threats, mainly human disturbance and food reduction, and is regionally classified as ‘Rare’ (Red Data Book of Russian Federation, 2008).

Bearded Vulture is a large cliff-nesting scavenger that occurs in Africa, West and Central Asia and Europe (Brink et al., 2020). Although globally ‘Near Threatened’ (BirdLife International, 2022), this raptor is assessed as ‘Critically Endangered’ in southern Africa (Brink et al., 2020), ‘Vulnerable’ in Europe (Vignali et al., 2021) and ‘Rare’ in Russia, including the Caucasus breeding range (Red Data Book of Russian Federation, 2008).

2.2   Nest data and environmental predictors

We used geographic records of nests and nesting site locations obtained during expedition surveys (expedition records), from the Global Biodiversity Information Facility (GBIF, 2022) (GBIF records) and reference sources (References data) (Table 1; Appendix Fig. S1). During expedition surveys in the central, western and eastern Caucasus, as well as in western, central and eastern Transcaucasia in 2003–2022, the localisation of nests or breeding colonies (in the case of Griffon Vulture) was recorded (Appendix Figs. S2 and S3). For the Lesser Caucasus and the Transcaucasian Highland, we used available data from GBIF and reference sources. The reference sources were current literature sources, supported by coordinates and photographs: Belik (2004, 2012); Komarov (2017); Ilyukh (2017); Karavaev and Potapenko (2018); Karyakin et al. (2018); Belik and Nasrulaev (2019, 2021); Belik et al. (2019); Perevozov (2020); Ilyukh and Shevtsov (2021); Karavaev and Vitovich (2021); Mnatsekanov (2022); Shevtsov and Ilyukh (2022) (Table 1). Geographic records from GBIF were checked for incorrect and shifted coordinates, inclusion in the study area. Records outside the natural range (zoos etc.) and observations outside the breeding period were also excluded from the analysis. As a result, the GBIF data was analysed only for Egyptian Vulture.

Table  1.  Nest data and species occurrence records used in the study.

Species	DOI from the GBIF	GBIF records	Expedition records	References data	Records in the analysis
Cinereous Vulture	–	0	7	7a–h	14
Egyptian Vulture	10.15468/dl.qg2j9c	152	8	12d, h, i, m	31
Griffon Vulture	–	0	10	10e, d, j–l	20
Bearded Vulture	–	0	4	6e, l–n	10
Scots Pine	10.15468/dl.ymbrx9	147	15	–	123
Alpine Chamois	10.15468/dl.4w7a8t	35	118	–	33
Caucasian Goat	10.15468/dl.72dcc4	187	64	–	123
Total		521	226	35	354
Notes: superscript letters denote the reference data, specifically: a–Perevozov (2020); b, c–Belik (2004, 2012); d, e–Belik and Nasrulaev (2019, 2021); f–Mnatsekanov (2022); g–Karavaev and Potapenko (2018); h–Karyakin et al. (2018); i–Shevtsov and Ilyukh (2022); j–Belik et al. (2019); k–Ilyukh (2017); l–Ilyukh and Shevtsov (2021); m–Karavaev and Vitovich (2021); n–Komarov (2017).

 | Show Table

DownLoad: CSV

We also collected occurrence data on wild ungulates, specifically Alpine Chamois (Rupicapra rupicapra) and Caucasian Goat (Capra caucasica), as the main natural food resources for the scavengers in the Caucasus (Belik, 2014; Karimov and Mammadov, 2019; Ilyukh and Shevtsov, 2021; Karavaev and Vitovich, 2021). Geographic records of Scots Pine (Pinus sylvestris), the preferred nesting tree species, were also used to identify Cinereous Vulture nesting sites (Table 1).

Using the ‘Clean duplicate’ function of the ntbox library in R (Osorio-Olvera et al., 2020), the occurrence records were checked for duplicates and sparse to one point per 1 ​km² grid cell. As a result, 354 geographic records were retained in the analysis (Records in the analysis) (Table 1).

Environmental predictors were sourced from the ENVIronmental Rasters for Ecological Modeling (Title and Bemmels, 2018; ENVIREM, 2022), one of the most effective datasets for studying biological objects (Adhikari et al., 2019; Tytar, 2021). ENVIREM comprises 16 climatic and two orographic variables, including the terrain roughness index (TRI), which were sourced as ready-to-use environmental layers from the dataset (Title and Bemmels, 2018). To minimise predictor collinearity and prevent model overfitting, we reduced the original set of 18 environmental variables by removing highly correlated predictors using the VIF (Variance Inflation Factor) test in R with a VIF threshold ≤3 (Table 2). The resolution of the resulting environmental layers was 1 ​km per pixel.

Table  2.  Noncorrelated ENVIREM variables in the analysis selected by VIF test (VIF threshold ≤3).

Variable code	Complete name, units	VIF
PETWettestQuarter	Mean monthly PET of wettest quarter, mm/month	1.78
embergerQ	Emberger's pluviothermic quotient	1.91
PETColdestQuarter	Mean monthly PET of coldest quarter, mm/month	2.06
TRI	Terrain roughness index	2.24
PETDriestQuarter	Mean monthly PET of driest quarter, mm/month	2.37
Notes: variable codes, complete names, and units followed Title and Bemmels (2018). PET denotes potential evapotranspiration.

 | Show Table

DownLoad: CSV

2.3   Modelling approach

In order to account for biotic and anthropogenic factors in SDM models, we considered building two types of models in Maxent, specifically A-models (abiotic) and BA-models (biotic-abiotic) according to the BAM concept of J. Soberón and A. Peterson (Soberón and Peterson, 2005; Peterson and Soberón, 2012). In the A-models we considered only five abiotic variables selected by the VIF test. In the BA-models we followed the same steps, but also used previously obtained SDM models of Scots Pine, Alpine Chamois and Caucasian Goat as biotic layers. We have previously used this approach in modelling the main forest-forming species of the Caucasus (Pshegusov et al., 2022), and more recently hierarchical modelling has been used to define the relative importance of prey resource distribution and abiotic factors on Harpy Eagle (Harpia harpyja) range limits in Central and South America (Sutton et al., 2023).

The availability of anthropogenic food resources (grazing livestock, landfills, dumps) was also considered in the BA-models. For this purpose, path_landuse and path_setlment factors were estimated via the Path Distance measure. Path Distance was calculated for each grid cell as the distance in metres to the nearest object (farms, pastures, landfills, settlements), taking into account elevation gradient (McCoy et al., 2001). The input data was represented by a spatial feature class from the NextGis vector map sets (NextGis, 2022) and the digital elevation model GMTED2010 (Amatulli et al., 2018). The resolution of biotic and anthropogenic variables was 1 ​km per pixel.

2.4   Model development and evaluation

In this study, we used Maxent (Baldwin, 2009; Phillips et al., 2017) as one of the most reliable, efficient and straightforward modelling methods based on presence-only data (Elith et al., 2006; Phillips and Dudík, 2008; Iverson et al., 2019; Komori and Eguchi, 2019). It proved robust with small sample sizes (Elith et al., 2006, 2011; Qin et al., 2017), up to three records for rare species (van Proosdij et al., 2016). Maxent provides reliable predictions when working with rare species whose absence data unreliably reflect habitat unsuitability (Vignali et al., 2021). Maxent software ver. 3.4.3 (Phillips et al., 2017) was applied in the R dismo package (Hijmans et al., 2017).

The best set of Maxent model parameters were defined in the R ENMeval package (Muscarella et al., 2014; Kass et al., 2021) using: (1) five replicates and 20% of occurrence data for model validating; (2) a default value of 10,000 background points; (3) features L, Q, H, LQ, and LQH; (4) regularization multiplier from 0.5 to 5 with 0.5 increments; (5) a maximum of 500 iterations. The application of five-fold cross-validation allows to obtain independent data sets, when all 100% of occurrence data are alternately included in 20% of test points (Phillips and Dudík, 2008). Splitting the original dataset into 5- or 10-fold is most common, since these values best fit the trade-off for obtaining not too high-variance nor too biased error rates (Guisan et al., 2017). The method avoids overfitting and underfitting of the model. The sample of 1000–10,000 background points should be sufficient to represent the full variety of ecological conditions in the study area (Ferrier et al., 2002; Phillips et al., 2006; Raes and ter Steege, 2007; Lissovsky and Dudov, 2020). Continuous Boyce index (CBItrain) (Boyce et al., 2002), Akaike's information criterion (AIC) (Akaike, 1974), the difference between AICc and its minimum value (deltaAICc), the area under receiver operating curve AUC from the training data (AUCtrain) (Fielding and Bell, 1997) and True Skill Statistics (TSStrain) (Allouche et al., 2006) were used to select the optimal models among the 50 models generated for each vulture species (ENMeval package in R (Muscarella et al., 2014; Kass et al., 2021)). To assess the Maxent model performance, we also calculated and evaluated the null models (Raes and ter Steege, 2007; Sillero et al., 2021) in the R ENMeval package (Kass et al., 2021). The calculation of null models was based on the Bohl et al. (2019) technique, which follows the original method (Raes and ter Steege, 2007) but evaluates null models on the same validation data that were used to estimate the species models. We used z-scores and p-values of differences to compare null and species models.

We calculated the percentage contribution and permutation importance of variables (Phillips et al., 2017) to analyse their significance in Maxent models. Optimal predictor values were sourced from the response curves illustrating the relationship between each variable and the probability of most suitable habitats. Since there is no uniform method for defining the habitat suitability threshold (Glover-Kapfer, 2015), in this study we used a fixed high-level threshold of 0.8 for optimal habitats. Such a threshold reduces the possibility of false-positives (Buhl-Mortensen et al., 2019). We used a threshold of 0.5 for potentially suitable locations.

For each model, distribution maps were created with a continuous scale of species occurrence probability from 0 (unsuitable habitats) to 1 (optimal habitats) in the Maxent default palette colour gradations. Distribution maps were generated by converting the Maxent output files into a netCDF files with visualisation in PanoplyWin software (PanoplyWin, 2021).

We used the Kernel density estimation (KDE) method (Blonder et al., 2014) to estimate the overlap of vulture ecological niches in environmental space. The KDE method was proposed to visualise ecological niches as point clusters in an n-dimensional space of environmental variables (biologically relevant independent axes), where the points represent acceptable values of these variables (Blonder et al., 2014). The method is conceptually simple, requires no absence data and suitable for SDM (Blonder et al., 2014). We applied Principal Component Analysis (PCA) axes, which combined the five ENVIREM variables previously selected by VIF test, as independent axes (complex environmental variables). PCA was run in the FactoMineR package in R (Le et al., 2008). The factoextra (Kassambara and Mundt, 2019) and ggplot2 (Wickham, 2009) packages in R were applied for axis extraction and visualisation.

3.   Results

3.1   Model selection and performance

Among the A- and BA-models created, 50 for each vulture species, the resulting models were presented in Table 3. Selection was based on the lowest AICc, deltaAICc and highest CBItrain, AUCtrain values, which identify models with the best balance between accuracy and complexity, and the best sensitivity and specificity in discriminating occurrence data from random data (Fielding and Bell, 1997; Boyce et al., 2002; Burnham and Anderson, 2002). The CBItrain values of the resulting models were 1 and the AUCtrain values ranged from 0.97 to 0.9, indicating the high predictive accuracy of the models. TSStrain did not decrease below 0.81 and, as reported by Leroy et al. (2018) for general situations, correlated with AUCtrain (Pearson's r for A- and BA-models were 0.98 and 0.91, respectively). Species (empirical) models differed from the null models (p ​< ​0.05) (Appendix Tables S1 and S2). Performance metrics (AUCtrain, AUCtest, AUCdiff, CBItrain and CBItest) significantly better for the species models, thus they demonstrates better performance than null models (Appendix Tables S1 and S2).

Table  3.  Predictive performance and settings of Maxent models with optimum complexity.

Species	Model	AICc	deltaAICc	CBItrain	AUCtrain	TSStrain	Feature type	RM
Cinereous Vulture	A-model	869	0.00	1.00	0.98	0.85	LQН	0.5
	ВА-model	120	0.00	1.00	0.99	0.97	LQН	0.5
Egyptian Vulture	A-model	1541	782	1.00	0.97	0.81	Н	0.5
	ВА-model	1724	49.5	1.00	0.97	0.93	LQН	0.5
Griffon Vulture	A-model	453	0.00	1.00	0.98	0.97	LQ	0.5
	ВА-model	358	0.00	1.00	0.99	0.99	LQН	0.5
Bearded Vulture	A-model	1059	23.7	1.00	0.98	0.97	LQН	0.5
	ВА-model	861	268	1.00	0.98	0.98	LQН	0.5
Note: Model performance was assessed by four indicators: AICс (Akaike's information criterion corrected), deltaAICc (the difference between AICc and its minimum value), CBItrain (continuous Boyce index from training data), and AUCtrain (AUC (area under the curve) values from training data); RM is regularization multiplier.

 | Show Table

DownLoad: CSV

3.2   Potential nesting areas and distribution maps

3.2.1   A-models

Modelling demonstrated that the nest and nesting site distribution of the studied vultures were influenced by orographic and humidity factors (Table 4). The most important variable in nesting-habitat suitability for all four species was the Terrain roughness index (TRI), with optimal values ranging from moderately (240–497) to highly (498–958) rugged slopes according to S. Riley classification (Riley et al., 1999). Among the precipitation variables, the mean monthly potential evapotranspiration of wettest quarter PETWettestQuarter contributed the most (20.5–34.1% contribution). The contribution of other abiotic predictors to the Maxent models did not exceed 10–13% including Emberger's pluviothermic quotient (embergerQ) with optimal values characterising subhumid and humid climates (Daget et al., 1988).

Table  4.  Variable importance represented as percentage contribution (PC) and permutation importance (PI) in the Maxent models.

Species	Environmental variable	A-Model				–	BA-Model
		PC (%)	PI (%)	OV	RT		PC (%)	PI (%)	OV		RT
Cinereous Vulture	TRI	53.7	29.8	300–650	/		13.7	12.9	300–650		/
	PETWettestQuarter (mm/month)	20.1	24.9	120–130	⁀		2.5	9.8	120–130		⁀
	PETColdestQuarter (mm/month)	12.3	15.5	15–20	⁀		0	0	15–20		⁀
	embergerQ	10.2	9.6	80–120	⁀		0.2	2.1	80–120		⁀
	PETDriestQuarter (mm/month)	3.6	10.3	18–20	⁀		0.4	0.5	18–20		⁀
	path_landuse (km)	–	–	–	–		10.2	9.2	20–30		⁀
	path_setlment (km)	–	–	–	–		9	6.3	10–15		⁀
	Scots Pine occurrence	–	–	–	–		29.6	30.3	0.8–1		/
	Alpine Chamois occurrence	–	–	–	–		20.4	14.9	0.6–1		/
	Caucasian Goat occurrence	–	–	–	–		13.9	14	0.6–1		/
	AUC ± SD	0.94 ± 0.03					0.92 ± 0.02
Egyptian Vulture	TRI	43.5	46.9	500–650	/		15.7	24.3	500–650		/
	PETWettestQuarter (mm/month)	34.1	27.2	125–135	⁀		6.4	4.2	125–135		⁀
	PETColdestQuarter (mm/month)	10.8	19.6	17–20	⁀		5.3	0.5	17–20		⁀
	embergerQ	5.9	2	100–140	⁀		1.4	2.5	100–140		⁀
	PETDriestQuarter (mm/month)	5.7	3.7	17–20	⁀		2.5	3.9	17–20		⁀
	path_landuse (km)	–	–	–	–		37.4	24.7	5–15		⁀
	path_setlment (km)	–	–	–	–		19.3	22.3	3–10		\
	Alpine Chamois occurrence	–	–	–	–		6	4.3	0–0.1		–
	Caucasian Goat occurrence	–	–	–	–		6	3.4	0–0.1		\
	AUC ± SD	0.93 ± 0.03					0.88 ± 0.11
Griffon Vulture	TRI	55.4	42.2	250–650	/		17.8	15.5	250–650		/
	PETWettestQuarter (mm/month)	24.2	31.1	120–130	⁀		6	14.1	120–130		⁀
	embergerQ	13.2	7.2	80–100	⁀		1.1	2.5	80–100		⁀
	PETColdestQuarter (mm/month)	6.1	15.2	15–20	⁀		0.1	1.1	15–20		⁀
	PETDriestQuarter (mm/month)	1.2	4.4	10–20	⁀		3.4	2.2	10–20		⁀
	path_landuse (km)	–	–	–	–		24.5	22.4	20–35		⁀
	path_setlment (km)	–	–	–	–		20.3	19.3	5–10		⁀
	Alpine Chamois occurrence	–	–	–	–		10.6	10.9	0.3–1		/
	Caucasian Goat occurrence	–	–	–	–		20.2	10.1	0.3–1		/
	AUC ± SD	0.96 ± 0.02					0.96 ± 0.03
Bearded Vulture	TRI	52.4	34.2	200–650	/		10	11.8	200–650		/
	PETWettestQuarter (mm/month)	20.5	33.4	120–135	⁀		5	13.4	120–135		⁀
	embergerQ	10.1	2.2	75–140	⁀		0.2	0.2	75–140		⁀
	PETDriestQuarter (mm/month)	10.3	16.9	20–22	⁀		1.4	5.6	20–22		⁀
	PETColdestQuarter (mm/month)	6.7	13.3	15–17	⁀		0	0	15–17		⁀
	path_landuse (km)	–	–	–	–		11.7	10.4	10–25		⁀
	path_setlment (km)	–	–	–	–		5.9	4.9	5–10		⁀
	Alpine Chamois occurrence	–	–	–	–		34.6	30.7	0.6–1		/
	Caucasian Goat occurrence	–	–	–	–		21.3	23	0.6–1		/
	AUC ± SD	0.92 ± 0.02					0.92 ± 0.02
Notes: Optimal values (OV) of variables (0.8 thresholds of habitat suitability) and Response types (RT) were sourced from the response curves. Response types are: “/” = increasing, “⁀” = hump-shaped, “\” = decreasing, “−” = flat. The occurrence of Scots Pine and ungulates represented the species occurrence probability in optimal habitats for vultures.

 | Show Table

DownLoad: CSV

3.2.2   BA-models

In the ВA-models, the contribution of biotic and anthropogenic factors exceeded that of TRI and especially humidity variables (Table 4). The occurrence probability of Alpine Chamois and Caucasian Goat was the main contributor to the nesting-habitat distribution models of Cinereous Vulture (cumulative contribution of about 34%) and Bearded Vulture (46%). Scots Pine localisation was also an important nesting-habitat determinant for Cinereous Vulture. Optimal distances from this vulture's nesting sites to agricultural facilities and settlements were at least 10–20 ​km, with the combined contribution of anthropogenic factors to the ВA-model of about 19%. The cumulative contribution of path_landuse and path_setlment to the ВA-model of Bearded Vulture was about 18%, with a minimum optimal distance of 5–10 ​km.

Modelling demonstrated that livestock carcasses and food waste from dumps play a greater role in the distribution of Griffon Vulture colonies than wild ungulate availability (Table 4). The accessibility of anthropogenic food resources, formalised through path_landuse and path_setlment, was the factor that presented the highest contribution (45% in total) to the ВA-model of this scavenger. Although the cumulative contribution of wild ungulate availability reached 31%, there was too high a variation in the occurrence probabilities of Alpine Chamois and Caucasian Goat in the predicted vulture nesting sites (0.3–1).

Egyptian Vulture was the most synanthropic of the studied species, with a percentage contribution of path_landuse and path_setlment to the ВA-model of about 57%. The optimal distance from the scavenger nesting sites to settlements and livestock facilities was only 3–10 ​km and 15 ​km, respectively (Table 4). The role of wild ungulate availability in the localisation of Egyptian Vulture nest and nesting sites was insignificant.

3.2.3   Distribution maps

According to the A-models, Bearded Vulture had the largest area of potentially suitable and optimal nesting territories in the Caucasus (Table 5). Nesting habitats of the species were predicted throughout the mountain systems of the Caucasus, with a core range in the western Caucasus (Fig. 2A). The second range centre was located in the eastern Caucasus (the Intra-Mountain Dagestan). In the western and central Transcaucasia, potential nesting sites were concentrated in the middle mountains and highlands, narrowing to the highland line in the east.

Table  5.  Predicted nesting areas of vultures in the Caucasus from results of Maxent modelling.

Species	Acceptable areas (103 km2)			Optimal areas (103 km2)
	A-Model	BA-Model		A-Model	BA-Model
Cinereous Vulture	35.9	13.2		10.7	3.0
Egyptian Vulture	16.7	2.5		4.6	0.8
Griffon Vulture	22.6	11.4		9.5	4.4
Bearded Vulture	48.8	8.2		16.9	2.3

 | Show Table

DownLoad: CSV

Figure  2.  Predictive maps of Bearded Vulture nesting site distribution in the Caucasus by A-model (A) and BA-model (B). 0–1 is the probability of nesting site/nest occurrence.

DownLoad: Full-Size Img PowerPoint

The second most common scavenger in the Caucasus, according to the A-models (Table 5), was Cinereous Vulture. Optimal nesting territories for this vulture were predicted in the middle mountains and highlands (along large river valleys) of the North Caucasus, with a large range centre in Dagestan (Fig. 3A). Relatively small areas of nesting habitats were localised in the Transcaucasia and the western ridges of the Lesser Caucasus.

Figure  3.  Predictive maps of Cinereous Vulture nesting site distribution in the Caucasus by A-model (A) and BA-model (B).

DownLoad: Full-Size Img PowerPoint

The distribution of optimal nesting sites for Egyptian Vulture was focal, with three potential core ranges from foothills to middle mountains in the western, central and eastern Caucasus (Fig. 4A). The A-model predicted the minimum nesting area for Griffon Vulture (Table 5). The scavenger nesting habitats were concentrated from foothills to middle mountains and highlands (along large river valleys) in the central and eastern Caucasus, as well as in the western ridges of the Lesser Caucasus and in a narrow mountain belt of the eastern Transcaucasia (Fig. 5A). BA-models predicted a significant reduction in the area of suitable and optimal vulture nesting sites almost uniformly across the range predicted by A-models (Table 5; Fig. 2B–5B).

Figure  4.  Predictive maps of Egyptian Vulture nesting site distribution in the Caucasus by A-model (A) and BA-model (B).

DownLoad: Full-Size Img PowerPoint

Figure  5.  Predictive maps of Griffon Vulture nesting site distribution in the Caucasus by A-model (A) and BA-model (B).

DownLoad: Full-Size Img PowerPoint

3.3   Differentiation of the species ecological niches

PCA identified four complex factors (components) with a cumulative variation of about 87% in the variables (Table 6). There was some differentiation in the ecological niches of the studied species by Factor 1 (Fig. 6), which included TRI, the distribution of Alpine Chamois and Caucasian Goat and the climate variable PETColdestQuarter. The correlation between these parameters was probably due to their coordinated variability on the altitude gradient, whereby the probability of wild ungulate occurrence increases and potential evapotranspiration decreases with increasing terrain roughness. The climate variables PETWettestQuarter, PETDriestQuarter and embergerQ formed the second and third PCA components, respectively. By these factors, the vulture ecological niches almost completely overlapped. Significant overlap of ecological niches was also noted by Factor 4, which characterised the distance of nesting sites from livestock farms, pastures and settlements (Table 6; Fig. 6).

Table  6.  Principal component analysis (Varimax normalized) of variables in vulture nesting sites.

Variables	Factor 1	Factor 2	Factor 3	Factor 4
Caucasian Goat occurrence	0.92	−0.01	0.18	0.15
embergerQ	0.06	0.17	0.94	0.10
path_landuse (km)	0.20	0.09	−0.06	0.75
path_setlment (km)	0.17	0.06	−0.23	0.77
PETColdestQuarter (mm/month)	−0.81	0.30	0.16	−0.28
PETDriestQuarter (mm/month)	−0.40	0.79	0.29	−0.19
PETWettestQuarter (mm/month)	−0.20	−0.90	−0.15	−0.13
Alpine Chamois occurrence	0.90	−0.04	0.20	0.19
TRI	0.89	−0.03	0.29	0.20
Eigenvalue	4.58	2.15	0.62	0.59
Explained variance (%)	53.57	21.54	6.23	5.81
Note: factor loadings >0.7 were indicated in bold.

 | Show Table

DownLoad: CSV

Figure  6.  Differentiation in ecological niches of Caucasian vultures by the PCA factors.

DownLoad: Full-Size Img PowerPoint

4.   Discussion

Previous studies covered different aspects of the vulture species ecology, such as orographic, climatic and vegetation requirements for nesting (Mateo-Tomás and Olea, 2010; Zuberogoitia et al., 2014; Rodríguez et al., 2018; Belik and Nasrulaev, 2021; Karimov and Mammadov, 2019; Ortiz-Urbina et al., 2020; Jha et al., 2021; Arslan and Kirazlı, 2022; Mnatsekanov, 2022; Shevtsov and Ilyukh, 2022), main food resources (Donázar et al., 2010; Tauler-Ametller et al., 2017; Stoynov et al., 2019; Ilyukh and Shevtsov, 2021; Karavaev and Vitovich, 2021; Xirouchakis et al., 2021), scavenger relationships with anthropogenic factors (Xirouchakis and Mylonas, 2005; Karyakin et al., 2018; Mishra et al., 2018; Ortiz-Urbina et al., 2020; Jha et al., 2021; Arslan and Kirazlı, 2022; Shevtsov and Ilyukh, 2022). Although most of these studies, particularly in the Caucasus, were conducted using ground-based observations, our research can be considered in the context of previous works. We assessed the potential distribution of vulture nesting sites in the Caucasus through the ecosystem approach to formalising environmental factors in species distribution models. This provided new insights into relative importance of the abiotic, biotic and anthropogenic requirements of species in the region.

4.1   Nesting-habitat suitability prediction

4.1.1   Abiotic predictors

In our study, the main abiotic predictors in the vulture A-models were the terrain roughness TRI (moderately to highly rugged slopes) and the potential evapotranspiration PETWettQuarter (120–135 ​mm/month) (Table 4). Many published sources reported that the Caucasian vultures inhabit mostly severely rugged mountains. The vulture usually nests in the highest sections of forested rocky ridges or cliff tops (Belik, 2004; Gavashelishvili and McGrady, 2006; Karavaev and Potapenko, 2018; Belik and Nasrulaev, 2021; Karimov and Mammadov, 2019; Mnatsekanov, 2022). Griffon Vulture, Egyptian Vulture and Bearded Vulture, typical cliff-nesting scavengers, inhabit steep rock ledges, niches and grottoes (Karyakin et al., 2018; Belik et al., 2019; Parfenov, 2019; Belik and Nasrulaev, 2021; Ilyukh and Shevtsov, 2021; Shevtsov and Ilyukh, 2022). Slope steepness (altitude range) mainly determined the distribution of Cinereous Vulture nesting sites in the Iberian Peninsula and Turkey (Morán-López et al., 2005; Yamac, 2007; Moreno-Opo et al., 2012; Kirazlı, 2016; Ortiz-Urbina et al., 2020). For Griffon Vulture, it was one of the most important predictors of nesting site selection in the eastern Iberian Peninsula, northern Sardinia and Spain (García-Ripollés et al., 2005; Mateo-Tomás and Olea, 2010; Aresu et al., 2021). Topographic features of the mountains (cliff availability) contributed most to the nesting-habitat suitability models of Egyptian Vultures in Nepal, Iraq and Spain (Zuberogoitia et al., 2014; Khwarahm et al., 2021; Sharma et al., 2023). According to researchers, steeply-sloping areas provide vultures with low-energy take-offs and landings, good visibility for detecting victims at a distance and unimpeded movement across open spaces due to prevailing updraft of air (Yamac, 2007; Kirazlı, 2016; Ortiz-Urbina et al., 2020; Aresu et al., 2021; Jha et al., 2021). Rugged terrain is also supposed to offer the nests better protection from human disturbance and predators (Morán-López et al., 2005; Yamac, 2007; Ortiz-Urbina et al., 2020; Jha et al., 2021). Accordingly, for all four vultures, habitat suitability increased with increasing terrain roughness TRI (see response types in Table 4).

In our opinion, the revealed dependence of nesting-habitat suitability on water regime in the wettest (summer) period is of interest for the breeding biology of the Caucasian vultures. Summer precipitation can affect the success of rearing young and their departure from the nest. Although some studies concluded that no climatic variables influence vulture nesting site selection (Moreno-Opo et al., 2012), other research revealed a significant contribution of precipitation variables to habitat models of these scavengers (Vignali et al., 2021; Jha et al., 2022; Sharma et al., 2023). Morán-López et al. (2005) reported that rainfall can limit foraging time and affect the reproductive efficiency of Cinereous Vulture. This was confirmed in Georgia by Gavashelishvili and McGrady (2006) who stated that low annual precipitation resulted in better conditions for the vulture soaring and breeding. In northern Sardinia, the breeding success of Griffon Vultures was also negatively correlated to mean annual precipitation because of low hunting efficiency during heavy rainfall (Aresu et al., 2021). On the other hand, precipitation was shown to play an important role in determining plant biomass and drought-related livestock mortality in Iraq and northern India, correlating positively with the occurrence probability of Cinereous Vulture, Griffon Vulture and Egyptian Vulture (Khwarahm et al., 2021; Jha et al., 2022). In Nepal, Egyptian Vulture was also more common in regions with high annual precipitation, which delayed rapid carcass depletion and favoured scavenging (Sharma et al., 2023). In the Caucasus, PETWettQuarter values of 120–135 ​mm/month correspond to humid and subhumid areas where precipitation is sufficient for vegetation development and probably does not hinder successful breeding of vultures. Accordingly, for all four vultures, response types of water regime parameters, including PETWettQuarter, were hump-shaped (Table 4).

4.1.2   Biotic and anthropogenic predictors

Vegetation cover, prey accessibility and anthropogenic activities are important factors in the nesting habitat selection by vultures (Jha et al., 2021). Human disturbance (grazing, persecution, landscape modification, forestry practices, etc.) may cause vultures to leave their nesting territories (Jha et al., 2021) and offspring are thrown from the nests (Arslan and Kirazlı, 2022). On the other hand, vultures show plasticity in their ecological requirements, including anthropogenic influences (Moreno-Opo et al., 2012). In the Caucasus, livestock carcasses usually attract large concentrations of various avian scavengers near settlements (Gavashelishvili and McGrady, 2006). Therefore, settlement availability may correlate positively with the probability of vulture occurrence.

The accessibility of pine forests and wild ungulates contributed most to the BA-model of Cinereous Vulture (Table 4). Habitat suitability of this vulture increased with increasing occurrence probability of Scots Pine forests, Alpine Chamois and Caucasian Goat (Table 4). In the study area, this scavenger prefers to nest on broken and flat tops of old freestanding pines or trees within clustered pine forests (Belik, 2004; Karavaev and Potapenko, 2018; Belik and Nasrulaev, 2021; Mnatsekanov, 2022). The nests on birches (Karavaev and Potapenko, 2018), junipers (Karimov and Mammadov, 2019), hawthorn (Appendix Fig. S3), cliffs and ground (Karimov and Mammadov, 2019) were rather special cases. Cinereous Vulture tended to prefer Pinus spp. as nesting trees also in Turkey and the Iberian Peninsula (Morán-López et al., 2005; Yamac, 2007; Ortiz-Urbina et al., 2020; Arslan and Kirazlı, 2022). Mature pine trees have many thick branches and relatively flatter tops, providing optimal nesting sites for the vulture (Yamac, 2007). However, Cinereous Vulture is known to select old (large and tall) trees that can support a heavy nest and facilitate landing and take-off, regardless of tree type or density (Morán-López et al., 2005; Moreno-Opo et al., 2012; Guerrero-Casado et al., 2013; Jha et al., 2021; Arslan and Kirazlı, 2022).

In the Caucasus, Cinereous Vulture mainly forages in the subalpine landscapes of the Lateral and Rocky Ranges (Belik, 2014; Ilyukh, 2017), where its food resources include Alpine Chamois and Caucasian Goat. Thus, widespread poaching of wild ungulates at the turn of the 19th and 20th centuries significantly reduced the scavenger population in the North Caucasus (Mnatsekanov and Tilba, 1998). On the contrary, the establishment of high-mountain nature reserves in the second half of the 20th century contributed to a considerable increase in the population of wild ungulates and Cinereous Vulture in the region (Belik, 2014). In principle, the vulture can use different food sources depending on the availability of wild animals or livestock (Costillo Borrego et al., 2011). In the Iberian Peninsula, Cinereous Vulture diet included deer, lagomorphs, rabbits and wild birds, and when wild animals were scarce, sheep, goats, domestic pigs and birds formed the main diet of the vulture (Morán-López et al., 2005; Costillo Borrego et al., 2011). Livestock also play a role in the scavenger diet in the Caucasian mountain pastures (Belik, 2004, 2014; Ilyukh, 2017; Mnatsekanov, 2022). However, the cumulative contribution of path_landuse and path_setlment factors to the BA-model of Cinereous Vulture did not exceed 19%, with optimal distance values of at least 10–20 ​km (hump-shaped response type) (Table 4). The vulnerability of tree-located nests forces the scavenger to inhabit sites away from densely populated areas of the Caucasus (Belik, 2004). Sensitive to disturbance, Cinereous Vulture usually chooses nesting sites away from human presence throughout its range, which ensures greater breeding success (Morán-López et al., 2005; Yamac, 2007; Kirazlı, 2016; Ortiz-Urbina et al., 2020; Arslan and Kirazlı, 2022).

The wild ungulate distribution contributed significantly (about 46%, increasing response type) to the BA-model of Bearded Vulture (Table 4). Preferring fresh carcasses, this scavenger usually occurs in the highlands of the Caucasus, where the degradation proceeds relatively slow (Karimov and Mammadov, 2019). Here, carcasses of wild animals killed in avalanches represent an important food source for the vulture. For example, in the western Caucasus, bones of Alpine Chamois and Caucasian Goat were found in the nests of Bearded Vulture, and the birds flew to the shots of ungulate hunters (Karavaev and Vitovich, 2021). In the central Caucasus, wild ungulate carcasses also formed the basis of the nestling diet at the upper border of the middle mountains (Komarov, 2017). Finds of young goat skulls in the nests indicated a significant role of active hunting in Bearded Vulture feeding (Komarov, 2017). The results of the study by Vignali et al. (2021) also showed that in the Swiss Alps the probability of chamois and ibex occurrence contributed most to the habitat models of Bearded Vulture during the warm season. In the lowland steppes of the Caucasus, however, livestock and food waste constitute the main food sources for the vulture and its nest location depends largely on the presence of livestock farms, landfills, picnic grounds and catering facilities near settlements (Karavaev and Vitovich, 2021). Several pairs nested for a long time within 5 ​km of the resort town of Kislovodsk (Ilyukh and Shevtsov, 2021). According to Margalida (2008) and Mishra et al. (2021), Bearded Vultures preferred bones of medium-sized ungulates (sheep) and birds, which are considered easier to transport. On the other hand, this cautious bird usually nests in areas away from humans (Ilyukh and Shevtsov, 2021). Consequently, the combined contribution of path_landuse and path_setlment factors to the BA-model of Bearded Vulture was only 18%, with a minimum optimal distance of 5–10 ​km and hump-shaped response type.

The great importance of livestock as a food source was revealed for Griffon Vulture, with a combined contribution of path_landuse and path_setlment factors to the BA-model of 45%. In the Caucasus, Griffon Vulture forages in subalpine and alpine meadows as well as in foothill steppes with intensive pastoralism (Belik, 2014; Ilyukh, 2017; Parfenov, 2019). At the end of the 20th century, a dramatic reduction in livestock numbers in highlands of the North Caucasus led to the disintegration of large scavenger colonies. Vultures relocated to the lowlands and made forced trophic migrations to the plains (Ilyukh, 2017; Parfenov, 2019). New colonies appeared in the vicinity of lowland villages and towns (Kislovodsk, Pyatigorsk, Khabaz, Lashkuta), where Griffon Vulture currently feeds at dumps, landfills and pastures (Belik, 2014; Ilyukh and Shevtsov, 2021; Mnatsekanov, 2022). In Georgia, this vulture, previously feeding on wild animal carcasses, is now attracted by more reliable food sources, such as waste dumps near settlements where food is constantly supplied by humans (Gavashelishvili and McGrady, 2006). The BA-model of Griffon Vulture accordingly predicted an optimum distance of 5–10 ​km from nesting sites to settlements (hump-shaped response type) (Table 4). Similar observations were obtained in the Mediterranean region, northern Spain and Greece where Griffon Vulture depended largely on livestock and landfills (Xirouchakis and Mylonas, 2005; Donázar et al., 2010; Xirouchakis et al., 2021).

Egyptian Vulture appeared to be the most tolerant of human activity, which is consistent with previous studies in the Caucasus (Shevtsov and Ilyukh, 2022). The combined contribution of the path_landuse and path_setlment factors to the vulture BA-model was about 57% with increasing response type (Table 4). In the study area, the scavenger forages in the foothills and plains, usually close to nesting sites (Belik, 2014). It nests on cliffs around villages and towns, near picnic grounds and livestock farms (Belik, 2014; Karyakin et al., 2018; Shevtsov and Ilyukh, 2022). A significant part of Egyptian Vulture diet consists of landfill waste and road-killed animals (Aghababyan and Ananyan, 2011). Habitat suitability of Egyptian Vulture increased with decreasing distances from the nests to human settlements and livestock facilities with optimum values of only 3–10 ​km and 5–15 ​km, respectively (Table 4). The scavenger nesting territories do not overlap spatially with the range of highland ungulates (Alpine Chamois and Caucasian Goat) (flat response type). In Nepal, Turkey, Kurdistan Region of Iraq, and north-eastern Iberian Peninsula anthropogenic trophic resources were also a major factor in nesting site selection for Egyptian Vultures near human settlements (Tauler-Ametller et al., 2017; Khwarahm et al., 2021; Sandesh et al., 2022; Sharma et al., 2023). In India, it was observed that this vulture often nested near sparsely populated settlements and used waste and anthropogenic materials for feeding and nest building (Mishra et al., 2018). In contrast, in southern Spain Egyptian Vulture predominantly fed on small- and medium-sized wild vertebrates (reptiles, rabbits and other), despite the abundance of livestock in the region (Margalida et al., 2011). In Iran, it also nested in sites most distant from populated areas and transformed lands (Farashi and Alizadeh-Noughani, 2018). This discrepancy in survey results may indicate a high plasticity of Egyptian Vultures in relation to biotic and anthropogenic factors.

4.1.3   Distribution maps

Optimal nesting sites of Bearded Vulture were predicted mainly in the middle mountains and highlands throughout the Greater Caucasus and in the western ridges of the Lesser Caucasus (Fig. 2), which is consistent with previous ground observations. In fact, in the Greater Caucasus the scavenger nesting sites were observed at altitudes of 550–2200 ​m on the Main, Lateral, Rocky, Front and rarely Pasture Ranges (Edisherashvili, 2011; Komarov, 2017; Belik and Nasrulaev, 2021; Ilyukh and Shevtsov, 2021; Karavaev and Vitovich, 2021). In the Lesser Caucasus in Azerbaijan, vulture nests were found at an altitude of 1200–2570 ​m (Karimov and Mammadov, 2019). As reported by Karavaev and Vitovich (2021), the species distribution in the Caucasian foothills, despite suitable rocky areas, is limited by human disturbance. The predicted core range in the western Caucasus, according to the existing literature (Tilba, 2014; Ilyukh and Shevtsov, 2021), concentrated the main nesting group of Bearded Vulture in optimal habitats. The second predicted range centre in the Intra-Mountain Dagestan, with its extensive network of rocky ridges, was the only Caucasian region where the Bearded Vulture population remained relatively stable during the 19th–20th centuries (Belik, 2014). In western Transcaucasia, according to Tilba (2014), the vulture preferred the highlands and probably nested on the southern slopes of Mount Chugush. In central Transcaucasia, the Bearded Vulture nesting sites occurred also in the highlands near the southern portal of the Roki tunnel (Edisherashvili, 2011). In the middle mountains and highlands of the Lesser Caucasus, six Bearded Vulture nesting sites were observed, while in the Talysh Mountains this scavenger did not actually occur (Karimov and Mammadov, 2019), which also agreed with our models.

Optimal nesting territories for Cinereous Vulture were predicted in the middle mountains and highlands of the North Caucasus (Fig. 3), which is in line with ground observations. Vulture nests were mainly observed in the middle mountains of the Rocky Range and in the highlands of the Main and Lateral Ranges at an altitude of 1300–2200 ​m (Belik, 2004; Karavaev and Potapenko, 2018; Perevozov, 2020; Mnatsekanov, 2022). The core range predicted in the Dagestan mountains (Avar Koisu and Karakoisu river valleys, forested mountains of the coastal region) is considered an optimal area for Cinereous Vulture nesting in the Caucasus (Belik, 2004). Small areas of nesting habitats were predicted in the Transcaucasia and the western ridges of the Lesser Caucasus. Only two nests of Cinereous Vulture were found here in Azerbaijan (Karimov and Mammadov, 2019) and two nesting sites in Georgia (Abuladze et al., 2011). Predicted nesting sites of Egyptian Vultures were distributed unevenly in the Caucasus (Fig. 4), as in the Nepal mountains (Sandesh et al., 2022). As in the Canary Islands (Rodríguez et al., 2018), in the Greater Caucasus this vulture mainly inhabits low altitudes (in the foothills and mid-mountains). In the predicted eastern core range, where Egyptian Vulture was considered a common scavenger during the 19th–20th centuries (Belik, 2014), it currently nests at altitudes below 1000 ​m (Karyakin et al., 2018; Belik and Nasrulaev, 2021). Rarely the species was recorded above 2000 m (upper reaches of the Avar Koisu and Karakasu rivers) (Belik et al., 2011; Karyakin et al., 2018), which is also consistent with our models. In the central Caucasus (central core range), Egyptian Vulture nesting sites were restricted to the Pasture and Rocky Ranges at altitudes up to 1800 ​m. In the Ciscaucasian part of the Pasture Range, the vulture nests in the vicinity of the resort town of Kislovodsk at an altitude of about 1200 ​m (Shevtsov and Ilyukh, 2022). The third scavenger core range was in the western Caucasus (Fig. 4), where Egyptian Vulture used to inhabit Mount Sauber-Bash and around the town of Maikop, but has now disappeared (Belik, 2014). Due to permanent human disturbance, the species almost disappeared from the Black Sea coast in western Transcaucasia (Belik, 2014), where it now occurs only on flyways and during summer migrations (Tilba, 2014). Few nesting sites were predicted in the Lesser Caucasus and Transcaucasia (Fig. 4). Here, 15 nesting sites at altitude of 160–1850 ​m were recorded in Azerbaijan (Karimov and Mammadov, 2019), 11 nests in Georgia (Abuladze et al., 2011), and 52–56 breeding pairs in Armenia (Aghababyan and Ananyan, 2011).

In fact, the main nesting sites of Griffon Vulture occur in the middle mountains throughout the Greater Caucasus, with the highest density in the central Caucasus (Belik, 2014; Belik et al., 2019; Parfenov, 2019), which agrees with our models (Fig. 5). At least 13 vulture colonies are known here on the Rocky and Pasture Ranges (Belik et al., 2019). In Intra-Mountain Dagestan, where the eastern scavenger core range was predicted, Belik et al. (2011) counted five vulture colonies along the Avar Koisu, Karakoisu, Sulak and Shuraozen rivers. In the Russian Black Sea coast, with a low probability of colony occurrence, Griffon Vulture is considered a rare migratory species (Tilba, 2014). In Transcaucasia, nine vulture colonies were observed in Georgia (Abuladze et al., 2011) and five nesting sites in Azerbaijan at altitudes up to 1645 ​m (Karimov and Mammadov, 2019), which also confirms the model results.

4.2   Overlapping ecological niches of vulture species

Some differentiation in the ecological niches of the studied vultures by the first PCA factor (TRI, wild ungulate distribution and PETColdestQuarter) (Fig. 6) was probably due to different food preferences of the scavengers. As shown above, the distribution of nesting sites was strongly dependent on wild ungulate availability in Cinereous Vulture, to a lesser extent in Griffon Vulture and Bearded Vulture, and weakly in Egyptian Vulture (Table 4). Significant overlap of ecological niches by the second, third and fourth PCA factors (Fig. 6) indicated similar abiotic conditions in vulture nesting sites and similar nest localisation relative to anthropogenic facilities. Thus, except for the wild ungulate availability, environmental conditions in the vulture nesting sites were generally comparable. Belik (2014) also concluded that in the Caucasus Cinereous Vulture often nest close to Griffon Vulture colonies. Joint nesting of all four studied vultures was observed on the outskirts of Orota village in the Intra-Mountain Dagestan (Belik and Nasrulaev, 2021). Moreover, in the Caucasus, early-nesting Bearded Vulture may occupy old Griffon Vulture nests, while Griffon Vulture sometimes uses vacant Bearded Vulture nests (Akbayev and Tkachenko, 2001). Such similarity in nesting sites may result in competition for nests between vultures, which is not an unusual phenomenon (Mishra et al., 2018). For example, in the southern Iberian Peninsula, Cinereous Vulture, Griffon Vulture and Egyptian Vulture can compete for the same resources, despite different adaptations (Morán-López et al., 2005). In the eastern Pyrenees, Griffon Vultures were responsible for the majority of Bearded Vulture nest usurpation, and the competition level for breeding space was positively correlated with the colony size (Bertran and Margalida, 2002).

At the same time, no conflicts between cliff-nesting Griffon Vulture, Egyptian Vulture and Bearded Vulture were recorded in the Caucasus due to the abundance of suitable nesting sites (Karavaev and Vitovich, 2021; Shevtsov and Ilyukh, 2022). Furthermore, the studied vultures often forage together on carcasses or dumps, forming large scavenger congregations and eating different parts of carcasses (Gavashelishvili and McGrady, 2006; Belik, 2014; Ilyukh, 2017; Ilyukh and Shevtsov, 2021; Mnatsekanov, 2022). Bearded Vulture and Egyptian Vulture often prefer smaller carcasses (Gavashelishvili and McGrady, 2006; Margalida, 2008; Mishra et al., 2021), which also reduces intra-specific competition within the avian scavenger guild. Vultures can benefit from joint nesting by using the behaviour of other scavengers as a signal of disturbance/safety and food availability (Gavashelishvili and McGrady, 2006). Feeding on dumps with an abundance of anthropogenic food resources enables more vultures to forage simultaneously and also reduces scavenger competition (Arkumarev et al., 2020). Thus, in the Caucasus, the interaction between Cinereous Vulture, Griffon Vulture, Egyptian Vulture and Bearded Vulture, which nest and forage together, is a good example of successful niche sharing among species with similar environmental requirements. Despite similar foraging and nesting requirements, the studied vultures are not pronounced nesting and trophic competitors.

5.   Management implications

Our findings and predictive maps may prove valuable for regional management of the four studied vultures. The most suitable nesting sites in the species core ranges in the Caucasus should be prioritised to establish buffer zones where human disturbance is restricted or prohibited. Conservation of current Cinereous Vulture nesting sites must be targeted at pine-covered moderately to highly rugged mountain slopes with preserved wild ungulate populations within 20–30 ​km from settlements and 10–15 ​km from agricultural facilities. These sites are most common in the middle mountains and highlands of the North Caucasus, the Transcaucasian highlands and the western ridges of the Lesser Caucasus. Buffer zones around Bearded Vulture nesting sites should include rocky highland areas with wild ungulate concentrations or middle mountain cliffs within 5–10 ​km from livestock farms and human settlements. The highest density of such areas was predicted in the middle mountains and highlands throughout the North Caucasus and Transcaucasia, as well as in the western part of the Lesser Caucasus. Conservation of rock formations in the foothills, middle mountains and upper river valleys at a distance of 5–10 ​km from settlements and 20–35 ​km from livestock farms may help to preserve Griffon Vulture nesting habitats in the North Caucasus and East Transcaucasia. Wildlife managers should use a buffer zone of 3–10 ​km from settlements and 5–15 ​km from agricultural facilities in the foothills and middle mountains of the North Caucasus to protect Egyptian Vulture. Wild ungulate availability is not an important element in the conservation of this vulture.

As the ecological niches of the four scavengers overlap significantly in terms of abiotic factors and nest localisation relative to anthropogenic facilities, each of the vultures can be considered an umbrella species for the protection of other vultures. Territories of high conservation status for one vulture may therefore be buffer zones for other species. Vultures in buffer zones should be protected by controlling human access to these sites (hunting and forestry, sport, tourism, etc.) especially during the breeding season.

Authors’ contributions

RP and VC conceived and performed this study, analysed the data and wrote and edited the paper. Both authors read and approved the final manuscript.

Ethics statement

The study was based on field observation without any experimentation or harm to the studied birds.

Declaration of competing interest

The authors declare that they have no competing interests.

Acknowledgements

We are grateful to R.A. Mnatsekanov for his help in organising and conducting a series of field studies in the central Caucasus.

Appendix. Supplementary data

Supplementary data to this article can be found online at https://doi.org/10.1016/j.avrs.2023.100131.