Ecological observations and paleontological data show that communities of organisms recur in space and time. Various observations suggest that communities largely disappear in extinction events and appear during radiations. This hypothesis, however, has not been tested on a large scale due to a lack of methods for analyzing fossil data, identifying communities, and quantifying their turnover. We demonstrate an approach for quantifying turnover of communities over the Phanerozoic Eon. Using network analysis of fossil occurrence data, we provide the first estimates of appearance and disappearance rates for marine animal paleocommunities in the 100 stages of the Phanerozoic record. Our analysis of 124,605 fossil collections (representing 25,749 living and extinct marine animal genera) shows that paleocommunity disappearance and appearance rates are generally highest in mass extinctions and recovery intervals, respectively, with rates three times greater than background levels. Although taxonomic change is, in general, a fair predictor of ecologic reorganization, the variance is high, and ecologic and taxonomic changes were episodically decoupled at times in the past. Extinction rate, therefore, is an imperfect proxy for ecologic change. The paleocommunity turnover rates suggest that efforts to assess the ecological consequences of the present-day biodiversity crisis should focus on the selectivity of extinctions and changes in the prevalence of biological interactions.

In modern ecosystems, species associate in recurrent communities that reflect biological interactions, overlapping environmental tolerances, or both (Liautaud et al., 2019). The fate of a community hinges on the persistence of some fraction of its species, but the strength of the relationship between extinction of taxa and disappearance of communities remains unclear, even though such knowledge may help in predicting the ecological consequences of taxonomic losses in the present day. The fossil record indicates that extinction (Bambach et al., 2004), morphological change (Brett et al., 1996), and ecological reorganization (Sheehan, 1996) tend to happen in discrete episodes. These observations imply that communities largely appear and disappear during relatively abrupt events and that their turnover rates have varied over the past 541 m.y. To date, however, no studies have directly estimated these rates. We fill this gap by conducting a network analysis of the entire Phanerozoic fossil record of marine animals, identifying paleocommunities, and quantifying their turnover through stratigraphy.

A paleocommunity is a unique group of fossil taxa that occur together at multiple locations arrayed, regionally or globally, in geographic and stratigraphic space. Paleontologists recognize such groups using a variety of methods (Boucot and Lawson, 1999). In quantitative analysis, a paleocommunity represents an aggregate of mathematically indistinguishable “local paleocommunities:” fossil assemblages at specific geospatial locations (Bambach and Bennington, 1996). Ideally, a local paleocommunity corresponds to all of the species observed in one bed at a single outcrop, but because assigning fossils to species can be challenging, and boundaries among beds and outcrops are unclear in some cases, paleocommunity analyses commonly focus on broadly defined samples of fossils classified by genus. Researchers have traditionally identified paleocommunities using multivariate clustering and ordination (Shi, 1993); however, the computational run times of these methods generally limit their application to data sets with no more than a few thousand samples. We instead use network-based methods (Muscente et al., 2018, 2019), which include a variety of “community detection algorithms” for identifying aggregates in large data sets (Yang et al., 2016).

Using the Paleobiology Database (https://paleobiodb.org; Alroy et al., 2001; Peters and McClennen, 2015), we compiled a global data set of collections of marine animal fossils of Phanerozoic age (Data Set S1, Fig. S1, and Tables S1 and S2 in the Supplemental Material1). Each collection in our data set represents one or more genera at a location (i.e., a local paleocommunity). We used genera instead of species, given that related species commonly resemble each other in terms of functional ecology due to phylogenetic niche conservatism (Wiens et al., 2010). In addition, data on genera are typically more accurate and complete than species data because species are more numerous and difficult to sample. Simulations suggest that when sampling is poor, genera may best represent evolutionary history (Sepkoski and Kendrick, 1993).

We structured the collections into a unipartite network where each node (n = 124,605) is a collection (Fig. S2) and each link (n = 47,294,900) connects two collections that share one or more genera (Fig. S3). The strength (weight) of each node equals its number of genera, and the weight of each connection equals its collections’ similarity, which we calculated using metrics (e.g., the Jaccard index) based on presence or absence of taxa (Shi, 1993). We identified modules, or clusters of nodes, in our network by applying community-detection algorithms. Because these modules represent aggregates of local paleocommunities, we treat them as paleocommunities. To find the best procedure, we repeated the analysis for 20 combinations of 5 algorithms and 4 similarity indices and evaluated the results based on the numbers, sizes, and modularity (Q) scores of the modules. The best result has the highest Q score, which is the fraction of links that connect nodes of the same module minus the corresponding fraction expected in an equivalent network with a random distribution of connections (Clauset et al., 2004). We assessed stochastic uncertainty by repeating the best procedure 200 times and ascertaining confidence intervals. To confirm that the modules represent paleocommunities, we examined the network for homophily, or the tendency of links to connect nodes with similar properties, by calculating “assortativity coefficients” for various variables (e.g., age, location, environment, and lithology). These coefficients range from 0 to 1, with homophily having high values.

To quantify the number and turnover of genera and paleocommunities in each Phanerozoic stage, we adapted common metrics of taxonomic change. We treat paleocommunities like taxa, using the oldest and youngest collections in the aggregates as analogs for their first and last appearance datums. Community diversity equals the corrected sampled-in-bin (CSIB) number (Alroy, 2008) after correction for unequal sampling, and appearance and disappearance rates are equivalent to “second-for-third proportions” of “origination” and “extinction,” respectively; simulations indicate these proportions are more accurate and precise than alternative metrics (Alroy, 2015). Prior to calculating metrics, we adjusted for unequal sampling by subsampling the data with the shareholder quorum method (Alroy, 2010).

Network Analysis

Assortativity coefficients (analogous to the Pearson correlation coefficient) indicate that absolute age (coefficient of 0.99) and chronostratigraphic system (0.60) best predict whether two nodes are connected; all other variables (Table S3) are associated with low coefficients (<0.25). We applied five algorithms and four similarity indices to the data, and consistently identified the most modules and those with the highest Q scores (Tables S4 and S5) with the Infomap algorithm (Rosvall and Bergstrom, 2008) using connection weights equal to Jaccard similarity (Shi, 1993). We completed 200 runs of this non-deterministic algorithm, but results (i.e., modules) typically varied by only 1% of information (Table S6), returning 3937 modules on average. The properties of these modules confirm that they represent paleocommunities with wide stratigraphic but limited geographic and environmental ranges (Supplemental Material text; Figs. S3–S5; Table S7).

Paleocommunity Turnover Rates

The number of paleocommunities rises and falls across the record, with peaks in the Ordovician, Permian, Cretaceous, and Paleogene (Fig. 1A; Fig. S6). Their CSIB number broadly covaries with the total CSIB number of genera, or gamma diversity (Whittaker, 1960); however, the number of paleocommunities significantly departs from gamma diversity in the middle Paleozoic (Silurian and Devonian) and Neogene. These departures may be related to changes in alpha diversity or beta diversity over time—the average number of genera per community (Fig. 1B) and the degree of community differentiation (Fig. 1C), respectively (Whittaker, 1960). Our work demonstrates that alpha diversity is highest in the Cenozoic, as generally thought (Bambach, 1977), but also shows that comparable levels exist in the mid-Paleozoic. The number of paleocommunities per genus is low in both intervals (Fig. 1D).

Our best estimates indicate that the appearance and disappearance rates of paleocommunities are both, on average, 20% per stage. However, paleocommunities disappear at nearly three times this rate (55%) in the “big five” extinction intervals, and rates of appearance generally rise to 42% in post-extinction recovery intervals. These differences are statistically significant, and the stochastic error of these estimates is low relative to the binomial uncertainty. The patterns closely resemble trends in the extinction and origination of animal taxa (Fig. 2; Figs. S6–S12). Indeed, reduced major axis regressions (Fig. 3; Figs. S11–S12) reveal positive correlations between genus extinction and paleocommunity disappearance (correlation coefficient, 0.83; p-value <0.001) and between genus origination and paleocommunity appearance (correlation coefficient, 0.86; p-value <0.001). Sensitivity analyses show that these results are robust (Figs. S13–S14; Table S8). Intervals of great taxonomic change are generally associated with high levels of paleocommunity reorganization and vice versa.

Our study provides the first estimates of paleocommunity appearance and disappearance rates across the Phanerozoic record (Fig. 2). These rates provide empirical evidence of ecological evolutionary units (EEUs), subdivisions of the record representing long (30–140 m.y.) intervals of relative stasis in communities at the genus level (Boucot, 1983). Boundaries between EEUs correspond to mass extinctions and major radiations, when interrelated communities disappeared and were replaced by new ones (Sheehan, 1996). During these transitions, paleocommunity turnover rates were roughly three times greater than normal. Although some community changes developed gradually over time, most reorganization happened during extinction and radiation events.

The correlation between community and genus turnovers indicates that ecologic and taxonomic change are not completely decoupled. For example, like taxonomic turnover (Bambach et al., 2004), community turnover was high from the Cambrian to the mid-Ordovician (Figs. 1 and 2), perhaps as a consequence of environmental dynamics (Gill et al., 2011; Saltzman et al., 2015). Taxonomic change, thus, provides a fair proxy for ecologic reorganization. Nonetheless, the variance is high (Fig. 3). The Capitanian (P7 in Fig. 3), Givetian (D5), and Induan (Tr1) biodiversity crises (McGhee et al., 2013) were times of significantly higher rates of taxonomic than community change (Fig. 3; Figs. S7, S8, S11). Conversely, the Chibanian (Q3)—the age that ended with the last interglacial period prior the current one (Holocene)—was not a time of great taxonomic extinction, but according to our results, it did see substantial community losses. Thus, our results corroborate studies (Droser et al., 2000; McGhee et al., 2013; Muscente et al., 2018) suggesting that there is limited (but not absolute) decoupling of taxonomic and ecologic turnover during some events.

According to our estimates, the “big five” mass extinctions (Bambach et al., 2004) may have markedly different rank orders in terms of community disappearance versus genus extinction (Table 1). The end-Permian mass extinction stands out as the highest-magnitude event. The end-Ordovician and end-Triassic events saw high losses of taxa but limited community change, whereas the end-Devonian and end-Cretaceous crises witnessed the demise of many communities but relatively few taxa. Whereas the end-Ordovician and end-Triassic extinctions involved declines in both alpha and beta diversity (Figs. 1B and 1C), the end-Permian and end-Cretaceous events primarily entailed losses in beta diversity; the Devonian event saw the largest drop in alpha diversity of the Phanerozoic Eon. Although there is uncertainty in these estimates, they suggest that marine ecologic structure has primarily been shaped not by the magnitudes of mass extinctions, but rather by the selectivity of such events (Clapham and Payne, 2011; Jablonski and Raup, 1995; Kitchell et al., 1986). A small number of extinctions can drive paleocommunity turnover on a large scale (Roopnarine, 2006; Vermeij, 2004).

We hypothesize that intervals of high taxonomic and low community diversity, like the middle Paleozoic and Cenozoic (Fig. 1A), represent times of ecospace-occupation expansion, when many niches were filled by specialist taxa. During the Devonian “mass depletion,” which involved a severe and protracted reduction in taxonomic origination (Bambach et al., 2004), marine transgressions may have limited speciation by vicariance and allowed invasive taxa to expand their geographic and environmental ranges (Stigall, 2012). Consistent with this hypothesis, the crisis involved a steep reduction in alpha diversity but only minor changes in beta diversity and the number of communities per genus (Figs. 1B–1D); the pre-crisis communities, which contained many specialist taxa with narrow ranges, were replaced with post-crisis communities dominated by generalist taxa with wide ranges (Stigall, 2012). The drop in alpha diversity, which continued through the late Paleozoic, may have also been related to extinctions that reduced epifaunal tiering (Ausich and Bottjer, 1985), thereby eliminating niches and limiting local diversity.

The cause of an event does not represent a good predictor of turnover. Ocean warming, acidification, and deoxygenation drove a number of biological crises over Earth history (Clapham and Renne, 2019). Some events involved more community change than taxonomic loss (e.g., end-Permian), while others saw the reverse (e.g., Toarcian and Capitanian). Although the magnitude of an event is related to its scale (e.g., regional or global), its impact on ecologic structure hinges on the prevalence of vulnerable taxa and the specificity of their interactions, which are functions of geography (Stigall, 2012) and evolutionary history (Foster et al., 2020). In this context, the paleocommunity record suggests that predictions about the current biodiversity crisis should not focus solely on numbers of extinctions, but instead on the selectivity of those extinctions and their impacts on biological interactions. Our understanding of the tempo and mode of ecologic change will ultimately improve with the amount and quality of data on local paleocommunities and their functional traits.

1Supplemental Material. Additional information on methods and results including supplemental figures and tables. Please visit https://doi.org/10.1130/GEOL.S.16985347 to access the supplemental material, and contact [email protected] with any questions.

This paper is dedicated to our late friend and colleague, Peter Fox. We thank Gene Hunt, John Huntley, Peter Wagner, and an anonymous reviewer for helpful comments. This work was supported by grants from the W.M. Keck Foundation (Los Angeles), the Deep Carbon Observatory, Alfred P. Sloan Foundation (New York), the Carnegie Institution for Science (Washington, D.C.), and the U.S. National Science Foundation (NSF grants EAR-1660005 and OAC-1835717). This is Paleobiology Database official publication number 413.