Global patterns and predictors of avian population density

Aim: How population density varies across animal species in the context of environmental gradients, and associated migratory strategies, remains poorly understood. The recent influx of avian trait data and population density estimates allows these patterns to be described and explored in unprecedented detail. This study aims to identify the main macroecological drivers of population density in birds.


| INTRODUC TI ON
Population density is a key parameter in ecology and conservation, and understanding how density varies across species and under various environmental conditions has long been a central focus of macroecology (Currie & Fritz, 1993;Damuth, 1981Damuth, , 1987;;Gaston et al., 2000;Santini, Isaac, Maiorano, et al., 2018;Silva et al., 1997).However, estimating population density is challenging and time-consuming, limiting the sample of species for which density estimates are available (Santini, Isaac, & Ficetola, 2018).At wider scales, studies applying a macroecological lens to the determinants of population density across species have largely focused on the relationship between body size and the average density of populations (Damuth, 1981;Juanes, 1986;Peters & Raelson, 1984;Peters & Wassenberg, 1983), dedicating relatively less attention to other potential biological and environmental drivers, and their relative influence.Only recently has it been demonstrated that population density may be influenced by other traits (Santini et al., 2022;Silva et al., 1997) and environmental conditions (e.g., climate and land cover, Currie & Fritz, 1993;Howard et al., 2015;Santini, Isaac, Maiorano, et al., 2018).Clarifying the drivers of population density across species, time and space, remains a fundamental question in ecology and conservation.
To date, several studies have attempted to disentangle the factors determining population density in birds.Besides body mass, trophic niche (i.e., diet) is an important predictor of bird population density, with the lowest densities found in higher trophic levels, including carnivorous (Santini, Isaac, Maiorano, et al., 2018) and insectivorous birds (Silva et al., 1997).This pattern is driven by diet-related differences in assimilation efficiency (Silva et al., 1997), as well as differences in resource availability and distribution, with higher trophic levels being on average less densely populated than lower trophic levels (Schoener, 1968).Although maximum densities are theoretically constrained by resource or energy availability (Lawton, 1990;Stephens et al., 2019), empirical evidence suggests that bird population densities peak at intermediate levels of net primary productivity (Santini, Isaac, Maiorano, et al., 2018).This may reflect a trade-off between increases in productivity and an increased number of competitors in more productive areas (Davies et al., 2007;Santini, Isaac, Maiorano, et al., 2018).
Other environmental factors potentially influencing local population density are climatic conditions.Climate change has been recently considered a more important driver of change in breeding bird abundance in Europe than land cover change (Howard et al., 2020).
At a global scale, bird populations living in areas with high precipitation seasonality usually have lower densities (Santini, Isaac, Maiorano, et al., 2018), perhaps because the availability of resources fluctuates seasonally and less predictably.A less explored angle is the influence of species traits other than body mass on population density.Theoretically, behaviours that have evolved to increase fitness by securing reproductive success should give rise to higher population densities.For example, ground-nesting birds may suffer high predation rates, hence their productivity (number of chicks per pair) may be lower on average than in species nesting in raised nests or cavities, whose clutches are more hidden from predators (Roos et al., 2018;Slagsvold, 1982).The outcome may depend critically on habitat type, with ground-nesting strategies being most successful in forests and raised or cavity nests being most successful in open habitats (Martin, 1993).Additionally, bird species with cooperative breeding behaviour are expected to live at higher densities as breeders benefit from helpers by enhancing survival of offspring (Cockburn, 1998) and the breeders (Crick, 1992), or by enabling breeders to increase their number of reproductive attempts within a season (Russell & Rowley, 1988).Territorial behaviour is also likely to limit density irrespective of trophic resources, thus territorial species are expected to occur at lower densities than non-territorial species (Errington, 1956).Understanding and quantifying how the environment and species' traits modify local population density is important for assessing which are the main drivers behind bird densities, helping to anticipate further biodiversity losses and improving conservation actions.
Recently there has been a surge of studies claiming that ecology should shift from merely descriptive and explanatory, to a more predictive science (Currie, 2019;Mouquet et al., 2015;Travers et al., 2019), in line with the progression of other traditional disciplines (e.g., climatology).Quantifying the relationships between species traits, environmental conditions, and population-specific densities, and testing the model's ability to predict independent samples allows for a shift to a more predictive science (Currie, 2019).
Improving the accuracy of ecological predictions is thus imperative to ensure a robust understanding of ecological mechanisms under scrutiny.Additionally, taxonomic biases in available density estimates (i.e., few well-known species; Santini, Isaac, & Ficetola, 2018) can bias our understanding of macroevolutionary and macroecological patterns and processes (e.g., Nakagawa & Freckleton, 2008).
Filling data gaps is thus essential for reconstructing large-scale diversity patterns and for unbiased conservation assessment (Callaghan et al., 2021;Cazalis et al., 2022;Santini et al., 2022).Ecological models able to predict species population densities and how they respond to drivers of global change are also fundamental to anticipate further biodiversity loss and its implications for ecosystem function.
The ongoing mobilization of large quantities of ecological and phenotypic data (e.g., Tobias et al., 2022) allows previously overlooked relationships between population density and a variety of traits to be investigated, potentially improving the predictive power of models.
to generate population density predictions for all mammal species for their application in macroecological and conservation analyses.
In a similar endeavour, Callaghan et al. (2021) recently attempted to quantify abundance for 92% of all extant bird species using citizen science data, expert-derived density estimates and missing data theory (but see Callaghan et al., 2022;Robinson et al., 2022).While the latter can be informative for global assessments, the results are not intended for use at smaller scales.
Here, we aim to identify the main macroecological drivers of population density in birds.Birds offer a useful system for further exploration of the macroecological drivers of population density because they are relatively well known, with an array of traits now quantified for the vast majority of species (Tobias et al., 2022), and often the subject of local-scale density estimation (Santini, Isaac, & Ficetola, 2018).We combine multiple trait and environment datasets to test whether species-specific traits or local environmental conditions are better predictors of population density in birds while controlling for phylogenetic and spatial autocorrelation.We apply a Bayesian modelling framework and focus on previously unexplored traits such as habitat preferences, nesting behaviour, primary lifestyle, cooperative breeding and territorial behaviour, along with body mass and diet.In addition, we also quantify the effect of primary productivity, precipitation of the driest month, and minimum annual temperature (Table 1), and predict how the average population density of birds varies along these gradients globally.Finally, we provide average population density predictions and their variation for 9089 (>80%) bird species for their use in future macroecological and conservation analyses.

| Data collection
We extracted avian population density estimates from an updated version of the TetraDENSITY database (Santini, Isaac, & Ficetola, 2018), the most comprehensive database on population densities of tetrapods to date.These estimates are all derived from local-scale field studies, therefore representing 'ecological densities' (i.e., densities in the species habitat, not at a regional scale such as in Callaghan et al., 2021).This first sample included >18K estimates.From these, we retained only estimates whose geographic coordinates were less uncertain than 50 km and collected after 1970 to match the spatial resolution and the period represented by the environmental covariates (see below).We further excluded non-native species, and seabird families (i.e., Alcidae, Diomedeidae, Fregatidae, Hydrobatidae, Laridae, Pelecanidae, Pelecanoididae, Phaethontidae, Phalacrocoracidae, Procellariidae, Rhynchopidae, Spheniscidae, Stercorariidae, Sternidae, Sulidae) because their densities on land are only reported as densities in the nesting colonies.
Finally, in order to avoid spatial and temporal pseudo-replicates, we averaged all estimates of the same species collected using the same sampling method and in the same location (exact set of coordinates) across years.This reduced the dataset to 5072 estimates for 1853 non-pelagic species (~19.2% of the total non-pelagic species) in 158 families (87.3%) and 35 orders (94.4%; Figure 1).Density estimates referring to pairs (i.e., pair density) were multiplied by 2.
We subdivided population density estimates into migratory (n = 715) and non-migratory populations (n = 4357) using the geographic range polygons from BirdLife International (BirdLife International and Handbook of the Birds of the World, 2017).Hence, the classification was conducted at the level of populations, not species, meaning that density estimates for resident populations of migratory species were considered non-migratory, that is, the resident portion of the range of a migratory species only includes density estimates for the resident population, hence we classified those population density estimates of such species as non-migratory.Other density estimates from the non-resident portion of such species were classified as migratory.Moreover, density estimates for migratory species are more commonly conducted within the breeding range compared to the non-breeding range because estimates are generally based on territory-holding (singing) individuals (Gregory et al., 2004).Therefore, among migratory species, we only retained estimates collected within their resident (classified as non-migratory) and breeding range (classified as migratory).All data sources for the density estimates included in the dataset are listed in the Supporting Information (Appendix S1).

| Trait and environmental covariates
We identified five traits and four environmental variables that can explain the variation in population density across species and space (Table 1).Environmental variables were all calculated year-round for non-migratory populations (sedentary species and resident portions of the range of migratory species), and for the breeding season only for migratory populations.The breeding period was defined as the months between April and September for the Northern hemisphere and, between October-March for the Southern hemisphere.This is clearly a simplification as some species breed outside these periods, but the method provides an estimate of environmental conditions at the time of reproduction for the vast majority of species breeding at higher latitudes.We used Google Earth Engine (Gorelick et al., 2017) to produce environmental variables at 0.5-degree resolution (in WGS84 lat/long), calculated as averages of the annual layers between 1981 and 2019.This period was constrained by the temporal window available for climatic and productivity data and is slightly narrower than the period from which density estimates are calculated .However, at the scale of our analyses, we do not expect this difference to have any effect on our results (i.e., if data from 1970 to 1980 were added, it would have a negligible effect on estimates of environmental variables or their geographic variation).
Finally, we extracted the values from the environmental variables using the location of the density estimates for resident populations (using year-round layers) and migratory populations (using breeding period layers) separately.Trophic niche/ level Dietary differences are linked to resource distribution in space and time, as well as differences in resource accessibility and assimilation efficiency (Silva et al., 1997) Primary  , exposed;Martin, 1993;Roos et al., 2018;Slagsvold, 1982), increase the likelihood of failure (e.g., brood parasites), or be a limiting resource constraining population density (e.g., cavities) Lower densities in species nesting on the ground due to higher mortality, lower densities in species nesting in cavities as cavities could be limited, or lower density in brood parasites due to high likelihood of failure.

Tobias and Pigot (2019) with updates and corrections
Cooperative breeding Cooperation over breeding is promoted in highdensity populations with saturated ownership of high-quality territories, suggesting that cooperative breeding is inherently correlated with population density (Komdeur, 1992).In addition, helpers can enhance the survival of offspring (Cockburn, 1998) and breeders (Crick, 1992), or enable breeders to increase their number of reproductive attempts within a season (Russell & Rowley, 1988).
Higher density in cooperative breeders due to increased survival and reproductive success, as well as increased spatial packing of individuals (i.e., larger numbers of individuals per unit area in comparison with pair-territorial species).
Territorial behaviour Poses a behavioural constraint to maximum density (Errington, 1956) Negative

| Model fitting
We modelled population density as a function of trait and environmental predictor variables (Table 1) using Bayesian inference and a mixed effect Generalized Additive Model (GAM) with a Gaussian family error distribution (Wood & Augustin, 2002).Population density, body mass and precipitation of the driest quarter were log 10transformed to meet normality and homoskedacity assumptions in the models' residuals (obtained using average posterior draws).All predictor variables were scaled to mean equal to 0 and standard deviation equal to 1 prior to modelling.To account for phylogenetic relatedness we used a phylogenetic variance-covariance matrix based on the consensus phylogenetic tree used by Stewart et al. (2022) derived from the phylogenies in Jetz et al. (2012).The consensus tree was obtained by calculating the maximum clade credibility tree topology and applying a 50% majority rule using median heights on 150 random trees out of the 10K provided tree topologies provided on www.birdt ree.org, based on the Hackett taxonomic backbone (Hackett et al., 2008).All predictors were modelled as smooth factors to account for non-linearity, but limiting k (dimension of the basis used to represent the smooth term) to 3 to avoid overfitting.
Since we expected a difference in the effect of environmental drivers for migratory and non-migratory species, we also included an interaction between environmental drivers and migratory status.
We included Species as a random intercept to account for pseudo-replicates at the level of species, as some species had multiple density estimates incorporated in the model.Similarly, we included Location (exact set of coordinates) as a random intercept to account for multiple estimates collected in the same location.Finally, we also included the Sampling method category as a random intercept to account for systematic differences in population density estimates due to different methods.Sampling methods were broadly classified into incomplete counts (e.g., strip transect counts, point counts within pre-defined buffers), census (e.g., territory mapping), distance sampling (from transect or point counts), mark-recapture and mist netting (i.e., minimum number known to be alive) (Santini, Isaac, & Ficetola, 2018).
To avoid overfitting, we applied weakly informative priors using a normal distribution with a standard deviation of 10 for the intercept, and a standard deviation of 1 for all slope coefficients, thereby limiting the range to a plausible gradient of variation considering the scaled coefficients (Lemoine, 2019).We ran 10 MCMC chains with 5000 iterations each, using the first 1000 as warmups.To limit the storage of MCMC chains we applied a thinning of 5. We checked models for chain convergence and parameter identifiability both visually and using the R-hat diagnostic.We summarized the posterior distributions of coefficient estimates using 95% credible intervals.
For all categorical variable levels, we also reported the probability of direction, a measure of evidence that varies from 50% to 100% and that indicates the probability of a coefficient being different from zero (Makowski et al., 2019).We ran the model using the 'brms' package v. 2.17.0 (Bürkner, 2017) in R v. 4.0.3(R Core Team, 2018).Since we had multiple residuals per location and species, residuals were aggregated per location and species by sum prior to the autocorrelation tests.

| Variable importance and population density predictions
We estimated variable importance by partitioning the variance explained across predictors (Santini et al., 2022).We first estimated a pseudo-R 2 by taking the square of the correlation coefficient between the median posterior prediction and the observed population densities used to train the model.Then, we repeated this for each predictor variable keeping it at its mean value for continuous variables.For categorical variables, we iteratively replaced all levels with the same level.For random effect variables, we replaced all random effect levels with a new level (not present in the training dataset), thus ignoring random intercepts of known levels by the model.We then subtracted the resulting pseudo-R 2 of predictors from the first pseudo-R 2 to obtain a measure of variable importance (average pseudo-R 2 for the level tested in categorical variables).We expressed the relative importance of predictors as a percentage contribution to the pseudo-R 2 .
To map areas of environmental extrapolation (predictions be- Note that the choice of the categorical level was irrelevant, since the purpose of this analysis was inspecting the geographic pattern of migratory and non-migratory populations, and the level of the intercept level does not influence this comparison.Random level effects (Location and Species) were not considered in the prediction.

| Model validation
We validated model predictions against taxonomically and spatially independent sets of data (Roberts et al., 2017).For taxonomically independent samples, we tested against independent orders, families and species.We divided the dataset into 10-folds, each including different taxonomic orders, families or species.Then fitted the model on all but onefold (training dataset), and predicted using the left-out fold (test dataset).Hence, onefold approximately equated to 10% of the training dataset, which corresponds to ~3-4 orders, 16 families, and ~507.We then calculated the average Root Mean Square Error (RMSE) for all 10-folds.To obtain spatially independent samples, we divided the dataset into 10 × 10-degree spatial blocks with a checkerboard pattern, allocated them into 10-folds and repeated the validation procedure (130 blocks with data and 13 iteratively removed).We compared prediction (block-validation) errors expressed as RMSE (Equation 1) with the standard deviation (Equation 2

| RE SULTS
The R-hat convergence diagnostic for all coefficients was equal or very close to 1 (Table S1), indicating that the MCMC chains suc-
As expected, population density in birds decreased with body mass (Figure 3a) and is higher in species eating vegetable items than those with a meat-based diet (Figure 3b).In particular, species feeding on energy-rich vegetable items such as fruit, nectar and seeds showed the highest densities, invertivore and folivore diet showed intermediate densities, followed by vertivores and scavengers, which exhibited the lowest density (Figure 3b).Wetland species and species present in human-modified habitats showed the highest densities, whereas species in rocky substrates showed the lowest (Figure 3c).Primary lifestyle did not show clear differences (Figure 3d).Non-cooperative breeders showed lower density than cooperative breeders.Brood parasites had lower densities than nesting birds, while we did not find differences among different types of  3f).Territorial species had lower densities than nonterritorial species (Figure 3g), and migratory populations had lower densities than non-migratory populations (Figure 3h; Table S1).
Environmental variables yielded different effects on migratory and non-migratory populations.Primary productivity showed a weak and uncertain relationship with density (Figure 3i).Precipitation of the driest month was positively related to population density, with a stronger effect in migratory populations, whereas it was almost flat in non-migratory populations (Figure 3j).Minimum annual temperature only showed a strong positive effect in migratory populations, and a weaker positive effect in non-migratory populations (Figure 3k).

| Population density predictions
We found no extrapolation in body mass values (Figure S3), and very limited extrapolation in environmental values (Figures S4- S6).Geographically, the training dataset does not capture well very dry (part of Sahara), extremely cold and low productive (part of Siberia and Greenland) considering the year-round conditions (Figure S4).In the case of seasonal conditions of the breeding period for migrant species, the dataset does not capture well very dry and warm areas in the southern hemisphere (e.g., Namibia; Figures S5 and S6).All cells falling into these regions were excluded from predictions.

| Predictive accuracy
The intraspecific variability in log10-transformed population density is substantial, it has a median of 0.4 corresponding to 2.5fold dispersion.The interspecific variability in log10-transformed population density estimates is 0.86, corresponding to 7.2-fold variation.
The posterior uncertainty around median predictions is higher than the observed intraspecific variability in the density estimates (standard deviation of log 10 density estimates per species) and more similar to the overall variability across all species (Figure S8).
Prediction errors per species highly overlap, but the extent of the interval changes substantially across species (Figure S9).
The validation of the model predictions against taxonomically and spatially independent sets of data yielded prediction errors that are comparable to the posterior errors for species, families and order level extrapolation (Figure S8), indicating that the posterior uncertainty adequately captures the extrapolation errors, thereby accurately reflecting our level of confidence in the predictions.

| DISCUSS ION
Our analyses provide the most in-depth assessment to date of the potential determinants of avian population density, including biological traits that reflect lifestyle, reproductive behaviour and ecology, as well as environmental predictors that directly and indirectly affect resource availability and mortality rates throughout the year.Our findings show that average population density is highly influenced by both species' traits and environmental conditions.Specifically, we found evidence for patterns described previously, including the predicted negative relationship between population density and body mass, along with substantial effects of diet (Silva et al., 1997) and climate (Forsman & Mönkkönen, 2003;Howard et al., 2015;Santini, Isaac, Maiorano, et al., 2018).In addition, our results reveal previously unexplored relationships between avian population density and habitat type, brood parasitism, territoriality, cooperative breeding and migratory behaviour.

| Biological traits determinants of population density
Our results support previous studies that identified body mass as an important predictor of avian population density (e.g., Juanes, 1986;Santini, Isaac, Maiorano, et al., 2018;Silva et al., 1997).Log 10transformed population density does not scale linearly with the log 10 -body mass, showing a minor relevance in discriminating population density of small species, and higher importance in explaining the low density of large species.Nonetheless, while body mass is an important predictor of bird population density, it seems to be less important than in mammals (Santini et al., 2022), perhaps reflecting the smaller range and more right-skewed distribution of body masses in birds (Currie, 1993).As such, the body mass variable in birds is less predictive of population density than in mammals.Interestingly, as previously shown in Greenwood et al. (1996) and Silva et al. (1997), at equal body mass, birds are substantially rarer than mammals, forming a relatively minor component of endotherm biomass in ecosystems, e.g., the average density of 500 g mammals ranges from 1 to 600 individuals per km 2 (Santini et al., 2022) whereas the density of 500 g birds ranges from 0.3 to 30 individuals per km 2 .This difference may be explained by considering the substantially faster metabolism of birds (McNab, 1988).Previous studies that have compared densities and energy use in birds and mammals suggest that granivorous birds account for a minor fraction of the energy transfer in grassland ecosystems (Pulliam & Brand, 1975;Wiens & Dyer, 1977), whereas small granivorous mammals are quantitatively much more important, using 10 times more energy than granivorous birds (Silva et al., 1997).Conversely, low avian densities compared to non-volant mammals suggests that birds compose a minor part of endotherm biomass in ecosystems (birds: ~0.002 gigatons of carbon, mammals: ~0.007 gigatons of carbon, Bar-On et al., 2018) even if they are able to maintain a high number of species in local communities (Silva et al., 1997).
Our results regarding trophic niche categories support previous studies showing clear differences between different feeding strategies (Silva et al., 1997).As expected, carnivores and scavengers have the lowest densities, with the latter partly reflecting their current critical conservation status (Buechley & Şekercioğlu, 2016).The other differences among other trophic niche categories broadly reflect what one would expect based on energy content and availability (Silva et al., 1997), with nectarivores, frugivores and granivores having the highest densities, followed by omnivores and invertivores, and those feeding on other vegetable items.Cooperative breeders have higher population densities than non-cooperative breeders, in line with studies concluding that helpers can enhance the survival of the nestlings, increase the survival of the breeders, and allow them to reproduce multiple times within the reproductive season, therefore, enhancing their reproductive success (Cockburn, 1998;Crick, 1992;Russell & Rowley, 1988).Territoriality yielded a negative effect on density, supporting the expectation of lower population density in territorial species (Errington, 1956), presumably due to a carrying capacity primarily determined by intra-specific competition for space, rather than resource availability.Contrary to our expectation, species with different nesting behaviours did not exhibit clear differences in population density, with the exception of the lower densities detected in brood parasites, which lay their eggs in the nests built by other species.While this behaviour confers an energetic advantage in terms of avoiding parental investment, it is also associated with a high chance of failure due to rejection of the eggs, limits posed by the available nests, and nest guarding behaviours (e.g., in cooperative breeders) (Feeney et al., 2013;Soler, 2009), all of them potentially resulting in lower productivity and recruitment.
In addition, brood parasites may have natural constraints on population density because their reproductive strategy is most successful when they are relatively rare compared to their hosts.

| Environmental determinants of population density
We found that precipitation in the driest month and minimum temperature influence avian population density is in line with previous evidence (Forsman & Mönkkönen, 2003;Howard et al., 2015;Santini, Isaac, Maiorano, et al., 2018).Increased precipitation and temperature are linked to higher plant productivity and thus greater resource availability which in turn is able to support denser populations of consumers.Populations of endothermic organisms, including birds, are also potentially limited in colder environments by the increased energetic cost of thermoregulation.The effect of temperature was particularly important in resident populations that do not escape from the harsher season.Contrary to previous findings (Santini, Isaac, Maiorano, et al., 2018), there was no clear effect of primary productivity.Interestingly, the density of bird populations appears to be considerably more dependent on climatic conditions than mammal populations do (Santini, Isaac, Maiorano, et al., 2018;Santini et al., 2022).
Overall, environmental drivers delineate a strong latitudinal effect with decreasing population densities at high latitudes for resident species (Figure 5), which was also observed in Forsman and Mönkkönen for Europe (Forsman & Mönkkönen, 2003), and a peak of population density at intermediate latitudes for migratory species.This pattern may have important implications for conservation, because protected areas, if aimed at ensuring long-term population persistence, should be larger in low-density areas in order to protect approximately the same number of individuals per species (Clements et al., 2018;Santini et al., 2014;Williams et al., 2022).However, we caution against the conclusion that protected areas should, therefore, be larger at higher latitudes because the opposite pattern appears to occur in forested habitats, with most tropical forest bird species occurring at much lower population densities than the equivalent temperate or boreal forest species, thus requiring larger reserves to protect viable populations (Tobias et al., 2013).Further studies are needed to explore latitudinal gradients of population density in the context of species ecology and habitat preferences.
One factor to consider that is particularly relevant to birds is migration.The average population density of migratory species is lower than non-migrants, a difference perhaps accentuated by high mortality during migration, with individuals exposed to multiple hazards along migration routes (Newton, 2008).In addition, the tendency for migratory species to have lower population density may arise partly from anthropogenic factors, since hunting and land-use change are currently driving rapid declines of migratory populations in temperate regions (Burns et al., 2021;Rosenberg et al., 2019).
The importance of our location-level random effect indicates that a number of drivers acting at smaller scales than the one used in our analyses explain important deviations from the global pattern.
Population densities of individual species may change substantially across habitat types, different habitat structures, and landscape configurations (Dolman, 2012).For example, some species tend to be more abundant in the core habitat while others prefer the habitat edges.Equally, some species may show local adaptations to anthropogenic drivers producing higher densities in one part of a species' range compared with another.Compensation effects due to the presence or absence of competitor species, as well as top-down (and release) effects due to presence (or absence) of predators may explain important deviations from patterns predicted on the basis of climate and habitat variables alone (e.g., Crooks & Soulé, 1999).
Similarly, spatial population dynamics such as source-sink, mainlandisland or classical metapopulation dynamics can explain important temporary fluctuations in species' population density that cannot be captured by environmental variables (e.g., Brawn & Robinson, 1996).Populations can also fluctuate substantially from year to year, especially in unstable environments (e.g., semi-deserts; Jordan et al., 2017).Finally, humans can directly (e.g., Benítez-López et al., 2017) or indirectly (e.g., Benítez-López et al., 2010;de Jonge et al., 2022) influence bird population densities, which in turn can cause trophic cascade effects (Estes et al., 2011).Therefore, while our model depicts a global pattern, this should not be used to make conclusions on local abundances since a multitude of specific drivers can introduce substantial uncertainty at the local scale.

| Density predictions
We were able to produce species-specific predictions of average population densities for 9089 species (>80% of the species).Even

F
I G U R E 1 Taxonomic and geographic coverage of population density data.(a) Geographic distribution of study localities from which published density estimates (N = 4488) are drawn, covering a total of 726 species; (b) percentage species coverage per order (with raw sample size [number of species sampled/total number of species] given for each order).datausing the 'bayesplot' package v. 1.7.2(Gabry & Mahr, 2022).We tested for spatial and phylogenetic autocorrelation in the residuals using the Moran I with DHARMa package v. 0.3.3.0(Hartig, 2018) and Pagel's Lambda with 'phytools' package v. 0.7.70 (Revell, 2012).
yond the environmental domain of the training dataset) we ran a Multivariate Environmental Similarity Surface (MESS) analysis(Elith et al., 2010).Here, positive values represent interpolation areas with similar environmental variation (NDVI, minimum temperature and precipitation of the driest quarter) than those in the training dataset.Negative values indicate localities where at least one environmental variable is outside the range of environmental variation in our data set.Further, we tested the overlap of the distribution of the only continuous trait variable (i.e., body mass) with that of the dataset used for predictions.We predicted average population densities while considering variation due to environmental factors for a set of 9108 non-pelagic bird species with information on traits, geographical distribution and phylogenetic relationships (BirdLife International and Handbook of the Birds of the World, 2017;Jetz et al., 2012;Tobias et al., 2022).For every species, we calculated a grid-level estimate across all 50 × 50 km cells within their resident and breeding portion of the geographic range using BirdLife range polygons (BirdLife International and Handbook of the Birds of the World, 2017).Subsequently, we excluded all cells with environmental values beyond the environmental domain covered in the training dataset and averaged the predictions across all other cells per species to obtain an average median population density and average 95% prediction intervals for each species.From the 9108, we excluded 19 species whose range falls entirely in environmental extrapolation areas, yielding a total of 9089 species with predicted population density estimates and prediction intervals (Appendix S2).Finally, to assess how environmental gradients influence the abundance of birds, we predicted the change in average population density across geographic space due to environmental conditions only.To that end, we kept all trait predictor variables at their median value except for climatic and NDVI variables and predicted population density for all cells globally.We repeated the procedure for the levels 'migratory' and 'non-migratory', and kept all other categorical variables at one level (i.e., Trophic niche = Nectarivore, Territorial behaviour = none, Primary habitat = Forest, Primary lifestyle = Aerial, Cooperative breeding = 'Cooperative', Nest placement = 'Cavity').
) of the posterior distribution of the predictions, and with the observed intra-specific and inter-specific variability in the density estimates (mean standard deviation of empirical density estimates per species and standard deviation across all species, respectively).Although used for different purposes, RMSE and the standard deviation are both measures of dispersion.The standard deviation measures the dispersion posterior predictions from the mean prediction, whereas the RMSE measures the dispersion of residuals around the predicted value.Therefore, they are directly comparable as they are effectively based on the same formula: Where p i is the predicted population density of a species i, y i is the observed, empirical density estimate, and N is the number of empirical estimates used to test model accuracy.Downloaded from https://onlinelibrary.wiley.com/doi/10.1111/geb.13688by Readcube (Labtiva Inc.), Wiley Online Library on [18/05/2023].See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions)on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License Where μ is the mean population density of a species based on N posterior predictions, x i is a single prediction of density for the species i.This allowed us to assess the accuracy of our model predictions for extrapolations (i.e., predictions beyond the range of independent variables in the training dataset) and interpolations (i.e., predictions within the range of independent variables in the training dataset).It is worth noting that the taxonomic block validation does not account for spatial extrapolation, and vice versa, spatial block validation does not account for taxonomic extrapolation, so maximum errors might be under-estimated.On the contrary, only 8 and 224 species belong to the missing orders and families, respectively.A very few areas of the globe are environmental extrapolations (see Section 3) and the exclusion of cells outside the environmental domain mitigates this problem.Additionally, a block validation tends to overestimate the error since the reduced training dataset is necessarily less representative than the full dataset.Hence, uncertainty around the estimate of the predictive errors remains.
cessfully converged.Posterior predictions closely matched well the observed log-transformed density estimates indicating a good fit of the model to the data (Figure S1).Model residuals did not show phylogenetic (Pagels' Lambda < 0.0001, p-value = 1) or spatial autocorrelation (observed Moran I = 0.0053, expected Moran I = −0.0025,p-value = 0.35; Figure S2).
Migratory and non-migratory populations showed similar but, nonetheless, different geographic patterns originating from different F I G U R E 3 Partial response plots and coefficient distributions.Levels in brackets for the categorical variables are used as intercepts (dashed line) and other levels displayed as contrasts.Population density is displayed in log 10 scale, continuous predictor variables are displayed with standardized values.The black dot in the probability density plots indicates the median of the distribution, the thick line is with 66% interval, and the thin line is with 95% interval.Interval bands for continuous predictors are set at 95% interval.responsesto environmental gradients (Figures3i-k and 5).Densities tended to be higher in tropical areas and lower in temperate areas.In temperate regions, maximum densities are found at intermediate latitudes and lowest densities are in the arctic tundra.In both tropical and temperate regions, densities are highest in wet areas and lowest in arid areas (Figure5).

F
Population density predictions for 9089 bird species.(a) Shows population density predictions displayed on the bird phylogeny to highlight phylogenetic clustering.(b) Violin plot showing the distribution of density predictions per order, and average estimates per species used to train the model on top displayed as dots.

F
Geographic pattern of population density in response to environmental conditions in (a) non-migratory and (b) migratory bird populations.The predictions are obtained by fixing all trait predictors at their mean value, so it represents how the average density of birds varies with environmental gradients.Yellow-orange indicates high population density and purple-black low population density.Dark grey land masses indicate areas outside the environmental domain of the dataset.The right margins indicate the median latitudinal value.
though our dataset only covered about 20% of the species, it represented well the full variation in species body mass, and the vast majority of the environmental variation globally (FiguresS4-S6), with exceptions limited to areas that support very few residents and migratory species, respectively(Somveille et al., 2018).While the predicted mean population density varies substantially across species (0.1-116 individuals per km 2 ), 95% predictive intervals span between 0.03 and 3819 individuals per km 2 and overlap widely across most species.This partly reflects statistical uncertainty around fixed effects, but mostly depends on factors that are not captured by 14668238, 0, Downloaded from https://onlinelibrary.wiley.com/doi/10.1111/geb.13688by Readcube (Labtiva Inc.), Wiley Online Library on [18/05/2023].See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions)on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License

Body mass Minimum Temperature Precipitation of the driest month Migratory Diet Primary habitat Territoriality
Variable importance score expressed as percentage contribution of each predictor variable on the total variance explained.Trait predictors are represented in blue colours and environmental predictors in green colours.14668238, 0, Downloaded from https://onlinelibrary.wiley.com/doi/10.1111/geb.13688by Readcube (Labtiva Inc.), Wiley Online Library on [18/05/2023].See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions)on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License nests (Figure