-
PDF
- Split View
-
Views
-
Cite
Cite
Dirley Cortés, Hans C E Larsson, Top of the food chains: an ecological network of the marine Paja Formation biota from the Early Cretaceous of Colombia reveals the highest trophic levels ever estimated, Zoological Journal of the Linnean Society, Volume 202, Issue 1, September 2024, zlad092, https://doi.org/10.1093/zoolinnean/zlad092
- Share Icon Share
Abstract
The Mesozoic Marine Revolution restructured the world’s ocean biodiversity into the complex marine ecosystems of today. This revolution began during the Triassic but the origin of this complexity is poorly understood due to a lack of detailed ecosystem reconstructions throughout time. We present the first site-specific ecological network for a marine Mesozoic fauna based on the Early Cretaceous Paja Formation biota of Colombia that preserves numerous, large-bodied, predatory marine reptiles. The trophic food-web was quantitatively reconstructed based on inferred trophic interactions of marine producers, consumers, and large apex predators. Compared to well-studied Caribbean reef ecosystem networks, the Paja biota network is missing a great proportion of benthic invertebrates and fishes, despite its rich higher trophic levels. We hypothesize that the ammonites from the Paja biota either mirrored the diversity represented by some fishes today or established a novel trophic unit with no living analogue. Recalibrating the Paja biota network to trophic analogues in the Caribbean, such as sea turtles, estimates that the largest Paja marine reptile hyper-apex predators occupied trophic levels a full tier higher than any extant marine apex predator. The Paja biota network is a starting point to tracing the evolution of marine ecosystems across the Mesozoic Marine Revolution.
INTRODUCTION
Extant biodiversity exhibits several patterns including latitudinal diversity gradients and regional hotspots (Pianka 1966, Myers et al. 2000, Hillebrand 2004, Willis et al. 2006, Cermeño et al. 2022). Yet the origin of these patterns is not obvious and palaeontological investigations have attempted to trace the deep time history of current biodiversity phenomena. The evolution of Phanerozoic biodiversity through time reveals several episodes of species originations and extinctions, and a generally increasing trend of biodiversity throughout the Mesozoic and Cenozoic (Sepkoski 1981, Benton 1995). Sampling and geological biases may account for some of these trends (Alroy et al. 2001). However, an upward trend in global marine biodiversity is still recovered from the Jurassic to the Present (Alroy et al. 2008). This pattern has been called the Mesozoic Marine Revolution (MMR) and its conceptual history began with the observation that this was a time of origination and diversification of several marine prey and predators and their assumed interactions (Vermeij 1977, 1987).
Marine food-webs in the Mesozoic oceans developed low trophic plankton systems and several higher trophic predators (Vermeij 1977, Walker and Brett 2002, Bardhan and Chattopadhyay 2003). The achievement of durophagy, in particular, is postulated to have enabled an arms’ race of marine fauna that may have shaped the marine faunas during the Mesozoic (van Valen 1973, Vermeij 1977, Walker and Brett 2002, Harper 2003, Buatois et al. 2016). Vermeij (2008) found evidence of the MMR beginning by the Late Triassic and recent studies have even suggested the MMR began as early as the early Late Triassic, following the Carnian Pluvial Episode (Benton and Wu 2022). It may be also plausible that these ecosystem changes increased predator abundance and, in turn, augmented primary production, termed the theory of escalation (Vermeij 1987, 2013, Knoll and Follows 2016). It should be noted that most of what has been studied comes from invertebrate faunas that dominate lower trophic levels, and thus the evolution of higher trophic levels has not yet been explored. During the Mesozoic, marine tetrapods, for example, found an opportunity to explore a wide range of morphological innovations, involving a diversity of bauplans, body sizes, ecological niches, appearance of ecological guilds not seen since the Cambrian, and also the origin of higher trophic levels (Kelley and Pyenson 2015, Stubbs and Benton 2016).
Biodiversity was driven in part by Continental Drift (involving tropical settings), sea-level change, and consequent biotic factors such as predation, competition, reproductive isolation, and ecospace evolution, among others (Hallam 1967, Heard and Hauser 1995, Barash 2008, Bush et al. 2016, Reeves et al. 2021). Temperature, precipitation, and area inversely correlate with latitude on land and are considered primary drivers of the terrestrial latitudinal biodiversity gradients prevalent today (Gaston 2000, Saupe et al. 2019). Temperature may be the most important driver in extant marine latitudinal biodiversity gradients, a gradient that has been present since the Early Paleozoic (Alroy et al. 2008). Sea levels, coastline extent, and temperatures have varied over time and their rises are associated with increased marine biodiversity (Littler et al. 2011, Föllmi 2012, Benson and Druckenmiller 2014, Haq 2014, Tennant et al. 2016, Vavrek 2016, Muscente et al. 2018). The fragmentation of Pangaea, for example, during the Mesozoic created novel continental and shoreline habitats that are expected to have concurrent increases in biodiversity (Vavrek 2016).
The physical and biotic factors listed above all took place during, and may have facilitated, the MMR. The MMR is assumed to have increased the energy budgets of marine systems (Finnegan et al. 2011) and might have also prompted great diversification (and rapid origination) of several marine prey and predatory taxa (Vermeij 1977, 1987). These prior studies have relied on diversity estimates and appearance of novel feeding guilds. However, none have examined how ecosystem-level patterns may have changed throughout the MMR. Until more site-specific ecological networks are produced, we can only speculate as to the tempo and complexity of the MMR through the Mesozoic.
A promising approach to examining ancient biodiversity is the use of quantitative ecological networks. These networks attempt to reconstruct energy flow throughout an ecosystem by approximating feeding interactions. Thus, ecosystem networks capture more than raw species richness (Whittaker 1972) but may also approximate many of their interactions, and ultimately energy distribution within ecosystems (Damuth 1981, Cohen et al. 2003, Hägen et al. 2012, Hatton et al. 2015, 2019, 2021, Gravel et al. 2019). They may also be used to reconstruct site-specific ecologies rather than time-averaging large deposits over broad geographic scales.
Quantitative ecological networks can also be used to calculate several parameters specific to each ecosystem. These include degrees of network connectance, modularity, node vulnerability, and trophic levels (Williams and Martinez 2000). Trophic levels, for instance, are important ecosystem properties because their numbers are direct indicators of ecosystem complexity (Dunne et al. 2004). Some of the highest trophic levels present in today’s marine ecosystems are occupied by large apex predators, such as orcas, sperm whales, and white sharks. These species are estimated at trophic levels of up to 5.6 (see below). Interestingly, some predatory invertebrates in Caribbean reefs even achieve this trophic level (Roopnarine and Dineen 2018).
Ecosystem networks describe the trophic interactions within a community, in which the taxa of the food web are represented as nodes, and feeding interactions are represented by links that connect resources (prey) to consumers (predators) (Roopnarine 2010, Delmas et al. 2019). Food webs reflect fluxes of biomass along with species compositions of communities, and provide useful data on ecological patterns over evolutionary timescales (Dunne et al. 2008, 2014, O’Connor et al. 2020).
Food webs have been primarily used to study extant communities (e.g. Yodzis 1998, Poisot et al. 2016, Muscente et al. 2018, Gravel et al. 2019, O’Connor et al. 2020, Tanentzap et al. 2020, Caron et al. 2022). However, recent work has applied these approaches to palaeontological data (Roopnarine et al. 2007, Caron and Jackson 2008, Dunne et al. 2008, 2014, Roopnarine and Angielczyk 2015, Kempf et al. 2020, Fricke et al. 2022). This is particularly relevant because the fossil record is an incomparable source for ecosystem-level food-web evolution across large timescales (Yeakel and Dunne 2015, Shaw et al. 2021).
Fossil food-web networks are useful in that: (i) since the fossil record is incomplete, a network may assist in ‘filling in the gaps’; (ii) energetic models can be estimated for these ancient systems and require a network model in place first; (iii) they offer a novel higher level (big picture) overview of palaeoecological systems; (iv) network evolution has never been addressed and offers a novel avenue of research comparing particularly well-sampled fossil beds; and (v) can be used to explicitly test ecosystem changes throughout the MMR, including large-scale biodiversity trends and ecosystem complexity.
In extant systems, a food web can be built by direct and indirect evidence, although even in these cases details of diets are difficult to assess (Nielsen et al. 2018). Taxon interactions in an ecosystem are made by a resource–consumer pairwise list (e.g. Dunne et al. 2008, 2014). Direct interaction evidence between taxa is rare in the fossil record. Indirect evidence of ancient trophic interactions can be inferred based on attributes such as feeding morphologies and body size. Some inferred trophic interactions may be based on evidence such as predation/scavenging marks, gut contents, fossil traces, among others, or by using extant analogues (e.g. Kelley et al. 2003, Fröbisch et al. 2013, Roopnarine and Dineen 2018; Cortés et al. 2019a; Klompmaker et al. 2019, Jiang et al. 2020, and references therein). Building ancient food-webs is challenging since the process relies primarily on indirect evidence of trophic interactions, biases in taxonomic representation in the fossil record, time-averaging, as well as a lack of a standardized framework for building and comparing different food-webs, since many incorporate varied methods (Kelley et al. 2003; Poisot et al. 2016, Delmas et al. 2019).
Ecological processes and ecological patterns can be studied throughout quantitative analyses of food webs (Thompson et al. 2012, Delmas et al. 2019). These can be assessed by considering multiple ecological ‘rules’ first. For example, generally, a predator consumes smaller prey (Cohen et al. 2003, Caron et al. 2022), abundances and taxonomic diversity of prey are greater than those of predators, and abundances of smaller-bodied taxa are higher than those of larger-bodied taxa (e.g. Briand and Cohen 1984, Cohen et al. 1993).
Quantitative ecological networks have had some success when applied to paleocommunities. Food webs of ancient systems have been analysed by empirical data (niche model) (e.g. Dunne et al. 2008, 2014), theoretical modelling (e.g. Solé et al. 2002), and data and modelling combined (Roopnarine and Dineen 2018).
Moderately complete marine networks have been reconstructed from several konservat-lagerstatten. These include the Cambrian Burgess Shale and Chengjiang Biotas, Ordovician Cincinnati Basin Biota, the Permo-Triassic Karoo Basin sequence, and Eocene Messel Biota (Roopnarine et al. 2007, Dunne et al. 2008, 2014, Roopnarine and Angielczyk 2015, Kempf et al. 2020). Some attempts have been made to examine limitations of fossil food-web models. Roopnarine and Dineen (2018) compared Caribbean coral reef communities from Jamaica to simulate fossil reefs to test the reliability of fossil food-web reconstructions and found that overall network structure, modularity, and guild diversity remained stable. Others have addressed sampling efforts and taphonomic biases while attempting to reconstruct food network models (Martinez et al. 1999, Shaw et al. 2021). A few of these fossil ecosystem networks have recovered some community network changes around mass extinctions, such as changes in ecosystem dynamics, disruptions of primary productivity, or network stability that assisted their survival across the mass extinction (Roopnarine et al. 2007, Roopnarine and Angielczyk 2015, Kempf et al. 2020). Additional Pleistocene terrestrial networks have been reconstructed and compared to contemporaneous networks for several fossiliferous localities demonstrating the impact of human-mediated megafaunal extinctions on mammalian ecosystems (Fricke et al. 2022). Banker and colleagues (Banker et al. 2022) recently presented a series of marine food-webs constructed from stage-wide ensembles of Western Tethys fossils in Carnian, Anisian, Bathonian, and Aptian deposits of Europe. Their results showed that trophic space generally expanded throughout the Mesozoic and that the Bathonian networks were dominated by longer food chains. The large temporal and spatial binning of Banker and colleagues (2022) may make direct inferences about their ecological networks difficult to interpret; however, the general results are promising. Qualitative food-webs have also been reconstructed for several Mesozoic marine deposits, including the Triassic (e.g. Hu et al. 2011, Lukeneder and Lukeneder 2021). Here we aim to apply current methods of quantitative ecosystem network modelling to the Early Cretaceous Paja Formation biota of Colombia. This is the first attempt at reconstructing a quantitative ecological network for any site-specific Mesozoic marine community. Our aim is to begin this food web to initiate discussions on how quantitative food-webs can be used to infer potential missing taxa, missing trophic levels, evolution of trophic complexity, and quantifying ecological changes throughout the MMR.
Geological setting
Geology
The Paja Formation is an Early Cretaceous sequence that outcrops extensively around Villa de Leyva (Boyacá, central Colombia) (Fig. 1). It consists of black to variegated shales, calcareous concretions, biomicrites, and abundant marine, and a few terrestrial, fossils (Etayo-Serna 1968, 1979, Gaona-Narvaez, Maurrasse, and Etayo-Serna 2013, Carballido et al. 2015, Cortés et al., 2023). The Paja Formation in the region of Villa de Leyva ranges from Hauterivian to the Aptian in age, overlying the Ritoque Formation (Valanginian–Hauterivian) and underlying the San Gil Inferior Formation (Aptian–Albian) (Etayo-Serna 1968, 1979, Patarroyo 2020) (Fig. 2A). The Paja Formation consists of three members. The lower Lutitas Negras Inferiores member is Hauterivian in age and consists of black shale that weathers to a reddish colour, with sandy interbeds representing the feather edge of a thick sandstone body. The Arcillolitas Abigarradas member, which is the middle of the three members, is Lower Barremian–Upper Aptian in age, thus spanning nearly 10 million years. It has five units (A–E) (Fig. 2A), and consists of light-coloured variegated claystone, which in its lower portion contains numerous thin beds of sandy clay, and in its upper portion contains intercalations of calcareous claystone and gypsum. The Arcillolitas Abigarradas member is the most fossiliferous of the three members and includes the highest taxonomic diversity. All of the vertebrate fossils of the Paja Formation biota are found in this member. Unfortunately, these fossils have little stratigraphic control within this member, and future work is needed for tighter stratigraphic controls of these localities. The upper member, the Arcillolitas con Nódulos Huecos member, Aptian in age, consists of claystone with hollow nodules (Etayo-Serna 1968, 1979).

Paja Formation geography and geological setting. A, location of Villa de Leyva within the Departamento de Boyacá and Colombia. B, location of Colombia within South America. C, lithological map of the Paja Formation, indicating the localities of the taxa found in the sequence (approximate location) across the Paja Formation’s members: Arcillolitas con Nódulos Huecos (An), Arcillolitas Abigarradas (Aa), and Lutitas Negras Inferiores (Ln). Lithological map based on Etayo-Serna (1968, 1979). Silhouettes of the taxa modified from Cortés et al. (2019b) and phylopic.org.

A, stratigraphic column of the Paja Formation. The Arcillolitas Abigarradas member, and its A–E units are highlighted (modified from Etayo-Serna 1968, 1979, Páramo-Fonseca et al. 2019). B, histogram plots of all ammonites known from these A–E units showing ammonite richness and body size distribution (log shell diameter against frequency).
The Paja Formation and its faunas
The Paja Formation was deposited over 15 million years throughout a sequence spanning the Hauterivian through to the Aptian in central Colombia (Etayo-Serna 1968, 1979). The Early Cretaceous was a transitional period in Earth’s history and of significant importance in terms of faunistic turnovers (Benson et al. 2013, Benson and Druckenmiller 2014, Tennant et al. 2017). This time represents the recovery at the end Jurassic mass extinction, high eustatic sea-levels, continental fragmentation, and high global temperatures—all factors that are expected to facilitate high levels of phyletic radiations (Littler et al. 2011, Föllmi 2012, Benson and Druckenmiller 2014, Haq 2014, Tennant et al. 2016, Vavrek 2016, Muscente et al. 2018). The Paja Formation preserves a rich marine fauna and outcrops extensively around Villa de Leyva, Boyacá. This formation is subdivided into several members, including the highly fossiliferous Arcillolitas Abigarradas member (Etayo-Serna 1968, 1979). This member preserves several species of marine invertebrates and vertebrates, which are here compiled in a single dataset for the first time. The invertebrate fauna consists of over 100 species of ammonites and several bivalves and crabs (Etayo-Serna 1968, 1979, Patarroyo 2000, 2009, 2020, Hoedemaeker 2004, Kakabadze et al. 2004, Luque et al. 2020). The vertebrate fauna includes massive short-necked plesiosaurs, large elasmosaurs, a teleosauroid crocodylomorph, several species of ichthyosaurs, several sea turtles, and at least four morphotypes of actinopterygian fish and one shark (Hampe 1992, 2005, Schultze and Stöhr 1996, Cadena 2015, Cadena and Parham 2015, Maxwell et al. 2016, Páramo-Fonseca et al. 2016, 2018, 2019, Cortés et al. 2019b, 2021, Carrillo-Briceño et al. 2019, this study). The ecological complexity of this palaeoecosystem hints at what this ecological network may look like. Several apex predators with ~10 m body lengths include the pliosaurs Monquirasaurus (=Kronosaurus) boyacensis and Sachicasaurus vitae, and a teleosaur (Hampe 1992, Páramo-Fonseca et al. 2018, Cortés et al. 2019b, Noè and Gómez-Pérez 2022). Mid-sized predators of >4 m include the pliosaurs Acostasaurus pavachoquensis and Stenorhynchosaurus munozi (=Brachauchenius sp.), and the ichthyosaur Kyhytysuka sachicarum (Hampe 2005, Páramo-Fonseca et al. 2016, Gómez-Pérez and Noè 2017, Cortés et al. 2021). Other predators include the large-bodied but small-headed elasmosaurs and stem sea turtles Desmatochelys padillai, and Leyvachelys cipadi (Cadena 2015, Cadena and Parham 2015).
Although the Paja Formation biota is most famous for its large marine reptiles, relatively little is known about the broad palaeoecological structure of the ecosystem, including its fishes and ammonites. Our goal with this ecological network is to add it to the short list of site-specific Phanerozoic ecosystem networks to help assess the evolution of marine ecosystems, highlight network properties of a marine ecosystem with exceptionally large predators, predict what trophic players are missing, and, ultimately, generate energy-flow models.
MATERIAL AND METHODS
Master taxon list:
To assemble the Paja Formation food-web datasets, we first compiled a list of all described and undescribed specimens known from the Paja Formation biota (Supporting Information, Table S1). We focused on specimens recovered in the region of Villa de Leyva because of the lengthy history of fossil exploration in the region (see major fossil localities in Fig. 1). A complete list of the specimens of the Centro de Investigaciones Paleontológicas (CIP) was provided by the institute (last update March 2022), and the data for collections of smaller museums were taken during visits and/or from the literature. We used the taxonomic referrals of the studied taxa from published literature and assessed the rest of the specimens individually for taxonomic identification. This is the first and most complete taxonomic approximation ever made for the Paja Formation biota.
The master list first included taxonomic referral, age, formation information (when available), and reference. We then estimated anatomical attributes and measurements either from first-hand data or from the literature. These were used to estimate body size, body mass, feeding attributes such as tooth guild (based on Massare 1987), tooth sizes, and jaw length, and postcranial attributes (if applicable). With these as evidence for trophic roles, we hypothesized diets, trophic roles, and trophic positions for all taxa. The vertebrates are composed of sauropterygians, ichthyosaurians, thalattosuchians, testudines, pycnodontiformes, aspidorhynchiformes, and lamniformes, and the invertebrates are composed of ammonites, belemnitides, classiduloids, decapods, littorinimorphs, lunicids, ostreids, pteriids, and trigonnids. The presence of echinoids, corals, crinoids, bryozoans, and forams is uncertain or completely unexplored and/or undersampled in the Paja Formation biota to date. The list of taxa is available in the Supporting Information, Table S2.
Food-web data matrix:
The Paja Formation food-web data (edgelist) was built following the approach by Dunne et al. (2008) that included consumer ID and resource ID (Supporting Information, Table S3). The base of the food web was taken from Dunne et al. (2008) hypothesizing the structure of the base of a healthy marine ecosystem, including phytoplankton, bacterioplankton, suspended organic matter, benthic detritus, algae, and zooplankton. We do not consider jellyfish in the analysis, although they were probably present. Taxon interactions are based on trait-based inferences with supporting comments for all. Trait-based inferences of taxon interactions are perhaps the most defensible for reconstructing palaeocommunity interactions. Ecological network analyses on extant communities have confirmed that feeding traits and body sizes are the best predictors of food webs (Caron et al. 2022). Exceptional konservat-lagerstatten yield numerous trait-based and some direct evidence, including gut contents, bite marks, and tooth wear (e.g. Carpenter 1998, Kelley et al. 2003, Fröbisch et al. 2013, Hone et al. 2015, Cortés et al. 2019a, Jiang et al. 2020).
A trait and body mass approach offers less reliance on preservation, and more reliance on anatomy and body size estimated from fossils, thus applying to a greater diversity of fossil localities and communities, and allows trophic stability analyses to inform how realistic body-size estimates are. A disadvantage of this approach is that more assumptions are needed about trophic compositions. Another limitation is that estimating body sizes and masses generally requires relatively complete specimens or preserved elements that correlate well to body sizes of related taxa.
We decided not to use the certainty scale proposed by Dunne et al. (2008) given the ambiguity of this number in terms of missing direct evidence in the Paja Formation biota. Instead, we established feeding rules among the taxa. For example, the only known feeding apparatus for ammonites and belemnites suggest that they were filter feeders of zooplankton (Kruta et al. 2011). The Colombian teleosauroid, on the other hand, of about 9.6 m in length, was assumed to be a generalist feeder, feeding on benthic animals larger than 10 cm, pelagic animals larger than 10 cm, but nothing larger than 1 m, since the longirostrine skull and piercing teeth in this clade are associated with piscivorous diets (Pierce et al. 2009, Ballell et al. 2019).
Ammonites are the most abundant fossil in the Paja Formation biota. In total, they represent 112 species distributed throughout the formation. The marine vertebrate fossils, however, are recovered only from the Arcillolitas Abigarradas member. Ammonites from this member have been identified from five discrete zones (A–E) (Etayo-Serna 1968, 1979). Using only ammonites from this member, we plotted the diversity and size variation to choose the zone that might best represent the ammonite community with as little time-averaging as possible. Unit A ammonites present the most diverse record with a relatively normal distribution of body sizes (Fig. 2B). We assume the ammonites from unit A represent the most completely sampled unit and the other units are either undersampled or have a greater taphonomic bias against ammonites. However, differences between units could also be due to lower and fluctuating ammonite biodiversity throughout this interval. Unfortunately, other workers have described several other ammonites present during this time interval but did not locate them to the specific stratigraphic units proposed by Etayo-Serna (1968, 1979), thus potentially under-representing those units.
The marine vertebrate diversity is restricted to the Arcillolitas Abigarradas member but has not been tied to any particular unit within the member. Thus, in an attempt to not over-represent ammonite diversity, we used only the ammonites known from unit A in comparison to all vertebrates recovered from the Arcillolitas Abigarradas member.
Extant marine system dataset:
We used, with some corrections, an edgelist of three large and well-studied extant marine ecosystems in northern Caribbean coral reefs: Cayman Islands, Cuba, and Jamaica (Roopnarine and Hertog 2012; Supporting Information, Table S4). These datasets provide the most complete tropical reef marine food-webs and were used as comparisons to the Paja Formation food-web. We also compiled the published data of Burgess Shale, Chengjian Shale, Messel Shale (lake), and Benguela systems, and re-ran all these matrices (Table 1) as controls for the estimated network properties stats (Yodzis 1998, Dunne et al. 2008, 2014).
Network values from the Paja Formation, Cayman Islands, Cuba, and Jamaica marine food web networks. S, species, the total number of nodes (species); L, total number of trophic links, as the size of the network; C, connectance; maxPATL, maximum prey-averaged trophic level (PATL), which assumes a predator’s diet fractions are equal across the prey (Williams and Martinez 2008); TL, trophic level; SWTL, short-weighted trophic level; meanSWTL, mean short-weighted trophic level (Dunne et al. 2008, Delmas et al. 2019, and references therein). The network statistics were calculated by using Network3D software (Yoon et al. 2004, Williams 2010). The properties output for the Paja Formation biota and the Caribbean systems are available in the Supporting Information, Table S5
Network . | S . | L . | LD . | C . | meanSWTL . | maxSWTL . | maxPATL . |
---|---|---|---|---|---|---|---|
Paja Formation biota | 74 | 977 | 13.2 | 0.177 | 2.95 | 5.12 Brachypterygiidae morphospecies 2 | 5.57 Monquirasaurus boyacensis, Sachicasaurus vitae |
Cayman Islands | 243 | 3905 | 16.1 | 0.066 | 3.93 | 5.69 blenny | 6.38 blenny |
Cuba | 233 | 3751 | 16.1 | 0.069 | 3.89 | 5.66 Gastropods with polychaete prey | 6.32 Gastropods with polychaete prey |
Jamaica | 248 | 4165 | 16.8 | 0.068 | 3.91 | 5.66 Gastropods with polychaete prey | 6.32 Gastropods with polychaete prey |
Network . | S . | L . | LD . | C . | meanSWTL . | maxSWTL . | maxPATL . |
---|---|---|---|---|---|---|---|
Paja Formation biota | 74 | 977 | 13.2 | 0.177 | 2.95 | 5.12 Brachypterygiidae morphospecies 2 | 5.57 Monquirasaurus boyacensis, Sachicasaurus vitae |
Cayman Islands | 243 | 3905 | 16.1 | 0.066 | 3.93 | 5.69 blenny | 6.38 blenny |
Cuba | 233 | 3751 | 16.1 | 0.069 | 3.89 | 5.66 Gastropods with polychaete prey | 6.32 Gastropods with polychaete prey |
Jamaica | 248 | 4165 | 16.8 | 0.068 | 3.91 | 5.66 Gastropods with polychaete prey | 6.32 Gastropods with polychaete prey |
Network values from the Paja Formation, Cayman Islands, Cuba, and Jamaica marine food web networks. S, species, the total number of nodes (species); L, total number of trophic links, as the size of the network; C, connectance; maxPATL, maximum prey-averaged trophic level (PATL), which assumes a predator’s diet fractions are equal across the prey (Williams and Martinez 2008); TL, trophic level; SWTL, short-weighted trophic level; meanSWTL, mean short-weighted trophic level (Dunne et al. 2008, Delmas et al. 2019, and references therein). The network statistics were calculated by using Network3D software (Yoon et al. 2004, Williams 2010). The properties output for the Paja Formation biota and the Caribbean systems are available in the Supporting Information, Table S5
Network . | S . | L . | LD . | C . | meanSWTL . | maxSWTL . | maxPATL . |
---|---|---|---|---|---|---|---|
Paja Formation biota | 74 | 977 | 13.2 | 0.177 | 2.95 | 5.12 Brachypterygiidae morphospecies 2 | 5.57 Monquirasaurus boyacensis, Sachicasaurus vitae |
Cayman Islands | 243 | 3905 | 16.1 | 0.066 | 3.93 | 5.69 blenny | 6.38 blenny |
Cuba | 233 | 3751 | 16.1 | 0.069 | 3.89 | 5.66 Gastropods with polychaete prey | 6.32 Gastropods with polychaete prey |
Jamaica | 248 | 4165 | 16.8 | 0.068 | 3.91 | 5.66 Gastropods with polychaete prey | 6.32 Gastropods with polychaete prey |
Network . | S . | L . | LD . | C . | meanSWTL . | maxSWTL . | maxPATL . |
---|---|---|---|---|---|---|---|
Paja Formation biota | 74 | 977 | 13.2 | 0.177 | 2.95 | 5.12 Brachypterygiidae morphospecies 2 | 5.57 Monquirasaurus boyacensis, Sachicasaurus vitae |
Cayman Islands | 243 | 3905 | 16.1 | 0.066 | 3.93 | 5.69 blenny | 6.38 blenny |
Cuba | 233 | 3751 | 16.1 | 0.069 | 3.89 | 5.66 Gastropods with polychaete prey | 6.32 Gastropods with polychaete prey |
Jamaica | 248 | 4165 | 16.8 | 0.068 | 3.91 | 5.66 Gastropods with polychaete prey | 6.32 Gastropods with polychaete prey |
Building the food-web network:
We used Network3D software that has visualization tools for ecological networks, modelling capabilities, and basic analyses for food webs (Yoon et al. 2004, Williams 2010). Network3D software was written by R.J. Williams and provided by J. Dunne (Santa Fe Institute).
Food-web properties:
Trophic network models are quantitatively interpreted following a series of properties or measurements, including: S, species, as the total number of nodes (species); L, total number of trophic links, as the size of the network; C, connectance; PATL, prey-averaged trophic level, which assumes a predator’s diet fractions are equal across the prey (Williams and Martinez 2008); TL, trophic level; SWTL, short-weighted trophic level; meanSWTL, mean short-weighted trophic level (Dunne et al. 2008, Delmas et al. 2019, and references therein). The network statistics were calculated by using Network3D software (Yoon et al. 2004, Williams 2010). The properties’ output for the Paja Formation biota and the Caribbean systems are available in the Supporting Information, Table S5.
RESULTS
Diversity:
In total, the Paja Formation biota has 157 species, including the six hypothetical base trophic levels of the chain (Supporting Information, Table S2). Of the 151 preserved fossil taxa, 73% (115) are identified to species level, 10% (16) only to genus, and the remaining 17% (26) only to family or any other higher taxonomic level. No marine plants or soft-bodied invertebrates are known from this formation. From the total list, ~75% of species correspond to cephalopods (including 65 ammonites from the Arcillolitas Abigarradas member, 29 only from unit A), over 7% species correspond to the rest of marine invertebrates, and ~15% species to marine reptiles, and almost 3% correspond to ray-finned fish (Fig. 3). The ammonite histograms and species diversity plots (Fig. 2B) served as a proxy for the ammonite community used to construct the ecological network. Using only ammonites from unit A brings the relative diversities to 45% of segment A, 26% of segment C, 6% of segment D, and 23% of segment E, respectively.

Pie charts of species diversity of the Paja Formation biota. A, Class–Order distribution. B, Order–Family distribution. Colours in (A) do not relate with colours in (B). Plots made in Rstudio (2022). Silhouettes of the taxa modified from Cortés et al. (2019b) and phylopic.org.
From this working taxonomic list, trophic levels are strongly biased towards predatory guilds (Fig. 4A). Apex predators make up about 3% of the biodiversity, whereas the first-level carnivores make up about 80% of total trophic roles; 16% of the taxa make up three other carnivorous consumer levels, spread across mid-level carnivores, and second- and third-level carnivorous consumers. Only about 1.3% of taxa comprise herbivorous consumers (Fig. 4).

Pie charts representing: A, trophic role distribution. B, trophic guild distribution by role. C, tooth guild distribution. Colours in (A) do not relate with colours in (B) and (C). Plots made in Rstudio (2022).
Tooth and trophic guild distributions of the reptilian carnivores are presented in Figure 4 (B, C). These distributions indicate a complex ecosystem dominated by piercing, clutching, and cutting dentitions. Notably absent are specialized crushing dentitions. The only taxon with this type of dentition is a rare, undescribed pycnodontiform fish (Supporting Information, Tables S1, S2). The full range of reptilian carnivore tooth guilds is present in the mid-level carnivores, which make up about half of the reptiles in the Paja Formation biota. Second-level reptile carnivores have only piercing dentition and third-level reptile carnivores have piercing, cutting, and clutch/tearing dentitions. The apex predators comprise about 10% of the total reptile diversity and present piercing, cutting, and crushing dentitions.
Food-web network:
The model representing the Paja Formation biota’s food-web network is shown in Figure 5. Figure 5A presents the network as an evenly distributed 3D network and Figure 5B as a flattened network. The network is composed of the same 74 taxa as above with 977 links. This yields a linkage density (links per species) of 13.2 and overall network connectance (C) of 0.177 (Table 1). The mean short-weighted trophic level (meanSWTL) for the Paja network is 2.95. The maximum short-weighted trophic level is 5.12 at Brachypterygiidae morphospecies 2 and the longest prey-averaged trophic level (PATL) is 5.57 at taxa Monquirasaurus (=Kronosaurus) boyacensis and Sachicasaurus vitae. We found that 2.7% of the taxa are apex predators consumed by no other taxa, 90.5% of taxa are both consumers and prey, and only 6.8% of taxa are basal trophic producers.

Paja Formation network. A, cylindrical projection, and B, rectangular projection. Taxa are represented by spheres that grade from red at the bottom to yellow at the top. The links among the taxa are represented by connecting lines in green. The vertical position of the taxa represents the approximate trophic levels. The numbering over the spheres corresponds to the taxa IDs (see Supporting Information, Tables S1, S2). The networks were computed in Network3D software written by R.J. Williams. See Yoon et al. (2004), Williams (2010), and Dunne et al. (2008) for details about the software’s algorithm, interpretations, and outputs.
To make comparisons to extant tropical marine ecosystems, three well-resolved marine food-webs were used from Roopnarine and Hertog (2012). These ecosystems were used because they are some of the most complete marine systems documented to construct food webs. They are also shallow and tropical, and expected to be good extant analogues to the Paja network. These systems reveal much more resolved food-webs (Fig. 6). Tables 1 and 2 present network values across the Paja and extant Caribbean networks. Some general results from the extant Caribbean networks are the higher number of taxa and links. They average about 241 species, 3940 links, a linkage density of 16.3, C of 0.068, meanSWTL of 3.91, and PATL of 5.67. About 4% of taxa are apex predators, 0.8% are basal producers, and the remaining are both consumers and prey. In the Cayman Islands system, the blenny reveals the highest SWTL of 5.69 and the highest PATL of 6.38. The oceanic whitetip shark has a lower SWTL of 4.95 and a PATL of 5.92. The Caribbean Reef shark has a SWTL of 3.78 and a PATL of 5.57. The barracuda has a SWTL of 5.52 and PATL of 6.05. Although the sharks and barracuda are apex predators, the blenny’s food-web is more complex with several benthic invertebrate interactions, leading it to rank at higher trophic levels. In the Cuba and Jamaica networks, the maxSWTL and maxPATL are both held by gastropods with polychaete prey, with the same values of 5.66 and 6.32, respectively. This is due to the diverse benthic diet of their prey. The maxSWTL and maxPATL for the apex vertebrate predators of these systems are 5.52 and 6.05, respectively, in the Cayman and Cuba networks, and 5.50 and 6.32 in the Jamaica network. The maximum SWTL and PATL values correspond to the barracuda in all the Caribbean networks (Table 2; Supporting Information, Table S5).
Comparisons of trophic level values across analogous taxa from the marine food-webs. Trophic levels are from the predator averaged trophic level (PATL) values. Mean values are calculated when multiple taxa are from the same network
Network . | Sea turtles . | Ray-finned fish . | Chondrichthyes . | Crustaceans . | Cephalopods . |
---|---|---|---|---|---|
Paja Formation biota | 4.76 Desmatochelys padillai, Leyvachelys cipadi | 5.00 Actinopterygii morphospecies 2 | 4.94 Protolamna ricaurtensis | 3.5 Bellcarsinus aptiensis, Planocarnicus olssoni | 4.00 Pulchellia leivaensis, Karsteniceras beyrichii |
Cayman Islands | 5.97 Loggerhead sea turtle | 6.05 Barracuda | 5.92 Oceanic whitetip shark | 4.54 Omnivorous crustacea | 5.13 Caribbean reef squid |
Cuba | 5.97 Loggerhead sea turtle | 6.05 Barracuda | 5.65 Blacktip shark | 4.54 Omnivorous crustacea | 5.15 Caribbean reef squid |
Jamaica | 5.96 Loggerhead sea turtle | 6.00 Barracuda | 5.88 Oceanic whitetip shark | 4.54 Omnivorous crustacea | 5.17 Caribbean reef squid |
Network . | Sea turtles . | Ray-finned fish . | Chondrichthyes . | Crustaceans . | Cephalopods . |
---|---|---|---|---|---|
Paja Formation biota | 4.76 Desmatochelys padillai, Leyvachelys cipadi | 5.00 Actinopterygii morphospecies 2 | 4.94 Protolamna ricaurtensis | 3.5 Bellcarsinus aptiensis, Planocarnicus olssoni | 4.00 Pulchellia leivaensis, Karsteniceras beyrichii |
Cayman Islands | 5.97 Loggerhead sea turtle | 6.05 Barracuda | 5.92 Oceanic whitetip shark | 4.54 Omnivorous crustacea | 5.13 Caribbean reef squid |
Cuba | 5.97 Loggerhead sea turtle | 6.05 Barracuda | 5.65 Blacktip shark | 4.54 Omnivorous crustacea | 5.15 Caribbean reef squid |
Jamaica | 5.96 Loggerhead sea turtle | 6.00 Barracuda | 5.88 Oceanic whitetip shark | 4.54 Omnivorous crustacea | 5.17 Caribbean reef squid |
Comparisons of trophic level values across analogous taxa from the marine food-webs. Trophic levels are from the predator averaged trophic level (PATL) values. Mean values are calculated when multiple taxa are from the same network
Network . | Sea turtles . | Ray-finned fish . | Chondrichthyes . | Crustaceans . | Cephalopods . |
---|---|---|---|---|---|
Paja Formation biota | 4.76 Desmatochelys padillai, Leyvachelys cipadi | 5.00 Actinopterygii morphospecies 2 | 4.94 Protolamna ricaurtensis | 3.5 Bellcarsinus aptiensis, Planocarnicus olssoni | 4.00 Pulchellia leivaensis, Karsteniceras beyrichii |
Cayman Islands | 5.97 Loggerhead sea turtle | 6.05 Barracuda | 5.92 Oceanic whitetip shark | 4.54 Omnivorous crustacea | 5.13 Caribbean reef squid |
Cuba | 5.97 Loggerhead sea turtle | 6.05 Barracuda | 5.65 Blacktip shark | 4.54 Omnivorous crustacea | 5.15 Caribbean reef squid |
Jamaica | 5.96 Loggerhead sea turtle | 6.00 Barracuda | 5.88 Oceanic whitetip shark | 4.54 Omnivorous crustacea | 5.17 Caribbean reef squid |
Network . | Sea turtles . | Ray-finned fish . | Chondrichthyes . | Crustaceans . | Cephalopods . |
---|---|---|---|---|---|
Paja Formation biota | 4.76 Desmatochelys padillai, Leyvachelys cipadi | 5.00 Actinopterygii morphospecies 2 | 4.94 Protolamna ricaurtensis | 3.5 Bellcarsinus aptiensis, Planocarnicus olssoni | 4.00 Pulchellia leivaensis, Karsteniceras beyrichii |
Cayman Islands | 5.97 Loggerhead sea turtle | 6.05 Barracuda | 5.92 Oceanic whitetip shark | 4.54 Omnivorous crustacea | 5.13 Caribbean reef squid |
Cuba | 5.97 Loggerhead sea turtle | 6.05 Barracuda | 5.65 Blacktip shark | 4.54 Omnivorous crustacea | 5.15 Caribbean reef squid |
Jamaica | 5.96 Loggerhead sea turtle | 6.00 Barracuda | 5.88 Oceanic whitetip shark | 4.54 Omnivorous crustacea | 5.17 Caribbean reef squid |

A, original trophic species network models of the Paja Formation biota, Cayman Islands, Cuba, and Jamaica systems plotted along (B) frequency of short-weighted trophic level (SWTL). The ammonites in the Paja Formation network have the highest value and it is highlighted in the histogram in yellow to show this distribution.
The trophic analogue for the Caribbean barracuda in the Paja network is the Actinopterygii ms 2, which has a SWTL and a PATL of 5.00. The Cayman network has SWTL and PATL for the oceanic whitetip shark of 4.96 and 5.92, respectively, and the Jamaica network of 4.94 and 5.88, respectively. The oceanic whitetip shark is absent in the Cuba network but the blacktip shark with values of 4.83 and 5.65, respectively, is the highest in these predators for Cuba. The values for the trophic analogue in the Paja network correspond to the only known shark for this system, Protolamna ricaurtensis with a SWTL of 4.47 and a PATL of 4.94. The marine turtles are placed higher in the Caribbean networks due to their diverse benthic prey. The loggerhead sea turtles in the Cayman and Cuba networks have a SWTL of 5.49 and a PATL of 5.97, and the Jamaica network has a SWTL of 5.48 and PATL of 5.96. These values contrast with the known marine turtles in the Paja network, Desmatochelys padillai and Leyvachelys cipadi, which both have a maxSWTL of 3.88 and a PATL of 4.76, indicating almost two trophic levels of difference between these taxa in the Paja and the Caribbean networks (Table 2; Supporting Information, Table S5).
Gastropods, crustaceans, and cephalopods are shared in the Caribbean and Paja networks. The only known gastropod in the Paja network has a SWTL 3.08 and a PATL 3.17. Extant Caribbean gastropods range from SWTL and PATL of 3.0 for herbivorous gastropods to 3.6 and 4.2, respectively, for omnivorous gastropods, and up to 5.7 and 6.3, respectively, for gastropods that are polychaete predators. The latter is, in fact, the invertebrate with the highest trophic levels in these Caribbean systems (Tables 1, 2; Supporting Information, Table S5).
The known crustacean taxa of the Paja network, Bellcarsinus aptiensis and Planocarcinus olssoni, have a SWTL of 3.25 and a PATL of 3.5. Crustaceans in all the Caribbean networks range between 3 and 4.60 of SWTL and 3 to 5.20 of PATL from herbivorous crustaceans at the bottom and carnivorous crustaceans at the top. All ammonites in the Paja network have SWTL and PATL values of 3.0. Although probably not good analogues, due to their range of diet, the highest values for Caribbean cephalopods correspond to the Caribbean reef squid with a maxSWTL of 4.6 and a PATL of 5.15 on average, and the lowest correspond to the slender inshore squid (SWTL of 4.45, PATL of 4.95) and common octopus (SWTL of 4.2, PATL of 4.5). These taxa are the closest analogues for the ammonites in the Paja network, which have a SWTL and a PATL of 4.0 (Table 2; Supporting Information, Table S5).
DISCUSSION
Comparisons of connectance to extinct and extant networks:
With the premise that ecological escalation results in greater taxonomic, trophic, and niche diversity, one network value that may capture these trends is connectance. Network connectance is simply the number of links per species2. Few extinct marine ecological networks have been published. We compared broad network results with them and extant marine ecosystems. The Chengjian Shale web (L = 99; C = 0.091) and the Burgess Shale web (L = 249; C = 0.108) are Cambrian marine networks whose connectance fall in the range of extant networks (0.071–0.315) (Dunne et al. 2008). The Caribbean networks average 3940 links with a connectance of 0.068. The Paja network has a connectance of 0.177. Examples of ancient marine ecosystem connectance values are too few to draw conclusions but may offer novel insights with the addition of more palaeoecosystems. Limited taxonomic sampling may also influence connectance. The Paja network has only about a third of the number of taxa and interactions of the Caribbean networks. The high numbers of links in the Caribbean networks are largely attributed to the rich benthic invertebrate and fish diversity. These are not recorded in more generally sampled extant systems (e.g. Benguela marine system) and certainly not present in the Paja faunal list. Excluding the benthos may elevate connectance values, especially in the presence of relatively high numbers of predatory links.
Connectance indices are single values of a higher dimensional network and may not be as relevant or comparable across systems, especially ones with taphonomic biases. Instead, we prefer a more trophic level scaled approach, similar to that presented in Figure 7. Vulnerability indices are direct downward predatory links to a specific taxon normalized by L/S. These may prove more useful when comparing ecosystem networks across such large temporal and geographic ranges.

Vulnerability (number of species it is consumed by normalized by links/species) plotted against short-weighted trophic level (SWTL) for: A, Paja Formation; B, Cayman Islands. Graphs made in Rstudio (2022).
Ammonite dominance in the Paja network
: In the Jurassic, there was an explosive adaptive radiation of ammonites, followed by a global diversification later in the Cretaceous (Ward and Signor 1983, Ward et al. 1985, Cecca et al. 2005). This diversification may have been facilitated by ecological escalation during the MMR (Vermeij 1977, Neige and Rouget 2015). The Early Cretaceous Paja Formation biota presents a stunning diversity of ammonites (Figs 2–4). Using even only a stratigraphically constrained subset assemblage of ammonites, the Paja Formation’s ammonites were recovered as important mid-trophic level taxa. They were all recovered at a fourth trophic level. This trophic level is discrete in the Paja network in that few other taxa are at or near this trophic level, forming a relatively isolated trophic level (Fig. 6). Histograms of trophic levels of the Paja and extant Caribbean networks show differences in both the median trophic levels and distributions about them. The Paja and extant networks are both centred about similar mean trophic levels of 3.95 and 3.94, respectively. However, the median trophic level of the Paja network is 4.0, whereas it is 4.27 for the Caribbean systems. The distribution of Paja trophic levels is biased toward ammonites, whereas several taxa fill the median trophic levels in the Caribbean systems and are relatively normally distributed about this median.
Paja Formation ammonites also reveal a unique connectivity within the ecological network. Examining the vulnerability of each taxon (number of species it is consumed by normalized by L/S) plotted against SWTL reveals a stark contrast between the Paja and extant networks (Fig. 7). In addition to the discrete trophic separation of ammonites from the remaining Paja network, ammonites also dominate high vulnerability indices of 2.0, compared to the SWTL. Only zooplankton (vul = 3.0) are higher. The contrasting zooplankton vulnerabilities in the extant networks are only 0.5 and occupy the same third trophic level. Although several taxa have higher vulnerabilities in the extant systems, the majority of taxa are concentrated at relatively low vulnerability values. The high vulnerability and trophic separation of ammonites imply these taxa either function as their own discrete trophic level without extant analogues, trophically similar taxa are unrepresented in the Paja Formation fossil record, or a combination of both. A near lack of understanding of the ecology of ammonites is also a factor. Currently, only a single feeding apparatus is known for an ammonite (Kruta et al. 2011) and this may not reflect the entire feeding variation for the clade.
The diversity of ammonites may also serve as an indicator of missing diversity of smaller-bodied animals. The near absence of smaller-bodied invertebrates and vertebrates is probably due to a combination of taphonomic and taxonomic biases (see below). However, future work could consider what diversity may be surrounding the diversity of ammonites in a stable marine ecosystem.
Taphonomic biases:
It is not surprising that ammonites dominate the Paja faunal record. Centuries of local collecting have contributed to an abundance of specimens in museum collections. Similar collecting biases are present for the large marine reptiles. However, notable depauperate taxa include benthic invertebrates and small-bodied vertebrates. Absence of soft-bodied invertebrates is not surprising; however, the near absence of non-ammonite shelled molluscs and arthropods is. The extant Caribbean systems have at least seven gastropods, over 30 crustaceans, and several infaunal suspension feeders and soft-bodied invertebrates. These taxa compose a rich network of lower trophic levels in the extant Caribbean systems that are not represented in the Paja biota. Only six bivalves, one gastropod, and two crustaceans are represented in the Paja biota. The few arthropods described are well aligned with extant lineages' diversity trends (Sepkoski 2000), suggesting there probably is a rich benthic invertebrate fauna in the Paja Formation still to be discovered. Paja fish diversity is also notably low. The Paja fauna has only one described shark and two described and two undescribed actinopterygian fishes. In comparison, the Caribbean systems have one ray, three to four sharks, and over 70 actinopterygian fishes. In spite of the much-rarefied samples of benthic invertebrates and fishes, turtles have comparable diversities in the Paja and extant Caribbean systems. Four turtles are present in the Paja fauna, whereas up to five are present in the Caribbean faunas.
Environmental biases:
One potential bias between the Paja and extant Caribbean marine faunas is the depositional environment. The Paja Formation was deposited in a shallow-sea setting, similar to the Caribbean systems. However, whereas the Caribbean systems are reef ecosystems, the Paja Formation has no indication of reef-forming taxa and was thus a soft-bottomed sea. Lack of hard substrates for the Paja biota may have also been a factor in its depauperate benthic invertebrate and fish faunas.
Evolutionary biases:
A further bias that requires consideration is that many extant taxa had not evolved by the Early Cretaceous. Although the Paja crustaceans and turtles share stem relations with extant lineages (Cadena 2015, Cadena and Parham 2015, Evers and Benson 2018, Luque et al. 2020), none of the ichthyosaur nor plesiosaur taxa do. Although actinopterygian fishes originated at least by the Devonian, common reef fishes, such as blennies, flatfishes, wrasses, and tetraodontiforms, did not evolve until the Late Cretaceous (Near et al. 2012). The remarkable diversification of extant reef fishes is dominated by biting feeders and has been shown to be a relatively recent event during the Eocene (Corn et al. 2022). By the time the Paja Formation was deposited, most lineages of extant teleost fishes were simply not present. The radiation of reef-dwelling teleost fish involves 24 lineages from which, apparently, based on ancestral-state reconstruction of feeding mode phylogeny, only two stem-based lineages could have been present in the Paja biota (see Corn et al. 2022: fig. 1).
Time-averaging:
We constrained the ammonite fauna to the Lower Barremian unit A within the Arcillolitas Abigarradas member to avoid time-averaging the entire ammonite fauna of this 10-million-year-long member. This unit is poorly constrained with units A–C in the Barremian and D and E in the Lower Aptian (Etayo-Serna 1979). The vertebrate fauna, however, has no stratigraphic control within this member, yet. One confounding factor in our study was the potential time-averaging of vertebrate taxa over the 10-million-year deposition of the Arcillolitas Abigarradas member. We argue that although this would inflate our taxonomic diversity, the trophic diversity is less prone to time-averaging. The presence of several large pliosaur and elasmosaur plesiosaurs, large ichthyosaurs, and turtles may be due to time-averaging. Several Late Jurassic and Early Cretaceous marine deposits have similar vertebrate taxa, including abundant large plesiosaurs and ichthyosaurs. However, even if several of the large marine reptiles in the Paja biota are artificially replicated due to time-averaging, the trophic disparity of the network would remain, retaining the high trophic levels and complex food chains of the large apex predators.
Comparisons with the Caribbean networks may be biased from the paleontological time-averaging of our Paja Biota. The Caribbean networks are derived from contemporaneous biodiversity surveys during the past 50 years (Roopnarine and Hertog 2012). If the Caribbean taxonomic lists were expanded to include all taxa present during the past 10 million years to establish similar possible time-averaged lists, they might include more large-bodied fish species such as megatooth sharks (Pimiento and Balk 2015) and extinct large marine mammals such as marine iniid whales (Pyenson et al. 2015) and manatees (Suarez et al. 2021). The list would exclude some recent human-mediated invasive taxa such as the lionfish. Human activity has also negatively impacted abundances of several Caribbean coral, macroalgal, and fish species (Gardner et al. 2003, Mora 2008, Hardt 2009, Dillon et al. 2021); however, pre-human Caribbean marine ecosystems have not yet been assembled. Relatively short time-intervals spanning only from the Late Pleistocene to the Present reveal significant ecosystem network changes in terrestrial mammal communities (Fricke et al. 2022). These terrestrial community changes are considered to have been largely human-mediated (Fricke et al. 2022, and references therein). Although no study has specifically assessed community changes in the Caribbean during the last several million years, we would expect it to have followed the patterns observed elsewhere with at least some extinctions and probably differences in particular species abundances during Holocene times. Quantifying these time-averaging biases with the Caribbean networks is clearly beyond the scope of this paper but worth noting for future work.
Time-averaging will be an important issue going forward with palaeoecosystem network modelling. Other published fossil systems also have this issue. For example, the Burgess Shale fauna has been suggested to range from 6 million years to as little as 200 000 years (Collom et al. 2009, Nanglu et al. 2020), the Chengjiang fauna may span 7 million years (Zhao et al. 2019), and the Messel biota spans 1 million years (Mertz and Renne 2005). However, when 10 discrete faunal horizons were examined in the Chengjiang biota, little taxonomic variation was observed, suggesting relative stability of the ecosystem (Zhao et al. 2014).
Recalibrating the Paja network:
In Figure 8, we aligned the Paja and Cayman networks using ecologically analogous taxa to calibrate missing trophic levels in the Paja network. PATL was used instead of SWTL because we wanted to maximize differentiation among apex predators. PATL values slightly augment trophic levels because the index assumes a predator’s diet fractions are equal across the prey. Five fossil taxa were chosen as ecological analogues to the Cayman network (Table 2; Supporting Information, Table S5). The lowest trophic level ecological analogue pairing was determined to be between the fossil decapods Bellcarsinus aptiensis and Planocarcinus olssoni with the extant omnivorous crustacean III. The latter, with a PATL of 3.9, was chosen because it represents the average for the 11 Cayman omnivorous crustaceans that range from PATLs 3.4 to 4.6. The Paja crustaceans were originally estimated at a PATL of 3.5, and their recalibrated position thus requires an addition of 0.4.

Recalibrated Paja Formation network (A) computed to Cayman network (B). The revised Paja Formation network is relatively linear. We calibrated it to the recovered predator-averaged trophic level (PATL) for the Cayman taxa. Trophic level (TL) scale is shown on each side of the networks, ranging from TL1 to TL8, and from red to yellow, meaning bottom to top.
All four Paja marine turtles had a PATL of 4.8. The range of marine turtle PATLs in the Cayman reef ecosystem ranged from 3.0 for the green sea turtle to 6.0 for the loggerhead sea turtle. We chose to align the Paja turtles with the hawksbill turtle with a PATL of 5.1 given the generalist diet of the taxon and relatively more conservative trophic position. The minor differences between the sea turtle trophic levels in the Paja and Cayman networks suggest that there is good correspondence between these putative trophic analogues of these very different marine ecosystems.
The fossil ecological analogues with the highest trophic levels are the shark Protolamna ricaurtensis and the undescribed Actinopterygii morphospecies 2 with PATLs of 4.9 and 5.0, respectively. These taxa were aligned within the ranges of the oceanic whitetip shark and barracuda, respectively, with extant SWTLs of 5.9 and 5.5, respectively. Increasing the trophic levels for these Paja taxa by a half trophic level brings them approximately in line with the Cayman taxa.
Recalibrating the Paja ecosystem network to the extant Cayman network elevated several trophic positions considerably by 0.3 to 0.5 levels. This recalibration has several assumptions, including the assumption that some taxonomically similar species are trophic analogues and that both ecosystems have similar topologies. The lower trophic level recalibration assumes a rich but undiscovered benthic invertebrate assemblage in the Paja biota. Another assumes a large missing fish fauna. Both these assumptions are supported by overlaying the body size distribution of taxa in the Paja and extant Caribbean networks (Fig. 9). Although the Paja biota occupies a large range of body sizes, the differences in the abundance of small invertebrates and fishes between the Paja and Caribbean ecosystems are clear.

Estimated body size distribution of the Paja Formation marine ecosystem (coloured) and the extant Caribbean ecosystems (grey). The diagram shows the Paja biota, abundant in large apex predators, occupies a large range of body sizes although is lacking a significant portion of invertebrates and fishes that are instead abundant in the Caribbean system, but missing the large apex predators.
Although the depauperate fossil record makes it difficult to further explore these assumptions, the Paja biota still has an abundance of apex predators that are not present in any of the extant Caribbean systems. Fourteen large predatory marine reptiles occupy trophic positions above the highest trophic positions of fossil fishes. The size disparity of the largest and top predators from the largest fishes suggests a large amount of missing trophic information. If we assume there are no missing fossil taxa, the apex predators would still be pushed up approximately a full trophic level due to the recalibration of the top fish predators. However, we can also hypothesize that there is at least one trophic level separating these fishes from the largest predatory reptiles. This expansion is based on the large size disparity within the large apex predators and between them and taxa at lower trophic levels (Fig. 9). This would push the largest apex predators, the pliosaurs Monquirasaurus (=Kronosaurus) boyacensis and Sachicasaurusvitae, to a trophic level of at least seven! The presence of such high trophic levels has no extant analogues.
The only extant large-bodied extant marine animals comparable in size to the large Paja marine reptiles are the predatory white sharks and suspension-feeding mobulid rays, basking sharks, megamouth sharks, whale sharks, and baleen whales (Stiefel 2021). None of these taxa interact directly with reef communities. No direct estimates of Caribbean whale shark and baleen whale trophic levels have been published. However, a global study based on diet presents a range of 3.2–3.4 for baleen whales, 3.8–4.6 for odontocete whales, with sperm whales and killer whales at the upper range (Pauly et al. 1998). Similar values are estimated for whale sharks and white sharks, respectively (Cortés 1999). Although methods to estimate trophic levels differ, plants were assumed to be at the first trophic level and estimated trophic levels of 2.2 for benthic invertebrates and large zooplankton, 2.7–3.3 for teleost fishes, and 4.0 for sea turtles and marine mammals that did not feed on other mammals (Pauly and Christensen 1993). These are approximately one trophic level less than those estimated for similar taxa in the Caribbean reef ecosystems because those systems, and all quantified palaeontological systems, included several unicellular and detrital nodes. The addition of a trophic level to recalibrate to the present study suggests that whale sharks and baleen whales may reach a trophic level of 4.2–4.4. The same calibration would have the largest predatory odontocetes and white sharks reaching a trophic level of 5.6. Such a level is not so different from the apex predators in the Caribbean ecosystems.
No large-bodied suspension feeder is yet known from the Paja Formation biota. The only Mesozoic taxa with this ecology are a Late Cretaceous shark (Shimada et al. 2015) and the Middle Jurassic to Late Cretaceous pachycormid fishes (Friedman et al. 2010). None of these taxa are yet reported from the Paja Formation biota, although two smaller pachycormid fishes are described in the Turonian Hondita Formation in Colombia (Páramo-Fonseca 2001).
The Eocene to Pliocene megatooth sharks (Otodus) may be more comparable to the largest Paja predators. Although absolute trophic levels have not been estimated for these sharks, nitrogen isotopes suggest that they occupied higher trophic levels than extant white sharks (Kast et al. 2022), whereas zinc isotopes suggest that these taxa occupied similar trophic levels (McCormack et al. 2022). In spite of these differences, the enormous size and predatory ecology of Otodus have been aptly described as ‘hyper-apex predators’. The great size of the largest predators in the Paja biota is also deserving of this title.
Future directions:
Simple food-webs tend to be less stable and more prone to extinctions than complex webs (Elton 1927, Roopnarine 2010, Marshall 2015). Our study is the first approximation of a marine Mesozoic quantitative food-web. In spite of probable taphonomic and artificial sampling biases, the Paja Formation biota yields a highly diverse and complex system. Many questions remain unexplored in terms of palaeoecosystems. Here, we highlight those that have the potential to be the starting point to explore ancient ecosystem complexity and the impacts of biotic and abiotic factors.
Extant food-webs outline the flow of energy from lower to higher trophic levels. The degree of basal productivity dictates the biomass of higher trophic levels that it can support. The greater the complexity and biomass of higher trophic levels, the greater the amount of basal productivity required (Cohen et al. 2003, Thompson et al. 2012, Guiet et al. 2016). Extant communities are generally composed of a high diversity of small-sized taxa that combine to a high biomass at lower trophic levels and less diversity and biomass at higher trophic levels. Presenting only the diversity of taxa present in the Cayman Islands reef ecosystem highlights this pattern (Fig. 9). Energy-flow modelling of predator–prey interactions has suggested stable scaling regressions of biomass change through aquatic and terrestrial ecosystems (Hatton et al. 2015). Metabolism and abundance of taxa also scale well with body size across a wide range of eukaryotes (Hatton et al. 2019). Extending these ecological models in conjunction with ecological networks to fossil biotas may help estimate both the extent of missing taxa within the network and potential abundances and biomasses of all taxa in the network. Ultimately, these values may work toward the total productivity of ecosystems.
The enormous body size of the hyper-apex predators in the Paja biota suggests a correspondingly complex ecological network supported by a large biomass at lower trophic levels (Young et al. 2013). No extant ecosystem may serve as a direct analogue. However, using extant ecosystem energetic models, we expect that the Paja biota and several other Mesozoic ecosystems could set upper limits to both ecosystem complexity and productivity.
The Paja Formation preserves several biotic and environmental changes that include a biotic recovery from the end of the Jurassic mass extinction, rapid rises in global sea-levels, continental fragmentation, and high global temperatures. All are expected to facilitate high rates of speciation; however, no ecosystem-level approach has assessed these influences. The ecological network introduced here for the Paja biota may be used as a starting point to examine how these abiotic and biotic factors may have influenced marine ecosystems. Several well-preserved and well-studied fossil marine assemblages are known from Late Jurassic to Late Cretaceous deposits. Additionally, trends in the diversity of specific lineages have documented rapid rises in diversity during the Early Cretaceous of marine invertebrates (Sepkoski 1981, Foote 2007, Alroy et al. 2008, Cermeño et al. 2022), fishes (Guinot and Cavin 2016), and reptiles (Benson et al. 2010). Integrating these lineage radiations within ecological networks may shed light on how ecosystems responded to and potentially drove these elevations in biodiversity.
Epicontinental seas were prevalent during the Cretaceous while Pangaea fragmented, creating new coastlines, and global sea-levels rose, rapidly flooding continents (Miller et al. 2005). Although open oceans have been estimated to have had higher origination rates during the Early Cretaceous than epicontinental seas, the latter had lower extinction rates (Miller and Foote 2009). These environmental and biotic patterns, in combination to pelagic successions, may have favoured ecologically dynamic and productive shallow marine ecosystems during the Cretaceous (Schlanger and Jenkyns 1976, Barron 1983, Takashima et al. 2006, Hou and Li 2018). Epicontinental seas flooded about 20% of continental areas during the Cretaceous (Barron 1983). A variety of habitats and ecological niches (habitat heterogeneity) could have been available for marine organisms to produce high levels of endemism (Hou and Li 2018). Thus, during the time in which the Paja Formation was deposited (Hauterivian–Aptian), a combination of geological and climatic factors may have promoted endemism in these assemblages at this time (Schlanger and Jenkyns 1976, Gasparini and Fernández, 1997, Takashima et al. 2006, Bardet et al. 2014). The fragmentation of the Tethys Ocean, for example, has been suggested to have augmented aquatic diversification via vicariant speciation, which in turn increased endemism, provincialism, and overall biotic diversity during the Mesozoic (Valentine and Moores 1970, Hou and Li 2018). The warmer tropical waters, along with an increased volume of new marine habitats resulting from continental fragmentation and marine transgressions, could have developed ideal ecological conditions for the marine fauna to radiate and explore new trophic levels. We suggest these factors contributed to the high diversity of the Paja Formation biota.
CONCLUSION
This is the first reconstructed food-web of the Early Cretaceous Paja Formation and the first quantitative site-specific ecological network for a Mesozoic marine fauna. The Paja network suggests this trophic network was diverse, complex, and possibly presented more trophic levels than any extant marine system. The Paja network yielded producers, consumers, and a significant proportion of large apex predators. Although the Paja Formation was deposited during the span of about 15 million years (Hauterivian–Aptian), the fossils used here are from a well-constrained stratigraphic member within the formation. Within this member, however, the possible range is still up to 10 million years of accumulation and presents several important biases of time-averaging in any fossil-based ecosystem network. The Paja Formation biota provides new insights into the Mesozoic trophic interactions and the food chains of the Early Cretaceous marine systems. These include highlighting the absence of the expected high diversity of small benthic invertebrates and fishes. The network also revealed surprisingly high trophic levels for the hyper-apex predators of the system with no extant trophic analogues. The results of this study provide useful data for understanding large-scale community dynamics and the stability of marine ecosystems during the Mesozoic Marine Revolution (MMR). Although this is the first taxonomic and palaeoecological compilation of the Paja Formation’s marine fauna, there is still need for additional basic taxonomic work on all vertebrate, invertebrate, and plant assemblages from this unique marine sequence. These datasets are meant to be starting points for work in progress and for collaborative initiatives to start reconstructing energy-flow models of marine ecosystems.
SUPPLEMENTARY DATA
Supplementary data are available at Zoological Journal of the Linnean Society Journal online.
Table S1. Dataset with master taxon list for the Paja Formation biota.
Table S2. Dataset with taxon list used to reconstruct the Paja Formation network.
Table S3. Dataset with Paja Formation biota edgelist.
Table S4. Dataset with Caribbean network edgelists.
Table S5. Dataset with network3D properties output. SWTL, short-weighted trophic level; PATL, prey-averaged trophic level, which assumes a predator’s diet fractions are equal across the prey.
ACKNOWLEDGEMENTS
This work would not have been possible without the immense contribution of the scientific and local people who have gathered data and fossils from the Cretaceous of Colombia, and who have strengthened extensively the palaeontological collections in the Departamento de Boyacá over the last 45 years. We appreciate the work done by those who collected, prepared, and/or described the fossil material known from the Paja Formation biota to date. We are grateful to M. Parra, J. Parra, and F.H. Parra from the Centro de Investigaciones Paleontológicas, Villa de Leyva, Boyacá (CIP) for allowing us access to their collections and facilities, and for making available a large amount of fossils and data to perform this study. We would like to express our deepest appreciation for the work and support of C.B. Padilla† and S. Padilla† (CIP) who dedicated a great part of their lives to the rescue and conservation of the fossils from Villa de Leyva and to developing the CIP. They will be missed. We thank J. Dunne (Santa Fe Institute) and R.J. Williams (Vibrant Data Labs) who gave us access to Network3D software and kindly clarified questions about the software. Relevant discussions made this work more focused and highly improved the main ideas. For these, we appreciate the important feedback given by C. Jaramillo (STRI) and A. Demers-Potvin (McGill). D.C. thanks the Paleontrópica members (CTPA), Larsson Lab members (Redpath Museum), Erin E. Maxwell (SMNS), and D. Mera-Rodríguez (UH) for their support and scientific input. Funding to D.C. was provided by BESS-NEO program, NSERC CREATE 46283-2015, the Smithsonian Tropical Research Institute, the Anders Foundation, the 1923 Fund, and Gregory D. and Jennifer Walston Johnson, the Fonds de recherche Nature et technologies Québec (FRQNT), the Redpath Museum’s Delise Alison Award-2019, the Sigma Xi Grant-in-aid-of-Research, Canada-2019 (GIAR), the Quebec Center for Biodiversity Science excellence award-2019 (QCBS), and the Society of Systematic Biologists Graduate Student Research Award, 2021 (GSRA). This work is part of DC.’s Ph.D. dissertation at McGill University. Funding for HL was supported by NSERC RGPIN-2022-04370. Relevant comments by the editorial staff and an anonymous reviewer highly improved this manuscript.
CONFLICT OF INTEREST
None declared.