Assessing the equilibrium between assemblage composition and climate: A directional distance-decay approach

J Anim Ecol. 2021;00:1–13. | 1wileyonlinelibrary.com/journal/jane Received: 3 February 2021 | Accepted: 23 April 2021 DOI: 10.1111/1365-2656.13509 RESEARCH ARTICLE Assessing the equilibrium between assemblage composition and climate: A directional distancedecay approach Andrés Baselga1,2 | Carola Gómez- Rodríguez3,2 This is an open access article under the terms of the Creat ive Commo ns Attri butio n- NonCo mmercial License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited and is not used for commercial purposes. © 2021 The Authors. Journal of Animal Ecology published by John Wiley & Sons Ltd on behalf of British Ecological Society 1Department of Zoology, Genetics and Physical Anthropology, Universidade de Santiago de Compostela, Santiago de Compostela, Spain 2CRETUS, Universidade de Santiago de Compostela, Santiago de Compostela, Spain 3Department of Functional Biology (Area of Ecology), Universidade de Santiago de Compostela, Santiago de Compostela, Spain Correspondence Andrés Baselga Email: [email protected] Funding information Ministerio de Ciencia, Innovación y Universidades, Grant/Award Number: C G L 2 0 1 6 - 7 6 6 3 7 - P Handling Editor: Jean- Philippe Lessard Abstract 1. The variation of assemblage composition in space is characterised by the decrease in assemblage similarity with spatial distance. Climatic constraint and dispersal limitation are major drivers of distancedecay of similarity. Distancedecay of similarity is usually conceptualised and modelled as an isotropic pattern, that is, assuming that similarity decays with the same rate in all directions. 2. Because climatic gradients are markedly anisotropic, that is, they have different strength in different directions, if species distributions were in equilibrium with climate, the decay of assemblage similarity should be anisotropic in the same direction as the climatic gradient, that is, faster turnover in the direction that maximises the climatic gradient. Thus, deviations from equilibrium between assemblage composition and climatic conditions would result in differences in anisotropy between distancedecay of similarity and climatic gradients. 3. We assessed anisotropy in distancedecay patterns in marine plankton assemblages, terrestrial vertebrates and European beetles, using two procedures: (a) measuring the correlation between the residuals of a distancedecay model and the angle in which pairs of sites are separated and (b) computing two separate distancedecay models for each dataset, one using only pairwise cases that are separated on North- South direction and another one using pairwise cases separated on East- West direction. We also analysed whether the degree of anisotropy in distancedecay is related to dispersal ability (proportion of wingless species and body size) and ecological niche characteristics (main habitat and trophic position) by assessing these relationships among beetle taxonomic groups (n = 21). 4. Anisotropy varied markedly across realms and biological groups. Despite climatic gradients being steeper in North- South direction than in East- West direction in all datasets, North- South distancedecays tended to be steeper than East- West distancedecays in plankton and most vertebrate assemblages, but flatter in European amphibians and most beetle groups. 5. Anisotropy also markedly varied across beetle groups depending on their dispersal ability, as the proportion of wingless species explained 60% of the variance in the difference between North- South and East- West distancedecay slopes. 2 | Journal of Animal Ecology BASELGA And GÓMEZ- ROdRÍGUEZ 1 | INTRODUCTION The similarity between two local species assemblages tends to be higher when they are close together, and decrease with spatial distance. This is because the closer two sites are, the higher the probability that a given species is present in both assemblages (Nekola & McGill, 2014). The decay of assemblage similarity with spatial distance is a widespread pattern across taxa and regions (Soininen et al., 2007). It has a central role in macroecology (Morlon et al., 2008), as this pattern is among the major predictions derived from any of the alternative biodiversity theories (McGill, 2010) that emphasise different causal processes (i.e. dispersal limitation vs. niche processes). This is because the decay of assemblage similarity can be driven either by dispersal limitation or by species sorting across environmental gradients via niche processes (Nekola & White, 1999). Distancedecay models are thus basic tools to understand how biological assemblages vary in space, and to infer the processes underpinning variation in community composition. With this aim, multiple studies have compared the relative explanatory power of climatic and spatial distances for specific biological groups, as dispersal limitation is predicted to produce higher correlations of assemblage similarity with spatial than with climatic distance, while species sorting via ecological niches is predicted to produce higher correlations with climatic than with spatial distance (e.g. Condit et al., 2002; Tuomisto et al., 2003). Another approach to infer the processes behind assemblage variation is to compare distancedecay patterns across biological groups, and using the crosstaxon variation in the relationship between assemblage similarity and climatic distance to infer the degree of equilibrium between assemblages and climatic conditions (Araújo & Pearson, 2005; Lenoir et al., 2012; Thuiller et al., 2019). These studies assume that the higher the degree of equilibrium between assemblage composition and climatic conditions, the higher the correlation between assemblage similarity and climatic distance. This relationship at the assemblage level derives from the fact that species ranges are either (a) constrained by their climatic niches and hence at equilibrium with their suitable climatic conditions or, on the contrary, (b) constrained by dispersal limitation and therefore not in equilibrium with their suitable climatic conditions. As a result, when species are fundamentally constrained by their climatic niches, assemblage composition is controlled by deterministic, nichebased processes, while stochastic processes gain relevance otherwise (Chase & Myers, 2011). In real assemblages, the relative relevance of both type of processes lies in a continuum between both extremes (Leibold & McPeek, 2006). This is a central question to understand how species will be able to track their suitable climatic conditions in the future (Svenning & Skov, 2004), and thus how climate change will affect the distribution of biodiversity (Faurby & Araújo, 2018; Lenoir et al., 2020; Thuiller et al., 2019). Previous research has shown that the degree of equilibrium between species ranges and climatic conditions varies across biological groups (Araújo & Pearson, 2005) and realms, being higher in marine than in terrestrial realms (Burrows et al., 2011; Poloczanska et al., 2013; Sunday et al., 2012). The decay of assemblage similarity with spatial distance is usually modelled as an isotropic pattern, even though the assumption is rarely made explicit (but see Baselga, 2007; Tovo & Favretti, 2018; Vellend, 2001) or tested (Burley et al., 2012; Svenning et al., 2011). Isotropic spatial patterns are spatially symmetric, and the magnitude and range of spatial autocorrelation do not vary with direction (Dale & Fortin, 2014). That is, from any focal point, the measured attribute varies with the same strength in any direction. In contrast, anisotropic spatial patterns are those in which the variation is more marked in one particular direction than in others. This is the case, for example, of environmental space (Soberón & Peterson, 2011), as environmental gradients are usually stronger in particular directions. This is especially true for climatic gradients at large spatial scales, because climate is mostly structured in latitudinal belts (Peel et al., 2007). This leads to the prediction that, if climatic niches were the major drivers of species distributions, then species ranges would tend to be constrained within these bands via their climatic niches, and spatial assemblage turnover (i.e. species replacement) should be faster in North- South direction (following the climatic gradients) than in East- West direction (within climatic belts). On the contrary, in an idealised scenario with a continuous land mass and no effect of climatic niches, species assemblages separated by a given spatial distance are expected to be equally similar in the East- West and North- South directions (isotropic distancedecay of similarity). Additionally, besides the general tendency of macroclimatic structure to favour anisotropy by leading to steeper decays of similarity in the North- South than the East- West direction, the specific spatial disposition of dispersal barriers (e.g. landmasses shape) can produce idiosyncratic anisotropic pattern in any direction. Therefore, we here introduce a novel directional distancedecay approach. We propose that assessing how the relationship between assemblage similarity and spatial distance varies across directions can help inferring the degree of equilibrium between assemblage variation and climatic gradients at macroecological scales (Figure 1). 6. Our results suggest that the degree of equilibrium decreases from marine to terrestrial realms, and is markedly different between vertebrates and beetles. This has profound implications on the expected ability of different groups to track their suitable climates, and thus on the impact of climate change on biodiversity. KEYWORDS anisotropy, beta diversity, climatic niches, dispersal limitation, equilibrium, isotropy, turnover | 3 Journal of Animal Ecolog y BASELGA And GÓMEZ- ROdRÍGUEZ Taking this into account, measuring anisotropy in distancedecay patterns and comparing it with the one observed along climatic gradients would allow inferring the relative importance of climatic gradients or dispersal limitation as drivers of assemblage variation. If a similar anisotropic pattern is observed in the climatic gradient and the decay of assemblage similarity, a major role of climatic gradients would be inferred. On the contrary, if the decay of assemblage similarity is isotropic or anisotropic in a different orientation than the climatic gradients, a major role of dispersal limitation would be inferred. Therefore, across biological groups, we expect that the degree of anisotropy will vary in accordance to their taxonspecific niche characteristics and dispersal abilities. To test this, we introduce two methods for assessing anisotropy in distancedecay models and address the following questions: (Q1) whether anisotropy in distancedecay models varies across biological groups, realms and geographical regions, (Q2) whether it differs from the anisotropy of climatic gradients in these regions, pointing to disequilibrium between assemblage composition and climatic conditions and (Q3) whether the degree of anisotropy in distancedecay patterns is related to morphological characteristics of the biological group associated with its dispersal ability (flight capacity and body size) or their ecological strategies (trophic level and habitat preference). FIGURE 1 Schematic representation of isotropic and anisotropic patterns of decay of assemblage similarity with spatial distance. For simplicity in the visual representation, only some of all possible pairs are considered, identifying the pairs separated in angles of 0° (blue), 45° (orange) and 90° (red) (a). If distancedecay of similarity is isotropic (b– d), similarity with a focal site would decrease at the same rate in all directions (b), so residuals of a distancedecay model would not be correlated with the angle (c) and the distancedecay curves fitted only with pairs separated in a specific direction (e.g. 0° or 90°) would be similar. However, if distancedecay of similarity is anisotropic (e– g), similarity with the focal site would decrease faster in some directions than in others (e), so residuals of a distancedecay model would be correlated with the angle (f), and the distancedecay curves fitted only with pairs separated in a specific direction (e.g. 0° or 90°) would be different. In this hypothetic example, residuals would be positive for site pairs separated in East- West direction and negative for pairs separated in North- South direction, and the slope of a distancedecay model fitted for pairs separated in North- South direction would by steeper than the slope in East- West direction (g) 4 | Journal of Animal Ecology BASELGA And GÓMEZ- ROdRÍGUEZ 2 | MATERIALS AND METHODS 2.1 | Data Two different types of data are used in this study. Differences in anisotropy among biological systems, taxonomic groups and geographical regions were first assessed with presence/absence data of marine plankton (Australian phytoplankton and copepods) and terrestrial vertebrates (amphibians, mammals and birds in Europe, North America and South America). The regular grid and high spatial resolution of this data (1° squares) provide high statistical power to detect and quantify anisotropy. However, because few biological groups are considered and they have disparate morphological traits, these data cannot be used to assess the relationship between the degree of anisotropy and the biological characteristics of the groups. With this aim, additional analyses were conducted based on presence/absence data of 21 beetle groups in continental European territories. These data allow comparing the degree of anisotropy among biological groups with different ecological requirements and dispersal ability (Baselga et al., 2012; Gómez- Rodríguez & Baselga, 2018; Gómez- Rodríguez et al., 2015), although the distancedecay patterns in each group are estimated with smaller sample size, lower spatial resolution and irregular size of sampling units, compared to vertebrate data. The marine plankton data were downloaded from the Australian Ocean Data Network (https://portal.aodn.org.au/). We built two separate presence/absence tables for phytoplankton (Richardson, Coman, Davies, McLaughlin, et al., 2020) and copepods (Richardson, Coman, Davies, Slotwinski, et al., 2020), respectively. To do this, we used the subset of species in the cited databases for which real presence/absences were recorded (i.e. we removed all species with no observations in any of the sampling points). We then aggregated observations in 1° × 1° cells, so we got two presence/absence tables in these 1° cells, for phytoplankton and copepods, respectively (hereinafter referred as ‘marine plankton datasets’, n = 2). The vertebrate data derive from species distribution maps of amphibians and terrestrial mammals from the IUCN (2013), and the breeding distributions maps of birds from the BirdLife International and NatureServe (2013) database. For each biological group (amphibians, mammals, birds) and in each land mass (Europe, North and South America), polygonal range maps were converted into presence/absence tables of 1° × 1° cells (hereinafter referred as ‘vertebrate datasets’, n = 9). In both plankton and vertebrate datasets, cell geographical coordinates were converted into UTM coordinates to allow computing Euclidean spatial distances among cells. The beetle data include species presence/absence tables for 21 beetle supraspecific taxa (large genera, tribes, subfamilies or families; 4,648 species in total) in 34 continental European territories (mostly countries > 2000 km2; extent: 11°W– 60°E longitude and 36°N– 72°N latitude) and was compiled from Löbl and Smetana (2003, 2004, 2006, 2010) and the Fauna Europaea Web Service (Fauna Europaea version 1.1, available online at https://fauna - eu.org/). The 21 supraspecific taxa (see Table 1) are representative of the major beetle lineages (superfamilies Caraboidea, Staphylinoidea, Scarabaeoidea, Chrysomeloidea and Curculionoidea) and include a wide spectrum of ecological strategies and dispersal abilities. Detailed information about this data can be found in Baselga, Lobo, Svenning, Aragón, et al. (2012). Supraspecific taxa were categorised according to their main habitat (epiedaphic, edaphic, aquatic and epiphytic) and trophic specialisation (coprophagous, zoophagous, xylophagous, phytophagous), and the proportion of apterous species (i.e. without functional wings) and species mean size were extracted from taxonomic monographs, as described in Baselga, Lobo, Svenning, Aragón, et al. (2012), Gómez- Rodríguez et al. (2015) and Gómez- Rodríguez and Baselga (2018). 2.2 | Q1. Measuring anisotropy in distancedecay models across biological groups, realms and geographical regions Distancedecay patterns were estimated from the relationship between the Simpson index of similarity and spatial distance. Simpson similarity (i.e. 1— Simpson dissimilarity) is independent of richness differences, thus reflecting only assemblage turnover (Baselga, 2010; Koleff et al., 2003; Simpson, 1960). Simpson similarity was computed with function beta.pair() in R package betapart (Baselga & Orme, 2012). The relationship between pairwise similarities and Euclidean spatial distance (between centroids of countries or square cells, depending on the dataset) was modelled as a negative exponential function with function decay.model() in betapart, which fits a GLM with log link and Gaussian error to the similarity matrix (Gómez- Rodríguez & Baselga, 2018). All calculations were done in R (R Development Core Team, 2018). We introduce two procedures to assess anisotropy in distancedecay patterns. In the first approach (anglebased method), anisotropy is measured as the correlation between the residuals of a distancedecay model and the angle in which pairwise dissimilarities and spatial distances are measured, from 0° in pairs of localities separated in East- West direction to 90° in pairs of localities separated in North- South direction. Under the hypothesis of isotropy in distancedecay patterns, there should be no relationship between residuals and the angle in which localities are separated (Figure 1b), while if the distancedecay is anisotropic, residuals would tend to be larger in a particular direction (angle), making residuals and angles to be correlated (Figure 1d). The robustness of this approach derives from the large sample size used in distancedecay computation (i.e. all pairwise distances are used). However, its power could be low if anisotropy occurs in a particular direction and the signal is lost among the noise. To account for this, a second approach (latitudeversus- longitude method) maximises the contrast between dissimilarities computed in bands following orthogonal directions (North- South and East- West in our case), as pairwise cases that are separated in intermediate angles (e.g. cases in the diagonal) are not considered. Thus, two different distancedecay models are built for each dataset, one using pairwise cases that are mostly separated on North- South direction and another one using pairwise cases that are mostly separated on East- West direction. North- South versus | 5 Journal of Animal Ecolog y BASELGA And GÓMEZ- ROdRÍGUEZ East- West directions were selected because climatic niche differences are expected to be stronger in the North- South direction due to latitudinal climatic gradients. In an isotropic distancedecay pattern, we would observe the same slopes in both models (Figure 1c) while in an anisotropic scenario, we would expect significantly different slopes in these distancedecay patterns (Figure 1e), given that higher turnover in species composition would be observed following the direction of stronger environmental variation (i.e. North- South). For the anglebased method, a distancedecay pattern was fitted using all pairwise similarities in each dataset. The degree of anisotropy was computed as the Spearman rank correlation between the decay model residuals and the angle in which site pairs were separated. Angles were computed as the arcsine of the ratio between North- South distance and the Euclidean distance (note that the former is the opposite cathetus and the latter the hypotenuse, so this ratio is the sine of the angle). The significance of the correlation was empirically assessed by comparing the observed correlation value with the distribution of nullmodel correlations resulting from a permutation test randomising distancedecay model residuals (n = 1,000 times for beetles, and 100 times for vertebrates and marine plankton, due to computational constraints derived from the high number of cases). Dataset rho p Difference in climatic slopes Difference in distancedecay slopes Marine plankton Australian phytoplankton −0.30 <0.01 0.00463 −0.00053 Australian copepods −0.41 <0.01 0.00468 −0.00075 Terrestrial vertebrates European amphibians 0.20 <0.01 0.00044 0.00015 European mammals −0.29 <0.01 0.00044 −0.00012 European birds −0.15 <0.01 0.00044 −0.00004 North American amphibians −0.20 <0.01 0.00014 −0.00019 North American mammals −0.43 <0.01 0.00014 −0.00033 North American birds −0.21 <0.01 0.00014 −0.00007 South American amphibians −0.21 <0.01 0.00008 −0.00025 South American mammals −0.22 <0.01 0.00008 −0.00015 South American birds −0.22 <0.01 0.00008 −0.00010 European beetles Amara 0.01 0.420 0.00074 −0.00001 Apionidae 0.04 0.179 0.00074 −0.00002 Carabus 0.21 <0.001 0.00074 0.00019 Cerambycinae 0.10 0.008 0.00074 −0.00001 Chaetocnema 0.04 0.191 0.00074 −0.00003 Chrysomelinae 0.16 0.001 0.00074 0.00010 Cryptocephalinae 0.13 <0.001 0.00074 0.00008 Hydraenini 0.15 0.001 0.00074 0.00014 Lamiinae 0.08 0.035 0.00074 0.00001 Lepturinae −0.16 <0.001 0.00074 −0.00005 Limnebiini −0.05 0.865 0.00074 0.00003 Ochtebiini −0.02 0.703 0.00074 0.00002 Otiorhynchus 0.17 <0.001 0.00074 0.00038 Phyllobiini 0.08 0.031 0.00074 0.00010 Polydrusini 0.10 0.009 0.00074 0.00006 Pselaphinae 0.21 <0.001 0.00074 0.00020 Pterostichus 0.19 <0.001 0.00074 0.00016 Scarabaeidae 0.20 <0.001 0.00074 0.00002 Scolytinae −0.22 <0.001 0.00074 −0.00013 Silphidae −0.20 <0.001 0.00074 0.00001 Trechus 0.33 <0.001 0.00074 0.00044 TABLE 1 Relationship between the distancedecay residuals and the angle in which the assemblages were separated (Spearman rho and p value), difference between N- S and E- W slopes of the relationship between climatic and spatial distance, and difference between N- S and E- W distancedecay slopes, for Australian marine plankton, American and European vertebrates and European beetles 6 | Journal of Animal Ecology BASELGA And GÓMEZ- ROdRÍGUEZ For the latitudeversus- longitude method, we fitted two separate distancedecay models for each dataset, one for the North- South distances (North- South distancedecay) and one for the East- West distances (East- West distancedecay). North- South distances and East- West distances were calculated as the difference between latitude coordinates (in km) and the difference in longitude coordinates (in km), respectively. To minimise the confounding effect of North- South distances on East- West distances, and vice versa, we restricted the East- West distance matrix to the subset of pairs with North- South distance <300 km and the North- South distance to the subset of pairs with East- West distance <300 km. This reduced the number of pairs considered, which still remained >100 in the beetle datasets (n = 102 in E- W bands, n = 111 in N- S bands) and >17,000 in all other datasets. The degree of anisotropy was computed as the difference between the slopes of the North- South and East- West decay models (Δb). 2.3 | Q2. Differences in anisotropy between distancedecay and climate Climatic gradients are expected to be steeper in North- South than in East- West direction. We tested this assumption by (a) assessing the relationship between climatic distance and North- South and East- West spatial distances, respectively and (b) assessing the relationship between assemblage similarity and climatic distance using the same subsets (North- South vs. East- West distancedecay patterns) as defined above. Climatic data for terrestrial vertebrates and beetles, originally at 2.5 arcminutes resolution, were downloaded from the WorldClim database (Hijmans et al., 2005) and rescaled to a 1° resolution grid, the same as the vertebrate data. To compute climatic distance, we submitted Mean Annual Temperature and Annual Precipitation to a PCA using command principal() in R package psych (Revelle, 2018), and computed climatic distance as the Euclidean distance between the scores of sites in the twodimensional climatic space. In the case of marine datasets, climate was characterised with monthly average sea surface temperature (downloaded from NOAA Global Surface Temperature database at https://psl.noaa.gov/). Then, climatic distance was regressed against North- South and East- West distances. In addition, assemblage similarity was modelled as a negative exponential function of climatic distance in North- South and East- West bands as explained above for spatial distances. Given the potential covariance between spatial and climatic distance, we also built a model with both predictors and partitioned the explained variation (pseudor2) into fractions uniquely explained by spatial or climatic distances, or jointly explained by both predictors. 2.4 | Q3. Relationship between anisotropy in distancedecay patterns and dispersal and ecological characteristics Finally, using the European beetle data, we assessed whether the anisotropy in distancedecay patterns could be predicted from dispersal and niche characteristics of taxa. To do that, we regressed the difference between the slopes of North- South and East- West distancedecay models (Δb) against two dispersalrelated characteristics (proportion of apterous species and mean body size) and two niche characteristics (trophic position and main habitat). Independent models were fitted within each set of predictors (dispersal or niche) using a forward stepwise procedure based on the F- statistic. If both dispersal and niche predictors explained a significant proportion of variance in anisotropy (Δb), variance partitioning was used to estimate the fractions of variance jointly or uniquely explained by each set of predictors. 3 | RESULTS 3.1 | Q1. Measuring anisotropy in distancedecay models across biological groups, realms and geographical regions The degree of anisotropy varied across biological groups and geographical regions. The relationship between distancedecay residuals and the angle in which the assemblages were separated (anglebased method, Table 1) was negative in the marine plankton datasets, and all American and European vertebrate groups, with the exception of European amphibians, for which the relationship was positive. A negative correlation implies that the similarity of sites separated in East- West direction tends to be larger than expected from the distancedecay model, but similarity of sites separated in North- South direction tends to be lower than expected. In other words, in the studied marine plankton datasets and most vertebrate datasets similarity decays faster in North- South than in East- West direction. In contrast, European beetles showed a significant, positive correlation between the residuals of the exponential decay model and the angle in 13 biological groups (out of 21). That is, sites separated in East- West direction tend to have lower similarity than predicted by the models, while sites separated in North- South direction tend to have higher similarity than predicted, that is, similarity decays faster in East- West than in North- South direction. In the remaining beetle groups, we found a significant negative correlation in three groups, and a negligible, not significant correlation (i.e. suggesting a high degree of isotropy) in five groups. The largest correlation values (i.e. stronger anisotropy) corresponded to North American mammals (Spearman rho = −0.43) and Australian marine copepods (rho = −0.41) among the negative values, and the European beetle genera Trechus (rho = 0.33) and Carabus (rho = 0.21), and the subfamily Pselaphinae (rho = 0.21) among the positive values. The latitudeversus- longitude method evidenced that in marine plankton and in terrestrial vertebrates the slope of the North- South distancedecay was significantly steeper (p < 0.01) than the slope of the East- West decay, again with the exception of European amphibians (Figures 2 and 3). The difference in slopes between North- South and East- West distancedecays was larger in the marine plankton | 7 Journal of Animal Ecolog y BASELGA And GÓMEZ- ROdRÍGUEZ FIGURE 2 North- South (red) versus East- West (blue) distancedecay models for Australian marine phytoplankton and copepods. Lines are negative exponential functions separately fitted for pairs of sites separated in North- South (red) or East- West (blue) direction. Goodness of fit was assessed with pseudor2 values (i.e. 1 – model deviance/null deviance). North- South decay slopes were significantly steeper in both phytoplankton and copepod datasets FIGURE 3 North- South (red) versus East- West (blue) distancedecay models for European, North and South American amphibians, terrestrial mammals and birds. Lines are negative exponential functions separately fitted for pairs of sites separated in North- South (red) or East- West (blue) direction. Goodness of fit was assessed with pseudor2 values (i.e. 1 – model deviance / null deviance). North- South decay slopes were significantly steeper in all cases, except European amphibians, for which the East- West decay slope was significantly steeper 8 | Journal of Animal Ecology BASELGA And GÓMEZ- ROdRÍGUEZ assemblages than in the terrestrial vertebrates (Table 1). In contrast, in nine groups (out of 21) of European beetles, the slope of the East- West distancedecay was significantly steeper than the slope of the North- South distancedecay (p < 0.01) while the slope of the North- South distancedecay was significantly stepper in only three groups and not significantly different in nine groups (Figure 4). In other words, as already indicated by the angledbased method, in plankton and vertebrate assemblages the decay in assemblage similarity tended to be faster in North- South direction, while in European beetles the decay of assemblage similarity tends to be faster in East- West direction, or very similar in East- West and North- South directions. 3.2 | Q2. Differences in anisotropy between distancedecay and climate Climatic data were highly anisotropic, as climatic distance was more tightly related to North- South distance in all the regions: climatic distance was well explained by North- South distances in both terrestrial (r2 = 0.42, 0.62 and 0.37 in Europe, North America and South America, respectively) and marine systems (r2 = 0.35 and 0.47 in Australian phytoplankton and copepod datasets, respectively), but not by East- West distances either in terrestrial (r2 = 0.13, 0.20 and 0.15 in Europe, North America and South America) or, particularly, in marine systems (r2 = 0.03 and r2 = 0.05 in Australian phytoplankton and copepod datasets, respectively). The tighter relationship between climatic and North- South distances reflects the existence of marked latitudinal climatic gradient in all the datasets. Contrary to what would be expected if anisotropy of distancedecay patterns was the result of anisotropy in climatic gradients, the difference between climatic North- South and East- West slopes was larger in Europe than in North and South America (Table 1; Figure S1), but the difference between North- South and East- West slopes of the distancedecay in assemblage similarity were larger in North and South America than in Europe. The decays of assemblage similarity with climatic distance within North- South and East- West bands showed steeper slopes and larger pseudor2 values in North- South than in East- West direction in all datasets (Figures S2– S4). In the marine plankton datasets, climatic and spatial distancedecay patterns showed virtually identical pseudor2 values, and the variation partitioning procedure revealed that assemblage similarity was jointly explained by spatial and climatic distances, with negligible unique contributions (Figure 5; Figure S5). In terrestrial vertebrates, pseudor2 values were generally larger for spatial than climatic distances, leading to unique contributions of spatial distances larger than those of climatic distances (Figure 5; Figure S5), both in North- South and East- West directions. In turn, in European beetles, climatic distance explained less variation in assemblage similarity than spatial distance, particularly in East- West direction. As a result, in many beetle groups, the unique contribution of spatial distance was the largest fraction, and the difference in pseudor2 was usually larger in East- West direction (Figure 5; Figure S5). 3.3 | Q3. Relationship between anisotropy in distancedecay patterns and dispersal and ecological characteristics When we assessed whether anisotropy (i.e. the difference between the slopes of North- South and East- West distancedecays, Δb) FIGURE 4 North- South (red) versus East- West (blue) distancedecay models for 21 European beetle groups. Lines are negative exponential functions separately fitted for pairs of sites separated in North- South (red) or East- West (blue) direction. Goodness of fit was assessed with pseudor2 values (i.e. 1 – model deviance/null deviance). East- West decay slopes were significantly steeper in nine cases (marked with *), North- South decay slopes were significantly steeper in three cases (marked with ^) and no significant differences were found in the remaining nine cases | 9 Journal of Animal Ecolog y BASELGA And GÓMEZ- ROdRÍGUEZ across beetle groups was related to dispersalrelated and niche characteristics of these taxa, it turned out that only the proportion of apterous species (r2 = 0.60, F1,19 = 29.06, p < 0.001) and main habitat (r2 = 0.23, F1,19 = 5.53, p = 0.030) were significantly related to Δb. A model including both variables (r2 = 0.61, F2,18 = 13.92, p < 0.001) explained almost the same amount of variance as the proportion of apterous species alone, implying that main habitat had a negligible unique contribution to the model (1% of the variance), while the proportion of apterous species uniquely explained 38% of the variance in Δb. 4 | DISCUSSION Our results show that distancedecay patterns in marine plankton, vertebrates and beetles deviate from isotropy in strikingly different ways. This has profound implications on our understanding of the processes driving global macroecological patterns, because distancedecay patterns are particularly informative about the spatial structure of ecological assemblages (Morlon et al., 2008). Moreover, we show that the variation in anisotropy of distancedecay patterns across datasets is not related to the anisotropy of climatic gradients, a relationship that would be expected to exist if species distributions were in equilibrium with climate (compare Figures 2– 4 to Figure S1). This suggests that the degree of equilibrium between biological assemblages and climate is not constant across realms (i.e. marine vs. terrestrial), biological groups and land masses. Deviations from equilibrium, as inferred from our analysis of anisotropy in distancedecay patterns, seem much higher in European beetles than in terrestrial vertebrates and marine plankton. Australian planktonic assemblages and American and European terrestrial vertebrates tend to have steeper distancedecay patterns in North- South than in East- West direction, and the difference is more marked in marine than terrestrial datasets. In contrast, most European beetles show the opposite trend, with steeper decays in East- West than in North- South direction. This suggests that ecological niches are a greater constraint in marine plankton and terrestrial vertebrates while dispersal limitation is probably more relevant in European beetles, which is further supported by the strong relationship between assemblage similarity and climatic distance in North- South bands, but not in East- West bands. This interpretation is also consistent with previous results showing (a) that the degree of equilibrium between species assemblages and climate depends on the relative dispersal abilities of the biological group (Araújo & Pearson, 2005) and (b) that marine species are able to track their suitable climatic conditions much faster than terrestrial species (Antão et al., 2020; Burrows et al., 2011; Lenoir et al., 2020; Sunday et al., 2012), which is probably also related to the greater dispersal potential of species in the ocean than on land (Kinlan & Gaines, 2003) and a smaller environmental heterogeneity in marine than in terrestrial systems (Antão et al., 2020). Therefore, our results show that anisotropy of assemblage turnover at large spatial scales is a relevant macroecological property that informs about the degree of equilibrium of biological assemblages with climate. This is a crucial parameter for predicting the effects of climate change on biodiversity (Bellard et al., 2012; Thuiller et al., 2019), and our results suggest that the degree of equilibrium of different biological groups can be estimated from the anisotropy of their distancedecay patterns. Within the terrestrial realm, amphibians, mammals and birds have steeper distancedecays in North- South than in East- West direction (except European amphibians), while most European beetles show similar slopes in both directions or even have steeper distancedecays in East- West than in North- South directions. Given that climatic gradients are more marked in North- South direction in all datasets, our results thus suggest that the replacement of vertebrate species with spatial distance is faster when climatic differences are larger. However, many European beetle groups have isotropic distancedecay patterns or show the opposite deviation from distancedecay isotropy, with faster species replacement in East- West than in North- South direction, even if climatic gradients are much more marked in North- South direction. This is particularly relevant because beetles are ectotherms so, if beetle species ranges were not dispersallimited, we would expect beetle assemblages to be more constrained by climatic gradients than endotherms. We observe the opposite result, which suggests a major effect of dispersal limitation on beetle assemblage similarity, as already pointed in previous studies (Gómez- Rodríguez & Baselga, 2018; Gómez- Rodríguez et al., 2015). The lack of marked turnover gradients in North- South direction is likely the result of incomplete postglacial recolonisation FIGURE 5 Percentage of explained variance accounted uniquely by spatial and climatic distance, or jointly by both, in North- South (left) or East- West bands (right). Each dot represents a dataset: colour denotes the biological system (blue for marine plankton, red for terrestrial vertebrates and orange for European beetles), and size is proportional to the total explained variance