The aim of our study was to identify the most important ecological factors influencing the breeding density of chukar partridges (Alectoris chukar) in Southeastern Bulgaria. Identified habitats were divided into two categories: 1) natural habitats with an extant population of the species, and 2) habitats in which it was absent. In each habitat, annual precipitation, lowest elevation in the habitat, Paliurus cover, grazing livestock per 100 ha, shrub height, shrub density, bedrock type, percentage of pastures, percentage of agriculture, percentage of wetlands, percentage of rocks, presence/absence of buildings, presence/absence of highways, plant species composition, plant layer coverage, and chukar breeding pairs were measured. We have found several habitat variables that probably affect the density of nesting chukar partridge. These are: 1) grazing intensity, 2) percentage of rock outcrops, 3) percentage of grasslands, 4) shrub cover, and 5) elevation. Grazing intensity and rock outcrop presence, together with the elevation gradient, influence the breeding pair distribution and abundance most significantly. Our study has confirmed only three types of habitats that the chukar currently prefers to nest in. Unknown factors probably changed the other sites in a way that they are no longer suitable for the species’ existence and reproduction.
Introduction
In Bulgaria, chukar partridge (Alectoris chukar) is found in Mediterranean forests, shrub lands and rocky habitat types. This bird also inhabits agricultural areas and grasslands (Tucker & Evans 1997). As a whole, 70 European priority bird species nest in these habitats in Bulgaria and the chukar partridge is one of them. It is classified by IUCN as a vulnerable species (Stoychev et al. 2007).
Chukar partridge has been numerous in Southeastern Bulgaria but in the late 40s of the last century, due to intensive hunting, its decline became obvious (Patev 1950). In 1959, its population reduction in Bulgaria coincided with the establishment of Austrian pine (Pinus nigra) plantations throughout most of its range (Ivanov 1985). Population number around 1989 was about 75000 individuals (Simeonov et al. 1990). In 1999, chukar partridge population was estimated at 39000 individuals (Yankulov & Irgeva 1999). During the same period, according to two independent ornithological estimates, population size was 1500– 3000 pairs (Birdlife International 2004) and 2000– 3000 pairs (Nankinov et al. 2004). Between 2006 and 2007 chukar partridge density was one pair per 100 ha in the breeding period and 2.2 individuals per 100 ha during the non-breeding period (Gruychev unpublished). During the 1996–2005 periods there was a significant reduction or complete extermination of the chukar partridge in some of its traditional habitats. The reason for its disappearance is unknown (Stojchev et al. 2007). Some studies (Tucker & Evans 1997) point that land abandonment and afforestation are crucial for shrub and rocky habitats of the chukar partridge. In other habitats as most important factors are considered the abandonment of agricultural land, pesticide use, intensity of agriculture and predator density (Tucker & Evans 1997). Our study aimed to identify the most important ecological factors influencing the breeding density of chukar partridges in Southeastern Bulgaria.
Material and Methods
Field techniques
During the 2007–2011 periods, a survey of potential suitable habitat of chukar partridge in Southeastern Bulgaria was conducted. Identified habitats were divided in two categories: 1) natural habitats with an extant population of the species, and 2) habitats in which it was absent. As suitable were considered all possible habitats matching the description of Patev (1950) and Simeonov et al. (1990). During the breeding season chukar presence in the surveyed areas was identified by transect with sound source and dogs. Density of the breeding populations was counted using transects (Bibby et al. 1992). Then it was equalized to 100 ha area and that figure was considered as maximum number of birds identified. The minimal altitude was determined as the lowest point found in each habitat. Means for temperature and precipitation were taken from the nearest weather stations. For each subject area bedrock type was described. Bedrock types were classified as follows: 1 = limestone, 2 = quartz, 3 = tarras, 4 = marl with branches, 5 = muscovite. Rocks, pastures, farmland and wetlands were expressed as percentage of the total habitat area. Plant species composition was identified and plant layer coverage was described in 35 (20 × 20 m) sampling plots, evenly distributed along the transects at intervals of 400 meters on average. Species coverage in each plot was determined visually using the Braun-Blanquet scale. Then it was transformed into 9 point cover-abundance scale of van der Maarel (van der Maarel 1979). Average shrubs height in the 20 × 20 m sampling plots was measured too. Shrub density was recorded in percentages at 1 m height from the ground at the base of the shrub in 30 × 50 cm subplots (Bibby et al. 1992). Identified grazing livestock during the field observation in each habitat was taken as the maximum value of this variable per unit area (100 ha).
Mathematical and statistical methods, classification
In the current study TWINSPAN classification (Two Way Indicator Species Analysis) (Hill & Šmilauer 2005) was used. The basic idea in TWINSPAN is that each group of samples can be identified on the basis of indicator species, i.e. such species that prevail at the one side of the dichotomy. TWINSPAN gives the opportunity of processing qualitative and quantitative data. The software TWINSPAN not only classify the samples but produces two-way ordered data table (samples × species). In construction of TWINSPAN table, two-way weighted average algorithm of Correspondence Analysis (CA) (Hill 1973) was used. Combination of the two has made the method one of the most popular among the ecologists nowadays (van Tongeren 2004).
Ordination
As ordination method we used Redundancy Analysis (RDA) which is the canonical form of Principal Component Analysis (PCA). RDA is technique selecting the linear combination of environmental variables that gives the smallest total residual sum of squares. In RDA site scores are restricted to a linear combination of the environmental variables. The species-environment correlation in RDA equals the correlation between the site scores that are weighed sums of the species scores and the site scores that are a linear combination of the environmental variables (ter Braak 2004). An advantage of RDA is that the method can handle variables that are measured in different units. It is a linear method where species and environmental variables are represented by arrows and samples are given by points on the so-called biplot or triplot (if the species are shown too). Employment of RDA in our analyses is justified by the short ecological gradient of our data which was measured using the first ordination axis of Detrended Correspondence Analysis (DCA).
Regression and correlation
Trying to reveal the detailed relationships between the chukar partridge distribution and abundance and the environment, we used multiple regression and correlation techniques. Because correlated variables were not normally distributed we used non-parametric correlation coefficient of Spearman, Rs (Spearman 1904).
We used two methods for multiple regressions — General Additive Models (GAM) and Locally Weighed Regression (LOESS). The aim of GAM model is to maximize the quality of the dependent variable description, which may have various distributions. GAM does so by developing unspecified nonparametric functions of the independent variables, which are “connected” with the depended variable by link function (Yee & Mitchell 1991).
LOESS is one of the contemporary methods of modeling intended to overcome some disadvantages of the classical ones. It combines the simplicity of the least square techniques with the flexibility of the non-parametric regression. This is a procedure for developing response surfaces by multivariable smoothing. The dependent variable is smoothed as a function of the independent variables in a moving fashion analogous to how a moving average is computed for a time series (Cleveland & Devlin 1988).
Obtained results were compared and tested for statistically significant differences with non-parametrical (Mann-Whitney Rank Sum Test), Monte Carlo Permutation test, and parametrical (F-test; t-test) tests. In all statistical tests the significance level was P < 0.05.
In our analyses we used the following specialized software products: STATISTICA, version 8.0 (StatSoft, 2007); TWINSPAN for Windows (Hill & Šmilauer 2005); CANOCO for Windows, version 4.51 (ter Braak & Šmilauer 2003); CanoDraw for Windows, version 4.1 (Šmilauer 2003); SigmaPlot for Windows, version 11.0 (Systat Software Inc. 2008).
Results
Classification of sampled vegetation is shown in Table 1. TWINSPAN clustering of all vegetation samples (N = 35) produced six groups or plant community types shown in the second top row of Table 1 with their names, composed of the latin name of the two or three dominant shrub and herb plant species. Two of the groups contained too few samples impeding further statistical analyses and were excluded. Table 1 is also showing the statistical comparison of the four left TWINSPAN groups by the most ecological variables that we consider as highly important for the distribution and abundance of chukar partridge. We found numerous statistically significant differences between these four plant community types shown with bold in the rows of Table 1.
Fig. 1.
RDA ordination all vegetation samples, breeding chukar partridge pairs and habitat variables. Cen dif= Centaurea diffusa. Bro tec = Bromus tectorum. Car aca = Carduus acanthoides, Pal spi = Paliurus spina-christi, Poa bul = Poa bulbosa, Ery cam = Eryngium campestre, Anth odo = Anthoxantum odoratum, Dich isch = Dichantium ischaemum, Chry gry = Chrysopogon gryllus. Ach com = Achillea compacta.
Breeding bird pairs were found only in three of the habitats and there were no significant difference between them. All four important plant communities were dominated by the same thorny shrub (Paliurus spina-christi) and different grass and herb species. As was already noted, significant differences are found in most variables mainly between the first two and the last habitat. However, most obvious they are for the precipitation, elevation, grazing, pasture, agriculture, and rock variables. For example, first two habitats, where breeding chukar pairs are most abundant, are found at significantly higher elevation compared to the last one. Similarly, first habitat receives significantly higher quantities of rain compared to the last one as well as the grazing pressure is significantly higher in the first two habitats in comparison with the last one. On the other hand, the first habitat has significantly lower pasture percentage than the last three. The same difference holds for the rock outcrops but the rock percentage in the first habitat is significantly higher compared to the last three. Regarding the agriculture and wetland variables, the significant differences are only between the first two habitats.
RDA ordination of the 35 vegetation samples, breeding chukar partridge pairs and measured habitat variables are presented in Fig. 1. The vegetation samples are shown with different symbols corresponding to the TWINSPAN group to which they belong. Habitat variables and bird breeding pairs are symbolized with arrows.
The arrow of the breeding pairs is pointing at the same direction in the same half (left one) of the biplot with the grazing, precipitation, elevation, rocks and wetland variables, which means that they are positively correlated (Fig. 1, Table 3). Symbols of two plant communities (Paliurus spina-christi-Poa bulbosa and Paliurus spina-christi-Eryngium campestre-Anthoxanthum odoratum) predominate there, meaning that these are the habitat types that the breeding chukar partridge prefers. On the other hand, in the opposite direction are pointed the arrows of variables negatively correlated with the distribution and abundance of the breeding birds (see also Table 3), i.e. these are the variables determining the ecological conditions in habitats which the partridge avoids. Statistical parameters of the ordination analysis show that the extracted axes are highly significant (Table 2).
Table 1.
Plant community types, breeding chukar partridge pairs and main ecological variables in the studied area.
Table 2.
Statistical parameters of the RDA ordination.
Table 3.
Spearman rank correlation between partridge breeding pairs, habitat variables and first two ordination axes.
Generalizing, the highest density of chukar partridge in Southeastern Bulgaria is found in Paliurus dominated plant communities, located in habitats with higher elevation and annual precipitation, with intensive grazing, with more wetlands and rock outcrops, with absent or minimal agriculture, and far away from human buildings and highways.
Trying to be more accurate about the power and direction of the relationships between the chukar partridge distribution and abundance we correlated most of the habitat variables recorded by us with the number of breeding pairs. All habitat variables were also correlated with the ordination axes aiming to reveal what latent gradient underlies the gathered data. Because all correlated variables had non-normal distribution we used Spearman rank correlation coefficient (Table 3).
First RDA axis correlates highly negatively with elevation, grazing, rock outcrops and annual precipitation and highly positively with buildings, highways, pastures and agriculture. Second ordination axis correlates highly positively with wetlands and annual temperature and highly negatively with pastures and shrub density. According to these results the highest variation in our data is explained by the elevation gradient, grazing intensity, rock outcrops and highway variables. On the other hand, breeding partridge pairs are highly positively correlated with elevation, grazing pressure and rock outcrops and highly negatively associated with the agricultural activities and proximity of highways (Table 3). Noted variables probably have the greatest importance on the birds' survival and reproduction and because of this in the following analyses we turn our attention mainly on them.
LOESS regression of breeding pairs on all other available data about the habitat variables and vegetation is illustrated on Fig. 2. The isolines represent the projected breeding pairs drawn by the LOESS model, which is highly accurate because the high regression coefficient (R 2= 0.788).
The graph shows that with the increasing of grazing intensity, elevation, annual precipitation, and rock outcrops the density of breeding partridge pairs will also increase. Conversely, with the increasing of highway proximity, agricultural activity, human building density and, probably, shrub height the likelihood for partridge breeding becomes close to zero. There is second probable pike associated with increasing pasture territory but this has not been supported by the correlation analysis (see Table 3).
Trying to be more precise in our claims we have undertaken further more detailed analyses of the spatial distribution of partridge breeding pairs according to some highly correlated with it variables (see Table 3) such as elevation, grazing intensity and rock outcrops (Fig. 3). Again, like the previous graph, here isolines represent the predicted partridge breeding pairs in the ecological space, defined by elevation and grazing pressure (left) and elevation and rock outcrops (right). The results of the two GAM models are highly significant which is confirmed by the F-test.
According to our model prediction, chukar breeding density will increase with increasing grazing intensity but will start to decrease after an optimum elevation of about 300–350 m. In other words the ecological response surface of the left graph is unimodal. Similar is the situation on the right graph where partridge distribution is defined by the elevation and rock outcrop gradients. According to the GAM model the density of breeding pairs should be maximal at elevation between 300–400 m with abundance of rock outcrops in the habitat.
Fig. 3.
GAM regression of partridge breeding pairs on the elevation, grazing (left) and rock (right) outcrop gradients.
Table 4.
Statistical juxtaposition of habitats where partridge is present and absent.
With the next analysis we tried to outline the most important variables for the presence of chukar partridge (Table 4). Doing so, we compared statistically two types of habitats: the ones where partridges were present with the ones where they were absent. These habitats where compared by all variables analyzed until this point, which we have measured at the two compared habitats.
There are only two habitat variables that have statistically significant difference between the two habitat types (Table 4). These are elevation and grazing gradients which must have principal influence on the survival and reproduction of chukar partridge in Southeastern Bulgaria. Alectoris chukar is supposed to be more abundant in habitats with elevation around 300 m and continuous grazing throughout the year.
Discussion
We have found several habitat variables that probably affect the density of nesting chukar partridge. These are: 1) grazing intensity, 2) percentage of rock outcrops, 3) percentage of grasslands, 4) shrub cover and 5) elevation. Grazing intensity and rock outcrop presence, together with the elevation gradient, influence the breeding pair distribution and abundance most significantly. Grazing intensity probably improves the habitat quality via herb composition maintenance, as well as providing low profile of vegetation cover. Intensive grazing prevents the shrub and tree species of assuming total dominance and turning the open habitat into woodland, shrub land or forest.
Grazing intensity, percentage of grasslands and shrub cover
With the increasing grazing intensity the number of chukar breeding pairs also increase. Sites with substantially smaller pasture area have significantly higher abundance of breeding chukars than those where the pastures cover larger areas. Lands where pastures cover smaller area require a higher grazing intensity but there is also higher number of nesting chukar partridges. However, the growth of breeding pair number influenced by the grazing intensity has its limit. For example, during the fieldwork in one of the habitats electric shepherd was built and consequently the number of grazing animals increased to 150 individuals per 100 ha. Despite the higher density accounted at first, thereafter chukar partridges have not been found in the area. Here, the higher number of grazing livestock most probably influenced chukar numbers negatively.
It is known that grazing and constant human presence in alpine habitats disturb gallinaceous birds and affect adversely their distribution (Bhattacharya et al. 2007). Overgrazing can change the structure and composition of the pastures followed by alteration in biodiversity and predator prey relationships (Rambo & Faeth 1999, Blaum et al. 2007). Insects that are valuable food for a number of animals also decrease their numbers with overgrazing (Rambo & Faeth 1999). Higher density of grazing animals leads to scarce grass cover, creating favorable conditions for accommodation of alien plant species (Chambers et al. 2007). In our case, the habitat with the highest nesting density of chukar partridge is occupied by 55–60 individuals of grazing livestock per 100 ha. The number of breeding birds per unit area decreases with increasing grazing livestock (1 breeding pair at density of 106.6 grazing individuals per 100 ha), but also decreases with decreased grazing (1 breeding pair at density of less than 10 grazing individuals per 100 ha).
Previous studies have found that intensive grazing affects adversely chukar density by changing composition of grass species, destroying shelters, nesting sites, and increasing predation (Gregg et al. 1994, Kirby & Grosz 1995). Several studies in the United States have confirmed that regular burning and intensive grazing reduced the density of nesting birds (Robbins et al. 2002, Svedarsky et al. 2003, Alexis & Powell 2008). Our results confirm these data, but also show that the lack of grazing can also be associated with the low density of chukars. Similar results were reported for other bird species in North Dakota (Salo et al. 2004).
Shrub cover is directly related to grazing intensity. Higher intensity leads to lower shrub cover and vice versa. However, the presence of shrubs determines the suitability of nesting sites, therefore highly disturbed habitats are positively correlated with lower nesting density. On the other hand, habitats with highest shrub cover are less accessible for the chukars leading again to low breeding density.
The density of grazing livestock in Bulgarian chukar habitats vary between 50–60 individuals per 100 ha and has positive impact on birds' nesting. However, maintaining favourable conditions for the breeding of chukar partridge is dependent on grazing livestock control and adequate management of rangelands.
Rock outcrops
Rock percentage in the habitat is the next factor that affects breeding density of the chukars. According to the model in Fig. 2, bird density increases with increasing rocks percentage in the habitat. Highest nesting density (two pairs) was found in the sites in which rocks cover 45.3 % of the area. Topography also plays an important role in the distribution of partridges. They prefer steep rocky slope because it offers easy escape from enemies flying down the slopes (Johnsgard 1973, Lee et al. 2003). In some parts of the United States, rock coverage in the chukar habitat reaches from 1/4 to 1/2 (Moreland 1950, Galbreath & Moreland 1953), which conform to our results.
Elevation
According to our model, chukar nesting density reaches its maximum at altitude between 300–400 m. In habitats from the low-elevation parts of species' range presence of partridges was not recorded. Chukar breeding at lower altitudes was reported in previous studies from Sakar mountain (Patev 1950), Strandja mountain (Milchev 1991), but later, Milchev (2010) lists the species as extinct from Strandja mountain. There are probably other environmental factors which led to the extinction of the species from these localities.
In our study wetlands had negligible influence on the nesting density (Figs. 1–2, Table 3). The presences of permanent wetlands have been pointed out as an important factor for the spreading of chukars in other parts of their range (Shaw 1971, Benolkin 1988, Lee at al. 2003). Some studies recognize the water presence as an important factor in the summer habitats of the species (Harper et al. 1958, Borralho et al. 1998, Larsen et al. 2010). In Bulgaria, however, permanent water sources are found across entire range of the chukar partridge. Most of them are small dams constructed by humans in the past for watering of livestock but today few of them are used for this purpose. They do not dry out through the year and provide enough water for the birds during the summer months. Some of these habitats are located along the shores of Kardzhali dam, Studen Kladenets dam and Ivaylovgrad dam. This ensures the year-round availability of water in them. For us, this is the reason why the presence of water does not play principal role for the chukar nesting density in the area.
Our study has confirmed only three types of habitats that the chukar currently prefers to nest in. Different factors changed the other sites in a way that they are no longer suitable for the species' existence and reproduction. Achievement of successful future recovery of the species is possible only if urgent measures for controlling the relevant factors are applied.
