On the mismatch in the strength of competition among fossil and modern species of planktonic Foraminifera

This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited. © 2019 The Authors. Global Ecology and Biogeography published by John Wiley & Sons Ltd 1Ocean and Earth Science, National Oceanography Centre Southampton, University of Southampton, Southampton, United Kingdom 2Center for Marine Environmental Sciences, University of Bremen, Bremen, Germany 3Department of Zoology, University of British Columbia, Vancouver, British Columbia, Canada 4Instituto de Física Teórica, Universidade Estadual de São Paulo, Sao Paulo, Brazil 5School of Geography, Earth and Environmental Science, University of Birmingham, Birmingham, United Kingdom


| INTRODUC TI ON
Ecological interactions play an important role in shaping the diversity of life across space and time. By interacting with each other, individuals change the dynamics of their populations and communities. Population dynamics, in turn, influence species' abundances and geographical ranges, which affect species' diversification rates (i.e., the balance between speciation and extinction rates) (Harnik, 2011;Jablonski, 2017). However, the resulting net effect of ecological interactions on species' diversification might not correspond to the fitness effects at the individual level. For example, positive interactions such as mutualism can have a negative effect on diversification, and negative interactions such as predation (for the prey) and parasitism (for the host) can enhance diversification (Jablonski, 2008a;Yoder & Nuismer, 2010). How individual-level interactions scale up to affect species' diversification remains a topic that has received scant attention (Jablonski, 2017), and we currently lack understanding of how congruent or discordant these effects are across ecological and geological time scales (Harmon et al. 2019;Weber, Wagner, Best, Harmon, & Matthews, 2017).
Competition among species has, by definition, a negative effect on the fitness of the interacting individuals. Yet, interspecific competition can affect species' diversification positively, by promoting niche partitioning and character displacement (Bailey, Dettman, Rainey, & Kassen, 2013;Meyer & Kassen, 2007) or negatively by driving species to extinction due to niche overlap (Bengtsson, 1989).
Within macroevolutionary research, interspecific competition has often been proposed to have a consistent negative effect across ecological and geological time scales, generating the pattern of a negative relationship between species diversity and diversification rates (Rabosky, 2013). This negative diversity-dependent diversification (DDD) pattern is thought to reflect how increases in species richness suppress diversification rates due to increased extinction or because opportunities for speciation become reduced once niches are filled, whereas decreases in species richness promote diversification because of weaker competitive interactions (Rabosky, 2013).
Many clades have been shown to display the negative DDD pattern Phillimore & Price, 2008;Pires, Silvestro, & Quental, 2017;Rabosky, 2013). Although most of these studies evoke interspecific competition to explain the observed macroevolutionary pattern, they do not explicitly test for it. Such an explicit test is difficult to perform, as it requires not only a sufficiently well-preserved fossil record to assess the deep-time individual-or population-level fitness responses to ecological interactions, but also an actual fossilized interaction proxy, both of which are rare (but see Liow, Martino, Voje, Rust, & Taylor, 2016). An alternative approach to test whether the negative DDD pattern is a result of upward causation from competition among individuals is to contrast macroevolutionary and present-day community ecology dynamics (Voje, Holen, Liow, & Stenseth, 2015;Weber et al., 2017). This approach involves investigating whether the ecological interactions among most (if not all) extant species of a clade (clade-wide) support interspecific competition, and are thus consistent with the clade's observed DDD pattern. To our knowledge, none of the clades that have shown the negative DDD pattern have been studied in a clade-wide community ecology framework.
Planktonic Foraminifera (PF) represent a useful model system for integrating macroevolution and community ecology (Yasuhara, Tittensor, Hillebrand, & Worm, 2017). They are marine unicellular zooplankton that produce calcium carbonate shells, which accumulate on the seafloor building up a remarkably complete fossil record Schiebel, 2002). The PF fossil record of the last 66 million years shows a negative DDD pattern (Figure 1), either by analysing the phylogenetic tree structure  considering biotic and abiotic correlates of diversification , or using nonlinear ecological models of alternative modes of competition (Ezard & Purvis, 2016). These ecological models suggest that although abiotic factors affect the dynamics of PF diversification, their diversification is driven by compensatory contest competition (Ezard & Purvis, 2016), which is analogous to a populational logistic growth where resource availability is limited and only successful individuals get the amount of resource they require (Hassell, 1975).
While competition among PF species likely affected their diversification dynamics Ezard et al., 2011;Ezard & Purvis, 2016), we do not know if (or how) species compete today. PF are difficult to culture because they do not reproduce under laboratory conditions (Hemleben, Spindler, & Anderson, 1989), hindering our understanding of their population dynamics and ecological interactions. Globally, there are 47 extant PF species (Siccha & Kucera, 2017), occurring in low densities when compared to other eukaryotic groups (Keeling & del Campo, 2017). PF prey on other plankton groups (e.g., phytoplankton, ciliates and copepods), and several species host photosymbionts (Hemleben et al., 1989;Takagi et al., 2019); however, we currently lack understanding of the role of ecological interactions in shaping extant PF diversity. Assessing if extant PF species compete is valuable not only to bridge macroevolution and community ecology, but also to improve ecological models used consider that diversification dynamics are influenced by groups that interact ecologically even when distantly related.
Here, we integrate PF abundance data globally in space and through time to investigate whether ecological signatures are expressed consistently across macroecological and macroevolutionary scales ( Figure 2). The temporal scale we study ranges from weeks (sediment trap data) to millennia (coretop data). We explore two community-level patterns expected under interspecific competition: local competitive exclusion (spatial pattern) and compensatory dynamics (temporal pattern). When resources are limiting, ecologically similar species (e.g., similar size or diet) are expected to compete more intensely, and ultimately exclude one another from the local community (Cavender-Bares, Kozak, Fine, & Kembel, 2009;Webb, Ackerly, McPeek, & Donoghue, 2002). Thus, we expect competitive exclusion to lead to the spatial segregation of ecologically similar species, resulting in overdispersed local communities (Cavender-Bares et al., 2009;Pearse, Purvis, Cavender-Bares, & Helmus, 2014;Webb et al., 2002). Yet if competing species coexist, population sizes of competitors are expected to show a negative relationship F I G U R E 1 Bounded diversification pattern observed in the planktonic Foraminifera fossil record. Number of species (black line) through geological time in millions of years (My) ago, zero represents the present. The coloured bars represent the diversification rate per million-year time bin, calculated by dividing the number of new and extinct species per bin . Blue bars represent positive diversification rates, while red bars represent negative rates. We removed the most recent time bin to avoid the pull of the present effect (i.e., species not having had time to go extinct yet). We used the Aze et al. (2011) lineage phylogeny of macroperforate Cenozoic planktonic Foraminifera for this figure F I G U R E 2 (a) Map of modern planktonic Foraminifera (PF) abundance data. Green dots, spatial data: 3,053 ocean-floor coretop samples, each with data on relative abundance of all the species present in the sample (usually larger than 150 μm; Siccha & Kucera, 2017). Orange triangles, temporal data: 35 sediment traps, each with data on abundance (shell flux) time series of PF species. Not all species were identified in the sediment trap samples. In total there are 370 population time series available (Supporting Information Table S2). (b) Latitudinal richness gradient of PF for each 5° bin, based on coretop samples with richness equal to or greater than two species through time (compensatory dynamics; Houlahan et al., 2007). Thus, we expect that, depending on the strength of interspecific competition, ecologically similar species will either show complete local competitive exclusion, or population abundances will covary negatively through time. If ecological interactions have a consistent effect across geological and current ecological time scales, we expect competition among PF species to play a key role in structuring their modern communities.

| MATERIAL AND ME THODS
We used two types of observational data: spatial data from 3,053 coretop samples and temporal data from 35 sediment traps (370 time series of population abundances in total) (Figures 2 and 3). We performed two distinct and independent analyses for each data type: spatial data were analysed using a community phylogenetics approach, considering random assembly null models. Temporal data were analysed by pairwise correlation of time series, considering a null model of species' annual seasonality. All analyses were conducted using R (version 3.3.3, R Core Team, 2017).

| Spatial data
Ocean-floor sediment samples are recovered using core sampling ( Figure 3); the coretop represents the most recent (surface) sediment. Because there is no precise estimate of local sedimentation rate in most sampled sites, core samples hold information on species relative abundances (relative to the sampled assemblage, as all individuals are identified) rather than absolute abundances. Coretop samples represent time-averaged assemblages on the order of hundreds to thousands of years depending on local sedimentation and bioturbation rates (Jonkers, Hillebrand, & Kucera, 2019). Such timeaveraged assemblages have the advantage of averaging out seasonal and inter-annual fluctuations in abundance that might blur macroecological patterns (Fenton, Pearson, Jones, & Purvis, 2016).
We used relative abundance data from the ForCenS database, which is a curated database, taxonomically harmonized and treated for sampling inconsistencies (Siccha & Kucera, 2017). Sites deeper than 3,500 m in the Pacific and Indian oceans and 4,500 m elsewhere were excluded to minimize the effects of calcium carbonate dissolution on assemblage composition (Tittensor et al., 2010). We also removed samples that only contained one species and those that did not differentiate between the species Globorotalia menardii and Globorotalia tumida; Globigerinoides ruber pink and white; or Turborotalita humilis and Berggrenia pumilio. With these restrictions, the data we analysed consisted of 3,053 unique sites of relative abundances of 41 PF species distributed across the globe (Figure 2a).  Table   S2) including 37 PF species, in total 370 time series of population F I G U R E 3 Schematic figure showing the two different types of data used in our analyses: sediment traps (temporal data) provide a virtually continuous time series of settling, dead planktonic Foraminifera shells. An upward-facing funnel (yellow triangle) directs sinking particulate material towards a sampling bottle. At pre-determined time intervals, a new sampling bottle is shifted into position. Core samples (spatial data) are a cylindrical section of the ocean floor sediment obtained by coring with a hollow tube. The top of the core (coretop) contains the most recently settled particles, and down core samples contain older, fossil material. The ship silhouette represents the German Research Vessel Meteor shell fluxes (11 species per trap on average). We used the sediment traps selected by Jonkers and Kucera (2015), who harmonized the taxonomy on these samples. The total duration of the sampling ranged from a minimum of 1 to 12 years (mean 3.4 years), and the mean length of the time series was 59 samples. Sampling intervals (resolution) varied between four and 58 days (mean 18 days). For more information, see Supporting Information Table S2 and Jonkers and Kucera (2015).

| Temporal data
Shell fluxes represent the settling of dead individuals and indirectly represent the population abundance in the water column (Jonkers & Kucera, 2015). PF have a short life span (typically 1 month; Hemleben et al., 1989), thus shell fluxes are likely to be a good proxy of population abundance (standing stock). The comparison of shell fluxes between different species makes the underlying assumption that species have similar life spans and settling time, which is reasonable given the average sampling resolution. As we expect competition to affect individuals' mortality, the dynamics of settling dead individuals are assumed to represent the dynamics of competition among living individuals.

| Ecological similarity data
We use three different measures to capture ecological similarity between species: phylogenetic relatedness, species' average size and presence or absence of photosymbionts. Phylogenetic relatedness can be used as a proxy for ecological similarity between species under the assumption of niche conservatism (Cavender-Bares et al., 2009;Wiens et al., 2010). PF are divided into two monophyletic groups, microperforates and macroperforates (Aurahs et al., 2009).  Figure S1). The final phylogeny includes 33 species (Supporting Information Figure S1); the 14 species not present in the phylogeny (e.g., microperforates) were excluded from the phylogenetic analysis.
We used similarity between species' average sizes as an alternative proxy for ecological similarity, because zooplankton body size is a functional trait related to physiological and ecological processes (Litchman, Ohman, & Kiørboe, 2013;Sauterey, Ward, Rault, Bowler, & Claessen, 2017). We estimated average size of 36 PF species, by compiling shell diameter data of 2,774 individuals from three different datasets (Baranowski, 2013;Rillo, Miller, Kucera, & Ezard, 2019;Weinkauf, Kunze, Waniek, & Kucera, 2016; Supporting Information Table S1). Although there is high intraspecific size variation within species, size differences among species are evident (Supporting Information Figure S2) and assumed to represent species-level niche differences. The 11 species for which no size data are available were excluded from the size analysis.
Thus, photosymbiosis seems to be a successful ecological strategy in nutrient-limited environments, possibly allowing for greater flexibility in resource acquisition for growth and reproduction (Takagi et al., 2019). We use presence or absence of symbionts (symbiotic strategy), recently measured for 30 species (Takagi et al., 2019), as a third proxy for ecological similarity between species. We considered species that display low to high intensity of photosymbiosis as phototrophs (19 species) and no photoactivity as heterotrophs (11   species

| Community phylogenetics analysis (spatial data)
To test for competitive exclusion in the spatial data, we use a community phylogenetics approach and the picante R package (version 1.6-2; Kembel et al., 2010). We test if species coexisting in a local community are more or less ecologically similar (clustered or overdispersed, respectively) than what would be expected given a random assembly process (null model). Overdispersed communities are predicted to emerge because of intense competition among ecologically similar species for limiting resources, whereas clustered communities contain species that share similar traits, hypothetically related to environmental tolerances (Cavender-Bares et al., 2009;Pearse et al., 2014;Webb et al., 2002).
We assessed the ecological similarity of species in the community given a species-by-species distance matrix (cophenetic distance, i.e., divergence time, for the phylogenetic tree and Euclidean distance for the species' average sizes). We assume that closely related species share more traits (niche conservatism) and therefore compete more strongly than distantly related species. Similarly, we also assume species of similar sizes compete more strongly.
For each of the 3,053 PF assemblages in the spatial data (Figure 2a), we calculated two community dispersion metrics: the mean pairwise distance (MPD) and the mean nearest (taxon) distance (MNTD), both calculated given the species-by-species distance matrices, and weighted by the abundance of each species in the local community. MPD is the average distance between all species in the assemblage and thus takes into account clade-wide patterns. MNTD is the average distance between each species and its nearest neighbour in the distance matrix, representing local patterns in the phylogeny or trait space.
To test if the observed MPD and MNTD values differ significantly from the values that would be expected under a random assembly process, we generated null models of community assembly.
If species are equivalent (neutral dynamics), we expect communities to be composed of a random subset of species capable of colonizing the local community (i.e., the source pool; Pearse et al., 2014;Webb et al., 2002). Our source pool includes all extant PF species, because the ranges of all but two species cover every ocean basin (Bé & Tolderlund, 1971). Thus, we assume in the null models that species can reach all oceanographic regions in the world. If there are local community assembly processes (competitive exclusion and environmental tolerances), the species' local establishment will be different from the random subset, because the local fitness of species differ (i.e., species are not equivalent).
We tested two random assembly models: 'taxa', which shuffles the species-by-species distance matrix, and tests if the observed communities show phylogenetic or size structure; 'richness', which shuffles the species' abundances in the sample (holding the richness constant), and tests if community structure could have been generated by randomly drawing species from the source pool (Kembel et al., 2010). We ran each null model 100 times, then standardized  (Locarnini et al., 2013) to assess the relationship between SST and community patterns.

| Time-series analysis (temporal data)
If two species compete in the local community, we expect that increases in the abundance of one species should have a negative impact on the abundance of the other species, generating a negative correlation of abundances through time. To test for such compensatory dynamics in the temporal data, we calculated the correlation between the differentiated time series (i.e., first differences) of co-occurring species pairs (sampled in the same sediment trap). We used the nonparametric Kendall rank correlation, which makes no assumption on the underlying pairwise distribution. In total, 2,303 correlations were calculated and analysed regarding ecological similarity between species pairs (i.e., phylogenetic relatedness, average sizes or symbiotic strategy).
PF species' abundances (fluxes) vary seasonally in response to SST and primary productivity (Jonkers & Kucera, 2015). When correlating these abundance time series, similar seasonal cycles can yield strong positive correlations even if species compete within their seasonal peaks. In order to correct for this bias, we built a null model assuming a seasonal cycle with a 1-year period, but allowing each species to have any number of abundance peaks within the year. This flexibility of species seasonality is consistent with observations that the same species can display unimodal or bimodal yearly abundance patterns depending on the regional oceanographic setting (Jonkers & Kucera, 2015). Thus, for each observed time series (i.e., each species in each sediment trap) we obtained a smoothing spline, which calculates the average abundance of the species for each day of the year, resulting in a yearly averaged series (spline series; Supporting Information Figure S7), without having to enforce an a priori fixed annual seasonality.
The smoothing spline parameter λ was estimated using the generalized cross-validation (GCV) method (Craven & Wahba, 1978), which finds the optimum parameter value by minimizing the mean squared error of the fitting (also used in generalized additive models). The mean error of the spline fitting is larger for time series of shorter length (number of samples) and shorter duration (number of years) (Supporting Information Figure S8), but not in a way that qualitatively affected our results.
To generate the seasonality null model, we calculated the residuals between the observed time series and its spline series (Supporting Information Figure S7). We then randomized these residuals and added them to the spline series. This procedure was repeated 100 times, to produce a set of surrogate time series with the same seasonality as the observed time series, but with random departures from it. For each species pair, we calculated the Kendall rank correlation coefficient between their surrogate series, resulting in a null distribution of 100 correlations for each co-occurring species pair (Supporting Information Figure S7). We report the standardized effect size (SES), which is the observed correlation standardized with respect to the null distribution. We compared the observed correlation with its null distribution for significance, and used a statistical significance level of 1% because of the large number of comparisons.
The residuals between the observed series and the spline series quantify the deviation of the observed abundances from the average annual abundance pattern. Thus, to calculate inter-annual variation in abundances, we correlated the residual series of all co-occurring species pairs within sediments traps that were sampled for longer than 2 years (21 traps in total, Supporting Information Table S2). A positive average correlation of residuals within the sediment trap indicates synchronous inter-annual variation in species abundances.

| Community phylogenetics (spatial data)
Out of the 3,053 PF assemblages studied, only two differ significantly from random assembly processes (Figure 4), that is, a negligible proportion of the studied assemblages. The majority of NRI and NTI values are positive, indicating a trend towards phylogenetic and size clustering in the PF assemblages. As there is no qualitative difference between the results of the 'richness' and 'taxa' null models (Supporting Information Figure S3), we discuss only the former. This single overdispersed assemblage has a remarkably low NTI value (around −4) and high relative abundance of the large species Globorotalia menardii (37%, Supporting Information Figure S6). In general, there is a tendency for negative, but not significant, NRI and NTI values at high SST (25-30 °C, Figure 4b,d), suggesting large size disparity. In fact, these tropical assemblages contain species of large average sizes (Supporting Information Figure S5), which is consistent with the known trend of increasing assemblage size towards the tropics (Schmidt, Renaud, Bollmann, Schiebel, & Thierstein, 2004).

| Time series (temporal data)
Out of the 2,303 standardized correlations, 634 differ significantly from the seasonality null model, of which 588 are positive and 46 negative. As such, even without correction for multiple testing, we find only 2% of the correlations to conform to the expected competi-

| Spatial patterns
Our community phylogenetics analysis of 3,053 assemblages revealed that PF co-occurrence patterns are indistinguishable from neutral assembly processes (Figure 4). The sediment samples used in this analysis are globally distributed (Figure 2a) and, importantly, contain data on relative abundance of species instead of presence and absence. Competitive exclusion over broad spatial extents is known to be a slow process (Yackulic, 2017) and particularly difficult to observe in the highly dynamical marine systems F I G U R E 4 Spatial data: most planktonic Foraminifera assemblages are indistinguishable from randomly assembled communities. Values show the net relatedness index (NRI) and nearest taxon index (NTI) calculated for 3,053 coretop assemblages, plotted against annual mean sea-surface temperature on each site, boxplots binned by 1 °C. (a) and (c) were calculated based on the phylogenetic distance among cooccurring species; (b) and (d) were calculated based on the shell size differences among co-occurring species. Grey dots show non-significant NRI or NTI values, the two red dots show significant negative values (overdispersed). No values were significantly positive (clustered). Significance level 1%, based on 100 runs of the null model 'richness' (Clayton, Dutkiewicz, Jahn, & Follows, 2013). However, by considering species local abundance, competitive exclusion need not have gone to completion. Overdispersed communities can emerge when ecologically similar species still coexist but have disparate abundances (Ulrich & Gotelli, 2010). Yet overdispersion was observed only in two assemblages (Figure 4), indicating that if competitive exclusion is happening in these communities, it leaves no phylogenetic or morphological signal in the seafloor sediments.
We observed no global significant relationship between the NTI and NRI values and SST (Figure 4), although PF show a clear latitudinal diversity gradient (Figure 2b). When considering a significance level of 5% (instead of 1%), NTI and NRI values regarding phylogenetic distance remain the same; however, regarding size differences, many significant clustered communities emerge (Supporting Information Figure S4). These clustered communities occur across the SST range, with no clear SST pattern; instead, they show a relationship with the local relative abundance of the species Globorotalia menardii and Globorotalia tumida. These disproportionately large species (Supporting Information Figure S2) bias the NTI and NRI metrics: when they are absent or in low abundances, communities are clustered (Supporting Information Figure S6). This bias in the PF co-occurrence patterns due to only two species suggests that our community phylogenetics analysis is limited by the low number of extant PF species.
The dichotomy of competition versus environmental filtering, pervasive in community phylogenetic approaches, is oversimplified (Cadotte & Tucker, 2017). Competition can instead reinforce co-occurrence of ecologically similar species, giving rise to clustered communities identical to those expected under environmental filtering (Mayfield & Levine, 2010). However, the fact that we observe no significant clustering or overdispersed patterns suggests that these plankton communities are influenced by neutral processes, such as dispersal and stochastic extinction (Rosindell, Hubbell, & Etienne, 2011). Indeed, dispersal assembly mechanisms have been shown to prevail over niche-assembly mechanisms for the past 30 million years of PF evolution (Allen & Gillooly, 2006), but the lack of support for competition found in our study does not necessarily indicate that PF dynamics are neutral. Other nonneutral mechanisms such as symbiosis, predation and competition with other groups can drive PF community assembly and need to be further investigated.

| Temporal patterns
The time-averaging of seafloor sediment samples (Jonkers et al., 2019) added to the complex spatio-temporal dynamics of the marine environment can mask local patterns that emerge from ecological interactions (Clayton et al., 2013;McGillicuddy, 2016). Therefore, we also investigated high-resolution time-series to test for competition among PF species. Previous sediment trap studies have shown that PF abundances (measured as shell flux) change synchronously through time largely in response to changes in temperature and food availability (Jonkers & Kucera, 2015;Marchant, Hebbeln, Giglio, Coloma, & Gonzalez, 2004;Marchant, Hebbeln, & Wefer, 1998;Northcote & Neil, 2005;Poore, Tedesco, & Spear, 2013;Rigual-Hernandez, Sierro, Barcena, Flores, & Heussner, 2012;Tedesco & Thunell, 2003;Wan, Jian, Cheng, Qiao, & Wang, 2010). However, other studies have shown that species abundance peaks can occur sequentially (Deuser & Ross, 1989;Jonkers et al., 2010;Jonkers, Heuven, Zahn, & Peeters, 2013) or be inversely related (Kincaid et al., 2000;Kuroyanagi, Kawahata, Nishi, & Honda, 2002;Mohiuddin, Nishimura, Tanaka, & Shimamoto, 2002;Rigual-Hernandez et al., 2012), suggesting that species have different ecological preferences. Species with opposite abundance patterns during the year are expected to exhibit low correlation values, whereas species with F I G U R E 5 Temporal data: most sediment traps show non-significant or positive correlations between planktonic Foraminifera population dynamics. Values show pairwise time-series correlations of species abundances (2,303 in total) plotted by sediment trap (see Supporting Information Table S2 for sediment trap name abbreviation). Grey dots show non-significant, blue significant positive, red significant negative standardized size effect (SES) values, based on the seasonality null model and significance level of 1% sequential abundance peaks (compensatory dynamics) should correlate negatively through time. Here, we show that the dynamics of PF abundances are largely explained by annual seasonality: 72% of the pairwise correlations did not differ significantly from our seasonality null model. The remaining significant correlations were predominantly positive (93%), indicating synchronous population dynamics.
This pattern was consistent across sediment traps with different sampling durations ( Figure 5 and Supporting Information Table S2) and across years (Supporting Information Figure S10). Moreover, species with similar ecologies did not correlate more negatively than other species pairs (Figure 6), as would be expected if they were partitioning the resource through time. Therefore, PF population dynamics seems to be influenced by external factors such as the abiotic environment and/or abundance of other distantly related plankton groups, rather than intra-clade density-dependent processes.
PF species might have different vertical niches and thus live in parapatry in the water column (Hemleben et al., 1989;Rebotim et al., 2017), yet they would still fall within the same sediment trap or coretop sample, and thus be observed as co-occurring species in our analyses. Species' vertical niches have been shown to cover a wide range of depths and to overlap considerably (Rebotim et al., 2017), making it difficult to define consistent depth habitats of species for all samples analysed here. Nevertheless, symbiont-bearing species are limited to the surface euphotic zone for the optimum maintenance of their photosynthetic algae (Bé, Spero, & Anderson, 1982;Takagi et al., 2019).
Assuming that photosymbiosis is a reasonable proxy for depth habitat, symbiont-bearing species did not show stronger compensatory dynamics than species with different symbiotic strategies (Figure 6c), supporting our results at an admittedly coarse depth habitat classification.

| Mismatch between macroevolutionary and community dynamics
The fossil record suggests that the diversification dynamics of PF are strongly regulated by competition among species Ezard & Purvis, 2016). However, we found no evidence of interspecific competition shaping spatial or temporal diversity patterns of extant PF. The global and clade-wide scales of our ecological analysis are comparable in scope to macroevolutionary analysis, yet these analyses differ markedly in temporal and taxonomic scales (Harmon et al. 2019). The PF fossil record showed evidence for competition when analysed at the species level (taxon counts) over millions of years (Ezard & Purvis, 2016). Here, we tested for competition among species at the population level and at shorter time scales, ranging from weeks up to a few thousands of years. While PF population dynamics could potentially show different emergent patterns when analysed over million-year time spans, interspecific competition operates at the individual level and requires population overlap in space and time, also as a diversity-dependent diversification process (Marshall & Quental, 2016). Moreover, the nested hierarchy of organizational levels (individuals, populations, species) implies F I G U R E 6 Temporal data: most planktonic Foraminifera species display seasonal and synchronized population dynamics. Values show pairwise time-series correlations of species abundances plotted against (a) the phylogenetic distance between species in millions of years (My) (1,843 pairwise correlations in total), boxplots represent divergence times with more than 40 species pairs; (b) the shell size difference between species (1,801 pairwise correlations), binned by 50 μm; (c) the photosymbiotic strategy of the species pair (1,721 pairwise correlations): both species do not host symbionts (heterotrophs), one species does not host symbionts and the other does (heterotroph & phototroph), or both species host symbionts (phototrophs). Grey dots show non-significant, blue significant positive, red significant negative standardized size effect (SES) values, based on the null seasonality model and significance level of 1% that dynamics at higher levels should always propagate downwards (Jablonski, 2008b). Thus, if a clade's diversification is limited by competition among species, this mechanism should be observable at the population level and ecological time scales.
It could be that interspecific competition limited PF diversification throughout the Cenozoic (Ezard & Purvis, 2016), but species are not competing at present because they are below the clade's carrying capacity (Rabosky, 2013). In fact, a total of seven PF species went extinct in the past 5 million years, and the clade is currently in a negative diversification phase (Figure 1; Aze et al., 2011). If the niche space was saturated, then the extinction of these species would have released the competitive pressure by making niche space available.
Nevertheless, we would still expect extant species to show niche partitioning and evidence for interspecific competition across space and/or ecological time, which we did not observe in our competition proxies. The macroevolutionary compensatory dynamics proposed by Ezard and Purvis (2016) leads to the hypothesis that empty niche space is filled by species occupying similar niches (incumbency advantages). Thus, we could expect that species that are ecologically similar to the recently extinct species expanded their niche by increasing in abundance and/or biogeographical range, compensating for the extinction of their competitors. This hypothesis remains to be tested in the PF fossil record.
Alternative mechanisms could generate the DDD pattern, without the need to evoke interspecific competition and niche partitioning. In a model where speciation and extinction rates emerge from population dynamics, Wang, Chen, Fang, and Pacala (2013) showed that the average speciation rate of a clade declines over time because species' abundances decrease at speciation. In this model, lineages with high speciation rates have the disadvantage of having high extinctions because of low population size (see also Harnik, 2011). Thus, a positive correlation between speciation and extinction rates naturally emerges (the macroevolutionary trade-off; Jablonski, 2017), which is also observed in the macroperforate PF fossil record (Supporting Information Figure S11).
The geography of diversification might also play a crucial role in the emergence of DDD patterns (Moen & Morlon, 2014). During allopatric speciation, cladogenesis not only reduces species' abundances but also species' biogeographical ranges. Larger ranges are more likely to experience the emergence of geographical or oceanographical barriers. Thus, as diversification proceeds, ranges are subdivided and per-species rates of speciation decline, generating the DDD pattern without competition among species (Moen & Morlon, 2014;Pigot, Phillimore, Owens, & Orme, 2010). PF have the ability to disperse long distances, evidenced by genetic and fossil data (Darling et al., 2000;Sexton & Norris, 2008). This high dispersal suggests that it may be difficult to geographically isolate pelagic populations for extended periods of time, a key component in allopatric speciation models. In fact, the PF fossil record has shown cases of sympatric speciation (Lazarus, Hilbrecht, Spencer-Cervato, & Thierstein, 1995;Pearson, Shackleton, & Hall, 1997), questioning the use of purely allopatric speciation models (e.g., Pigot & Etienne, 2015). Nevertheless, PF genetic studies revealed phylogeographical patterns of allopatry (Darling, Kucera, & Wade, 2007;De Vargas, Norris, Zaninetti, Gibb, & Pawlowski, 1999;Weiner et al., 2014), suggesting that ocean circulation and directional flow create dynamical barriers to populations. Moreover, species' abundances and geographical ranges might interact in a negative way, and rare but widespread species (e.g., globorotaliids) might experience oceanographical gene-flow barriers more often than abundant species. Genetic data provide new perspectives for speciation modes in the plankton, and could reveal PF ecological interactions at other organizational levels (Alizon, Kucera, & Jansen, 2008;Aurahs, Treis, Darling, & Kucera, 2011;Darling, Kucera, Wade, von Langen & Pak, 2003).
Overall, our results suggest that PF species are not competing in the modern oceans, which is consistent with the known low densities of PF in the water column (Keeling & del Campo, 2017;Schiebel & Hemleben, 2005). Abiotic factors such as SST have been shown to affect PF diversity patterns in space and time (Fenton et al., 2016;Jonkers & Kucera, 2015;Yasuhara, Hunt, Dowsett, Robinson, & Stoll, 2012), as well as their diversification dynamics Ezard & Purvis, 2016). However, phytoplankton abundance also plays an important role in PF population dynamics (Jonkers & Kucera, 2015), but their effect on PF diversification has not been explored. To further understand the processes governing biodiversity dynamics, we need to include species that are coupled ecologically, even when belonging to distantly related clades (Marshall & Quental, 2016). By defining the species pool of macroevolutionary dynamics ecologically rather than only phylogenetically, we can potentially elucidate the mismatch found in our study.

ACK N OWLED G M ENTS
We thank all the researchers who have made their data openly avail-

DATA AVA I L A B I L I T Y S TAT E M E N T
The data that support the findings of this study are openly available and previously published: refer to Siccha and Kucera (2017) for coretop abundance data; Supporting Information Table S2 for