Pierre Gaüzère

Louise O'connor

Christophe Botella

Giovanni Poggiato

Tamara Münkemüller

Laura J Pollock

Ulrich Brose

Luigi Maiorano

Michael Harfoot

Wilfried Thuiller

The diversity of biotic interactions complements functional and phylogenetic facets of biodiversity

Keywords: biodiversity, vertebrates, macroecology, biogeography, conservation

Taxonomic, functional and phylogenetic diversities are important facets of biodiversity. Studying them together has improved our understanding of community dynamics, ecosystem functioning and conservation values [START_REF] Faith | Conservation evaluation and phylogenetic diversity[END_REF][START_REF] Tilman | Functional diversity[END_REF][START_REF] Pollock | Protecting Biodiversity (in All Its Complexity): New Models and Methods[END_REF] . In contrast to species, traits, and phylogenies, the diversity of biotic interactions has so far been largely ignored as a biodiversity facet in large-scale studies. This neglect represents a crucial shortfall because biotic interactions shape community dynamics, drive important aspects of ecosystem functioning 4- 7 , provide services to humans, and have intrinsic conservation value [START_REF] Thompson | The Ecological Basis of Conservation: Heterogeneity, Ecosystems, and Biodiversity[END_REF][START_REF] Tylianakis | Conservation of species interaction networks[END_REF] . Hence the diversity of interactions can provide crucial and unique information with respect to other diversity facets. Here, we leveraged large datasets of trophic interactions, functional traits, phylogenies and spatial distributions of >1000 terrestrial vertebrate species across Europe at a 10km resolution. We computed the diversity of interactions (Interaction Diversity, ID) in addition to functional (FD) and phylogenetic diversities (PD). After controlling for species richness, surplus and deficits of ID were neither correlated with FD nor with PD, thus representing unique and complementary information to the commonly studied facets of diversity. A threedimensional mapping allowed for simultaneously visualizing different combinations of ID-FD-PD. Interestingly, the spatial distribution of these diversity combinations closely matched the boundaries between ten European biogeographic regions, and revealed new, interaction-rich areas in the European Boreal region and interaction-poor areas in central Europe. Our study demonstrates that the diversity of interactions adds new and ecologically relevant information to multi-facetted, large-scale diversity studies with implications for understanding eco-evolutionary processes and informing conservation planning.

Results and discussion

Figure 1. Conceptual workflow for a joint analysis of phylogenetic, functional, and interaction diversity. (a) Occurrences and probability of presence for 1149 terrestrial vertebrate species on 117,000 10×10km cells across Europe are combined with (b) the phylogenetic tree, a set of functional traits, and the trophic interactions of species. (c) We combined species distribution with phylogenetic, functional and trophic species attributes to compute local terrestrial vertebrate diversities using Hills numbers (q = 0, i.e "richness") and statistically corrected the diversity values by the local species richness. Note that the expected relationships (gray ellipses) are not necessarily linear. (d) We projected the diversity values in a 3dimensional space with each axis representing a diversity facet and a color in the Red-Blue-Green space (x = PD / blue, y = FD / green, z = ID / red), and discretized particular types of combinations based on surplus and deficits of each diversity. Red identifies surpluses of ID and FD associated with deficits in FD and PD (ID>0, FD<0, PD<0); Yellow identifies surpluses of ID and FD associated with deficits in PD (ID>0, FD>0, PD<0); Green identifies deficits in ID and PD associated with FD surpluses (ID<0, FD>0, PD<0); Black identifies deficits in ID, PD and FD, Pink identifies surpluses of ID and PD associated with FD deficits; Dark blue identifies surpluses of PD associated with ID and FD deficits; Light blue identifies surpluses of PD and FD associated with ID deficits; white identifies surpluses in ID, PD and FD.

Biodiversity -the diversity of life on Earth-was originally used to refer to species diversity, but it is now used to reflect a multi-faceted concept [START_REF] Pollock | Protecting Biodiversity (in All Its Complexity): New Models and Methods[END_REF] . Given the evidence that species diversity alone cannot appropriately describe community assembly, ecosystem functioning and variation in community composition [START_REF] Hooper | Species diversity, functional diversity and ecosystem functioning[END_REF] , several complementary measures of biodiversity have emerged in the last three decades [START_REF] Faith | Conservation evaluation and phylogenetic diversity[END_REF][START_REF] Tilman | Functional diversity[END_REF] . The most important are the diversity of species' evolutionary histories (i.e phylogenetic diversity, PD) and their ecological functions (i.e functional diversity, FD), but while PD and FD are becoming central to many studies [START_REF] De Palma | Dimensions of biodiversity loss: Spatial mismatch in land-use impacts on species, functional and phylogenetic diversity of European bees[END_REF][START_REF] Schweiger | Plant spectral diversity integrates functional and phylogenetic components of biodiversity and predicts ecosystem function[END_REF] , the diversity of biotic interactions (i.e interaction diversity, ID) has been poorly considered as a biodiversity facet in large-scale studies (but see refs [START_REF] Thompson | The Ecological Basis of Conservation: Heterogeneity, Ecosystems, and Biodiversity[END_REF][START_REF] Dyer | Diversity of Interactions: A Metric for Studies of Biodiversity: Interaction Diversity[END_REF] ). This is a major gap since biological interactions are tightly linked to species coexistence [START_REF] Jacquet | No complexity-stability relationship in empirical ecosystems[END_REF] , ecosystem productivity and functioning [START_REF] Grime | Biodiversity and Ecosystem Function: The Debate Deepens[END_REF][START_REF] Tilman | Biodiversity and Ecosystem Functioning[END_REF][START_REF] Brose | Predicting the consequences of species loss using sizestructured biodiversity approaches[END_REF][START_REF] Schneider | Animal diversity and ecosystem functioning in dynamic food webs[END_REF] . In its simplest form, ID is the total number of interactions shared by all species of a given assemblage [START_REF] Tylianakis | Conservation of species interaction networks[END_REF] . Interactions considered can be of different types and nature, e.g antagonistic (competition for resources), mutualistic (pollination [START_REF] Burkle | The future of plant-pollinator diversity: understanding interaction networks across time, space, and global change[END_REF][START_REF] Burkle | The beta-diversity of species interactions: Untangling the drivers of geographic variation in plant-pollinator diversity and function across scales[END_REF] ), or trophic (predation [START_REF] Sandom | Mammal predator and prey species richness are strongly linked at macroscales[END_REF] ). Although the concept of interaction diversity is not novel [START_REF] Thompson | The Ecological Basis of Conservation: Heterogeneity, Ecosystems, and Biodiversity[END_REF][START_REF] Dyer | Diversity of Interactions: A Metric for Studies of Biodiversity: Interaction Diversity[END_REF] and has its own methodological tools [START_REF] Ohlmann | Diversity indices for ecological networks: a unifying framework using Hill numbers[END_REF][START_REF] Poisot | The dissimilarity of species interaction networks[END_REF] , the lack of information available on biotic interactions [START_REF] Hortal | Seven Shortfalls that Beset Large-Scale Knowledge of Biodiversity[END_REF][START_REF] Poisot | Global knowledge gaps in species interaction networks data[END_REF] has limited its study across large taxonomical and spatial scales [START_REF] O'connor | Unveiling the food webs of tetrapods across Europe through the prism of the Eltonian niche[END_REF][START_REF] Braga | Spatial analyses of multi-trophic terrestrial vertebrate assemblages in Europe[END_REF][START_REF] Cumming | Network analysis in conservation biogeography: Challenges and opportunities[END_REF][START_REF] Gravel | Bringing Elton and Grinnell together: a quantitative framework to represent the biogeography of ecological interaction networks[END_REF][START_REF] Baiser | Ecogeographical rules and the macroecology of food webs[END_REF] . Here, we leveraged unique and valuable data combining spatial distributions [START_REF] Maiorano | Threats from climate change to terrestrial vertebrate hotspots in Europe[END_REF] (Figure 1.a), trophic interactions (Figure 1.b), functional traits [START_REF] Thuiller | Conserving the functional and phylogenetic trees of life of European tetrapods[END_REF] , and phylogenies [START_REF] Roquet | One tree to link them all: a phylogenetic dataset for the European tetrapoda[END_REF] of most terrestrial vertebrate species in Europe 30 at a 10 km resolution. Within each 10km cell, we computed interaction diversity (ID, as the number of trophic interactions), functional diversity (FD, as the sum of functional pairwise Gower distances between species in the cell), and phylogenetic diversity (PD, as the sum of the branch lengths of the phylogenetic tree containing all species present in the cell 1 ) using Hill numbers 31 18 . We statistically corrected each diversity by the local species richness [START_REF] Zupan | Spatial mismatch of phylogenetic diversity across three vertebrate groups and protected areas in Europe[END_REF][START_REF] Safi | Understanding global patterns of mammalian functional and phylogenetic diversity[END_REF] in order to measure and map deficits and surpluses of ID, FD, and PD (Figure 1.c). We also investigated the correlation and complementarity between the three facets, and created a 3-dimensional diversity space that reveals different local combinations of ID-FD-PD (Figure 1.d) and their distribution across biogeographical regions in Europe.

Surpluses and deficits of diversities

Trophic networks of terrestrial vertebrates found within 10km cells in Europe contained up to 4834 trophic interactions with an average of 1958 interactions across cells (Figure 2a). Once corrected for species richness, ID ranged from a deficit of -942 interactions (1667 observed interactions with 202 species involved) to a surplus of +968 interactions (3730 interactions with 210 species involved, see Supplemental information -Trophic network examples). Because highly connected assemblages are often considered as the signature of functional and resilient ecosystems [START_REF] Gilbert | Connectance indicates the robustness of food webs when subjected to species loss[END_REF][START_REF] Thébault | Stability of ecological communities and the architecture of mutualistic and trophic networks[END_REF] , areas with high ID are important from a conservation point of view [START_REF] Gilbert | Connectance indicates the robustness of food webs when subjected to species loss[END_REF][START_REF] Heleno | Connectance of species interaction networks and conservation value: Is it any good to be well connected?[END_REF] . Further, comparing spatial distributions of surplus and deficit IDs with those of FD or PD can complement our understanding of community dynamics and underlying processes. Because phylogenetic and trait data contain information about evolutionary history and species niches, the spatial distribution of their diversity (Figure 2b-c) is thought to hold the signature of the eco-evolutionary drivers that shape biodiversity patterns [START_REF] Mcgill | Rebuilding community ecology from functional traits[END_REF][START_REF] Webb | Phylogenies and Community Ecology[END_REF][START_REF] Münkemüller | Dos and don'ts when inferring assembly rules from diversity patterns[END_REF] . For example, for a given species richness, observed surplus of FD (Figure 1d, and Figure 2b, green color) could result from competitive exclusion between species with similar traits, while a deficit of FD might result from environmental filtering constraining the range of locally viable traits or hierarchical competition where a given set of traits is the best adapted locally [START_REF] Mayfield | Opposing effects of competitive exclusion on the phylogenetic structure of communities[END_REF] . PD surplus (Figure 1d and Figure 2c, dark blue color) could result from slow extinction rates of old and distant lineages (i.e., museums of biodiversity [START_REF] Stebbins | Flowering Plants: Evolution Above the Species Level[END_REF] ), and PD deficit from rapid recent speciation (i.e., cradles of biodiversity). ID surplus and deficit brings additional information, as observed ID surplus (Figure 1d and Figure 2a, red color) indicates particularly dense or long trophic networks, such as those emerging from high levels of omnivory and intraguild predation [START_REF] Duffy | The functional role of biodiversity in ecosystems: incorporating trophic complexity[END_REF] , or from bottom-up control when large amounts of basal resources sustain longer trophic chains and the presence of top predators. ID deficits can result from weakened topdown control when top predators are absent from local assemblages, for example following human-induced removal [START_REF] Estes | Trophic downgrading of planet Earth[END_REF][START_REF] De Visser | The Serengeti food web: empirical quantification and analysis of topological changes under increasing human impact[END_REF][START_REF] Pauly | Fishing down marine food webs[END_REF] . Overall, the different facets of diversity are shaped by eco-evolutionary drivers which are not mutually exclusive [START_REF] Tylianakis | Ecological Networks Across Environmental Gradients[END_REF] . Any combination of ID-FD-PD could potentially exist locally and bring complementary information to the others, although one can expect the facets of diversity to be (partly) correlated when similar drivers influence multiple diversity facets. We showed that FD and PD were clearly and positively correlated (Figure 3b). This correlation is due to the fact that species tend to retain their ancestral traits through evolution [START_REF] Webb | Phylogenies and Community Ecology[END_REF][START_REF] Wiens | Niche Conservatism: Integrating Evolution, Ecology, and Conservation Biology[END_REF][START_REF] Peterson | Conservatism of ecological niches in evolutionary time[END_REF][START_REF] Lavergne | Biodiversity and Climate Change: Integrating Evolutionary and Ecological Responses of Species and Communities[END_REF] , and suggests an important effect of evolution and phylogenetic niche conservatism on biodiversity patterns [START_REF] Wiens | Niche Conservatism: Integrating Evolution, Ecology, and Conservation Biology[END_REF][START_REF] Wiens | Niche conservatism as an emerging principle in ecology and conservation biology[END_REF][START_REF] Münkemüller | Common conceptual and methodological pitfalls in the analysis of phylogenetic niche conservatism[END_REF] . While one could expect ID to be related to FD (because of the link between trait similarity and competition for resources) or PD (because biotic interactions can drive the (co)evolutionary history of the species [START_REF] Tylianakis | Ecological Networks Across Environmental Gradients[END_REF][START_REF] Williams | Simple rules yield complex food webs[END_REF] ), this was not what we observed (Figure 3c-d). Instead, ID represented unique and complementary information to the commonly studied facets of diversity.

Bottom left: Distribution of deficits and surpluses of diversities, where model residuals correspond to "corrected diversity" values with deficits (dark shades) and surpluses (red for ID, green for FD, blue for PD). Right: spatial distribution of corrected values for each biodiversity facet, color corresponds to distributions on the left.

Distribution of diversity combinations

To investigate the congruence between the interaction, functional and phylogenetic facets of biodiversity, we created a 3-dimensional space where each dimension represents one diversity facet. We further attributed a color channel for each diversity facet (red = ID, green = FD, blue = PD) to visualize all possible combinations of biodiversity facets (Figure 1d). Each combination of three color channels (Red, Blue, Green) resulted in a particular color in the RGB color space that corresponds to a given combination of three diversity facets, and allowed us to identify a continuum of ID-FD-PD combinations (Figure 1d). We also interpreted particular types of combinations by discretizing colors based on the combinations of surplus and deficits of each diversity facet (Figure 1d). This joint analysis of diversity facets highlighted various local combinations of ID-FD-PD, with all kinds of combinations being observed in different proportions (Figure 3a). The most commonly observed combinations were ID surpluses with FD and PD deficits (covering 21.8% of the total study area); surpluses in ID, FD, and PD (white, 21.6%); surpluses of FD and PD with deficits in ID (light blue 19.6%); and deficits in ID, FD, and PD (black, 17.3%), which is consistent with the positive correlation observed between FD and PD (Figure 3b). The spatial structure of diversity combinations aligned well with many boundaries of European biogeographical regions (Figure 3a), a striking spatial congruency considering that the identification and delimitation of bioregions are based on geographic distribution of vegetation types [START_REF] Noirfalise | Map of the Natural Vegetation: Of the Member Countries of the European Economic Community and the Council of Europe[END_REF] . Beyond species distribution, biodiversity facets such as phylogenetic diversity already have been shown to match some ecological regions across the globe [START_REF] Holt | An update of Wallace's zoogeographic regions of the world[END_REF] . ID strongly varies between different regions (e.g between the Mediterranean region and the Alps, or between the Continental region and the Carpathian mountains) and thus further refines boundaries between them. These results suggest that species interactions (along with species co-occurrences and phylogeny) could have a strong structuring effect on (bio)regional species pools. Such a question, however, would require a deeper analysis based on the turnover of interactions within and between regions as regional diversity is connected to local diversity by the turnover in composition between locations. Interestingly, the mapping of diversity combinations also revealed the specificity of several sub-regions within their biogeographical region, e.g the Balkan peninsula sub-region in the Mediterranean region, or the Carpathian mountains in the Alpine region. These results further highlight that biotic interaction diversity adds new and independent information and that a dense network of trophic interactions can occur in areas of poor functional and phylogenetic diversity.

black shows lowest ID-FD-PD values, white shows highest ID-FD-PD, and so on for each combination. Black lines show the boundaries of the European biogeographical regions. (b) Pair plot of corrected FD (y-axis) VS corrected PD (x-axis), (c) corrected FD (y-axis) VS corrected ID (x-axis), (d) corrected ID (y-axis) VS corrected PD (x-axis). In top right, r is the value of Pearson's product-moment correlation between y and x. Points colors correspond to colors in the map.

Southern Europe showed strong diversity surpluses in all diversity facets (white / light color shades in Figure 3a), which confirms the Mediterranean bioregion as a multifaceted biodiversity hotspot [START_REF] Myers | Biodiversity hotspots for conservation priorities[END_REF][START_REF] Mittermeier | Global Biodiversity Conservation: The Critical Role of Hotspots[END_REF] . This result shows that, for a given number of species, local assemblages of Mediterranean terrestrial vertebrate species were particularly rich in terms of ecological strategies, contained long evolutionary history, and had particularly dense trophic networks. In the Mediterranean basin, the warm climate and the geographical proximity with Africa and Asia explains the high diversity of amphibians and reptiles, as well as the presence of unique evolutionary lineages, leading to high functional and phylogenetic diversities compared to the rest of Europe (Supplemental information -Groups and species distributions across Europe). In addition to these high levels of functional and phylogenetic diversities, the Mediterranean region showed surpluses in interaction diversity, in particular in the subregion of the Balkan peninsula. The densely connected trophic networks observed in the Mediterranean region resulted from (i) numerous top predators in this region (Supplemental information -Maps of relevant network properties) previously identified as birds, felids, and snakes predating upon small reptiles and rodents [START_REF] O'connor | Unveiling the food webs of tetrapods across Europe through the prism of the Eltonian niche[END_REF] ; and (ii) to a lesser extent from a high degree of omnivory in the Iberian peninsula [START_REF] Braga | Spatial analyses of multi-trophic terrestrial vertebrate assemblages in Europe[END_REF] .

Conversely, the northernmost areas tended to show low levels of diversities (black areas in North of Scandinavia and Iceland, Figure 3a). The Boreal and Arctic bioregions showed deficits in functional and phylogenetic diversities, but tended to sustain surpluses in interaction diversity (red areas, Figure 3a). In these regions, FD deficits were likely to be driven by the cold climate constraining the range of functional traits that can be found in these regions, and similarly for PD via trait conservatism. In particular, the fact that cold temperature limits the presence of ectotherms (amphibians and reptiles) in high latitudes reduces functional and phylogenetic diversities, in line with the expected effect of environmental filtering on these diversity facets [START_REF] Mcgill | Rebuilding community ecology from functional traits[END_REF][START_REF] Webb | Phylogenies and Community Ecology[END_REF] . The consideration of ID brings additional and complementary information since FD-PD deficits are associated with ID surpluses in Northern Europe. The presence of ID surpluses in the Boreal and Arctic bioregions likely resulted from a high degree of omnivory (Supplemental information -Maps of relevant network properties), which is known to increase trophic network connectance [START_REF] Baiser | Ecogeographical rules and the macroecology of food webs[END_REF][START_REF] Cirtwill | Latitudinal gradients in biotic niche breadth vary across ecosystem types[END_REF] . Species that live under high latitudes tend to be trophic generalists [START_REF] Braga | Spatial analyses of multi-trophic terrestrial vertebrate assemblages in Europe[END_REF] because the higher seasonality in high latitudes promotes the evolution of larger niche breadth, in accordance with the latitude-niche breadth hypothesis 57,58 .

Within the Alpine bioregion, different mountain ranges displayed contrasting diversity combinations. The marked differentiation between the Alps and the Carpathian mountains subregions is a striking example supporting the consideration of interaction diversity in biodiversity studies and conservation biogeography. These two mountain ranges located in Central Europe are part of the same Alpine bioregion, which partly explains their similarity in terms of functional surpluses and phylogenetic deficits (Figure 2b-c). Based on functional and phylogenetic diversities alone, these two mountain ranges would be considered as similarly diverse -but they are markedly different in terms of interaction diversity. The Carpathians displayed a clear ID surplus in (Figure 3, yellow), while the Alps were deficitary (Figure 3, green). The proximate cause of such difference was the rarity of top predators in the Alps compared to the Carpathians (see Supplemental information -Maps of relevant network properties). Human influence likely explains this discrepancy because many apex-predators (bears, wolves, lynx) that are often trophic generalists are still present in the Carpathians, while they were exterminated in the Alps [START_REF] Kaczensky | Status, Management and Distribution of Large Carnivores -Bear, Lynx, Wolf & Wolverine -in Europe[END_REF] .

Potential drivers of diversity facets

While environmental filtering is likely to drive the decrease of FD and PD observed in high latitudes, ID might be more influenced by human activities than climate. As such, local deficits of trophic interactions appeared as a marker of high human impact across Europe. This is in line with the negative correlation between connectance and human influence previously reported for the same study system [START_REF] Braga | Spatial analyses of multi-trophic terrestrial vertebrate assemblages in Europe[END_REF] , and suggests that the diversity of interactions is influenced by different drivers than functional and phylogenetic diversity. It is, however, noteworthy that other studies reported higher connectance in more human impacted systems [START_REF] Tylianakis | Conservation of species interaction networks[END_REF][START_REF] Tylianakis | Ecological Networks Across Environmental Gradients[END_REF] . Indeed, the human-induced relative increase of generalist intermediate predators could counterbalance the decrease in ID due to the loss of a few top predators. The human influence on large-scale diversity has been considered and studied in terms of phylogeny and traits [START_REF] Fritz | Geographical variation in predictors of mammalian extinction risk: big is bad, but only in the tropics[END_REF] . However, its consequences on large-scale patterns of interaction diversity have been largely overlooked, although they are probably stronger. Indeed, human activities have been (and still are) particularly detrimental to large-bodied species [START_REF] Fritz | Geographical variation in predictors of mammalian extinction risk: big is bad, but only in the tropics[END_REF][START_REF] Gaston | Birds, body size and the threat of extinction[END_REF][START_REF] Cardillo | Multiple causes of high extinction risk in large mammal species[END_REF][START_REF] Böhm | Correlates of extinction risk in squamate reptiles: the relative importance of biology, geography, threat and range size[END_REF] . While this observation is generally viewed as a trait-induced consequence (humans are more detrimental to larger animals), it might also be a trophic-induced consequence (humans are more detrimental to apex and generalist predators) [START_REF] Pauly | Fishing down marine food webs[END_REF][START_REF] Brose | Improving nature conservancy strategies by ecological network theory[END_REF] .

The importance of interaction diversity

A clear understanding of the impact of human activities on ID has yet to emerge. More generally, ID is likely to be highly context and taxa dependent, and the understanding of its multi-scale drivers represents a research agenda for the years to come. Among others, the Eltonian shortfall is one big challenge that currently limits the description of ID in many parts of the world where information on biotic interactions is lacking [START_REF] Poisot | Global knowledge gaps in species interaction networks data[END_REF] . Here, we overcame this challenge for trophic interactions by inferring local interactions from species distributions and their known potential trophic interactions from the literature and expert knowledge (as commonly done, see for example refs [START_REF] Poisot | The dissimilarity of species interaction networks[END_REF][START_REF] Gravel | Bringing Elton and Grinnell together: a quantitative framework to represent the biogeography of ecological interaction networks[END_REF] ). While this approach overestimates interactions at a given time, "realized" and "potential" number of interactions are very likely to converge in the long term. On the contrary, a field sampling approach would underestimate the realized ID. This underestimation can be quite severe and a massive sampling effort is required to detect most interactions [START_REF] Chacoff | Evaluating sampling completeness in a desert plant-pollinator network[END_REF] . Combining both approaches (inferring interactions from a metanetwork and species distribution, vs. observing interactions), and comparing their accuracy across a range of temporal and spatial scales will provide valuable insights in community ecology and biogeography [START_REF] Fortuna | Partner Fidelity and Asymmetric Specialization in Ecological Networks[END_REF] .

Although ID patterns appear robust to data depletion and spatial contexts (see Supplemental information -Robustness of diversity patterns), whether the patterns described in this study can be extrapolated to other biomes remains an open question. For example, our conclusions from European terrestrial vertebrates might not hold true for tropical rainforests which shelter many trophic specialist species with narrow ecological niches (but comprehensive data on traits and interactions are lacking). Nonetheless, we argue that interaction diversity is a particularly valuable facet for biogeography and conservation planning. Although this view has been empirically challenged [START_REF] Heleno | Connectance of species interaction networks and conservation value: Is it any good to be well connected?[END_REF] , more densely connected trophic networks are generally considered as desirable from a conservation point of view [START_REF] Tylianakis | Conservation of species interaction networks[END_REF] . Areas with surpluses of interactions represent interaction networks that are expected to be more robust to cascading species extinctions [START_REF] Gilbert | Connectance indicates the robustness of food webs when subjected to species loss[END_REF] , and consequently more resilient to perturbations. Coupled with its apparent sensitivity to human activities [START_REF] Braga | Spatial analyses of multi-trophic terrestrial vertebrate assemblages in Europe[END_REF] , interaction diversity might be viewed as a marker of both ecosystem degradation and resistance to future degradation. We argue that a general consideration of interaction diversity as an important and meaningful diversity facet alongside the functional and phylogenetic diversities should be a priority for macroecology and conservation biogeography.

STAR Methods

Study area and data

Study area. The study area, hereafter referred to as "Europe", included the entire European subcontinent (with Macaronesia and Iceland) plus Anatolia to include a complete picture of the North Mediterranean coast (Figure 1a). The study area was divided into 117,000 cells on a 10×10 km equal-size area grid (ETRS89). Within the study area, we considered ten biogeographical regions defined by the European Environment Agency 67 : Alpine, Anatolian, Arctic, Atlantic, Boreal, Black Sea, Continental, Macaronesia, Mediterranean, and Steppic. These bioregions are large scale ecological units based on an interpretation of geobotanical data [START_REF] Roekaerts | The Biogeographical Regions Map of Europe. Basic principles of its creation and overview of its development[END_REF] , and represent areas with homogeneous ecological context. Species distributions. We extracted the distributions for all terrestrial vertebrates naturally occurring within the study area from Maiorano et al. (2013). Species distributions for 509 bird, 288 mammal, 250 reptile and 104 amphibian species were mapped by combining the IUCN extent of occurrence for each species with their habitat requirements. A species was considered potentially present in a 10×10 km cell if the grid cell met the three following criteria: i) is within the species extent of occurrence, ii) contains at least one 300x300m area of primary habitat for the species, i.e. habitat where the species can persist (defined by experts and published literature) and iii) meets species requirements in terms of elevation and distance from water. A full description of species distribution data and definition of primary habitat can be found in Maiorano et al. 2013 [START_REF] Maiorano | Threats from climate change to terrestrial vertebrate hotspots in Europe[END_REF] . In addition, we used the percentage of primary habitat of the species in each cell as a proxy for the probability to find the species in a random locality within this cell. For example, we considered that if the primary habitat of a species covered 80% of the cell, the probability to find the species in a random locality of the cell was 0.8. As such, it represents a proxy for the probability of presence of the species within the cell and was used as a weight in the entropy-based diversity measures (i.e when q=1) provided as supplementary analyses.

Functional traits. We gathered biological trait data from Thuiller et al. 2015 [START_REF] Thuiller | Conserving the functional and phylogenetic trees of life of European tetrapods[END_REF] , excluding traits describing diet (and thus trophic interactions) and traits for larvae and juveniles. Our analysis was based on four lifehistory and ecological traits common to mammals, amphibians, birds and reptiles. The only quantitative trait was body mass [grams, log-transformed]. The three other traits were multichoice nominal variables coded by binary values. Feeding behavior was coded by four binary columns: opportunistic feeder, active hunter, browser, grazer. Nesting location was coded by eleven binary columns: tree/hole/fissure in bark, ground, rocks, building/artificial, underground water, cave/fissures/burrows, lodge, temporary water, brooks/springs/small rivers, puddles/ponds/pools/small lakes, brackish waters. Activity time was coded by four binary columns: nocturnal, crepuscular, diurnal, arrhythmic. These traits were selected because they represent informative niche dimensions linked to the use and acquisition of resources in space and time, and are related to ecosystem functioning [START_REF] Luck | Improving the application of vertebrate trait-based frameworks to the study of ecosystem services[END_REF][START_REF] Sekercioglu | Increasing awareness of avian ecological function[END_REF] . A thorough description of traits and the list of publications where the data were gathered is available in supplementary material from [START_REF] Thuiller | Conserving the functional and phylogenetic trees of life of European tetrapods[END_REF] available at https://royalsocietypublishing.org/doi/suppl/10.1098/rstb.2014.0005. We computed the pairwise dissimilarities (distances) of this trait matrix using a mixed variable coefficient of distance (using function dist.ktab in ade4) that generalizes Gower's general coefficient of distance to allow the treatment of various statistical types of variables when calculating distances [START_REF] Pavoine | On the challenge of treating various types of variables: application for improving the measurement of functional diversity[END_REF] . Euclidean distance was used for body mass, and Jaccard index was used for the four other multichoice nominal variables (S3 coefficient of in Gower and Legendre 1986 [START_REF] Gower | Metric and Euclidean properties of dissimilarity coefficients[END_REF] ).

Phylogenetic tree. We used the 100 phylogenetic trees for European terrestrial vertebrates assembled and published by Roquet at al. 2014 [START_REF] Roquet | One tree to link them all: a phylogenetic dataset for the European tetrapoda[END_REF] . We chose these phylogenetic trees as they are the only species-level phylogenies encompassing all european vertebrates, and have already been valuably used to depict phylogenetic diversity in this context vertebrates in the past [START_REF] Thuiller | Conserving the functional and phylogenetic trees of life of European tetrapods[END_REF][START_REF] Zupan | Spatial mismatch of phylogenetic diversity across three vertebrate groups and protected areas in Europe[END_REF][START_REF] Pearman | Phylogenetic patterns of climatic, habitat and trophic niches in a European avian assemblage[END_REF][START_REF] Saladin | Environment and evolutionary history shape phylogenetic turnover in European tetrapods[END_REF] .

Trophic networks. We used data on species trophic interactions from the metaweb of European terrestrial vertebrates, (Tetra-EU 1.0, [START_REF] Maiorano | TETRA-EU 1.0: A species-level trophic metaweb of European tetrapods[END_REF] ). This metaweb is based on expert knowledge, published information and field guides. Potential trophic links between a predator and a prey were identified from published accounts of their observation, morphological similarities between potential prey and literature-referenced prey or -in the absence of this information-the diet of the predator's sister species. The metaweb of European terrestrial vertebrates contained 1,164 species and a total of 50,408 potential trophic interactions. The full dataset and methods description can be found in ref [START_REF] Maiorano | TETRA-EU 1.0: A species-level trophic metaweb of European tetrapods[END_REF] .

In order to maximize the species coverage for each diversity, we allowed for different sets of species to be used to compute ID, FD and PD. For ID we retained 1149 species for which we had information on their European distribution range and trophic interactions; for FD we retained 1009 species for which we had information on their European distribution range and functional traits; for PD we retained 993 species for which we had information on their European distribution range and phylogeny. This varying set of species should have low impact on the assessment of diversities as ID, FD and PD were corrected by their corresponding taxonomic richness to compute surpluses and deficits. In order to investigate the potential bias resulting from the variation of species coverage across space and diversities, we computed diversities based on the same set of 884 species for which we had all shared information. The resulting diversity patterns were similar when considering the 884 species or varying set of species (see Supplemental information -Diversities based on the same set of 884 species ).

Diversity measures and the 3-dimensional diversity space

Within each 10x10km cell, we used Hill numbers [START_REF] Chao | Unifying Species Diversity, Phylogenetic Diversity, Functional Diversity, and Related Similarity and Differentiation Measures Through Hill Numbers[END_REF] to compute FD, PD, and ID. In this framework, diversity values are converted into effective numbers of species, the Hill numbers. When considering taxonomic diversity, the effective number of species is the number of equally abundant species necessary to produce the observed value of diversity (an analogue to the concept of effective population size in genetics). This approach has then been generalized to incorporate species phylogenetic relatedness and species functional distances. We used the framework from Chao et al. 2014 [START_REF] Chao | Unifying Species Diversity, Phylogenetic Diversity, Functional Diversity, and Related Similarity and Differentiation Measures Through Hill Numbers[END_REF] implemented in the R package hillR for phylogenetic and functional diversity and in the package econetwork [START_REF] Ohlmann | Diversity indices for ecological networks: a unifying framework using Hill numbers[END_REF] for interaction diversity. We computed each diversity as a Hill number analogous to a measure of richness by setting q=0 (ignoring abundance). The ID richness was the sum of trophic links formed by the species present in the cell, the FD richness was the sum of functional pairwise gower distances between species in the cell, and the PD richness was the mean sum of the branch lengths of the phylogenetic tree connecting all species present in the cell 1 across the 100 trees. We focused our study on richness-based results (q=0) as they are the easiest to interpret, but we also analyzed and showed results based on Shannon entropy in the Supplemental information -Results based on Shannon entropy. To compute the results as a Shannon entropy, we set q=1 and used the % of species' primary habitat within the cell as the probability to find the species in the cell. More precisely, when q=1 the ID entropy is the Shannon entropy over the interaction weights (product of the two species abundances), the FD is the Shannon entropy of effective number of species-pairs with unit-distance between species, and the PD is the mean Shannon entropy of the effective total branch length across the 100 trees. More details on the calculations of FD and PD can be found in Chao et al. 2014 [START_REF] Chao | Unifying Species Diversity, Phylogenetic Diversity, Functional Diversity, and Related Similarity and Differentiation Measures Through Hill Numbers[END_REF] , and in Ohlmann et al. 2019 [START_REF] Ohlmann | Diversity indices for ecological networks: a unifying framework using Hill numbers[END_REF] for ID.

We corrected FD, PD, ID richness and Shannon entropy for the number of species in the cell (i.e taxonomic richness) based on the set of species used to compute each diversity. We fitted a thin plate spline regression, a particular Generalized Additive Model (GAM), to predict each diversity measure from species richness. The residuals of each model (one for each diversity facet and order q) were retained as the species richness corrected value of the diversity, with positive residuals considered as surplus and negative residuals considered as deficits given the species richness [START_REF] Zupan | Spatial mismatch of phylogenetic diversity across three vertebrate groups and protected areas in Europe[END_REF][START_REF] Safi | Understanding global patterns of mammalian functional and phylogenetic diversity[END_REF] . In other words, a deficit (or surplus) indicates a lower (or higher, respectively) diversity value than expected given the local species richness (Figure 1) .

To investigate the congruence between the interaction, functional and phylogenetic facets of biodiversity, we created a 3-dimensional space where each dimension represents one diversity facet. In order to visualize all possible combinations of biodiversity facets, we attributed a color channel for each diversity facet (red = ID, green = FD, blue = PD) where the residual values for each diversity were rescaled to 0-255 value in the corresponding color channel (Figure 1). Hence, each combination of three color channels (Red, Blue, Green) results in a particular color in the RGB color space that corresponds to a given combination of three diversity facets, and allows us to identify a continuum of ID-FD-PD combinations depicted in figure 1.d. We can also interpret particular types of combinations by discretizing colors based on the combinations of surplus and deficits of each diversity. As shown in Figure 1d and Figure 3a, Red identifies surpluses of ID and FD associated with deficits in FD and PD (ID>0, FD<0, PD<0); Yellow identifies surpluses of ID and FD associated with deficits in PD (ID>0, FD>0, PD<0); Green identifies deficits in ID and PD associated with FD surpluses (ID<0, FD>0, PD<0); Black identifies deficits in ID, PD and FD, Pink identifies surpluses of ID and PD associated with FD deficits; Dark blue identifies surpluses of PD associated with ID and FD deficits; Light blue identifies surpluses of PD and FD associated with ID deficits; white identifies surpluses in ID, PD and FD.

Figure 2 .

 2 Figure 2. Patterns of Interaction diversity ID (a, in red), Functional diversity FD (b, in green), and Phylogenetic diversity PD (c, in blue). Top left: Relationship between each diversity facet and the species richness. Dotted lines show relationships as fitted by Generalized Additive Models. Bottom left: Distribution of deficits and surpluses of diversities, where model residuals correspond to "corrected diversity" values with deficits (dark shades) and surpluses (red for ID, green for FD, blue for PD). Right: spatial distribution of corrected values for each biodiversity facet, color corresponds to distributions on the left.

Figure 3 .

 3 Figure 3. (a) Spatial projection of the 3-dimensional diversity space. In the top left barplot, we created 8 discrete categories based on the combinations of deficits (-) and surpluses (+) of each diversity and reported the number of cells falling in each category. In the map, points are colored by their location in the Red-Green-Blue 3-dimensional color space, with each diversity facet corresponding to a distinct channel: Red channel = Interaction Diversity, Green Channel = Functional Diversity, Blue channel = Phylogenetic Diversity. black shows lowest ID-FD-PD values, white shows highest ID-FD-PD, and so on for each combination. Black lines show the boundaries of the European biogeographical regions. (b) Pair plot of corrected FD (y-axis) VS corrected PD (x-axis), (c) corrected FD (y-axis) VS corrected ID (x-axis), (d) corrected ID (y-axis) VS corrected PD (x-axis). In top right, r is the value of Pearson's product-moment correlation between y and x. Points colors correspond to colors in the map.

ACKNOWLEDGEMENTS

This study has received funding from the ERA-Net BiodivERsA -Belmont Forum, with the national funder Agence Nationale de la Recherche (ANR-18-EBI4-0009), part of the 2018 Joint call BiodivERsA-Belmont Forum call (project 'FutureWeb'). WT was also supported by the Agence Nationale pour la Recherche through the FORBIC (ANR-18-MPGA-0004) and EcoNet (ANR-18-CE02-0010) projects. We warmly thank Florian Maderspacher for handling the paper, Brittany Pugh, Richard Field, Kamran Safi and one anonymous referee for their valuable reviews and comments. Interaction data depletions

SUPPLEMENTAL INFORMATION

Supplemental Information includes six appendix and can be found at:

AUTHOR CONTRIBUTIONS

P.G and W.T. conceived the study, with early advice from L.O.C, C.B, G.P. and T.M. P.G. performed all analyses and wrote the first version of the manuscript with inputs from W.T, L.O.C, C.B, G.P. and T.M. All authors contributed substantially to the interpretation of the results and to the writing of the manuscripts and its revisions.