Functional habitat suitability and urban encroachment explain temporal and spatial variations in abundance of a declining farmland bird, the Little Bustard Tetrax tetrax

. Species response to land use can be examined under a functional perspective, where habitats are described according to species´ resource dependencies. Distribution or abundance models based on resource availability rather than land use types can be more informative about the ultimate processes behind observed population or distribution trends. Habitat use may depend on resources available, as well as disturbances that affect accessibility to such resources. Increasing human presence and urban encroachment may thus alter the relationships between habitat suitability and species abundance. Using 10 years of field data, we investigated whether variability in Little Bustard ( Tetrax tetrax ) abundance was explained by functional habitat suitability (assessed through resource-based models) and urban encroachment. We found that spatial and temporal variations in Little Bustard abundance were explained by functional habitat suitability and avoidance of urban areas, but that the significance of each variable varied with spatial scale. Little Bustard abundance at each observation point significantly increased with local nesting but not foraging habitat suitability, and decreased with increasing proportion of urban areas. At larger spatial scales, temporal changes in Little Bustard abundance were highly significantly related to changes in foraging habitat suitability. Moreover, the positive relationship between foraging habitat suitability and Little Bustard abundance weakened as the proportion of urban areas increased, and almost disappeared when the proportion of urban areas was more than 5%. Our results underline the benefits of using resource-based models to better understand processes that relate animal abundance and habitat suitability, while simultaneously considering avoided elements of the landscape. 10 ans de données prises sur le terrain, nous avons cherché à savoir si la variabilité de l’abondance de l’Outarde canepetière ( Tetrax tetrax ) s’expliquait par l’adéquation fonctionnelle de l’habitat (évaluée par des modèles basés sur les ressources) et l’empiètement urbain. Nous avons constaté que les variations spatiales et temporelles de l’abondance de l’Outarde canepetière s’expliquaient par l’adéquation fonctionnelle de l’habitat et l’évitement des zones urbaines, mais que l’importance de chaque variable variait selon l’échelle spatiale. L’abondance de l’outarde à chaque point d’observation a augmenté significativement avec l’adéquation de l’habitat local pour la nidification mais pas pour l’alimentation, et a diminué avec la proportion croissante de zones urbaines. À des échelles spatiales plus grandes, les changements temporels de l’abondance de l’outarde étaient fortement liés aux changements de l’adéquation de l’habitat pour l’alimentation. De plus, la relation positive entre l’adéquation de l’habitat pour l’alimentation et l’abondance de l’outarde a diminué à mesure que la proportion de zones urbaines augmentait, et l’espèce était presque totalement absente lorsque la proportion de zones urbaines était supérieure à 5 %. Les présents résultats soulignent les avantages de l’utilisation de modèles fondés sur les ressources pour qu’on comprenne mieux les processus qui relient l’abondance des animaux et l’adéquation de l’habitat, tout en considérant simultanément les éléments du paysage qui sont évités.


INTRODUCTION
Habitat change is one of the most important drivers of alterations in species distribution or abundance. The potential effects of changing landscapes on species' dynamics have been frequently addressed through correlative models whereby species-habitat associations are estimated by relating presence or abundance to the availability of particular land cover types (Dormann et al. 2012). While these models have provided valuable information for understanding species ecological requirements and conservation threats (e.g., Carone et al. 2014, Phifer et al. 2017, Tarjuelo et al. 2020, they commonly fail when extrapolated outside the area or habitat conditions for which the model has been calibrated (Graf et al. 2006). Furthermore, even when significant correlations between species abundance or trends and land use changes are found, those associations do not always allow for an understanding of the processes underlying them (Buckley et al. 2010, Dormann et al. 2012. As an alternative to using structural land cover types in modeling approaches, land use-species relationships might be better examined under a functional perspective, where habitats are described on the basis of resource dependencies of a species or a group of species (Fahrig et al. 2011, Butler andNorris 2013). Developing distribution or abundance models based on resource availability rather than land use types is likely to allow more robust predictions, even under changing environmental conditions, and can be more informative about the processes behind observed population or distribution trends (Butler and Norris 2013); e.g., the relative importance of modified nesting versus foraging resources.
Nonetheless, habitat use may not just depend on resources provided by different land covers but also on disturbance levels that might affect accessibility of individuals to such resources. In particular, increasing human presence and urban encroachment in many landscapes may alter the relationships between habitat suitability and species abundance. Exposure to human-related disturbances may increase stress in wildlife, alter behavior, and modify habitat use (Arlettaz et al. 2007, Carrete and Tella 2010, Tarjuelo et al. 2015, Casas et al. 2016, Arroyo et al. 2017, Rabdeau et al. 2021. The presence of human infrastructures may also lead to avoidance of certain areas, even if they are covered by land uses that provide valuable resources (e.g., Eldegard et al. 2012, Díaz-Ruiz et al. 2016, Filla et al. 2017. Avoidance of those areas may result in changes in distribution and, ultimately, population sizes (e.g., López-Jamar et al. 2011). In a time when the presence of humans and their infrastructures is increasing everywhere, this should be taken into account to better understand the relationships between habitat suitability and species abundance (Torres et al. 2016).
Farmland is a particularly dynamic social-ecological system that under certain conditions can support rich biodiversity (Pain and Pienkowski 1997). However, the process of agricultural intensification experienced from the second half of last century has led to much higher agricultural production (Tilman et al. 2011) but also to a decline in farmland biodiversity, including birds (Donald et al. 2001, Benton et al. 2002, Gregory et al. 2019 Vermouzek 2019), which it is ongoing. There is a crucial need to assess the effect of current and future potential changes in agricultural landscapes on farmland biodiversity to identify contexts where both agronomic and ecological sustainability may be simultaneously attained (Firbank 2005). With this aim, it is important to have a better understanding of the processes responsible for declines in farmland bird populations that are associated with changes in land use or management in order to design solutions that are compatible with present-day agriculture (Henle et al. 2008). Following this perspective, Cardador et al. (2014) developed a resource-based model framework to estimate habitat suitability for target farmland bird species, according to information on two of their key resource requirements: diet and preferred vegetation height in foraging or nesting habitats. This approach provided habitat suitability estimates that correlated with local species distribution (Cardador et al. 2014) and allowed the impacts of several alternative farming systems on farmland bird habitat suitability to be anticipated (Cardador et al. 2015). Such a framework could also be used to explain variations in abundance of target species. Moreover, buildings and human infrastructures are increasingly abundant within farmland areas, as in many other habitat types, and this may alter the relationship between resource availability and species abundance patterns (Rabdeau et al. 2021). Exploring whether species abundance is related to patterns of habitat suitability assessed through resource-based models, and whether this relationship is modulated by human infrastructures in the farmland matrix, may help identify the most relevant processes behind abundance changes (e.g., loss of foraging resources versus loss of available places for nesting), which is a crucial step for developing more effective policies to protect declining species.
We explored the associations between habitat suitability inferred from resource-based models and abundance of a farmland bird of conservation concern, the Little Bustard (Tetrax tetrax). This species has suffered strong declines throughout its breeding range in western Europe (Morales and Bretagnolle 2021), which are associated with changes in agricultural land uses and practices that are thought to reduce availability of suitable nesting areas or food supply (e.g., Wolff et al. 2001, Martínez and Tapia 2002, although the specific processes behind the declines or the relative contribution of different processes are still not clear in all areas. Additionally, this species is known to be sensitive to human disturbance (Casas et al. 2009, Tarjuelo et al. 2015, Estrada et al. 2016, and increasing pressure from infrastructures in farmland is considered a threat to the species (Silva et al. 2021). The proximity of houses and built-up areas, beyond the direct effects of disturbance, may enhance the negative effects of habitat degradation by reducing access to available resources.
Using 10 years of field data on this species in a study area in central Spain (Casas et al. 2019), we evaluated whether spatial or temporal variability in Little Bustard abundance was explained by functional habitat suitability (i.e., variations in nesting or foraging habitat suitability), and whether urban encroachment modified the relationships between habitat suitability and the species' abundance.

Study area and species
The Little Bustard is a charismatic species of European agroecosystems, and its most important Western European breeding populations are in Spain (García de la Morena et al. 2018). The species suffered strong population declines throughout its range during the last century (Morales and Bretagnolle 2021) and has disappeared from many European countries (Schulz 1985). It is classified as "Near Threatened" worldwide and "Vulnerable" in Europe (BirdLife International 2015). At the end of the last century, major declines in the Little Bustard population occurred in the highly intensified farmland of western France (Jolivet and Bretagnolle 2002); at the beginning of the 21st century, a decline was also first noticed in the cereal and long-term fallow mosaic areas of southwestern Iberia (De Juana 2009, Delgado andMoreira 2010). A recent national survey indicated that populations in Spain had declined approximately 50% in 10 years (García de la Morena et al. 2018).
Our study was conducted in Campo de Calatrava (Ciudad Real, central Spain, approximately 38º54' N, 3º55' W) (Fig. 1) within a Special Protection Area. The area, which covers 54.66 km 2 , is flat to slightly undulating (590-685 masl); it is used primarily for the cultivation of dry cereal (> 75%) and to a lesser extent includes olive groves (Olea europaea, 3-5% of the area), leguminous crops (1-5% of the area), and vineyards (Vitis vinifera, approximately 3% of the area). Cereal is grown in a traditional way, which creates a mosaic of sown, ploughed, stubble, and fallow fields of different ages. Grape vines are grown either in a traditional way-free standing without any attachment (goblet-shaped)-or as trellised vines with partial irrigation. The area has been highlighted as a hot spot for steppe birds (Traba et al. 2007) and supports a significant population of breeding Little Bustards, which has declined steeply in recent years (Casas et al. 2019).

Little Bustard abundance and habitat data
Field surveys were conducted at the peak of Little Bustard breeding display activity (April) from 2002 to 2011. Each year, a survey of the whole study area was conducted over 2-7 consecutive days. Surveys were based on point count methods and were carried out in the first 3 hours after sunrise and the last 3 hours before sunset; the hottest central hours of the day when bird activity is lowest were avoided, as were adverse weather conditions (Bibby et al. 1992). Using binoculars, observers counted Little Bustards for a period of 10 minutes from points situated along tracks. Six highly experienced observers participated over the years, four of them in multiple years. Observation points (totaling 120) were regularly distributed throughout the study area (Fig. 1). Each year, 110 ± 8 (mean ± SD) of those points were surveyed (range 91-116). Little Bustards observed from each point were plotted on a map, and double counting among consecutive points was avoided (if a bird was spotted in the same place, it was considered to be the same individual as in the nearby point and was not counted again [see Casas et al. 2019 for more details]). We detected both males and females, but most observations were of males; females accounted for about 13% of the total 1489 Little Bustard observations over the years. Given their different detectability, we used only males for analyses. However, it is known that males of this species display at locations that include resources that are useful to females for breeding , and that female observations are more common in areas with higher numbers of males (Jiguet and Bretagnolle 2006), so our assumption was that observed male abundance at different scales is also a proxy of female abundance.
Land use (the cover of each agricultural plot within the study area) was determined in the field each year (except in 2007, the year for which land use data were not available) when conducting the surveys, and was then entered to geographical information system software (QGIS 2.8.0). We considered the following land uses: cereal (rain fed), ploughed fallow field (a field left as annual fallow, which was ploughed to avoid the growing of weeds), unmanaged fallow (a field left as annual or pluriannual fallow, which was not ploughed and thus had vegetation), traditional vineyard, trellised vines, olive groves, pastures, alfalfa, other leguminous crops (Vicia spp. and Pisum sativum), cereal mixed with leguminous crops, beetroot, maize, riverine vegetation, buildings (isolated constructions, mainly secondary houses with a garden, all of them surrounded by a fence), and roads and tracks.

Resource-based modeling
We followed the framework described and validated in Cardador et al. (2014) to estimate functional habitat suitability in the study area for different periods. This approach first identifies the pool of potential land uses based on their agronomic characteristics in a particular study area (in our study case, those mentioned in the previous paragraph). The framework consists of four steps: (1) construction of a matrix to describe the species' resource requirements (r) for each vital activity (i.e., nesting and foraging) and for each categorized time period, (2) quantification of resource availability in each land use and for each time period, (3) calculation of habitat suitability indices for each vital activity in each land use and for each time period, and (4) temporal and/ or spatial integration of habitat suitability indices to encompass temporal and spatial variation at the scales at which the vital activities occur. We adapted the model developed by Cardador et al. (2014), which considered vegetation height and diet as variables for estimating resources, by integrating an additional variable, vegetation cover, which is known to be important in determining habitat use by the Little Bustard , Lapiedra et al. 2011, Morales et al. 2013, 2014. The following provides details of each step: Step 1. We categorized resource requirements (r) for Little Bustards on the basis of preferred vegetation height and cover for foraging or nesting, or presence of a given food resource in the diet, based on published literature (Jiguet 2002, Silva et al. 2007, Morales and Traba 2009, Martinez 2011, Ponjoan et al. 2012, Bravo et al. 2017, Cabodevilla et al. 2021. Four vegetation height categories (0-25 cm, 25-50 cm, 50-100 cm, >100 cm) and five vegetation cover categories (0-5%; 5-25%, 25-50%, 50-75%, 75-100%) were defined to describe foraging and nesting resources related to habitat characteristics. For each vegetation height or cover category, the r i value assigned reflected an assessment of the capability of Little Bustards to use vegetation height or cover i (0 -not preferred; i.e., vegetation heights/cover never used or avoided; 0.5 -vegetation heights/ covered used occasionally; 1 -preferred; i.e., vegetation heights/ cover where a high proportion of the individuals forage or nest) (see Table A1 for values used). We considered nesting resources to analyze male abundance because males display at locations that include resources that are useful to females for breeding . For dietary resources, we considered four food types (plants, seeds, invertebrates, and vertebrates), with values that reflected an ordinal measure of the degree of preference for each food type by Little Bustards (0 -not used or very rarely consumed; 0.5 -consumed secondarily, when usual food is not widely available; 1 -a primary and frequent food resource). Values were derived for two periods, spring (April-June) and summer (July-September), to reflect temporal changes in resource requirements through the breeding season; e.g., invertebrates are more important in summer because they are essential in chicks´ diet (see Table A2 for values used).
Step 2. Resource availability data for both habitat and dietary resources are described in a habitat type × resource table (A = [a kj ]) for each time period. Each a kj value is a measure of the availability of resource j in habitat unit k. We used available information on agricultural practices applied to different cover types in our study area (i.e., sowing and harvesting dates, fertilizers used, irrigation, and ploughing), in combination with our expert knowledge based on 10 years of field surveys, to qualitatively describe the probability of a given land use k having a given vegetation height or cover j in each time period (0 -not possible; i.e., vegetation height/cover category never or very rarely present; 0.5 -rare; i.e., infrequent or marginal vegetation height/cover category; 1 -usual; i.e., dominant vegetation height/cover category) (see Table A3 for values used). We then transformed these values to relative frequencies by dividing the score of each category by the sum of scores of all categories in a given time period and land cover type so that the sum of all categories was 1. Values were derived monthly according to vegetation growth patterns and land management.
For dietary resources, a kj values indicate the relative abundance of resource j in cover type k in a given time period. We assumed that the abundance of dietary resource j in cover type k was inversely related to both the number of agricultural practices that negatively affect that resource (n) and the intensity of the production system (f). Specifically, we calculated relative food abundance for each resource as a kj = 1/(n · f + 1), where n · f + 1 was used to avoid infinite a kj values. For these calculations, we considered the effect of the main practices known to be directly related to food abundance: agro-chemical use, irrigation, ploughing, and harvest (Table A4). These practices can lead to a reduction in food supply for Little Bustards directly (e.g., reduction in weed availability through the use of herbicides) or indirectly (e.g., elimination through competition of many broadleaved plant species and invertebrates associated with them by stimulation of crop growth due to crop irrigation or fertilizer use). We used a scaling factor for production system intensity (f), based on expert criteria, considering yield (tonnes/ha) and the intensity of implementation of management practices. This scaling factor was set to 1 for fallow lands and pastures (this, however, did not take into account whether fallow land plots were subject to intensive management when cultivated in previous years, so this may be an overestimation). Expected food abundance was calculated for spring (April-June) and summer (July-September) ( Table A4).
Step 3. For each ecological requirement considered (i.e., dietary resources, foraging vegetation, or nesting vegetation) and time period, suitability s i of habitat type i was defined as the scalar product of the corresponding vectors of matrices A (availability) and R (requirements). In the case of foraging suitability, we calculated the product of the suitability derived from habitat characteristics and that from food abundance. Suitability values derived from vegetation characteristics were, by definition, bounded between 0 and 1. Similarly, suitability values associated with food abundance (which represent the cumulative suitability of all types of food resources, and thus can sum to more than 1) were bounded between 0 and 1 (namely, values above 1 were set to 1); given that a species can use complementarily different food resources within the species' trophic niche, the maximum value would indicate that food requirements are fully provided in a given land use and month (regardless of by which food type/s).
We calculated suitability values monthly according to the temporal resolution of vegetation height and cover data. However, we assumed that (1) food preferences remained constant across months within both spring (April-June) and summer (July-September), and (2) expected food abundance was also constant within each period. Nesting habitat suitability was calculated only for those months when nesting activity occurs (April-June), whereas foraging habitat suitability was calculated for all breeding months (April-August). Table A5 presents the final suitability values for each land use/period. We attributed the foraging and nesting habitat suitability estimates for each plot in the study area in relation to the land use. Areas covered by non-farmland land uses (i.e., buildings, tracks or roads, riverine vegetation) were considered to have suitability 0.
Step 4. We calculated the suitability for each observation point and year as the average suitability in a radius of 400 m around each point (83% of Little Bustards were observed within that distance from the observation point). For this, information about suitability was rasterized (to pixels of 10 m). We generated two raster layers for each year: one for nesting habitat suitability, and the other for foraging habitat suitability. We then performed zonal statistics on each raster layer for each observation point buffer. Similarly, we calculated foraging and nesting suitability for the whole study area as the average value for all pixels within. Additionally, we divided the study area into four subareas of roughly the same size (Fig. 1), and we calculated foraging and nesting suitability for each subsector of the study area each year as the average value for all pixels within.
All indices were implemented in R 3.4.3 software (R Development Core Team 2018), and GIS manipulations were done with QGIS 2.8 (QGIS Association 2015).

Analyses of factors affecting variation in Little Bustard
Habitat suitability was calculated in relation to resources provided by each land use, and areas covered by land uses other than crops (buildings, tracks and roads, riverine vegetation) were considered to have suitability 0 for the Little Bustard. However, built-up areas, in addition to not providing resources for Little Bustards, may disrupt key behavioral processes for the species, thereby preventing the effective use of resources provided by nearby plots. We explored this by considering the variable "proportion of surface covered by buildings" (hereafter urban land) as an explanatory variable in models, in addition to nesting or foraging suitability. Tracks or roads may also locally influence Little Bustard occurrence and thus spatial distribution (Martínez-Marivela et al. 2018). Therefore, we also included the proportion of surface covered by roads and tracks (hereafter "tracks-roads") in the initial models for analyzing spatial variations in Little Bustard abundance.
To analyze factors that affected spatial variations in Little Bustard abundance (i.e., variations in abundance among observation points), we first used a generalized linear mixed model in R 3.4.3 (glmer in package lme4) (Bates et al. 2015), with the number of Little Bustard males counted at each observation point during yearly bird censuses as the dependent variable with a Poisson error distribution (log link function). The model included two random factors: "observation point identity", to account for repeated measures at the same point in different years, and "year", to account for variations among years in variables not measured in this study. Nesting and foraging suitability, and the proportion of built-up areas and tracks-roads were included as fixed explanatory variables. We assessed collinearity among explanatory variables with the variance inflation factor. All variance inflation factor values were lower than 1.5, indicating lack of collinearity between them (Zuur et al. 2009). Spatial distribution of Little Bustards may be aggregated due to conspecific attraction (Morales et al. 2001, Martínez-Marivela et al. 2018; therefore, there may be a potential lack of independence among sampling points. We conducted a Moran´s test of the generalized linear mixed model residuals (with the DHARMa package), and it confirmed a significant (P < 0.00001) spatial correlation of abundance observations up to 2.5 km (Fig. A1). Therefore, we fitted a spatial regression model (fitme function in spaMM package), where a spatial component was included as a random term including the X and Y coordinates in a Matern correlation function (Rousset and Ferdy 2014). The model also included the two random factors and the four fixed explanatory variables previously described. We checked that the model fit was correct by inspecting plots of residuals.
Spatial distribution of Little Bustards may change interannually in relation to changes in land use associated with crop rotation (Morales et al. 2005). To evaluate temporal trends in abundance, it is advisable to analyze numbers at a larger spatial scale than the observation point to avoid the noise created by distribution changes among years (Morales et al. 2005). We therefore assessed abundance each year as the total number of Little Bustards observed in the study area. This rendered a data set of 10 data points (one per year) for analyzing trends, which did not allow us to adequately explore additive effects of several explanatory variables and their interactions. In order to further explore results obtained with this limited data set, we divided the study area into four subareas of roughly the same size ( Fig. 1), each of which included a sizeable number of observation points and differed in land use characteristics, and assessed abundance each year in each subarea as the total number of Little Bustards observed there. To analyze trends in abundance, we used generalized linear models (GLMs), and fit the response variable (number of males observed) with a Poisson distribution, and used the number of observation points (log transformed) as an offset (to account for sampling effort variations because not all observation points could be monitored each year, and to have a response variable that was comparable across subareas and years). The abundance model at the study area level included year as a continuous explanatory variable. The abundance model at the subarea level included year, subarea identity, and their interaction as explanatory effects. Temporal trends in habitat suitability or urban land area were also analyzed with GLMs at the two spatial scales by fitting response variables to a normal error distribution and using the same explanatory variables as those specified in the models for temporal trends in abundance. Interactions were removed from the models when they were not significant in order to correctly assess the significance of direct effects (Zuur et al. 2009).
Finally, we tested whether abundance at the study area and subarea levels was related to habitat suitability. For this, we constructed GLMs with abundance as the response variable (Poisson distribution) and (log-transformed) number of observation points as an offset. When analyzing abundance at the study area level, we included foraging suitability, nesting suitability, and proportion of built-up areas as explanatory variables. When analyzing abundance at the subarea level, we included foraging suitability, nesting suitability, proportion of urban land, subarea identity, and the two-way interactions as explanatory variables. The latter model was over-dispersed (dispersion parameter > 2). We therefore fitted the response variable in this model with a negative binomial distribution (with function glm.nb). When interactions were not significant, they were removed from the models.
R 2 values of models were calculated with function rsq (Zhang 2020). The only exception was the spatial regression model, for which this function did not work and for which we could not estimate R 2 values.

Habitat suitability, urban land, and spatial variations in abundance
The abundance of Little Bustard males varied between 0 and 8 individuals per point (mean ± SD: 1.31 ± 1.54).
Little Bustard abundance at each observation point significantly increased with nesting habitat suitability, and significantly decreased with increasing proportion of urban land around each point, even when controlling for spatial autocorrelation (Table 1). No significant relationships were found with either foraging habitat suitability or the surface covered by tracks-roads (Table 1).

Temporal trends in abundance, habitat suitability, and urban land
Little Bustard abundance in the study area declined significantly over the study period (Chi 2 1 = 49.17, P = 0.0001), with numbers in 2011 approximately 50% lower than those in 2002 (Fig. 2a). During the same period, suitability of foraging habitats in the study area also declined significantly (Chi 2 1 = 24.6, P < 0.0001) (Fig. 2b), whereas no such decline was observed for nesting habitat suitability (Chi 2 1 = 0.65, P = 0.41). Values for nesting habitat suitability were, overall, lower than those for foraging habitat suitability (Fig. 2b). The proportion of urban land in the whole study area increased approximately 1% over the study period (from 2.28% to 3.34%; Chi 2 1 = 59.6, P < 0.0001) (Fig. 2c). The area and distribution of tracks and roads did not change with time during the study period (Table A5).
Generalized linear models showed significant differences in Little Bustard abundance and foraging and nesting habitat suitability among subareas. Abundance and foraging habitat suitability temporally declined in all subareas (although the rate of decline differed among subareas), whereas nesting habitat suitability did not (Table 2). Both average values and temporal trends in urban land differed significantly among subareas (Table 2).
Temporal changes in Little Bustard abundance in the whole study area were positively correlated with changes in foraging and nesting habitat suitability at that scale but not with changes in urban land (Table 3). At the subarea level, GLM results showed a significant positive relationship between Little Bustard abundance and foraging habitat suitability, even when taking into account the strong subarea differences in Little Bustard abundance, the slope of which (but not the intercept) was modulated by the proportion of urban land (Table 3); the positive relationship between foraging habitat suitability and Little Bustard abundance was marked when the proportion of urban land was small, but this relationship became much weaker as the proportion of urban land increased (Fig. 3). The additive effect of nesting habitat suitability was only marginally significant (Table 3).

DISCUSSION
Spatial and temporal variations in Little Bustard abundance were explained by functional habitat suitability and avoidance of urban land, but the relative importance of each variable changed as a function of spatial scale.
We observed significant spatial autocorrelation in the abundance of males at each observation point. This may be related to the reproductive system of the species, exploded leks (Morales et al. 2001), where social or behavioral aspects such as conspecific attraction may also be important in explaining the distribution of males. Nonetheless, even when taking this spatial autocorrelation into account, spatial variations in Little Bustard abundance were significantly related to variations in nesting habitat suitability and the occurrence of buildings in the farmland matrix. Male Little Bustards set up territories in areas containing resources that are potentially used by females to enhance the probability of encounter between the sexes , which may explain the influence of nesting habitat suitability on male abundance. It has been previously shown that male distribution relates to high food abundance ), but we did not find a relationship with foraging habitat suitability Table 3. Results of generalized linear mixed models explaining variations in Little Bustard annual abundance in relation to habitat suitability and proportion of urban areas at the study area and subarea levels. Results of type III Anova tests for each variable and model are presented. The subarea model initially included twoway interactions between variables, but they were removed when not significant for assessment of the direct effects. The interaction between urban land and suitability without the direct effect of urban land in Model 2 indicates a common intercept at zero suitability, irrespective of urban area, but differing slopes for the relationship between abundance and foraging habitat suitability under different proportions of urban area. (assessed for the whole of the breeding season). Also, and in contrast to another study that showed that the presence of tracks and roads was negatively related to the probability of Little Bustard male occurrence (Martínez-Marivela et al. 2018), we did not find a significant effect of tracks and roads on the spatial distribution of males. In contrast, our results indicate a strong avoidance of built-up areas at the observation point level. A negative effect of building density (but not track density) on nest distribution was also found for Montagu´s Harriers (Circus pygargus) nesting in farmland, modulated by female personality (Rabdeau et al. 2021). Negative effects of urban encroachment on Little Bustards have also been found in previous studies (Silva et al. 2004, Devoucoux 2014, although such effects have been absent or less marked in other studies (Martínez 1994, Faria andRabaça 2004), which suggests that the relationship may be context dependent. Our results showed that the negative effect of urban areas at this spatial scale was additive, which indicates that this effect goes beyond that of habitat loss occurring from urbanization. Little Bustards have been shown to be sensitive to disturbance caused by human activities (Casas et al. 2009, Tarjuelo et al. 2015, which altered their short-term physiological stress, as well as their behavior and habitat use. Therefore, an avoidance of built-up areas may be associated with the presence of humans, but built-up areas in farmland are also often associated with the presence of cats and dogs, the former behaving in semi-feral ways, even when they are domestic pets. Therefore, avoiding the proximity of built-up areas may also be a way of reducing predation risk. Variations in abundance at larger spatial scales through time were better explained by variations in habitat suitability than by urban encroachment. Annual male abundance trends at the study area level were explained by an additive effect of both nesting and foraging habitat suitability, although the effect size was stronger for foraging habitat suitability, which significantly declined through time, in contrast to nesting habitat suitability. A decline in nesting habitat suitability could have been expected due to the decline in fallow land after the abolishment of set-aside within the European Union´s Common Agricultural Policy in 2008 . However, in the study area, the decline affected mainly ploughed (managed) fallow, which is not suitable for nesting, whereas the availability of fallow with vegetation cover remained relatively stable (Table A5). The fact that no longterm decrease in nesting habitat suitability was detected could be related to the simplicity in describing this type of functional resource, as if a certain structural plant cover is maintained, the functionality would remain according to our model. On the other hand, and despite the simplicity in our approach, our results suggest that feeding resources have declined more over time, and support that this is an important process behind observed declines. The additive effect of both variables means that even if nesting habitat suitability is maintained through time (as happened in our study area during the study period), a decline in foraging habitat suitability may lead to temporal declines. In our case, a 15% reduction in estimated foraging habitat suitability during the study period coincided with an approximately 60% decline in Little Bustard abundance over the years. Although other factors not considered in our study may have also contributed to the decline of the Little Bustard population, the additive relationship we found means that measures aimed at enhancing foraging habitat suitability could potentially be very effective in helping Little Bustard populations recover (e.g., Bretagnolle et al. 2011) if nesting habitat suitability is at least maintained.
Results at the subarea scale, which allowed us to test for interactions between variables, also showed that the presence of urban land modulated the relationship between habitat suitability and species abundance: when the proportion of urban land was high, abundance was much lower than expected from foraging habitat suitability. This means that the local-scale avoidance of urban land may lead to an overall reduction in abundance, an effect particularly marked in good quality areas. The latter has marked relevance because those good-quality areas may be of utmost importance for maintaining the species (Inchausti and Bretagnolle 2005). Construction of built-up areas should therefore be regulated in the areas where Little Bustards are currently abundant in order to not accelerate the decline of this vulnerable species.
More generally, our results represent a validation for the use of resource-based models for inferring habitat suitability. Cardador et al. (2014) showed that describing the farmland landscape in the view of only two variables that described resources selected by a given species (vegetation height and food abundance) rendered habitat suitability values that correlated with occupancy probability for a suite of farmland bird species. We acknowledge that estimating availability of food resources based mainly on potential effects of management practices is a simplification, which should be re-assessed with information from future studies. For example, accumulated effects of agricultural intensification over the years may result in the homogenization of specific food resources across land uses (Tarjuelo et al. 2019). Regardless of this limitation, we showed that habitat suitability calculated with three simple variables (vegetation height, cover, and food abundance) allows for the production of habitat suitability indices that have biological meaning and explain variations in abundance of Little Bustards. Furthermore, this approach allows for the discussion of the relative importance of foraging versus nesting habitat suitability in explaining these variations in abundance, which will aid in the design of management measures that help reverse the species' population trend beyond the maintenance of favored specific land uses. Such an approach, however, may overshadow other important processes in the species' biology that are not related to resources used. Thus, our study shows the benefits of simultaneously assessing the influence of preferred resources with that of avoided elements of the landscape. Even if we have included the latter as an additive explanatory variable in our models, they could potentially be included in the resourcebased suitability model as elements that modulate the use of the available resources in each plot. This way, it would be possible to use those models to explore the effects of different land use change scenarios (including urban or infrastructure development) on a range of species in a manner similar to that used by Cardador et al. (2015).
Responses to this article can be read online at: https://www.ace-eco.org/issues/responses.php/2243   Table A3. Availability of vegetation height and cover categories in different land uses and months. A "0" indicates that a particular cover or height category is never present in that month for that land use. A "0.5" indicates that the category occurs sometimes. A "1" indicates that the category occurs frequently.  Table A4. Food availability values used in the resource-based model, and parameters used to calculate them. Specifically, we note for each land use and period (spring, Sp; and summer, Su) whether agricultural practices applied are likely to negatively influence food resources (seeds, plants, invertebrates and vertebrates) () or not (empty cells) (the latter may reflect the practice not being applied, or not leading to loss of food supplies). We then summarize the total number of practices that may negatively affect availability of different food types (n), and the intensification scale used (f). Availability of seeds, plants, invertebrates and vertebrates is calculated as 1 / (n • f + 1).