Introduction
Large-scale taxonomic turnover between modern and fossil faunal assemblages can make it difficult to recognize the ecological similarities and differences between otherwise comparable ecosystems. Ecologists could leverage the fossil record to address questions such as how communities change on long timescales in response to environmental changes, if they could accurately identify paleoecological analogues to modern systems (Williams and Jackson Reference Williams and Jackson2007; Rull Reference Rull2010; Kidwell Reference Kidwell2013; Seddon et al. Reference Seddon, Mackay, Baker, Birks, Breman, Buck and Ellis2014). These fossil–modern analogues can also help us overcome the preservational biases of the fossil record to better reconstruct past ecosystems (Overpeck et al. Reference Overpeck, Webb and Prentice1985; Williams and Jackson Reference Williams and Jackson2007). To accomplish both these goals, we must compare the ecological structure and composition of assemblages independent of taxonomic identity. Paleoecologists use guild-based analyses and theoretical ecospaces to reveal ecological changes through the Phanerozoic and across extinction events (Bambach Reference Bambach, Tevesz and McCall1983; Droser et al. Reference Droser, Bottjer, Sheehan and MeGhee2000; Bush et al. Reference Bush, Bambach and Daley2007; Novack-Gottshall Reference Novack-Gottshall2007; Bush and Novack-Gottshall Reference Bush and Novack-Gottshall2012; Foster and Twitchett Reference Foster and Twitchett2014); however, this method is rarely used to examine community-scale ecology between taxonomically disparate fossil assemblages. Here we use a guild-based ecospace to compare the functional ecology of modern cold-seep communities and Cretaceous Western Interior Seaway (WIS) cold-seep assemblages.
Cold-seep communities are characterized by invertebrates dependent upon symbiotic chemoautotrophs that use methane or hydrogen sulfide from the cold-seeps as their energy source. Recent studies have shown the importance of modern hydrocarbon seep communities to the surrounding ecosystems and to global biogeochemical cycles (Levin Reference Levin2005; Boetius and Wenzhöfer Reference Boetius and Wenzhöfer2013; Levin et al. Reference Levin, Baco, Bowden, Colaco, Cordes, Cunha and Demopoulos2016). Seeps on continental margins contribute to local fisheries’ productivity, and many predators acquire prey from both chemoautotrophic and non-chemoautotrophic environments (Levin et al. Reference Levin, Baco, Bowden, Colaco, Cordes, Cunha and Demopoulos2016). In addition, seep communities add habitat complexity to the seafloor through archeal bacterial–mediated precipitation of carbonate that remains after fluid flow ceases and through the deposition of organism hard parts such as shells after death (Levin et al. Reference Levin, Baco, Bowden, Colaco, Cordes, Cunha and Demopoulos2016). Thus, seep communities have lasting impacts on ecosystems even after the seeps are no longer active, and it is critical that we understand the long-term stability of these communities though deep time to preserve their ecosystem services in the modern.
Cretaceous fossil cold-seep faunas are often viewed as analogues for modern cold-seep communities despite significant differences in taxonomic composition and in environmental conditions, including water depth and sulfide concentrations (Beauchamp and Savard Reference Beauchamp and Savard1992; Powell et al. Reference Powell, Callender and Stanton1998; Callender and Powell Reference Callender and Powell1999; Kiel Reference Kiel2015; Meehan and Landman Reference Meehan and Landman2016). Some of the major clades of chemosymbiotic invertebrates that occupy cold-seeps today likely did not originate until the early Cenozoic, ecologically important Cretaceous taxa such as ammonites have gone extinct, and major changes in the composition of cold-seep communities occurred after the Cretaceous as seawater sulfide availability increased (Vrijenhoek Reference Vrijenhoek2013; Kiel Reference Kiel2015). Thus, there may be significant differences in ecological function between cold-seep communities in the Cretaceous and the modern, despite broad similarities between fossil and modern cold-seeps faunas, which must be explored before drawing a precise ecological analogue across deep time (Campbell Reference Campbell2006).
Modern seep communities are diverse, and assemblages within a region vary with differences in hydrocarbon fluid flow, seep age, water depth, and sedimentary characteristics (Sahling et al. Reference Sahling, Galkin, Salyuk, Greinert, Hilmar, Piepenburg and Suess2003; Campbell Reference Campbell2006; Cordes et al. Reference Cordes, Hourdez, Predmore, Redding and Fisher2005, Reference Cordes, Bergquist, Predmore, Jones, Deines, Telesnicki and Fisher2006, Reference Cordes, Bergquist and Fisher2009, Reference Cordes, Becker, Hourdez and Fisher2010; Dando Reference Dando and Kiel2010). The primary drivers of faunal variation in modern cold-seep communities are the rate of methane release, sulfide concentration, and oxygen availability (Bergquist et al. Reference Bergquist, Fleckenstein, Knisel, Begley, MacDonald and Fisher2005; Cordes et al. Reference Cordes, Hourdez, Predmore, Redding and Fisher2005, Reference Cordes, Bergquist, Predmore, Jones, Deines, Telesnicki and Fisher2006, Reference Cordes, Bergquist and Fisher2009). Because these properties change directionally over the lifetime of a seep, changes in the associated faunas are often described as successional stages (Cordes et al. Reference Cordes, Hourdez, Predmore, Redding and Fisher2005, Reference Cordes, Bergquist, Predmore, Jones, Deines, Telesnicki and Fisher2006, Reference Cordes, Bergquist and Fisher2009). In the early stages of seep development, when soft substrate is most common, communities are dominated by chemosymbiotic mussels that require minimal hard substrate, such as Bathymodiolus childressi. Vestimentiferan tubeworms colonize the site as carbonate precipitation establishes sufficient hard substrate, typically leading to more diverse communities than the mussel beds. As the seep further ages, decreased methane and sulfide concentrations and increased oxygenation permit colonization by organisms less tolerant of seep fluids (Cordes et al. Reference Cordes, Bergquist and Fisher2009).
Similar stages of seep development have been proposed for Cretaceous seeps in the WIS and Japan, as well as for Jurassic to Pliocene cold-seeps from western North America. In these studies, early stages of seep development are characterized by infaunal chemosymbiotic bivalves that used the high concentration of hydrogen sulfide in the local environment (Beauchamp and Savard Reference Beauchamp and Savard1992; Campbell and Bottjer Reference Campbell and Bottjer1993; Campbell et al. Reference Campbell, Farmer and Des Marais2002; Hikida et al. Reference Hikida, Suzuki, Togo and Ijiri2003). In later stages of seep development, lithification of the sediment hindered the settlement of infaunal seep taxa and facilitated colonization by epifaunal taxa such as serpulid tubeworms. In the Pierre Shale of the WIS, baculitid ammonites, inoceramids, and lucinids, which also commonly occur outside seep settings, dominate seep assemblages, but the seep assemblages vary in the relative abundance and diversity of secondary organisms that attach to exposed potions of the carbonate cap, reflective of changing substrate conditions at the seep over its lifetime (Meehan and Landman Reference Meehan and Landman2016).
Despite these potential parallels between the modern and fossil cold-seep systems, quantitative comparisons necessary to determine the degree of analogy in functional ecology, ecological structure, and patterns of community succession through geologic time are limited. Here we quantitatively compare 45 modern faunal assemblages representing five cold-seep biofacies in the Gulf of Mexico (GoM) and one cold-seep biofacies from the North Sea to 53 fossil assemblages collected in the Cretaceous Pierre Shale in the Black Hills region of South Dakota. We quantitatively compare the composition and structure of the fossil cold-seep assemblages to those of the modern biofacies using ecological guilds (after removing taxa from the modern faunas that are not readily fossilizable) and reconstruct parallel patterns of succession between these modern and fossil systems.
Seep Data and Collection Methods
Cretaceous Black Hills Cold-Seep Sites
Cold-seeps in the geologic record are recognized by deposits of authigenic carbonate with low δ13C values, microbial textures, and dense macrofossil accumulations (Beauchamp and Savard Reference Beauchamp and Savard1992; Kauffman et al. Reference Kauffman, Arthur, Howe and Scholle1996; Conti and Fontana Reference Conti and Fontana1999; Shapiro and Fricke Reference Shapiro and Fricke2002; Campbell Reference Campbell2006; Hammer et al. Reference Hammer, Nakrem, Little, Hryniewicz, Sandy, Hurum, Drunkenmiller, Knutsen and Høyberget2011; Larson et al. Reference Larson, Brezina, Landman, Garb and Handle2013). In the Cretaceous WIS Pierre Shale, the carbonate deposits are commonly expressed as a conical erosional feature in the surrounding shale and are referred to as “tepee buttes” (Kauffman et al. Reference Kauffman, Arthur, Howe and Scholle1996; Gilbert and Gulliver Reference Gilbert and Gulliver1894; Metz Reference Metz2010; Landman et al. Reference Landman, Cochran, Larson, Brezina, Garb and Harries2012b). These cold-seep systems formed along the Front Range of the present-day Rocky Mountains and around the flanks of the Black Hills during the Campanian (83.5–70.6 Ma; Fig. 1A; Cobban et al. Reference Cobban, Walaszczyk, Obradovich and McKinney2006; Metz Reference Metz2010) and may be associated with faulting within the Niobrara Formation or lower Pierre Shale during the early Larmide Orogeny (Kauffman et al. Reference Kauffman, Arthur, Howe and Scholle1996; Metz Reference Metz2010; Larson et al. Reference Larson, Brezina, Landman, Garb and Handle2013).
Distribution of hydrocarbon seep deposits near the Black Hills. A, The three geographic areas where tepee buttes outcrop in the Pierre Shale near the Black Hills (adapted from Metz Reference Metz2010; © 2010 by the University of Chicago). Ellipse shows the area encompassed by American Museum of Natural History (AMNH) samples included herein (Meehan and Landman Reference Meehan and Landman2016). B, South Dakota School of Mines and Technology Museum of Geology (MoG) collection locations. RR, Richards Ranch locality. Modified image ArcGIS® software © 2018 by ESRI. C, Sample locations at RR. Modified image ArcGIS® software © 2018 by ESRI.
Cold-seep faunal assemblages from the WIS are consistent with water depths of 30 to 100 m (Callender and Powell Reference Callender and Powell1999; Landman et al. Reference Landman, Cobban and Larson2012a), shallower than the maximum 300 m depth of the seaway (Wright Reference Wright1987; Jewell Reference Jewell1993). The carbonate core of a tepee butte is often concentrically zoned into vuggy carbonate, lucinid (Nymphalucina) bivalve-dominated coquina, and thrombolytic microbialite facies (Fig. 2; Kauffman et al. Reference Kauffman, Arthur, Howe and Scholle1996; Shapiro and Fricke Reference Shapiro and Fricke2002).
Zonal distribution of carbonate facies in a typical tepee butte. VC, vent carbonate; NLC, Nymphalucina coquina; and TM, thrombolitic microbialite. The infaunal chemosymbiotic bivalve Nymphalucina occidentalis dominates the Nymphalucina coquina facies (after Shapiro and Fricke Reference Shapiro and Fricke2002).
Original Field Collections at Richards Ranch (RR)
In situ abundance counts of taxa were collected from 14 Pierre Shale cold-seep deposits at RR northeast of Newell in Butte County, South Dakota (Fig. 1B,C, Table 1). To standardize sampling, a 0.5 m quadrat was placed once at each carbonate facies at each seep, and all individuals within the quadrat were identified to genus. This counting method was repeated until 300 individuals or the entire exposed surface was counted. We defined individual bivalves as either an articulated individual or the mold of a weathered-out individual; disarticulated valves were counted as 0.5 individuals to avoid double counting. We count individual ammonites as isolated shell fragments of at least 0.5 cm diameter. Individual gastropods are isolated internal molds, fragments, and complete shells of at least 0.5 cm diameter. Most genera are monospecific in these assemblages, thus genus-level identifications were sufficient for ecological comparisons. Two species of Hoploscaphites were present, however, H. gilli and H. gilberti have similar mouth parts and thus likely had a similar trophic ecology (Landman et al. Reference Landman, Kennedy, Cobban, Larson and Jorgensen2013). We used species-level identifications of Baculites and Didymoceras specimens to place assemblages within ammonite biostratigraphy zones where possible (Table 1). Representative specimens and specimens of interest are cataloged in the South Dakota School of Mines and Technology Museum of Geology (Supplementary Table 1).
Campanian ammonite zones of the Western Interior Seaway (WIS) represented by Black Hills fossil data sets herein following Cobban et al. (Reference Cobban, Walaszczyk, Obradovich and McKinney2006). Some associated species found in tepee buttes near the Black Hills during this time range are included in their appropriate ammonite zone(s) (Kauffman et al. Reference Kauffman, Sageman, Kirkland, Elder, Caldwell and Kauffman1993; Larson et al. Reference Larson, Jorgensen, Farrar and Larson1997; Landman et al. Reference Landman, Kennedy, Cobban, Larson and Jorgensen2013). Samples MoG 4–6, 10, and 16 did not contain taxa necessary for assignment to a specific ammonite zone. American Museum of Natural History (AMNH) sites are from Meehan and Landman (Reference Meehan and Landman2016). MoG, South Dakota School of Mines and Technology Museum of Geology; RR, Richards Ranch.
Museum Collections from the Museum of Geology (MoG)
Fossil assemblages from tepee buttes associated with the Gale A. Bishop Collection, held at the South Dakota School of Mines and Technology Museum of Geology (MoG), are from Fall River and Butte Counties in South Dakota and from Weston and Niobrara Counties in Wyoming (Fig. 1B, Supplementary Table 1). Samples were gathered through surface collection or breaking apart pieces of the exposed carbonate (Bishop and Williams Reference Bishop and Williams2000). A fossil assemblage from this collection was included in the analyses herein if the original collectors restricted sampling to a single tepee butte and retained individuals from all encountered taxa. The MoG collections yielded 16 cold-seep deposits from seven localities and two ammonite zones (Fig. 1B, Tables 1 and 2). At least one locality is present in each of the three exposed tepee butte regions flanking the Black Hills and spans two ammonite zones (Fig. 1A,B). We identified fossil specimens to the lowest taxonomic rank possible, but as in the new field collections, most genera were monospecific. The collection methods for the MoG samples did not ensure all individuals were retained, and thus we excluded these data from our analyses that require relative abundances.
Summary of each cold-seep faunal data set.
Published Data Sets from the American Museum of Natural History (AMNH)
Published faunal counts from 23 methane cold-seeps from two of the three Black Hills tepee butte regions across Custer, Pennington, and Butte Counties in South Dakota are also included (Fig. 1A, Tables 1 and 2, Supplementary Table 2); however, sampling methods differ significantly from our field collections (Meehan and Landman Reference Meehan and Landman2016). This data set includes counts based on “extensive sampling” of seeps in outcrop; some seeps were also exposed in cross section (Meehan and Landman Reference Meehan and Landman2016). Fossil collections were made either 1 m from the main conduit or across the carbonate cap within <1 m of the carbonate edge and restricted to 0.5 m2 of collected material (Handle Reference Handle2014). Disarticulated bivalves were counted as 0.5 individuals, and only whole or nearly whole ammonites were counted as individuals. Species-level identifications were provided for the majority of specimens. None of the AMNH samples are from the RR site; we verified this by comparing GPS coordinates published in Meehan and Landman (Reference Meehan and Landman2016).
Modern Cold Seeps
Extensive work in the GoM has made its cold-seep faunas some of the best described (Cordes et al. Reference Cordes, Bergquist and Fisher2009) and thus represent an excellent candidate for modern-fossil comparisons. In addition, successional changes in species composition with seep age are well quantified in modern GoM seeps (Cordes et al. Reference Cordes, Bergquist, Predmore, Jones, Deines, Telesnicki and Fisher2006, Reference Cordes, Bergquist and Fisher2009), allowing us to test the hypothesis that similar successional patterns are reflected in fossil seeps through quantitative comparisons. Previous workers report five different faunal assemblage types, or biofacies, from GoM cold-seeps: lucinid bivalve, thyasirid bivalve, vesicomyid bivalve, mussel-bed (dominated by bathymodiolins), and vestimentiferan tubeworm biofacies, each with different environmental associations (Callender Reference Callender1992; Callender and Powell Reference Callender and Powell1997; Cordes et al. Reference Cordes, Bergquist and Fisher2009). We compiled faunal assemblage data from each of these biofacies on the modern Louisiana slope from literature sources (Table 2, Supplementary Table 3). This modern cold-seep data set includes live-collected assemblages and faunal counts based upon the remains of dead organisms (death assemblages). The inclusion of death assemblages provides a link between the living seep fauna and its potential fossil record.
Detailed studies of modern shallow-water cold seeps are rare (Dando Reference Dando and Kiel2010), thus there is little quantitative ecological data for these ecosystems. In the GoM, cold-seeps at water depths <500 m are dominated by microbial mats (Sassen et al. Reference Sassen, Roberts, Aharon, Larkin, Chinn and Carney1993), and predation at shallow depths can prevent cold-seep macrofaunal communities from establishing (Sahling et al. Reference Sahling, Galkin, Salyuk, Greinert, Hilmar, Piepenburg and Suess2003). To broaden the environmental range of our comparisons, we include five shallow-water seeps (~150 m) from North Sea pockmarks (Dando et al. Reference Dando, Austen, Burke, Kendall, Kennicutt, Judd, Moore, O'Hara, Schmalijohann and Southward1991). These assemblages are thyasirid dominated, and we refer to them as the thyasirid biofacies in our analyses (Supplementary Table 3).
Live-collected Assemblages
The six included mussel-bed samples were collected from two localities about 29 km apart on the Louisiana slope (Bergquist et al. Reference Bergquist, Fleckenstein, Knisel, Begley, MacDonald and Fisher2005). The number of individuals >2 mm in size were reported, except for B. childressi, which was reported in terms of biomass. We estimate the number of B. childressi individuals using the area sampled and the number of individuals per square meter reported. Mussel-bed assemblages were collected from the Green Canyon (GC) petroleum lease area at 540–650 m water depth using a scoop and suction methodology, which may not capture the diversity of infauna and the pelagic community at these sites (Bergquist et al. Reference Bergquist, Fleckenstein, Knisel, Begley, MacDonald and Fisher2005).
The 23 included vestimentiferan tubeworm assemblages were collected from a 2532 km2 area across the GC and Garden Banks (GB) petroleum lease blocks at 528–571 m water depth (Bergquist et al. Reference Bergquist, Ward, Cordes, McNelis, Howlett, Kosoff, Hourdez, Carney and Fisher2003; Cordes et al. Reference Cordes, Hourdez, Predmore, Redding and Fisher2005). As for mussel beds, all individuals >2 mm are reported, while the tubeworm component is reported as biomass. We chose to exclude tubeworms from comparisons with the fossil record due to their low probability of preservation (see “Analytical Methods”) and thus do not calculate their absolute abundance. Tubeworm assemblages were collected with a mesh net attached to a ring that was used to scoop the aggregation, its associated epifauna, and some of its infauna; however, the infaunal and pelagic components of these communities may also be underrepresented (Cordes et al. Reference Cordes, Hourdez, Predmore, Redding and Fisher2005, Reference Cordes, Bergquist, Predmore, Jones, Deines, Telesnicki and Fisher2006, Reference Cordes, Bergquist and Fisher2009).
The five North Sea thyasirid assemblages were collected in 150–166 m water depth via box coring within a single pockmark 600 m by 300 m across (Dando et al. Reference Dando, Austen, Burke, Kendall, Kennicutt, Judd, Moore, O'Hara, Schmalijohann and Southward1991). Macrofauna were collected from the sediments using a 0.5 mm sieve and reported as the number of individuals per square meter.
Death Assemblages
Death assemblages used herein were collected from GC and GB petroleum lease blocks on the Louisiana continental slope at 513–754 m water depth using box cores and submersible-operated grab samples (Callender Reference Callender1992; MacDonald et al. Reference MacDonald, Schroeder and Brooks1995; Callender and Powell Reference Callender and Powell2000). Relative abundance data from five lucinid biofacies death assemblages were captured from published bar plots (MacDonald et al. Reference MacDonald, Schroeder and Brooks1995; Callender and Powell Reference Callender and Powell2000) using Image J (Rasband Reference Rasband2016). Relative abundance data for individual death assemblage samples were not available for tubeworm, mussel, and vesicomyid biofacies, but pooled data within each biofacies were provided (Callender Reference Callender1992); the thyasirid biofacies death assemblage was derived from a single box core (Callender Reference Callender1992).
Analytical Methods
Comparing Fossil and Modern Assemblages
Despite difficulties collecting from slope sites, modern sampling is still likely to record more individuals and more taxonomic diversity than fossil assemblages due to taphonomic loss of soft-bodied and small-shelled individuals in lithified fossil assemblages. Thus, we selectively removed individuals from modern assemblages so they more closely approximated the fossilizable assemblage for analysis. Only those taxa with fossilizable hard parts greater than 5 mm in length were included in the data sets, because lithified assemblages are biased against small taxa (Sessa et al. Reference Sessa, Patzkowsky and Bralower2009). Vestimentiferans and Serpula sp. were also removed from all data sets, because “tubeworms,” while potentially present at some WIS fossil seeps (Beauchamp and Savard Reference Beauchamp and Savard1992; Kauffman et al. Reference Kauffman, Arthur, Howe and Scholle1996; Larson et al. Reference Larson, Brezina, Landman, Garb and Handle2013; Meehan and Landman Reference Meehan and Landman2016), do not always preserve well at fossil seep sites (Shapiro and Fricke Reference Shapiro and Fricke2002; Hammer et al. Reference Hammer, Nakrem, Little, Hryniewicz, Sandy, Hurum, Drunkenmiller, Knutsen and Høyberget2011) and some reported “tubeworms” may be fluid conduits (Krauss et al. Reference Krauss, Clark, Sayegh and Perez2009). We did, however, retain the fossilizable component of the fauna of the tubeworm biofacies in our analyses. In addition, fossil taxa for which individuals were not well quantified (sponges, microbial mats, bryozoans) in our sampling of RR and MoG sites or in the modern cold-seep data (Cordes et al. Reference Cordes, Hourdez, Predmore, Redding and Fisher2005; Bergquist et al. Reference Bergquist, Fleckenstein, Knisel, Begley, MacDonald and Fisher2005) were removed from the AMNH data set.
For all data derived from the literature, we began with the original identifications of the authors, which we report in Supplementary Tables 2 and 3. The modern data sets often provide species-level, and at times subspecies-level, identifications, but this taxonomic specificity was not always feasible in the fossil assemblages due to poor preservation of key morphological features. Thus, we completed analyses using taxonomic identity at the generic level, except for taxa identified only to the family level (i.e. Vesicomyidae, Raninidae). We do not re-identify taxa, although we do reassign species to different genera if the World Register of Marine Species (Horton et al. Reference Horton, Kroh, Ahyong, Bailly, Boury-Esnault, Brandão and Costello2018) indicated changes in synonymization since publication of the original work (Supplementary Table 3).
Guild Assignment and Occupation
To facilitate direct comparison of modern and fossil assemblages despite little to no overlap in taxonomic composition, we assigned each taxon to a broad ecological guild based upon feeding mode, relationship to the substrate, and mobility following Bush et al. (Reference Bush, Bambach and Daley2007) (Table 3, Supplementary Tables 2 and 3). We assigned modern taxa to guilds using ecological information gathered from Neogene Marine Biota of Tropical America (Todd Reference Todd2001), the Biodiversity of the Gulf of Mexico Database (Moretzsohn et al. Reference Moretzsohn, Sánchez Chávez and Tunnell2015), the Marine Species Identification Portal (de Kluijver et al. Reference de Kluijver, Ingalsuo and de Bruyne2018; de Kluijver and Ingalsuo Reference de Kluijver and Ingalsuo2018), the European Register of Marine Species (Costello et al. Reference Costello, Emblow and White2001), and the literature (Dando et al. Reference Dando, Austen, Burke, Kendall, Kennicutt, Judd, Moore, O'Hara, Schmalijohann and Southward1991; Bergquist et al. Reference Bergquist, Ward, Cordes, McNelis, Howlett, Kosoff, Hourdez, Carney and Fisher2003, Reference Bergquist, Fleckenstein, Knisel, Begley, MacDonald and Fisher2005; Cordes et al. Reference Cordes, Hourdez, Predmore, Redding and Fisher2005, Reference Cordes, Bergquist, Predmore, Jones, Deines, Telesnicki and Fisher2006; Kraeuter Reference Kraeuter, Felder and Camp2009; McEachran Reference McEachran, Felder and Camp2009; Pawson et al. Reference Pawson, Vance, Messing, Solís-Marin, Mah, Felder and Camp2009; Rosenberg et al. Reference Rosenberg, Moretzsohn, García, Felder and Camp2009; Turgeon et al. Reference Turgeon, Lyons, Mikkelsen, Rosenberg, Moretzsohn, Felder and Camp2009)."
Principal coordinates analysis (PCO) scores for ordination analysis of guild abundances from fossil Black Hills and modern cold-seep assemblages. “G no.” refers to guild codes in Supplemental Tables 2 and 3. An “x” indicates guild presence in fossil or modern samples.
Fossil taxa were assigned to guilds using information gathered from the published literature (Kase and Ishikawa Reference Kase and Ishikawa2003; Lockwood Reference Lockwood2003; Kruta et al. Reference Kruta, Landman, Rouget, Cecca and Larson2010, Reference Kruta, Landman, Rouget, Cecca and Tafforeau2011; Landman et al. Reference Landman, Kennedy, Cobban, Larson and Jorgensen2013) and the Paleobiology Database. The ammonites Didymoceras nebrascense and Placenticeras spp. were assigned to a predatory guild based on their “cutting jaws,” which are consistent with a macrophagous diet (Landman et al. Reference Landman, Tsujita, Cobban, Larson, Tanabe and Flemming2006; Kruta et al. Reference Kruta, Landman, Rouget, Cecca and Larson2010). All other aptychophoran ammonites were assigned to a suspension-feeding guild based on the food contents found in the buccal mass of a baculite that were consistent with a suspension-feeding diet (Kruta et al. Reference Kruta, Landman, Rouget, Cecca and Tafforeau2011; Landman et al. Reference Landman, Cobban and Larson2012a). We assign Inoceramus spp. to an “unknown” guild, because while some authors suggest inoceramids harbored chemosymbionts (MacLeod and Hoppe Reference MacLeod and Hoppe1992; Kauffman et al. Reference Kauffman, Harries, Meyer, Villamil, Arango and Jaecks2007), others find no evidence for chemosymbiosis and instead suggest members of this genus were filter feeders (Berrocoso et al. Reference Berrocoso, MacLeod, Calvert and Elorza2008; Knight et al. Reference Knight, Morris, Todd, Howard and Ball2013). Our results do not change significantly across analyses in which we considered several possible guild assignments for Inoceramus spp.
Guild-Level Community Composition
Multivariate ordination of the proportional abundances of guilds from each assemblage can summarize similarities in community function between Cretaceous and modern seeps despite taxonomic dissimilarities. Faunal gradients in the live-collected tubeworm and mussel biofacies data sets were previously investigated using principal coordinates analysis (PCO) (Cordes et al. Reference Cordes, Bergquist, Predmore, Jones, Deines, Telesnicki and Fisher2006, Reference Cordes, Bergquist and Fisher2009); thus we focus on the PCO for comparison. All PCO analyses were performed on the Bray-Curtis dissimilarity matrix of the square-root-transformed proportional abundances using the cmdscale function in the ‘stats’ package of the R programming environment (R Core Team 2015). We used the same dissimilarity matrix in a cluster analysis based on unweighted pair-group method with arithmetic mean using hclust in the ‘stats’ package of the R programming environment (R Core Team 2015).
Assemblage Structure and Heterogeneity
To test whether genera were similarly partitioned into guilds within Cretaceous and modern seeps, we calculated the guild redundancy for each sample by dividing the number of genera in a sample by the number of guilds in the sample. The guild redundancy within modern and fossil data sets were compared with Mann-Whitney tests using the pairwise.wilcox.test function with Bonferroni correction in the ‘stats’ R package (R Core Team 2015).
We calculated the taxonomic and guild richness of individual seep assemblages using size-standardized subsamples of 47 individuals, which is the smallest abundance count of any one assemblage (we excluded samples with fewer than 47 individuals after culling for taxa unlikely to enter the fossil record). The rarefied faunal richness and the 95% confidence interval of the subsampling procedures were calculated using the rarefy function in the ‘vegan’ package for R (Oksanen et al. Reference Oksanen, Blanchet, Kindt, Legendre, Minchin, O'Hara, Simpson, Solymos, Stevens and Wagner2015). To determine whether richness differed among assemblage types, we applied the Mann-Whitney test to the rarefied richness using the wilcox.test function with Bonferroni correction in the base R package (R Core Team 2015).
We measured the taxonomic and guild evenness of abundance distributions from individual seep assemblages using the probability of interspecific encounter (PIE). PIE determines the probability that two random individuals drawn from an assemblage will be different species with lower values signifying a less even assemblage that is dominated by one or a few species (Hurlbert Reference Hurlbert1971). PIE is relatively unaffected by sample size as long as the sample size is larger than 20 individuals (Bulinski Reference Bulinski2007). The Mann-Whitney test was used to determine whether the evenness differed among assemblage types. The occurrence-based MoG assemblages were not included in PIE calculations.
We quantified the heterogeneity in generic and guild composition among samples within a biofacies using pairwise Simpson dissimilarities. Simpson dissimilarity was chosen over other metrics because it is relatively insensitive to differences in sample size and species richness (Koleff et al. Reference Koleff, Gaston and Lennon2003; Hammer and Harper Reference Hammer and Harper2006). Differences in dissimilarity can be driven by differences in sampling intensity (Beck et al. Reference Beck, Holloway and Schwanghart2013; Darroch and Wagner Reference Darroch and Wagner2015); thus we rarefied each abundance-based data set to 47 individuals and subsampled to five sites for 100 iterations using the rrarefy function in the ‘vegan’ package for R (Oksanen et al. Reference Oksanen, Blanchet, Kindt, Legendre, Minchin, O'Hara, Simpson, Solymos, Stevens and Wagner2015). This procedure randomly samples sites without replacement in each iteration and calculates an average dissimilarity matrix of a subsample of the data. The occurrence-based MoG data set was also subsampled for five sites through 100 iterations using the rrarefy function in the ‘vegan’ package for R (Oksanen et al. Reference Oksanen, Blanchet, Kindt, Legendre, Minchin, O'Hara, Simpson, Solymos, Stevens and Wagner2015).
Results
Cretaceous Black Hills Cold-Seep Assemblages
Richards Ranch Assemblages
The vuggy carbonate, Nymphalucina coquina, and thrombolytic microbialite facies were sampled at all 14 sites at RR. At 11 tepee buttes, 219 to 248 individuals were counted; however, site 8 had only 172 individuals and sites 10 and 14 had fewer than 100 individuals, even after exhausting the exposed surface of those tepee buttes. Didymoceras nebrascense occurred at sites 3 and 4. Baculites pseudovatus, which is only found in the Didymoceras nebrascense Zone (Kauffman et al. Reference Kauffman, Sageman, Kirkland, Elder, Caldwell and Kauffman1993; Larson et al. Reference Larson, Jorgensen, Farrar and Larson1997), occurred at sites 9 and 15. These ammonite identifications place the RR sites within the D. nebrascense Zone.
The taxonomic richness pooled across all RR sites is nine genera, representing seven guilds (Table 2). RR sites have a median of 1.1 genera per guild (Fig. 3). Nymphalucina occidentalis, an infaunal chemosymbiont, dominates all 14 tepee buttes and accounts for 82% to 99% of individuals exposed at each site (Fig. 4). Inoceramus spp. and Baculites spp. also frequently occur in the RR sites. Fossil specimens with ambiguous identifications are excluded from analyses. For example, cylindrical structures 4–7 mm in diameter occurred at site 5 (SDSM 109790) and at site 13 (not collected); these may be tubeworm trace fossils. In addition, site 12 contained a 7-mm-long fragment resembling a poorly preserved decapod cheliped (SDSM 109800).
Guild redundancy measured as the number of genera per guild. Boxes are defined by the interquartile range of the data, while whiskers are defined as 1.5× the interquartile range, and circles are data beyond 1.5× the interquartile range. Notches represent the 95% confidence interval of the true median's value. Fossil data sets are in gray. AMNH, American Museum of Natural History; MoG, South Dakota School of Mines and Technology Museum of Geology; RR, Richards Ranch.
Proportional abundance of ecological guilds for Black Hills fossil assemblages (Richards Ranch [RR] and American Museum of Natural History [AMNH]), the fossilizable component of modern live-collected assemblages (thyasirid, mussel, tubeworm), and death assemblages (lucinid). A, Feeding modes; B, substrate relationships. Inoceramus was not assigned to a guild. The deep infaunal guild is present in the lucinid biofacies, but is <0.5%. AMNH, American Museum of Natural History; RR, Richards Ranch.
Museum of Geology Collections
Baculites scotti occurs in samples from MoG 2, MoG 7, MoG 8, MoG 13, and MoG 15, which places these sites within the B. scotti Zone (dated to both 75.56 ± 0.11 Ma and 75.84 ± 0.26 Ma; Cobban et al. Reference Cobban, Walaszczyk, Obradovich and McKinney2006). Didymoceras nebrascense occurs at MoG 1 and MoG 3, which places these sites within the D. nebrascense Zone (dated to 75.19 ± 0.28 Ma; Cobban et al. Reference Cobban, Walaszczyk, Obradovich and McKinney2006). Hoploscaphites gilberti occurs in material from MoG 5, suggesting the sample is from either the B. scotti or D. nebrascense Zone (Landman et al. Reference Landman, Kennedy, Cobban, Larson and Jorgensen2013).
Pooling across all MoG collections, there are 19 genera representing 10 guilds. Site-level richness is 3–10 genera and 3–7 guilds. The median guild redundancy of MoG sites is 1.3 genera per guild (Fig. 3), which is indistinguishable from the RR and AMNH data sets (see Supplementary Table 4 for all Mann-Whitney test results). MoG 4 and MoG 14 contain fragments suspected as tubeworm traces, but we exclude them from analysis as we did for the RR samples. Decapods are well-preserved, and all could be confidently identified at least to the family Raninidae. Diagnostic features for genus-level identifications are indiscernible, with the exceptions of the type specimens SDSM 11016, 11017, 11018, 11019, 11021, 11022, and 11023 (Bishop and Williams Reference Bishop and Williams2000) from MoG 3, MoG 4, and MoG 14. To maintain consistency in identification, we lumped the types into Raninidae for analyses. Of the 19 taxa found across the MoG and RR collections, only Heus foersteri and other members of Raninidae are restricted to cold-seep localities (Bishop and Williams Reference Bishop and Williams2000; Meehan and Landman Reference Meehan and Landman2016).
American Museum of Natural History Collections
Published AMNH cold-seep deposits include samples from the Baculites compressus, Didymoceras cheyennense, and Baculites scotti–Didymoceras nebrascense ammonite zones (Meehan and Landman Reference Meehan and Landman2016). Those authors used the physical structure of the seeps to place them into three categories, reflecting different stages of seep development: (1) conduit with few secondary pipes and no carbonate cap implying strong advective flow, (2) conduit with moderate number of secondary pipes and a small carbonate cap implying diffusive flow, and (3) conduit with multiple secondary pipes and a broad carbonate cap implying continued flow but dampening at the sediment–water interface (Handle Reference Handle2014; Meehan and Landman Reference Meehan and Landman2016).
The taxonomic richness pooled across all AMNH sites is 28 genera, representing 13 guilds with a median guild redundancy of 1.3 genera per guild. Like the RR seeps, AMNH sites are often dominated by Nymphalucina spp., but also have high proportional contributions from Inoceramus spp., pelagic suspension feeders such as Baculites spp. and Hoploscaphites spp., and predators such as Didymoceras spp. and Placenticeras spp. (Fig. 4); all of these taxa are common across the WIS and are not restricted to seep sites (Meehan and Landman Reference Meehan and Landman2016).
Modern Cold-Seep Assemblages
Of the 22 guilds identified in this study, 10 occur in the modern data sets that do not appear in the fossil data sets (Table 3). The original taxonomic lists from modern sites included 238 taxa, but only 96 had potentially fossilizable hard parts greater than 5 mm in length (see Supplementary Table 3 for complete faunal lists).
The 23 vestimentiferan seep samples from modern GoM cold-seeps yielded a total of 37 fossilizable genera across 15 guilds, but are dominated by surficial guilds with only minor contributions from infaunal and pelagic guilds (Fig. 4). The median guild redundancy is 1.6 genera per guild, significantly higher than the RR and MoG collections but comparable to all other assemblages (Fig. 3). The tubeworm biofacies death assemblage is less diverse than the live-collected tubeworm biofacies (Table 2).
The six mussel biofacies samples from the GoM yielded 10 fossilizable genera across six guilds with a median guild redundancy of 1.3 genera per guild. In contrast to the tubeworm biofaces, the mussel biofacies death assemblage is richer at both the genus level and guild level than the live-collected mussel biofacies (Table 2). Like the tubeworm biofacies, the mussel biofacies is dominated by surficial guilds, but the mussel-bed assemblages have a higher proportion of fossilizable chemosymbionts (Fig. 4).
The five live-collected thyasirid biofacies from the North Sea yielded 10 genera in seven guilds with a median guild redundancy of 1.7 (Fig. 3, Table 2). These assemblages are dominated by the infaunal chemosymbiotic guild with contributions from other infaunal guilds and predators (Fig. 4). The GoM thyasirid biofacies death assemblage is monogeneric and contains only Thyasira spp. The five lucinid biofacies death assemblages from the GoM yielded a total of 21 fossilizable genera across 12 guilds, with a median guild redundancy of 1.4 (Fig. 3, Table 2). The lucinid biofacies has a similar guild composition as the thyasirid facies, but it has a greater contribution from surficial guilds and suspension feeders (Fig. 4).
Comparisons among Modern and Fossil Cold-Seep Assemblages
Guild Composition
Guild-level relative abundances from all modern and fossil samples (except the MoG data set) were incorporated into (1) a single PCO (k = 3) and (2) a single cluster analysis to examine similarities in guild composition among assemblages. PCO axis 1 summarizes 57.2% of the faunal variation and separates the live-collected tubeworm and mussel biofacies (negative axis 1 scores) from the modern lucinid, modern thyasirid, and fossil Black Hills assemblages (Fig. 5). Similarly, the first major division in the cluster analyses separates the live-collected tubeworm and mussel biofacies from all other samples (Fig. 6). Pelagic guilds, infaunal guilds that are present in both modern and fossil data sets, and all suspension feeders have positive PCO axis 1 scores. All surficial guilds that are not suspension feeders and the infaunal guilds that are absent from the fossil data sets have negative PCO axis 1 scores (Fig. 5, Table 3). Pooled modern death assemblages from the vesicomyid, mussel, and tubeworm biofacies have intermediate PCO axis 1 values (Fig. 5). The cluster analysis most closely affiliates these death assemblages and the lucinid biofacies death assemblages, which together form a distinct cluster (cluster 2.1), with the fossil assemblages and thyasirid biofacies (Fig. 6).
Principal coordinates analysis (PCO) of fossil Black Hills and modern cold-seep assemblages. A, Sample scores on PCO axes 1 and 2. For modern assemblages, symbol shapes represent the different biofacies, with open symbols indicating live-collected assemblages and gray-filled symbols indicating death assemblages. For fossil assemblages, symbols differentiate Richards Ranch (RR) samples (this study) from American Museum of Natural History (AMNH) samples; AMNH samples are further coded to indicate the seep categories assigned by Meehan and Landman (Reference Meehan and Landman2016). B, Guild scores on PCO axes 1 and 2. Symbol shape indicates substrate relationship, and symbol color indicates feeding mode. Black filled, predators; dark gray filled, chemosymbionts; light gray filled; suspension feeders; open symbols are other feeding modes, including grazing and mining.
Cluster analysis of fossil Black Hills and modern cold-seep assemblages. Modern assemblages include five biofacies: tubeworm (TW), mussel (MB), lucinid (LB), vesicomyid (VB), and thyasirid (TB). DA denotes modern samples collected as death assemblages. Fossil assemblages include the Richards Ranch (RR) samples (this study) and American Museum of Natural History (AMNH) samples (Meehan and Landman Reference Meehan and Landman2016). “Cat#” refers to the seep categories assigned by Meehan and Landman (Reference Meehan and Landman2016). Sample numbers otherwise as included in Supplementary Tables 2 and 3.
PCO axis 2 summarizes 16.9% of the faunal variation. Among samples with low PCO axis 1 scores, PCO axis 2 differentiates live-collected tubeworm biofacies assemblages (more positive PCO 2 scores) from live-collected mussel biofacies assemblages (more negative PCO 2 scores), although the groups overlap in PCO space (Fig. 5A). Similarly, the primary mussel-bed and tubeworm cluster is broken into two secondary clusters, with one composed entirely of tubeworm biofacies samples and the other a mix of both biofacies (Fig. 6). Within the tubeworm biofacies, PCO axis 2 scores are positively correlated with the age of the assemblage (Spearman's rho = 0.7263, p < 0.001; Supplementary Table 3). Among samples with high PCO axis 1 scores, PCO axis 2 differentiates the modern lucinid and thyasirid biofacies and RR fossil assemblages (more negative PCO axis 2 scores) from the majority of AMNH fossil assemblages (Fig. 5A). This division is also reflected in the cluster analysis, where all but four of the AMNH assemblages form a cluster distinct from the RR assemblages, thyasirid biofacies, and remaining AMNH samples (cluster 2.2 in Fig. 6). Higher relative abundances of predatory guilds occur in samples with high PCO axis 2 scores (Fig. 5B, Table 3). Inoceramus spp., which we assigned to an “unknown” guild, has a high PCO axis 2 score. Low PCO axis 2 scores are largely related to higher relative abundances of the chemosymbiotic guilds (Fig. 5B, Table 3). We do not interpret PCO axis 3, as it summarizes 10.9% of the variation in our data and collapses the variation among fossil assemblages.
The physical seep categories defined by Handle (Reference Handle2014) for the AMNH seeps do not appear related to gradients in faunal composition defined by our PCO analyses (Fig. 5A). One AMNH seep (no. 3343), listed as category 3 by Meehan and Landman (Reference Meehan and Landman2016), overlaps with RR seeps, but category 3 seeps occur throughout the ordination space occupied by AMNH seeps sites as a whole. In the cluster analysis, three of the four AMNH seeps that are grouped with the RR assemblages and thyasirid biofacies are category 3, but category 3 seeps are also closely affiliated with category 1 and 2 seeps in the AMNH exclusive cluster (cluster 2.2 in Fig. 6).
Assemblage Structure and Heterogeneity
After sample-size standardization, RR samples have significantly lower genus and guild richness than all other assemblage groups except for the thyasirid biofacies in terms of guild richness (Fig. 7). Similarly, RR samples have significantly lower evenness at the genus level than all other assemblages (Fig. 8), whereas at the guild level RR has significantly lower evenness than all other assemblages except the thyasirid biofacies after Bonferroni correction (see Supplementary Table 4 for all Mann-Whitney test results).
Rarefied richness (n = 47) of samples by data set. A, Genus level; B, guild level. Boxes are defined by the interquartile range of the data, while whiskers are defined as 1.5× the interquartile range, and circles are data beyond 1.5× the interquartile range. Notches represent the 95% confidence interval of the true median's value. Fossil data sets are in gray. *Lacks relative abundance data, so not rarefied. AMNH, American Museum of Natural History; MoG, South Dakota School of Mines and Technology Museum of Geology; RR, Richards Ranch.
Evenness calculated as the probability of interspecific encounter (PIE) for each site by data set. A, Genus level; B, guild level. Boxes are defined by the interquartile range of the data, while whiskers are defined as 1.5× the interquartile range, and circles are data beyond 1.5× the interquartile range. Notches represent the 95% confidence interval of the true median's value. Fossil data sets are in gray. AMNH, American Museum of Natural History; RR, Richards Ranch.
In contrast to the RR samples, AMNH samples have significantly greater genus and guild richness than the mussel and thyasirid biofacies, but have richness comparable to the lucinid and tubeworm biofacies (Fig. 7). AMNH evenness is also not significantly different from any of the modern assemblage groups (Fig. 8).
At the genus level, RR has significantly lower among-sample dissimilarity than all other assemblage groups except for the mussel biofacies (Fig. 9). At the guild level, RR, MoG, and the mussel biofacies all have low dissimilarity among constituent samples, and their within-group dissimilarity is not significantly different. AMNH samples have significantly greater among-sample dissimilarity than the RR and MoG assemblages but are significantly more homogeneous than the lucinid, thyasirid, and tubeworm biofacies (Fig. 9).
Simpson dissimilarity among samples from the same data set. A, Genus level; B, guild level. Boxes are defined by the interquartile range of the data, while whiskers are defined as 1.5× the interquartile range, and circles are data beyond 1.5× the interquartile range. Notches represent the 95% confidence interval of the true median's value. Fossil data sets are in gray. *Lacks relative abundance data, so not sample-size standardized. AMNH, American Museum of Natural History; MoG, South Dakota School of Mines and Technology Museum of Geology; RR, Richards Ranch.
Discussion
The Cretaceous Black Hills and modern cold-seep assemblages examined here share only the genus Pecten; however, by examining the assemblages with taxon-free ecological metrics and comparing assemblage composition using a guild-defined ecospace, we can assess the ecological similarity of these communities across deep time.
To lessen the impact of taphonomically driven differences between modern and fossil assemblages, we remove taxa from modern data sets that do not have fossilizable hard parts large enough to reliably recover from lithified material. Despite this, some guilds are present in the modern assemblages that are not in the fossil assemblages; our modern data sets span a much greater range of biofacies than our fossil material, and thus we would expect additional ecological modes. In addition, new lineages have evolved physiological tolerances for seep environments since the Cretaceous (Kiel Reference Kiel and Kiel2010a, Reference Kiel2015). Pelagic guilds are also more common at Cretaceous seeps than modern seeps, which is partly due to the extinction of ammonites and partly due to modern live-collection methods that can exclude fast-moving pelagic organisms.
The modern death assemblages we incorporate in our multivariate analyses come from the same GoM region as the live-collected assemblages. Their inclusion allows us to assess changes in faunal composition expected as the live assemblages transition into the sedimentary record, which may not be wholly accounted for by excluding soft-bodied taxa. We find that modern mussel and tubeworm biofacies death assemblages incorporated proportionally more infaunal taxa, which distanced them compositionally from their living counterparts along PCO axis 1 and in the cluster analysis, but both biofacies remained distinct from the lucinid death assemblages and thyasirid living assemblages in the PCO analysis. Thus, the hard-substrate and soft-substrate assemblages should be distinguishable in the fossil record even if the sometimes-dominant soft-bodied components are absent. The thyasirid GoM death assemblage is also compositionally similar to the living thyasirid North Sea assemblages despite these assemblages being from very different localities and water depths.
We additionally exclude ambiguous taxa from fossil data sets and remove taxonomic groups from the AMNH samples that we did not devote effort to observe in the RR and MoG assemblages in order to make the data sets more comparable. While we have removed any suspected “tubeworms,” the lack of limpets that used tubeworms as their habitat (e.g., Jenkins et al. Reference Jenkins, Kaim and Hikida2007a) further suggests tubeworms were indeed rare or absent from the Black Hills assemblages. Furthermore, we analyze fossil collections obtained by different workers separately so we can address potential differences due to collection strategies. For example, RR and MoG collections were obtained from surface samples, or broken pieces of surface samples, but AMNH assemblages were collected from rock volumes (Meehan and Landman Reference Meehan and Landman2016), which could enhance the richness of the AMNH samples. However, richness is comparable between the AMNH and MoG data sets, and thus the low diversity of the RR samples cannot simply be due to its restriction to surface samples or differences in taxonomic expertise, because the same workers evaluated both the RR and MoG collections.
Comparing the Guild Composition of Black Hills Cold-Seeps to Modern Biofacies
Shallow- and deep-water seeps from the Late Cretaceous are often different in their taxonomic composition (Kiel Reference Kiel2010b; Kiel et al. Reference Kiel, Wiese and Titus2012), and this may influence our guild-based analyses. For example, in the present analyses, all suspension-feeding guilds have positive PCO axis 1 scores, with their highest relative abundances occurring in some AMNH seep assemblages; the distribution of this feeding mode may reflect the greater influence of organic carbon production from photosynthesis at shallow-water seeps (Levin et al. Reference Levin, Mendoza, Gonzales, Thurber and Cordes2010; Dando Reference Dando and Kiel2010). However, the deeper-water modern lucinid biofacies has proportionally higher abundances of suspension feeders than the shallow-water RR seeps; thus suspension-feeding guilds do not entirely covary with water depth in this data set. In addition, modern assemblages from the same water depth occur across the full length of PCO axis 1, indicating that water depth does not entirely explain the primary variation in the relative abundances of guilds in our analyses.
Black Hills cold-seep assemblages, modern thyasirid assemblages from shelf depths, and modern lucinid and thyasirid assemblages from slope depths overlap on PCO axis 1 and form a coherent group in the cluster analysis. These assemblages are dominated by shallow-infaunal guilds, consistent with a high availability of soft substrate, and share abundant chemosymbiotic and suspension-feeding guilds. In contrast, the modern mussel and tubeworm biofacies samples, which are from the same water depth and geographic region as the lucinid and thyasirid slope assemblages, have much lower PCO axis 1 scores and form a distinct cluster. These two biofacies are dominated by surficial taxa and have abundant surface-grazing guilds, consistent with a high availability of hard substrates. Thus, by examining faunal composition using the relative abundances of guilds, we recognize an ecological gradient that reflects differences in substrate availability at cold-seeps despite water-depth differences that may drive differences in taxonomic composition.
Some RR samples overlap with modern thyasirid assemblages in both PCO axes 1 and 2, which combined summarize 74% of the variation in guild composition, indicating the RR collections are nearly identical in guild composition to the fossilizable component of these modern assemblages. Despite their similarity in guild composition, RR seeps have significantly lower evenness and lower richness at the genus level than the thyasirid biofacies, although they are more similar in richness at the guild level. The infaunal chemosymbiotic bivalves in the thyasirid biofacies were originally all assigned to the genus Thyasira (Dando et al.Reference Dando, Austen, Burke, Kendall, Kennicutt, Judd, Moore, O'Hara, Schmalijohann and Southward1991), but these taxa were later split into four different genera (Horton et al. Reference Horton, Kroh, Ahyong, Bailly, Boury-Esnault, Brandão and Costello2018); using the original taxonomy, the RR and thyasirid assemblages would be indistinguishable in genus richness. This taxonomic decision may also explain the higher guild redundancy in thyasirid biofacies than in the fossil assemblages. All RR and thyasirid assemblages are contained in a distinct cluster accompanied by the four AMNH assemblages with the lowest PCO axis 2 values. Overall, the AMNH assemblages have similar PCO axis 1 scores as the modern thyasirid and lucinid biofacies, but have, on average, higher PCO axis 2 scores and higher richness more comparable with the modern slope-depth lucinid assemblages from the GoM than with the shallow-water thyasirid assemblages from the North Sea.
Previous workers also concluded that Colorado tepee butte assemblages are most similar to lucinid and thyasirid biofacies (Callender and Powell Reference Callender and Powell1997; Powell et al. Reference Powell, Callender and Stanton1998) based upon the dominance of lucinids and low diversity. Other similar Cretaceous shallow-water cold-seep assemblages are also reported from Utah (Kiel et al. Reference Kiel, Wiese and Titus2012) and Antarctica (Little et al. Reference Little, Birgel, Boyce, Crame, Francis, Kiel, Peckman, Pirrie, Rollinson and Witts2015). While cold-seep assemblages similar to our Black Hills faunas are more common at depths <300 m (Majima et al. Reference Majima, Nobuhara and Kitazaki2005), some deep-water Campanian seeps from Hokkaido, Japan (Omagari site; Hikida et al. Reference Hikida, Suzuki, Togo and Ijiri2003; reviewed in Majima et al. Reference Majima, Nobuhara and Kitazaki2005) and some slope to bathyal Early Cretaceous and Cenozoic seeps from California and the Pacific Northwest (Hickman Reference Hickman1984; Campbell and Bottjer Reference Campbell and Bottjer1993) are similarly dominated by thyasirids and lucinids. These deeper-water instances of fossil cold-seeps dominated by infaunal chemosymbionts further suggest that the assemblage composition of WIS seeps is not purely a result of their shallow water depth. Instead, the dominance of infaunal chemosymbionts may relate to the successional stage of seep development in which seepage was stable but significant carbonate cements had not yet precipitated (early stage C of Beauchamp and Savard [1992]; earliest diagenetic stage of Campbell et al. [2002]; stage 2 of Hikida et al. [2003]; type III of Majima et al. [2005]).
Modern shelf-depth seeps are ephemeral and often do not persist long enough for communities dominated by guilds other than infauna to establish (Levin et al. Reference Levin, James, Martin, Rathburn, Harris and Michener2000b; Campbell Reference Campbell2006); if WIS seeps were similarly short-lived, this could explain the dominance of early successional stage guilds. The consistent dominance of infaunal taxa in Black Hills seeps further suggests that sediment covered the authigenic carbonates that formed beneath the sediment–water interface and created little hard substrate for epifaunal attachment. This hypothesis is reinforced by the rarity of taxa associated with hard substrates, such as some echinoderms (Landman et al. Reference Landman, Cochran, Larson, Brezina, Garb and Harries2012b; Meehan and Landman Reference Meehan and Landman2016; Thuy et al. Reference Thuy, Landman, Larson and Numberger-Thuy2018). Sedimentation rates of surrounding shale that kept up with carbonate production (Gilbert and Gulliver Reference Gilbert and Gulliver1894; Metz Reference Metz2010) may have prevented exposure of sizable amounts of hard substrate, despite a growing carbonate deposit during the seep's active period. High sedimentation rates coupled with a relatively short seep-active period could thus maintain the dominance of soft-sediment faunas at the Black Hills seep sites.
Comparing the Ecological Structure of Modern and Fossil Cold-Seeps
WIS marine faunas are generally more depauperate than other Cretaceous faunas due to the basin's restricted circulation, brackish surface waters, dysoxic benthos, and episodic connection between the northern and southern regions, which would restrict dispersal (Jeletzky Reference Jeletzky1970; Hay et al. Reference Hay, Eicher, Diner, Caldwell and Kauffman1993; Lowery et al. Reference Lowery, Leckie, Bryan, Elderbak, Parker, Polyak, Schmidt, Snoeyenbos-West and Sterzinar2018). An even smaller pool of species with the physiological tolerances necessary to colonize Black Hills seep habitats would have been present in the WIS; for example, tolerance of sulfide is important for colonization of shallow seeps (Dando Reference Dando and Kiel2010). In addition, seep-obligate taxa are rare to absent at WIS seeps (Bishop and Williams Reference Bishop and Williams2000; Kiel Reference Kiel2010b; Meehan and Landman Reference Meehan and Landman2016), similar to modern shallow-water seeps, which have few seep-obligate symbiont-bearing species due to the high flux of organic matter from photosynthetic production that can support organisms otherwise tolerant of sulfidic conditions (Sahling et al. Reference Sahling, Galkin, Salyuk, Greinert, Hilmar, Piepenburg and Suess2003; Levin et al. Reference Levin, Mendoza, Gonzales, Thurber and Cordes2010; Dando Reference Dando and Kiel2010). Thus, we might expect WIS seeps to be less diverse than the modern deep-water seeps we include.
While RR seeps do have lower richness and evenness than the modern seeps we examine, the other two Black Hills data sets (MoG and AMNH) have sample-level richnesses comparable to, or higher than, some of the modern seeps. In addition, evenness is similar among AMNH samples and the majority of modern biofacies. Thus, the lower diversity is not characteristic of WIS seeps in general. The number of genera per guild is also comparable among the three fossil data sets and the modern lucinid and mussel biofacies, suggesting modern genera are not more finely partitioning guild space. Furthermore, faunal lists from Late Cretaceous seep sites outside of the WIS are not richer than those in the Black Hills data set: a previous compilation reports genus richness of 1 to 12 (median = 5.5; Kiel et al. Reference Kiel, Wiese and Titus2012), well within the range of our WIS sample-level genus richness. Thus, it appears that richness at cold-seeps in the WIS is not generally lower than modern cold-seeps nor is it lower than richness of other seeps in the Late Cretaceous; however, sample-size standardized measures of richness are necessary to make these broader comparisons.
Among-Sample Heterogeneity
In addition to having lower richness and evenness, RR seep assemblages are also more homogeneous among samples than the compositionally similar shallow-water thyasirid biofacies. The spatial and temporal extent of collections within a data set may influence the faunal dissimilarity among samples with larger areas and time spans expected to produce greater dissimilarity (MacArthur and Wilson Reference MacArthur and Wilson1963; Qian and Ricklefs Reference Qian and Ricklefs2012). However, the low dissimilarity among seep assemblages is not simply due to the restricted geographic area of the RR seeps; MoG seeps, which encircle the Black Hills region, are similarly homogeneous at the guild level, although they are more heterogeneous than RR samples at the genus level. The geographically broader AMNH assemblages are also more similar to one another than are samples from all modern biofacies, except the mussel biofacies. This homogeneity exists despite an expectation of higher dissimilarity among the fossil samples based upon their temporal breadth, which samples multiple ammonite zones, whereas the modern samples were collected within years of one another. Thus, the fossil seeps are more homogeneous in their compositions than would be expected from these modern seeps.
Environmental conditions can contribute to compositional homogeneity of seep communities. For example, higher fluxes of organic carbon from the photic zone at shallow-water seeps can increase the heterogeneity of seep faunas, because the distribution of organic matter across seeps is itself heterogeneous and leads to microhabitat differences within the sediments (Dando Reference Dando and Kiel2010). In contrast, high predation pressure in modern shallow seeps can homogenize faunal composition. For example, infaunal chemosymbiotic bivalves often dominate at modern shallow seeps because their subsurface habitat protects them from pelagic predators (Sahling et al. Reference Sahling, Galkin, Salyuk, Greinert, Hilmar, Piepenburg and Suess2003; Dando Reference Dando and Kiel2010). The presence of ammonite and gastropod predators at the Black Hills seeps investigated here and the low δ13C values of ammonite shells found preserved in seep environments (Landman et al. Reference Landman, Cochran, Larson, Brezina, Garb and Harries2012b, Reference Landman, Cochran, Slovacek, Larson, Garb, Brezina and Witts2018; Cochran et al. Reference Cochran, Landman, Larson, Meehan, Garb and Brezina2015) suggest predation did occur at WIS seeps. Thus, predation pressure may contribute to homogeneity among the Black Hills seeps; however, heterogeneity is greater at AMNH seeps where predators are more frequently preserved, making predation a less likely driver.
Among the modern assemblages we analyze, only the slope-depth mussel biofacies, which persists in the highest sulfide and methane conditions of the modern seeps explored here (Cordes et al. Reference Cordes, Bergquist and Fisher2009), is as homogeneous as the fossil Black Hills samples. This structural similarity with the mussel biofacies may suggest that the RR seeps reflect a level of environmental stress from seep fluids similar to that of the mussel biofacies. The MoG and AMNH assemblages, which have greater dissimilarity among samples than the RR and mussel biofacies, may represent seeps with generally lower seep-fluid influence, although they are still more homogeneous than the other modern biofacies, including the compositionally similar lucinid and thyasirid biofacies.
In addition to stress from seep-fluid composition, the restricted circulation and low-oxygen conditions common in the WIS could be an additional stressor leading to more homogeneous assemblages (Jeletzky Reference Jeletzky1970; Wright Reference Wright1987; Hay et al. Reference Hay, Eicher, Diner, Caldwell and Kauffman1993; Lowery et al. Reference Lowery, Leckie, Bryan, Elderbak, Parker, Polyak, Schmidt, Snoeyenbos-West and Sterzinar2018). Well-developed oxygen minimum zones in the modern ocean tend to homogenize faunal assemblages (Levin et al. Reference Levin, Gage, Martin and Lamont2000a); thus the low-heterogeneity of the Black Hills seeps may instead reflect lower-oxygen conditions during their formation. However, this compositional heterogeneity may not be unique to WIS seeps but may instead reflect a difference between modern and Cretaceous cold-seep ecosystems. For example, two shelf-edge seep faunal lists reported from the Campanian of Japan (Hikida et al. Reference Hikida, Suzuki, Togo and Ijiri2003; Jenkins et al. Reference Jenkins, Kaim, Hikida and Tanabe2007b; analyzed in Kiel et al. Reference Kiel, Wiese and Titus2012) are also fairly homogeneous in their taxonomic composition, as they share 10 of 14 genera among the 12 genera reported from each site. Additional multi-seep comparisons within Cretaceous regions outside the WIS are necessary to test whether the low heterogeneity we find in the WIS is a general feature.
Evidence for Community Succession among Black Hills Seep Assemblages
Modern seep faunas in the GoM progress through a series of successional stages as exposed carbonate increases and the influence of seep fluids wanes (Cordes et al. Reference Cordes, Bergquist, Predmore, Jones, Deines, Telesnicki and Fisher2006, Reference Cordes, Bergquist and Fisher2009). Differences between two Cretaceous fossil seep sites in Hokkaido, Japan, have also been inferred to reflect the longevity of the seep and changes in seep discharge (Jenkins et al. Reference Jenkins, Kaim, Hikida and Tanabe2007b). The lack of an obligate seep fauna in the WIS makes it difficult to recognize a decrease in seep-fluid influence over time, which is recognizable in modern, deeper-water seeps from the decreasing dominance of seep-obligate taxa. A minor increase in the proportion of epifaunal taxa through a suspected seep succession among the AMNH samples used here suggests some compositional change related to seep complexity, but the dominance of infaunal, foundational, taxa persists, and the transition to seeps dominated by secondary taxa as observed at modern seeps does not occur (Meehan and Landman Reference Meehan and Landman2016). However, by combining modern seeps that display faunal succession with fossil seep assemblages in the same multivariate ordination, we are able to elucidate parallel ecological changes among modern and fossil seeps.
PCO axis 2 in our guild-based ordination reflects ecological changes expected with the aging of a seep. Among the modern assemblages, samples from the tubeworm biofacies have higher PCO axis 2 scores than those from the mussel biofacies. In addition, the tubeworm biofacies samples with the highest PCO axis 2 scores tend to be the oldest seeps, and there is a strong correlation between seep age and PCO axis 2 score (Supplementary Table 3). Cordes et al. (Reference Cordes, Bergquist, Predmore, Jones, Deines, Telesnicki and Fisher2006, Reference Cordes, Bergquist and Fisher2009) reported a similar gradient in these same live-collected assemblages that reflected a gradient in exposure to seep fluids. In the modern GoM, seep communities associated with the mussel biofacies are more tolerant of hypoxic conditions and elevated methane and sulfide concentrations. Communities associated with tubeworm biofacies experience lower epibenthic sulfide concentrations, and the tubeworms create habitat that permits colonization by nonendemic taxa and higher trophic levels such as predators (Cordes et al. Reference Cordes, Bergquist and Fisher2009). All chemosymbiotic guilds, whether infaunal or surficial have negative or weakly positive PCO axis 2 scores, while predatory taxa and the majority of suspension feeders have higher PCO axis 2 scores, consistent with a successional interpretation of PCO axis 2. Thus, we interpret our PCO axis 2 as a gradient in seep succession and suggest that the parallel gradient in fossil assemblages along PCO axis 2 similarly reflects a gradient from relatively higher sulfide at negative values to lower epibenthic sulfide and greater habitat complexity at positive values. Similar spatial clines in community composition from a Cretaceous seep in Hokkaido, Japan (Yasukawa site; Jenkins et al. Reference Jenkins, Kaim, Hikida and Tanabe2007b) and in modern seeps (Sahling et al. Reference Sahling, Rickert, Lee, Linke and Suess2002; Dando Reference Dando and Kiel2010; Decker et al. Reference Decker, Morineaux, Van Gaever, Caprais, Lichtschlag, Gautherier, Andersen and Olu2012; Levin et al. Reference Levin, Baco, Bowden, Colaco, Cordes, Cunha and Demopoulos2016) reflect the decreasing influence of seep fluids with distance, supporting our inference that the faunal differences we observe among Black Hills seeps reflect differences in seep-fluid exposure. The parallelism between the modern hard substrate–dominated seeps and the soft sediment–dominated Black Hills seeps demonstrates that these Cretaceous seeps had a similar trajectory of faunal succession as seen in the modern hard substrate–dominated seeps.
Inoceramus spp., which was not assigned to a guild in this study out of caution, does not covary with the chemosymbiotic guilds in ordination space and instead has one of the highest PCO axis 2 scores. This suggests Inoceramus spp. has environmental preferences quite different from the chemosymbiotic taxa and is not itself chemosymbiotic. This is consistent with recent evidence from gill morphology that inoceramids do not have chemosymbionts (Knight et al. Reference Knight, Morris, Todd, Howard and Ball2013).
Given this guild-based gradient along PCO axis 2, we interpret the RR samples as being from a younger successional stage than the majority of the AMNH samples. This is consistent with the lower evenness and higher homogeneity of the RR samples, characteristics expected of a more stressed community (Odum Reference Odum1985) and characteristics that the RR samples share with the early successional stage mussel biofacies in the modern. While the AMNH samples occupy a large part of the PCO axis 2 gradient, the physical seep characteristics described by Handle (Reference Handle2014; and related to faunal assemblages in Meehan and Landman [2016]) are unrelated to this axis; thus those physical categories are independent of the faunal gradient recognized here. Future geochemical or petrographic studies could be used to test these interpretations of seep successional stages made from faunal composition.
Conclusions
Translating taxonomic identities into ecological guilds allowed disparate assemblages separated by tens of millions of years to be included in a single multivariate ordination and compared using taxon-free ecological metrics. We found that guild composition of the Black Hills seeps most closely resembles the modern thyasirid and lucinid biofacies, primarily in the dominance of infaunal chemosymbiotic bivalves. The persistence of these soft-sediment communities without a transition to hard-substrate communities likely reflects high sedimentation rates that kept up with carbonate production over the life of the seep.
Faunal composition is more similar among fossil WIS seep assemblages than would be expected from the modern seeps we explored. This homogeneity may reflect environmental stresses such as low oxygen or high fluxes of seep fluids in the WIS assemblages, but this should be tested more broadly across cold-seep faunas to determine whether this is a general difference between modern and Cretaceous seep communities.
Guild-based analyses also allow us to apply ecological gradients emergent from modern faunal samples to the fossil assemblages, because the data can be analyzed in a common ordination space. Modern hard-substrate assemblages define a gradient from young seeps with high sulfide and methane and low oxygen to older seeps less influenced by seep fluids with more structural complexity. The Black Hills assemblages form a parallel gradient in ordination space, allowing us to recognize a similar trajectory of faunal succession despite differences in the primary substrate type, water depth, and taxonomic composition. Black Hills assemblages analogous to younger successional stage modern seeps on this gradient (mussel biofacies) also have similarly low among-sample dissimilarity, further suggesting these fossil assemblages developed in a more stressed environment than those analogous to older successional stage modern seeps. Thus, by standardizing modern and fossil data to account for differences in preservation and using taxon-free methods, we are able to make more detailed paleoecological interpretations and recognize a more precise analogy between modern and fossil systems. Broader application of these methods to other seep faunas and to other ecological settings may help strengthen the integration of modern ecology and paleoecology.
Acknowledgments
We thank S. Danise, K. Meehan, and an anonymous reviewer for their valuable suggestions and review of this publication. We also wish to thank K. Antonen, L. Anderson, and the rest of the staff at the South Dakota School of Mines and Technology Museum of Geology for their help and support. Special thanks go to J. Brezina and the Richards family for introductions and graciously allowing access to their lands, respectively. J.D.L. was supported by Institute of Museum and Library Services grant IMLS-MA-30-13-0500-13.
