Global distribution and conservation of avian diet specialization

Ecologically specialist species are more prone to extinction than generalist species, yet the global distribution and conservation of ecological specialism is poorly understood. Here, we show that the global distribution of avian dietary specialization is roughly congruent with overall bird species richness for resident and breeding species, as well as for non‐breeding species. However, some areas harbour a higher number of diet specialist birds than expected given overall species richness (e.g. the Amazon, Gabon and Cameroon in Central Africa, extensive parts of Indonesia and some parts of northern Eurasia, Baltic cost and Mediterranean areas for resident and breeding birds, and tropical zone and south part of subtropical zone in South America for non‐breeding birds). These areas represent hotspots of avian specialization that need to carefully be considered in conservation strategies. We found that overall, 49.6% of resident and breeding species and 45.5% of non‐breeding diet specialist species are adequately represented by the global protected area system, but that this percentage is lower for the most threatened species. Policies that modify conservation planning approaches to include measures of specialization alongside other more traditional metrics of biodiversity could improve the protection of biodiversity in the face of rapidly accelerating anthropogenic threats.

In the last few decades, policies to identify priority locations for protected areas have focused mainly on locations characterized by high species richness or that contribute to representing a target proportion of species' geographic distributions (Maes et al., 2005;Fleishman et al., 2006). Yet species differ enormously in the extent to which they are vulnerable to threats, and full incorporation of this issue into conservation planning is still in its infancy. A useful starting point might be consideration of species' functional characteristics (Kosman et al., 2019). For example the degree of specialization of the species in a community could profoundly affect their resilience to threats, and thus influence the return on investment of protecting them. A species-trait approach focused on functional aspects of biodiversity could be adapted to evaluate specialization in communities for this purpose, and a conservation policy based on a multi-dimensional assessment of species assemblages, considering not only the number of species but also their functional traits, could be much more reliable than using simple metrics of species occurrences alone (Kosman et al., 2019).
Ecologically specialist species can be defined as those with a narrow niche, or using a restricted range of resources, whereas ecologically generalist species can in contrast exploit a broad variety of available resources (Irschick et al., 2005). Data on the breadth of resource exploitation in different species have become increasingly available over the last few decades, especially for birds (Reif et al., 2010;Clavero et al., 2011), and new tools for measuring specialization have been developed that can integrate multiple dimensions of species' traits, yielding a continuous index of specialization .
Specialist species are more prone to extinction than generalist species (Sekercioğlu et al., 2004;Colles et al., 2009;Clavel et al., 2011). Generalism can predict a species' capacity to respond to rapid environmental change, and adapt to disturbance (Hammond et al., 2018). Conversely, specialization leads to an increased dependence on a specific and limited range of resources (Begon et al., 2006), making the species more prone to conservation risks should there be natural or anthropogenic variability in the availability of those resources (Balisi et al., 2018;Chichorro et al., 2019). In the specific, some studies in mammals suggested that dietary specialization is associated with a higher risk of extinction, reducing species durations (Balisi et al., 2018). This negative relationship was attributed to a lower dispersion capacity of the most dietary specialized species (Balisi et al., 2018). Additionally, studies in bats and insects also found a higher risk of local extinction for dietary special-ists species (Biesmeijer et al., 2006;Boyles & Storm, 2007). Dietary specialists should be more sensible than generalists species to the loss of prey (Sierro & Arlettaz, 1997). For the same reasons, in birds, we can expect that dietary specialization and ecological plasticity could be associated to the extinction risk of species (Owens & Bennett, 2000;Sekercioglu, 2007;Ducatez et al., 2020).
Incorporating ecological specialization into IUCN Red List assessments has been proposed , but not yet implemented. However, consistent estimates of specialization for all species in any taxon are not yet available, and -most important -very little is known about the spatial distribution of specialization and how well existing protected areas are representing specialist species.
Here, we map avian diet specialization globally and assess the extent to which diet specialists are represented in the global protected area system. We (i) develop an index of diet specialization for the world's bird species, (ii) map the global distribution of avian diet specialization and (iii) assess how representation in protected areas varies along the gradient of diet specialization.

Global distribution of birds and protected areas
We created a 30 × 30 km 2 grid over the Earth (World Cylindrical Equal Area coordinate system, SR-ORG:8287), and overlaid polygons delineating the global distribution of 11,120 bird species (BirdLife International & Handbook of the Birds of the World, 2016). We removed uncertain, possibly extinct and extinct distributions, and additionally excluded passage, uncertain seasonal distributions from analyses, resulting in a set of 10,933 species for analysis. We created separate maps of the distribution of avian species using locations where species are (i) present year-round as resident and/or during the breeding season, and (ii) during the non-breeding season (the species is known or thought very likely to occur regularly during the non-breeding season; BirdLife International & Handbook of the Birds of the World, 2016).
Data on protected area boundaries were obtained from the World Database on Protected Areas (https://www. protectedplanet.net; Figure S1). Following standard procedures for cleaning the protected area dataset (e.g. Butchart et al., 2015), (i) sites with unknown or proposed designations were omitted, (ii) UNESCO Biosphere Reserves were omitted (Coetzer et al., 2014), (iii) sites represented as point localities were buffered to their reported area, (iv) boundaries were dissolved to prevent issues due to overlapping areas, and (v) slivers were removed (code available at https://github.com/jeffreyhanson/globalprotected-areas). Protected area data cleaning was completed using the World Behrman Equal Area coordinate system (ESRI:54017). Protected area data were then overlaid with a 30 × 30 km 2 grid covering the Earth. All spatial data processing was completed using ArcMap (version 10.3.1;ESRI, 2012), Python (version 2.7.8) and the R statistical computing environment (version 3.5.3; R Development Core Team, 2019) with the following R packages 'sf', 'raster', and 'fasterize'.

Mapping avian diet specialization
We derived a set of species traits describing diet specialization for each bird species (Table S1), expressing the composition of diet using percentages of 10 major food items (invertebrates, vertebrates (endotherm), vertebrates (ectotherm), vertebrates (fish), vertebrates (unknown), scavenger, frugivore, nectarivore, granivore, folivore), founded on data provided in a relatively recent publication (Wilman et al., 2014). For determining the diet of each bird species, the diet proportions were scored from secondary literature based on word order in sentences describing the diet. Thus, the diet data are based on semi-quantitative information assessing the relative importance of each item for the whole life history, characterizing a large portion of the 'Eltonian' niches of species (see more details in Wilman et al. (2014)). A similar type of data were used in a recent study focusing on bird trophic niche (Pigot et al., 2020). Following Morelli et al. (2019), we estimated the degree of diet specialization of each species using the Gini index of inequality. This index is a measure of statistical dispersion that ranges between 0 and 1, representing low to high specialization (respectively). Developed by Corrado Gini in 1921, it is probably the best single measure of inequality (Gastwirh, 1972), commonly used in studies of economic health of nations and recently adopted also in ecology and conservation studies (Barr et al., 2011;Morelli et al., 2019;Roeland et al., 2019). The Gini coefficient is estimated with the following formula: where 'x' is an observed value, 'n' is the number of values observed and ' ′ ' is the mean value.
In the specific case of avian diet traits used in this study, if every diet component has the same value, the index is 0, indicating the maximum degree of generalization possible. At the other extreme, a Gini coefficient of 1 indicates perfect inequality in the use of resources (i.e. high diet special-ization), in which a species has a diet represented by just a single component.
The geographical data and the available trait data used different underlying taxonomies, and to minimise the loss of species in the merged dataset, we sourced up to three synonyms for each species to facilitate matching. Using the diet specialization index for each bird species, we estimated two quantities to describe the avian community present in each 30 × 30 km 2 . Firstly, we estimated the richness of diet specialist species ('specialist richness') in the community by counting the number of bird species with a Gini diet specialization of 1 (species exploiting only a single type of item in the diet). Secondly, we calculated the 'community diet specialization index' (DSI) or mean diet specialization of the community as the sum of all diet specialization indices divided by species richness.
For the visualisation of each variable in the maps of resident and breeding, and non-breeding birds, we used the Jenks optimisation method based on natural breaks classification criteria.

Statistical analyses
We explored associations among diet specialist richness, DSI and overall species richness of resident and breeding, and non-breeding bird communities using spatially explicit Mantel tests (Valcu & Kempenaers, 2010). The Mantel test compares distance matrices among variables, using permutations of N rows and columns of each dissimilarity matrix. We compared distance matrices with differences between each variable among 30 × 30 km squares, by applying Monte Carlo permutations with 999 randomizations. Mantel tests were performed separately for each continent, using the package 'ade4' in R. Because the number of diet specialist species in a community is correlated with overall species richness, to reduce any effect associated with this spatial association, we used the residuals of a regression of specialist richness on overall species richness, a tool commonly used to highlight spatial patterns in ecology (Jetz & Rahbek, 2002;Caro, 2010). We ran a generalized linear mixed model between both variables, introducing continent as a random effect to account for possible consistent regional differences. Models were fitted by maximum likelihood using the package 'lme4' in R. The residuals of the regression indicate the areas where there are more or fewer specialists than expected by chance, considering the total number of species. Positive residuals indicate the number of diet specialist species above the expected, whereas negative values indicate areas where specialist richness is lower than expected considering overall bird species richness.
We assessed how well diet specialist bird species are represented by existing protected areas, considering the distribution of resident and breeding, and non-breeding birds. We compared the percentage of each species' distribution currently covered by protected areas to a target threshold percentage. First, we calculated the proportion of distribution range of each species for each category of threat, overlapped with the worldwide protected area network. Then, following standard practices (Venter et al., 2014;Butchart et al., 2015), targets for species with distributions smaller than 1,000 km 2 were set as 100% of their geographic range size, targets for species with ranges larger than 250,000 km 2 were set as 10% of their geographic range size and targets for species with intermediate range sizes were log-linearly interpolated between these two thresholds. Additionally, targets for species with geographic range size larger than 10,000,000 km 2 were capped at 1,000,000 km 2 . After calculating the targets, we identified species that are inadequately represented by existing protected areas by comparing the percentage of their distribution covered by protected areas to their species-specific coverage target. In other words, if the proportion of a species' distribution covered by protected areas was smaller than its representation target, it was considered inadequately represented by existing protected areas. Antarctica was excluded from all coverage statistics following standard practices. Because the percentage of conservation target met did not fulfil the normality and homoscedasticity requirements for parametric tests, we used a Wilcoxon test to compare whether the level of coverage by the network of protected areas is similar for diet specialist and non-specialist species.
All statistical tests were performed with R statistical programming environment version 3.6.0 (R Development Core Team, 2019).  Figure 1a). On the other hand, richness of non-breeding bird species peaked in Central America (mean = 78 species), Southeast Asia and Southern Asia (mean = 61 species for both areas), followed by Eastern, Western and Middle Africa (mean = 58, 46 and 43 species, respectively; Table 2 and Figure 2a).

Resident and breeding bird species richness peaked in
Of the 9,993 bird species assessed (those species that matched between maps of the global distribution and the dataset of traits), approximately one-third (34.16%) were classified as diet specialists, consuming only a single type of item. Diet specialization among the non-specialist species showed a wide range of variation, with a central tendency around 0.47, and SD = 0.19 (Table S1 and Figure S2). The spatial distribution of diet specialization at the global scale broadly corresponds with that of overall bird species richness (compare for example resident and breeding birds in Figures 1a and S3a), and diet specialist richness was spatially correlated with the overall species richness in all continents, for both resident and breeding, and non-breeding species (Table S2). Yet the pattern of residuals of the regression of diet specialist species on overall species richness highlighted some areas where there are more diet specialist species than expected from the total number of species, for both resident and breeding, and non-breeding species (Figures 1b and 2b). Areas with elevated diet specialist species richness for resident and breeding birds included the Amazon Basin, parts of Central Africa (e.g. Gabon and Cameroon), Southeast Asia (e.g. extensive parts of Indonesia) and some parts of the Northern Hemisphere (e.g. northern Eurasia, Baltic cost and some portion of the Mediterranean areas; Figure 1b). For non-breeding birds, the areas with a higher concentration of diet specialists than expected from overall species richness were located mainly in the tropics and the austral subtropics (e.g. in South America; Figure 2b).
DSI was spatially congruent with overall species richness and number of diet specialist species in the resident and breeding communities in all continents, with the exception of Oceania (Table S2 and Figure 1a-c). A different pattern of spatial correlation was observed for nonbreeding birds, with spatial incongruence mainly between the metrics of richness (overall and diet specialist species) and DSI in Europe and North America (Table S2 and Figure 2a-c). For resident and breeding birds, there was particularly high diet specialization (community DSI) in the far north (Arctic), in central America, New Guinea, Australia, some areas of sub-Saharan Africa and more comprehensively across South America than expected from specialist richness (Figure 1c). The percentage of diet specialist species in the resident and breeding avian communities showed in Figure S3 is complementary to the information about community diet specialization and is independent of the total richness of species in each spatial unit mapped. For non-breeding birds, the hotspots of high diet specialization are more clustered in some parts of northern Eurasia and North America (Figure 2c).
The proportion of the geographic range overlapping with protected areas, for diet specialist and non-specialist species, is showed separately for resident and breeding, and non-breeding birds in Figures S4 and S5 (Table S3). For non-breeding birds, only one of the seven diet specialist species classified as CR resulted adequately covered by the global network of protected areas: Gyps africanus (Table S3). On the other hand, overall, the F I G U R E 2 Global distribution of non-breeding bird species richness (a), residuals of the regression model between diet specialist richness versus overall bird species richness (b) and community 'diet specialization index' (DSI) (c). The values are presented in a coloured gradient from lower (dark blue), moderate (light yellow) to higher (red), with the ranges for each category based on a Jenks classification. The spatial unit is a grid of 30 × 30 km squares geographic distribution of 20 diet specialist species classified as CR remains completely outside the global network of protected areas (Table S3). Additionally, our results showed that resident and breeding specialist bird species most efficiently covered by the current global network of protected areas are the species classified as a lower threat (e.g. Near Threatened and Least Concern; IUCN & BirdLife International, 2017; Figure 3), with a similar pattern found for non-specialist bird species (Wilcoxon test statistic = 14, P = 0.834; Figure 3). Among non-breeding specialist birds, a high proportion of those classified as Near Threatened and Least Concern met their conservation target (Figure 3), and diet specialist birds classified as Endangered were actually more effectively covered by F I G U R E 3 Coverage of protected areas for resident and breeding bird species (above) and non-breeding bird species (below) for diet specialists and non-specialists, broken down by IUCN Red List Category. The bar-plots show the percentage of species which met the conservation target, for diet specialist (red) and non-specialists (dark blue). Critically Endangered (CR), Endangered (EN) and Vulnerable (VU) species are considered to be threatened with imminent global extinction. Other IUCN categories are as follows: Near Threatened (NY), Least Concern (LC), and Data Deficient (DD) the global network of protected areas than non-specialist species in their non-breeding distribution (Figure 3).

DISCUSSION
Based on a diet-trait characterization of species, we found that approximately one-third of the world's bird species can be classified as diet specialist, and thus potentially at higher risk facing an environmental change. Functional traits are especially useful for arranging species and communities along a specialization gradient, because such traits are the consequence of a series of diversification process of adaptive radiation (Castiglione et al., 2017). For this reason, the level of specialization of single species and overall species assemblages can provide an indication of the capacity of communities to deal with sudden climatic or ecological changes (Clavero et al., 2011;Cooke et al., 2019). Additionally, the degree of diet specialization is a measure that could be positively associated with other levels of ecological specialization in birds .
The main novelty of this study is the map of the global richness of diet specialist birds, highlighting a number of areas of high avian diet specialization, which could merit further consideration in the development of conservation policies. Globally, areas with high overall bird species richness of resident and breeding birds were in South America, some regions of Africa, Southeast Asia, Australia and New Zealand and Melanesia. However, the areas characterized by relatively high number of non-breeding bird species were in Central America, Southeast Asia and Southern Asia, followed by Eastern, Western and Middle Africa. This fact highlights the importance of considering separately the different stages of the annual cycle of migratory birds (Runge et al., 2015).
Regarding the spatial distribution of avian diet specialist species, our results agreed with previous studies despite a different strategy for classifying diet specialization (Belmaker et al., 2012). Diet specialist species richness is rather congruent with overall species richness, and this spatial correlation is higher for non-breeding species. However, exploring cases where the number of specialists is higher than expected from overall species richness, we were able to highlight some areas characterized by unusually high avian diet specialization (Figures 1b and 2b). Some areas (e.g. the Amazon, parts of Central Africa and Southeast Asia) are already recognised areas of conservation importance, but other regions identified in the analysis that harbour relatively high diet specialist species richness, especially those with relatively low overall species richness, suggest some important additional conservation priorities (e.g. some areas of the Northern Hemisphere in Europe, and Australasia for resident and breeding birds, and some areas in the tropics and south of the subtropical zone for non-breeding birds). Ecological specialization could be integrated directly into approaches for identifying additions to protected area systems, for example by modifying estimates of the benefit of protecting candidate planning units to include ecological specialization as well as more traditional biodiversity metrics based on the occurrence or abundance of species.
Areas with high diet specialization in avian communities were not necessarily the same as areas characterized by a high species richness (see Table 1 and Figures 1a-c  and 2a-c). This fact implies that the number of species in a community is not necessarily a good proxy for the level of specialization of a given community, and again highlights the need to explicitly incorporate specialization into conservation planning policy. For example the mismatch between overall species richness and DSI is relatively high in some north parts of the Northern Hemisphere, meridional areas of the Southern Cone, Madagascar and Australia for resident and breeding birds (Figure 1a-c), a result confirmed by the spatial correlations (Table S3). Similar discrepancies between overall species richness and DSI were found also for non-breeding species (Figure 2a-c  and Table S3). All these differences in the spatial distribution of avian specialization need to carefully be assessed when considering the establishment of a natural sanctuary or a protected area. In Figure S8, we provide two examples of areas characterized by a high mismatch between the overall species richness and the total number of diet specialist species in the bird communities, overlapped with the network of protected areas. Clearly, decisions around protected area designation could usefully incorporate information on hotspots of specialization.
A second novelty of this study is the assessment of whether diet specialist species are better covered by the protected area network than non-specialist birds. Specialist birds have fewer ecological options, so they arguably need greater protection. The fact that diet specialist richness was correlated with overall bird species richness across all continents is a first suggestion that the network of protected areas, which was built mainly based on the distribution of species richness or simply situated in remote areas, incidentally offers substantial protection to avian specialists. However, the most threatened species have relatively little protected area coverage, independent of whether they are specialist or non-specialist, and independently of whether they are resident, breeding or nonbreeding species. Only six Critically Endangered diet specialist species were adequately protected by the global network of protected areas, among the resident and breeding bird species. The Magenta petrel Pterodroma magentae is one of the world's rarest seabirds, with a diet mainly of squid, and threatened from predation by feral cats (Johnston et al., 2003). The Santa Marta wren Troglodytes monticola is a small passeriform inhabiting the upper elevations of the Santa Marta Massif in northern Colombia, with a diet characterized by invertebrates (Kroodsma et al., 2019). The last four diet specialist species are the Indigo-winged parrot Hapalopsittaca fuertesi, a frugivorous species from Colombia (Benavidez et al., 2018), the White-backed vulture Gyps africanus, which is an Old-World vulture with a diet of carrion and bone fragments of ungulate carcasses in African savannahs (Phipps et al., 2013), the Orangebellied parrot Neophema chrysogaster that is a small parrot endemic to southern Australia and the Sumatran ground cuckoo Carpococcyx viridis, a large and very elusive species of cuckoo endemic to Indonesia (Erritzøe et al., 2012). The remaining 85% of Critically Endangered diet specialist resident or breeding species were inadequately protected, highlighting the importance of integrating ecological specialization into conservation planning policy and practice. Interestingly, the White-backed vulture was the only Critically Endangered diet specialist to meet the conservation target among the non-breeding species, with other six species failing to meet such conservation target.
Inevitably, the quality of the diet specialization index is linked to the reliability and completeness of the underlying data, and because the data quality varies from species to species, the accuracy of the index is no doubt uneven across species. Furthermore, our classification of a species as specialist when associated with only a single major diet category could be rather too coarse to reflect variation in the types of food consumed within that broad category. Underlying diet data are limited mainly by missing information or by uncertainty in the expert assessments about some species' diet used in Wilman et al. (2014). However, missing information for a small number of can be estimated using interpolated values based on taxonomic relations among species (Molina-Venegas et al., 2018). In addition, some species show strong seasonal variation in diet, for example the blackcap Sylvia atricapilla is primarily insectivorous during the breeding season and frugivorous outside this period.
We have shown that high latitude and equatorial regions are characterized by avian communities with high diet specialization, perhaps representing hitherto unrecognised conservation priorities. We suggest paying more attention to the spatial distribution of the areas where the number of specialists is above that expected from overall species richness, and where overlap with the current network of protected areas is inadequate. Designation of new or expanded protected area might be necessary to cover some of these core areas of avian diet specialization. Conserving ecological specialization is potentially more efficient than purely species-based approaches because it captures the attributes of communities that could be more sensitive to environmental change. Additionally, we highlight the need to better understand the local and regional distribution of avian specialization, as well as the ecological mechanisms determining where specialists occur.

A C K N O W L E D G E M E N T
F.M. and Y.B. were financially supported by the Czech Science Foundation GAČR (project number 18-16738S).

A U T H O R C O N T R I B U T I O N S
F.M., Y.B. and R.A.F. conceived the idea and designed methodology. F.M., Y.B. and J.O.H. prepared the data and performed data analyses. All authors contributed critically to the drafts and gave final approval for publication.

D ATA A N D C O D E AVA I L A B I L I T Y
All data generated or analysed during this study are included in this published article (and its Supporting Information files). Code available at https://github. com/jeffreyhanson/global-protected-areas. Other codes are available from the corresponding author on reasonable request. https://github.com/jeffreyhanson/globalprotected-areas

C O N F L I C T O F I N T E R E S T
The authors declare no conflict of interest.