Global climate model comparisons of niche evolution in Turritelline gastropods across the end-Cretaceous mass extinction

Paleo-Ecological Niche Modeling (PaleoENM) aims to map the distributions of extinct species using paleo-coordinates of fossils and local environmental data. While General Circulation Models (GCMs) have been widely used to estimate climate conditions in deep time, they have primarily been applied to the terrestrial vertebrate record. Furthermore, variations in paleo-elevation models used in GCM construction can significantly influence the outcomes of PaleoENM. This study addresses two main objectives: (1) to investigate whether changing climatic factors drove niche shifts following the end-Cretaceous mass extinction in the shelly marine invertebrate group, the tower snails (Turritellidae: Turritellinae), and (2) to compare the effects of two different paleo-elevation models on the results of GCM-based predictions of species distribution. Fossil occurrence data from the Maastrichtian and Danian time periods were obtained from the Paleobiology Database, supplemented by museum collections and published literature. Environmental data were extracted from atmosphere-ocean General Circulation Model (GCM) simulations using the HadCM3L model, applying two different sets of paleogeographic and CO2 boundary conditions: Scotese-based and Getech-based. Additional sedimentology and depositional environment data were sourced from the Paleobiology Database (PBDB). We predicted the distributions of Turritellines using the maximum entropy (MaxentMaxEnt) algorithm and performed niche similarity analysis using principal component analysis and kernel density estimation. We found significant differences in the spatial arrangement of suitable habitats between the Maastrichtian and Danian time periods across GCMs. The results also showed that the Getech-based GCM outperformed the Scotese-based GCM in terms of model metrics. Niche overlap across both time periods was high, with niche similarity and equivalency being higher than expected by chance within both GCMs. Our results also suggest that differences in elevation model boundary conditions led to variations in the predicted distribution and niche patterns. This study provides a novel approach to understanding ecological resilience and niche change in invertebrate taxa after mass extinction events. It also explores the robustness of varying GCM boundary conditions on PaleoENM studies and offers a framework for future paleoecological research on fossil invertebrate taxa.

Methods
Taxonomy of Turritellinae—Although turritellids are well known from the fossil record, the intrafamilial relationships among turritellids remain an active area of taxonomic research. Turritellidae was divided into five subfamilies, Turritellinae, Protominae, Pareorinae, Turritellopsinae, and Vermiculariinae by Marwick (1957), with Turritellinae including most nominal genera and subgenera. Subsequent analyses, informed by both molecular and morphological data, have found that the Vermiculariinae are well nested within the Turritellinae, which may also be the position of the Protominae (Anderson and Allmon 2024). Many generic and subgeneric names exist in the literature, but very few are used consistently. Recent work, supported by new molecular analyses (Anderson 2018, Anderson and Allmon 2024, Lieberman et al. 1993, Sang et al. 2019) has supported Marwick’s (1957) prioritization of protoconch form, order of appearance of spiral ornamentation, and growth line characters as indicative of taxonomic affinity, in the absence of significant morphological differentiation in the teleoconch (e.g. Vermicularia, Caviturritella) (Friend et al. 2023). However, the comparative lack of well-preserved protoconch material available leaves many species broadly assigned to “Turritella sensu lato” on the basis of plesiomorphic or convergent teleoconch form, making “Turritella” a true “wastepaper basket” genus (Allmon 1996, 2011, Allmon and Cohen 2008, Hendricks et al. 2014, Plotnick and Wagner 2006). Yet, species assigned to Turritella s.l. (given their absence of teleoconch characters that assign them to other subfamilies) are all expected to be encompassed by a monophyletic Turritellinae. Furthermore, extant species within Turritellinae all have similar ecologies as shallow infaunal, suspension (and occasionally deposit) feeding gastropods with the overwhelming majority also having highly similar life histories (Allmon 1988, 2011, Anderson and Allmon 2020). Therefore, we consider the Maastrichtian and Danian members of Turritellinae, a clade first appearing in the Jurassic (Das et al. 2018) and diversifying in the Cretaceous, to be broadly like a modern gastropod genus in terms of species richness, clade age, and ecological similarity among constituent species (Allmon 2011).

Unfortunately, confusion in the literature and large databases like the PBDB is not limited to the proper generic assignments of species. Large numbers of fossil occurrences of “Turritella sp.” do not represent members of monophyletic Turritellinae (or even Turritellidae) but may represent other high-spired gastropods mistakenly identified as “Turritella.”  Furthermore, the name “Turritella" has been erroneously applied to entire rock layers such as ‘turritella agate’ when such layers belong to the freshwater gastropod Elimia tenera Hall, 1845 (Allmon and Knight 1993). Therefore, we chose to restrict our analyses to only those occurrences which had species-level identification.

Because of these taxonomic and temporal record uncertainties, we have decided to analyze the subfamily Turritelline as a single taxon. We recognize that ENM techniques are best interpreted at the species-level, however, taxonomic reality and spatiotemporal resolution support a more robust analysis at the clade level (Hendricks et al. 2014).

Occurrence Records—We acquired occurrence records of species belonging to Turritellinae within the Maastrichtian and Danian periods from the Paleobiology Database (PBDB; paleodb.org), Global Biodiversity Information Facility (GBIF; gbif.org, downloaded July 2022), additional occurrences from C.M. (Myers et al. 2013), as well as an exhaustive literature search for the U.S. coastal plain occurrences by W.A. and K.C. We selected genera belonging to Turritellinae based on taxonomic, morphological (Allmon 1988, 2011, Friend et al. 2023, Harzhauser and Landau 2019) and recent molecular phylogenetic evidence from Anderson (2018) and Anderson and Allmon (2024) (Table 1). Phylogenetic evidence suggests a deep split between the majority of turritellid gastropod species, assigned to the subfamily Turritellinae, and the subfamily Pareorinae (Marwick, 1957) which includes the genera Mesalia and Sigmesalia (Anderson 2018, Anderson and Allmon 2024); we omitted both genera from the analysis (Table 1). Molecular evidence also suggests that two other traditional subfamilies, Vermiculariinae Dall 1913 and Protominae Marwick 1957, are nested within the Turritellidae (Anderson and Allmon 2024), however, neither of these has Maastrichtian or Danian fossil records and therefore this inclusion does not affect the analyses. Since the genus Turritella s.s., is nested within Turritellinae, and fossils are commonly misidentified as “Turritella sp.,” we only selected taxa possessing a species epithet (in the ‘accepted name’ category within PBDB), and that could be verified in the literature as being present within the selected time periods. Further occurrence filtering consisted of omitting specimens for which latitude and longitude could not be computed. We binned the time intervals to the Maastrichtian (72.1 - 66 mya) which possessed 2032 occurrences, and the Danian (66-61.6 mya) consisting of 1883 occurrences (Table 2, See Supplemental Information).

Paleogeographic maps of the Phanerozoic predominantly follow the paleo-digital elevation model PALEOMAP, which estimates Earth’s past paleoceanography, the changing area of land, mountains, shallow seas, and deep oceans through time (Scotese 2016, Scotese and Wright 2018). Fossil databases such as the Paleobiology Database (PBDB; paleodb.org) follow PALEOMAP for paleo-rotation of taxa to their estimated place of deposition within the Phanerozoic. The GCM simulations of past climate from (Valdes et al. 2021) have utilized Scotese’s model of paleo-elevation. Paleo-elevation models such as those from the Getech Plc (Getech.com) have also been used in GCM simulations from Lunt et al. (2016) and Farnsworth et al. (2019); these provide alternative reconstructions of Earth’s past climate. However, differences in the paleo-coastlines of each model (particularly epicontinental sea boundaries) may result in the reconstruction of coastal invertebrate occurrences on land instead of in the ocean (Scotese and Wright 2018).

Using the rgplates package v 0.3.2 (Kocsis et al. 2023) in R v. 4.0 (Team 2021), we paleo-rotated both datasets to their respective time periods using Scotese’s PALEOMAP (Scotese 2016) and Getech’s paleo-rotation model provided by Farnsworth et al. (2019). PALEOMAP consists of 120 unique paleo-digital elevation models (PaleoDEMs) representing three-million-year time slices roughly equating to different stratigraphic stages. We paleo-rotated the occurrences of the Maastrichtian to the 69mya time slice, and the Danian to the 66mya time slice. We repeated this process using Getech’s in-house model with the same time slices (Farnsworth et al. 2019) (See Supplemental Information for all paleorotated datasets).

Climate Data—We generated environmental layers using output from the HadCML3 climate model version 4.5 (Valdes et al. 2017), with inputs from the solar luminosity, and both the PALEOMAP paleogeographic atlas (Scotese 2016, Scotese and Wright 2018) (hereafter referred to as the Scotese simulations) and from Getech Plc (Getech.com) (hereafter referred to as the Getech simulations) for the Maastrichtian and Danian periods. For both sets of simulations, we use the GCM “HadCM3LB-M2.1aD,” described in detail in (Valdes et al. 2017), for which surface resolution is 3.75° longitude × 2.75° latitude (grid box size of ~420 × 220 km; at the equator, reducing to ~200 × 280 km at 45° latitude).

The Scotese simulations are similar to the latest Maastrichtian (Map number 22) and Danian (Map number 20) simulations described in Valdes et al. (2021). These are the simulations that prescribe a pCO2 concentration according to Foster et al. (2017). Compared to the Valdes et al. (2021) study, simulations were run for an additional 2,000 years, and with modified atmospheric and ocean physics by applying methods similar to those described in Sagoo et al. (2013), resulting in an improved representation of polar amplification in deep-time climates. In addition, they have islands defined correctly for the purposes of calculating the ocean barotropic stream function, according to Foreman (2005).

The Getech simulations are identical to those described in Farnsworth et al. (2019). The pCO2 concentration for both Maastrichtian and Danian Getech simulations were each 1,120 ppmv, which is within the typical range of uncertainty in Foster et al. (2017) pCO2 data, which itself approximates the actual evolution of pCO2 through time with some uncertainty (see Fig. 1 in (Foster et al. 2017)).

Climate model output variables chosen within the Maastrichtian and Danian for both model sets were: annual mixed layer depth (meters), annual maximum and minimum (warmest and coldest) sea-surface temperature within a 3-month interval (season; Celsius), and monsoon seasonality index (difference in precipitation between the three driest and three wettest months of the year); monsoon seasonality is a proxy for preference of seasonal variability around the tropics. Furthermore, we acquired annual potential temperature (Celsius), and annual salinity (PSU, practical salinity unit) averaged from 5 – 5800m depth (See Valdes et al. (2017) for in-depth definitions of variables). Although modern Turritellines are found predominantly between depths of 2 – 100m (Allmon 2011), occurrences affected by paleo-rotation and placed at deeper depths would be mistakenly removed if layers were restricted to shallow depths including epeiric seas found in North America and Eurasia. Finally, we incorporated bathymetry (meters), derived from the digital elevation models (DEMs) which underlie both the Getech and Scotese GCM simulations. These variables were chosen to encapsulate the abiotic restrictions of modern Turritelline species (Allmon 1988, 1996, 2011, Allmon and Cohen 2008, Allmon and Knight 1993). We rotated all variables from 0 - 360 longitudinal format to -180 - 180 longitudinal format to avoid cutoff at the international date line.

Sedimentology Data—We generated two distinct sedimentological, interpolated datasets derived from lithology and depositional environment proxies at fossil localities within the Maastrichtian and Danian. We downloaded all marine occurrences within the Cretaceous/Paleogene boundary stages from the PBDB and assigned each collection to either carbonate or siliciclastic using the primary lithologic data associated with the occurrence, recorded as ‘lithology1’ field in the collection record (see Hopkins (2014) and Foote (2006) for criterion). Within the PBDB, each lithological term possesses specific definitions so that the collection is consistently described in ways relative to all other collections in the database, allowing for consistent translation to the two lithological categories. We also assigned each collection to different depositional environments recorded from the ‘environment’ field to determine preferred depositional settings for Turritellines within the Maastrichtian and Danian. In total, we used 58 depositional categories. We supplemented the PBDB lithological and environmental datasets with localities from Paleo Reefs (https://www.paleo-reefs.pal.uni-erlangen.de/; (Kiessling et al. 1999)), and the Sedimentary Geochemistry and Paleoenvironments Project (SGP) (https://sgp.stanford.edu/; (Farrell et al. 2021)). We classified all paleo-reef collections as being carbonate, and only selected collections from the SGP belonging to one of the two lithological categories. We categorized paleo-reef collections as ‘reef’ while the SGP only possessed ‘basinal’, ‘fluvial’, ‘inner shelf’, ‘lacustrine’, and ‘outer shelf’ environmental categories, all of which are consistent with the PBDB ‘environment’. Definitions for lithologies and depositional environments can be found at https://paleobiodb.org/public/tips/lithtips.html, and https://paleobiodb.org/public/tips/environtips.html, respectively. It should be noted that depositional environment categories as defined and used in the PBDB are not strictly independent of one another and are defined in a hierarchical fashion. Thus, our analysis provides an example of a very simplified use of this data. In the Supplemental materials, we provide a detailed description of how these categories can be condensed into more directly comparable categories and return to this topic in the Discussion (Supplemental Table S1).

We paleo-rotated the lithological and depositional environment datasets to the Getech and Scotese GCM paleo-rotation models within the Maastrichtian and Danian. We assigned each lithological and environmental category an integer (ex: 1 = carbonate, 2 = siliciclastic) and performed a nearest-neighbor interpolation using the ‘interpNear’ function in the R package terra (Hijmans et al. 2022) to generate categorical raster layers. When lithological or depositional categories vary within a single pixel, the nearest-neighbor interpolation assigns the pixel’s value based on Euclidean distance—selecting the occurrence closest to the pixel’s center. Nearest-neighbor is an efficient method for interpolating categorical data since it preserves original categories, avoids blending (i.e., averaging), and is less computationally intensive (Johnson and Clarke 2021). We resampled the categorical rasters to the same resolution of the GCM layers using the ‘resample’ function in the raster package in R (Hijmans et al. 2021). Climate model, lithological, and depositional environment output variables can be obtained in this dataset.

Ecological Niche Modeling—Before modeling, we spatially thinned occurrences to the resolution of the raster layers to reduce sampling bias, artificial clustering, and subsequent spatial autocorrelation (Veloz 2009) using the spThin package (Aiello-Lammens et al. 2015), which resulted in more even sample sizes of occurrence points across time periods (see Table 2). ENM modeling is more accurate when occurrences are thinned to the resolution of pixels within our raster layers, such that each occurrence gets a singular value of each of the climate and sedimentological values used (Araújo et al. 2019). Similarly, after thinning, only a single species’ occurrence may exist in any single pixel, which is an important step to prevent model over-sampling of the same pixels within the study extent (Aiello-Lammens et al. 2015). Given that the collated occurrence data often contained already-assigned environmental conditions (e.g., lithology and depositional environment in PBDB-based records), we checked whether (1) multiple occurrences within the same pixel had the same environmental designations, and (2) thinning occurrences caused discrepancies between the environmental values derived from our lithology and depositional environment layers and the original PBDB-based record designations. We found no discrepancies between the PBDB environmental designations linked to occurrences within a pixel and the interpolated pixel values from our lithology and depositional environment layers.

Turritellines possess a global distribution (Allmon 2011), thus we chose a study extent to encapsulate the entire preserved geographic range of Turritellinae in both the Maastrichtian and Danian time periods. Furthermore, to include potentially undersampled areas despite dispersal limitations of individual species (Peterson and Soberón 2012), we defined a bounding box buffered by 8.12° (~100,000 m2) as this is the maximum dispersal capability of modern-day Turritellines (‘M’ training region as defined by  Soberón and Nakamura (2009)) (Allmon 1996, 2011). From the total region encompassing all thinned occurrences enveloped by the bounding box, we randomly sampled background points for modeling (n=888 for Scotese and n=1842 for Getech study extents); background points differ between niche models due to the number of occurrence points differing between paleorotation models which produced the resultant bounding box. Background points are a random sample of the environmental conditions across our entire study area which help define what is suitable and unsuitable habitat in model output. Our model compares environmental conditions at our presence (occurrence) points to that of our background points to identify which conditions are more strongly associated with the presence of the species (Elith et al. 2011, Phillips 2021, Phillips et al. 2017, Phillips and Dudík 2008). We extracted environmental values from our background points to calculate correlations between variables using the ‘vifcor’ function in the usdm package (Naimi 2017) and filtered out variables with correlation coefficients higher than 0.7. In ENM analyses, the balance of the correlation threshold against the number of environmental variables retained for analysis is an area of active study with coefficients ranging from 0.5 - 0.9 (Graham 2003, Yan et al. 2020). Since fewer variables create underfitting models, we chose a moderate correlation threshold (0.7), which acted as a natural break in which higher thresholds significantly reduced the number of environmental variables retained. We used the same set of variables for subsequent model building, retaining only those that remained after assessing correlations between GCM type and time period.

We used the machine learning algorithm MaxEnt v3.4.4. (Phillips et al. 2017), which remains one of the top-performing algorithms for fitting presence-background ENMs (Valavi et al. 2021). Even though a substantial number of occurrences were removed during thinning, the sample size was large enough to produce robust results (Table 2) (Shcheglovitova and Anderson 2013). At small spatial and taxonomic scales, fewer than 10 occurrences can be sufficient to create highly accurate models (Hernandez et al. 2006, Pearson et al. 2007, van Proosdij et al. 2016). A recent study by Qiao et al. (2017) suggested that the incorporation of occurrences of related species increases model fidelity on undersampled taxa. Nonetheless, the sample size should be considered carefully because the number of occurrence points varies greatly depending on the study system, and too few paleoENM studies have been published to generate a consensus (See Supplemental Information of Myers et al. (2015)). Within ENMs, partitioning data involves splitting the dataset into two subsets: the training dataset is the subset of data used to build the model, while the testing data is used to evaluate the performance of each model iteration (Elith et al. 2011, Phillips 2021, Phillips et al. 2017, Phillips and Dudík 2008). Partitioning reduces the risk of model overfitting and assesses its performance with unseen data. We used the ‘checkerboard2’ strategy of spatial partitioning which divides up our study extent into a grid of spatial cells whereby occurrence points are then alternatively assigned to different partitions (the default setting used here assigns four partitions).

All final models were fitted to the full datasets (training + testing) (Phillips et al. 2017). As the combination of two key complexity settings in MaxEnt models, feature classes and regularization multipliers, can strongly influence model outputs (Radosavljevic and Anderson 2014, Warren and Seifert 2011), we tuned model complexity to find optimal settings. For tuning, in order of increasing complexity, we chose the feature classes linear (L), quadratic (Q), and hinge (H), as well as regularization multipliers 1 through 5 (higher numbers penalize complexity more) (See Supplemental Information for a more detailed explanation of MaxEnt model building). In brief, feature classes determine the shape of the model fit, while regularization multipliers control how much complexity is penalized—this can result in predictor variable coefficients shrinking to 0 and thus dropping out of the model (Phillips and Dudík 2008).

We assessed optimal models using sequential criteria that included threshold-dependent (omission rate) and threshold-independent (AUC) performance metrics. We first filtered models that possessed a delta Akaike Information Criterion (with correction for small sample sizes) less than 2 to compare models between time periods and GCM type (AICc; Warren and Seifert (2011)). We then chose models that possessed the smallest 10-percentile omission rate. If a series of models possessed near-identical 10-percentile omission rates, we chose the simplest model or the one with the fewest non-zero lambda values (model coefficients). For each time and GCM model, we documented variable importance and plotted marginal response curves to better understand the modeled relationships between the predictor variables and the data. We recorded the permutation importance metric output by MaxEnt, which is calculated by randomly permuting the values of all environmental variables but one, building a new model, and then calculating the difference between each model’s training AUC and that of the empirical model (Phillips et al. 2017). Marginal response curves are generated by constraining all predictor variables to their means except for one, then making model predictions along the full range of the focal variable associated with the training data. These curves show the modeled relationship of each variable individually with the occurrence data when all other variables are held constant and are affected by the complexity of the model settings (Phillips et al. 2017).

We calculated a multivariate similarity surface (MESS; Elith et al. (2010)) to detect the degree of similarity in extracted empirical occurrence (presences) environmental values and extracted background environmental values for each time period and GCM type. Increasingly negative values in MESS scores suggest dissimilar environmental conditions (non-analogue environments resulting in model extrapolation), while positive values suggest more similar environmental conditions (interpolation). In the case of our data, we calculated MESS plots by extracting environmental values from occurrence points and projecting them to their respective time period. However, to prevent over-extrapolation, we performed an ‘informed’ MESS analysis, by clamping each environmental variable a priori based on the completeness of the response curve (ex: extrapolation can be used if the response curve for an environmental variable reaches 0 or 1 for full modeled behavior). Code for our informed MESS analysis is provided in our Supplemental Information, which was modified from the ‘mess’ function in the R package dismo (Hijmans et al. 2017).

We made habitat suitability predictions for Turritellinae using our environmental predictor variables. We reclassified continuous prediction output to binary (0 = unsuitable; 1 = suitable), determined by the maximum sum of sensitivity and specificity (MaxSSS) threshold calculated from model evaluation metrics. MaxSSS remains one of the best threshold selection methods for presence/absence data (Liu et al. 2016). We also compared our MaxSSS predictions to those calculated from the 10-percentile omission rate, which is one of the strictest binary thresholds (Phillips et al. 2017). Finally, we generated continuous predictions by transforming our raw MaxEnt predictions to a scale of 0 – 1 to approximate the probability that Turritellinae will be present at a particular location (‘cloglog’ transformation) (Phillips et al. 2017).

Niche Overlap—We compared niche overlap of the Maastrichtian and Danian occurrences using an ordination framework by first reducing dimensionality within the datasets via a Principal Component analysis (PCA). Using the ‘espace_pca’ function in the Wallace v2.0.5 and ade4 package v1.7 in R (Dray and Dufour 2007, Kass et al. 2018), we generated convex hulls of principal component space using the extracted environmental variables from the Turritellinae occurrences within both time periods and plotted with correlation loadings to infer the degree of influence particular environmental variables possess in the distribution of background points within niche space (see Supplemental Information for convex hull reconstructions). Using the ‘espace_occDens’ function in the package ecospat v3.2 (Di Cola et al. 2017) an occurrence density grid was estimated for both the environmental values at each occurrence point and background extent points using a kernel density estimation approach. Niche overlap between occurrence density grids of environmental values at each occurrence point and background points were compared using Schoener’s D (Schoener 1968) using the ‘espace_nicheOv’ function.

Using the ‘ecospat.plot.overlap.test’ function, we conducted niche similarity and equivalency tests to assess the degree of niche differentiation of Turritellinae between time periods. For the similarity test, a null distribution was generated by randomly shifting the occurrences of Turritellinae from one time period (Maastrichtian) within the combined background extent of both time periods, while keeping the niche from the other time period (Danian) fixed. This process is permuted 1000 times to assess whether the observed niche overlap was greater than expected by chance. This process was then repeated by swapping time periods, shifting occurrences from the Danian while keeping the Maastrichtian niche fixed. A p-value < 0.05 indicates that the niches from both time periods are more similar than expected. In contrast, the equivalency test pools the occurrences of Turritellinae from both time periods and randomly assigns them within their combined background extent to assess whether the niches are statistically indistinguishable; this process is also permuted 1000 times. A p-value < 0.05 suggests that the niches are significantly different from each other, indicating they are less equivalent than expected by chance.

References
Aiello-Lammens, M. E., R. A. Boria, A. Radosavljevic, B. Vilela, and R. P. Anderson. 2015. spThin: an R package for spatial thinning of species occurrence records for use in ecological niche models. Ecography 38(5):541-545.

Allmon, W. D. 1988. Ecology of Recent Turritelline gastropods (Prosobranchia, Turritellidae): current knowledge and paleontological implications. Palaios:259-284.

Allmon, W. D. 1996. Systematics and evolution of Cenozoic American Turritellidae (Mollusca: Gastropoda) I: Paleocene and Eocene coastal plain species related to" Turritella mortoni Conrad" and" Turritella humerosa Conrad". Paleontological Research Institution.

Allmon, W. D. 2011. Natural history of Turritelline gastropods (Cerithiodea: Turritellidae): a status report. Malacologia 54(1-2):159-202.

Allmon, W. D., and P. A. Cohen. 2008. Palaeoecological significance of Turritelline gastropod-dominated assemblages from the mid-Cretaceous (Albian-Cenomanian) of Texas and Oklahoma, USA. Cretaceous Research 29(1):65-77.

Allmon, W. D., and J. L. Knight. 1993. Paleoecological significance of a Turritelline gastropod-dominated assemblage in the Cretaceous of South Carolina. Journal of Paleontology 67(3):355-360.

Anderson, B. M. 2018. THE EVOLUTION OF UNUSUAL SHELL MORPHOLOGIES IN FOSSIL AND LIVING TURRITELLIDAE (GASTROPODA).

Anderson, B. M., and W. D. Allmon. 2020. High calcification rates and inferred metabolic trade-offs in the largest turritellid gastropod, Turritella abrupta (Neogene). Palaeogeography, Palaeoclimatology, Palaeoecology 544:109623.

Anderson, B. M., and W. D. Allmon. 2024. Phylogeny and systematics of fossil and recent Vermicularia (Caenogastropoda: Turritellidae). Malacologia 66(1-2):1-59.

Araújo, M. B., R. P. Anderson, A. M. Barbosa, C. M. Beale, C. F. Dormann, R. Early, R. A. Garcia, A. Guisan, L. Maiorano, B. Naimi, R. B. O’Hara, N. E. Zimmermann, and C. Rahbek. 2019. Standards for distribution models in biodiversity assessments. Science Advances 5:eaat4858.

Dall, W. H. 1913. New species of the genus Mohnia from the North Pacific. Proceedings of the Academy of Natural Sciences of Philadelphia: 501-504.

Das, S. S., S. Saha, S. Bardhan, S. Mallick, and W. D. Allmon. 2018. The oldest Turritelline gastropods: from the Oxfordian (Upper Jurassic) of Kutch, India. Journal of Paleontology 92(3):373-387.

Di Cola, V., O. Broennimann, B. Petitpierre, F. T. Breiner, M. d'Amen, C. Randin, R. Engler, J. Pottier, D. Pio, and A. Dubuis. 2017. ecospat: an R package to support spatial analyses and modeling of species niches and distributions. Ecography 40(6):774-787.

Dray, S., and A.-B. Dufour. 2007. The ade4 package: implementing the duality diagram for ecologists. Journal of Statistical Software 22(4):1-20.

Elith, J., M. Kearney, and S. Phillips. 2010. The art of modelling range‐shifting species. Methods in Ecology and Evolution 1(4):330-342.

Elith, J., S. J. Phillips, T. Hastie, M. Dudík, Y. E. Chee, and C. J. Yates. 2011. A statistical explanation of MaxEnt for ecologists. Diversity and Distributions 17(1):43-57.

Farnsworth, A., D. Lunt, C. O'Brien, G. Foster, G. Inglis, P. Markwick, R. Pancost, and S. A. Robinson. 2019. Climate sensitivity on geological timescales controlled by nonlinear feedbacks and ocean circulation. Geophysical Research Letters 46(16):9880-9889.

Farrell, Ú. C., R. Samawi, S. Anjanappa, R. Klykov, O. O. Adeboye, H. Agic, A. S. C. Ahm, T. H. Boag, F. Bowyer, and J. J. Brocks. 2021. The sedimentary geochemistry and paleoenvironments project. Geobiology 19(6):545-556.

Foote, M. 2006. Substrate affinity and diversity dynamics of Paleozoic marine animals. Paleobiology 32(3):345-366.

Foreman, S. J. 2005. Unified Model Documentation Paper Number

40, The Ocean Model, Report, The Met.

Foster, G. L., D. L. Royer, and D. J. Lunt. 2017. Future climate forcing potentially without precedent in the last 420 million years. Nature Communications 8(1):1-8.

Friend, D. S., B. M. Anderson, E. Altier, S. Sang, E. Petsios, R. W. Portell, and W. D. Allmon. 2023. SYSTEMATICS AND PHYLOGENY OF PLIO-PLEISTOCENE SPECIES OF TURRITELLIDAE (GASTROPODA) FROM FLORIDA AND THE ATLANTIC COASTAL PLAIN. Bulletins of American Paleontology (402).

Graham, M. H. 2003. Confronting multicollinearity in ecological multiple regression. Ecology 84(11):2809-2815.

HALL, J. 1845. Description of organic remains collected by Capt. JC Fremont in the Geographical Survey of Oregon and North California: Twenty-Eighth Cong., 2nd Sess., House Ex. Doc 166:304-7.

Harzhauser, M., and B. Landau. 2019. Turritellidae (Gastropoda) of the Miocene Paratethys Sea with considerations about turritellid genera. Zootaxa 4681(1):1–136-1–136.

Hendricks, J. R., E. E. Saupe, C. E. Myers, E. J. Hermsen, and W. D. Allmon. 2014. The generification of the fossil record. Paleobiology 40(4):511-528.

Hernandez, P. A., C. H. Graham, L. L. Master, and D. L. Albert. 2006. The effect of sample size and species characteristics on performance of different species distribution modeling methods. Ecography 29(5):773-785.

Hijmans, R. J., R. Bivand, K. Forner, J. Ooms, E. Pebesma, and M. D. Sumner. 2022. Package ‘terra’. Maintainer: Vienna, Austria.

Hijmans, R. J., S. Phillips, L. Leathwick, and J. Elith. 2017. Package ‘dismo’. Circles 9.1:1-68.

Hijmans, R. J., J. van Etten, M. Mattiuzzi, M. Sumner, J. Greenberg, O. Lamigueiro, A. Bevan, E. Racine, and A. Shortridge. 2021. Raster package in R. Version.

Hopkins, M. J. 2014. The environmental structure of trilobite morphological disparity. Paleobiology 40(3):352-373.

Johnson, J. M., and K. C. Clarke. 2021. An area preserving method for improved categorical raster resampling. Cartography and Geographic Information Science 48(4):292-304.

Kass, J. M., B. Vilela, M. E. Aiello‐Lammens, R. Muscarella, C. Merow, and R. P. Anderson. 2018. Wallace: A flexible platform for reproducible modeling of species niches and distributions built for community expansion. Methods in Ecology and Evolution 9(4):1151-1156.

Kiessling, W., E. Flügel, and J. Golonka. 1999. Paleoreef maps: evaluation of a comprehensive database on Phanerozoic reefs. AAPG bulletin 83(10):1552-1587.

Kocsis, Á., N. Raja, and S. Williams. 2023. rgplates: R interface for the GPlates Web service and desktop application.

Lieberman, B. S., W. D. Allmon, and N. Eldredge. 1993. Levels of selection and macroevolutionary patterns in the turritellid gastropods. Paleobiology 19(2):205-215.

Liu, C., G. Newell, and M. White. 2016. On the selection of thresholds for predicting species occurrence with presence‐only data. Ecology and Evolution 6(1):337-348.

Lunt, D. J., A. Farnsworth, C. Loptson, G. L. Foster, P. Markwick, C. L. O'Brien, R. D. Pancost, S. A. Robinson, and N. Wrobel. 2016. Palaeogeographic controls on climate and proxy interpretation. Climate of the Past 12(5):1181-1198.

Marwick, J. 1957. Generic revision of the Turritellidae. Journal of Molluscan Studies 32(4):144-166.

Myers, C. E., R. A. MacKenzie, and B. S. Lieberman. 2013. Greenhouse biogeography: the relationship of geographic range to invasion and extinction in the Cretaceous Western Interior Seaway. Paleobiology 39(1):135-148.

Myers, C. E., A. L. Stigall, and B. S. Lieberman. 2015. PaleoENM: applying ecological niche modeling to the fossil record. Paleobiology 41(2):226-244.

Naimi, B. 2017. Package ‘usdm’. Uncertainty Analysis for Species Distribution Models. Wien: www. cran. r-project. org.

Pearson, R. G., C. J. Raxworthy, M. Nakamura, and A. Townsend Peterson. 2007. Predicting species distributions from small numbers of occurrence records: a test case using cryptic geckos in Madagascar. Journal of Biogeography 34(1):102-117.

Peterson, A. T., and J. Soberón. 2012. Species Distribution Modeling and Ecological Niche Modeling: Getting the Concepts Right. Natureza & Conservação 10(2):102-107.

Phillips, S. J. 2021. A brief tutorial on Maxent. AT&T Research 190(4):231-259.

Phillips, S. J., R. P. Anderson, M. Dudík, R. E. Schapire, and M. E. Blair. 2017. Opening the black box: an open-source release of Maxent. Ecography 40(7):887-893.

Phillips, S. J., and M. Dudík. 2008. Modeling of species distributions with Maxent: new extensions and a comprehensive evaluation. Ecography 31(2):161-175.

Plotnick, R. E., and P. J. Wagner. 2006. Round up the usual suspects: common genera in the fossil record and the nature of wastebasket taxa. Paleobiology 32(1):126-146.

Qiao, H., A. T. Peterson, L. Ji, and J. Hu. 2017. Using data from related species to overcome spatial sampling bias and associated limitations in ecological niche modelling. Methods in Ecology and Evolution 8(12):1804-1812.

Radosavljevic, A., and R. P. Anderson. 2014. Making better Maxent models of species distributions: complexity, overfitting and evaluation. Journal of Biogeography 41(4):629-643.

Sagoo, N., P. Valdes, R. Flecker, and L. J. Gregoire. 2013. The Early Eocene equable climate problem: can perturbations of climate model parameters identify possible solutions? Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences 371(2001):20130123.

Sang, S., D. S. Friend, W. D. Allmon, and B. M. Anderson. 2019. Protoconch enlargement in Western Atlantic turritelline gastropod species following the closure of the Central American Seaway. Ecology and Evolution 9(9):5309-5323.

Schoener, T. W. 1968. The anolis lizards of bimini: resource partitioning in a complex fauna. Ecology 49(4):704-726.

Scotese, C. R. 2016. Tutorial: PALEOMAP paleoAtlas for GPlates and the paleoData plotter program. Technical Report, 56. Available at: https://www. earthbyte. org/paleomap ….

Scotese, C. R., and N. Wright. 2018. PALEOMAP paleodigital elevation models (PaleoDEMS) for the Phanerozoic. PALEOMAP Proj.

Shcheglovitova, M., and R. P. Anderson. 2013. Estimating optimal complexity for ecological niche models: a jackknife approach for species with small sample sizes. Ecological Modelling 269:9-17.

Soberón, J., and M. Nakamura. 2009. Niches and distributional areas: concepts, methods, and assumptions. Proceedings of the National Academy of Sciences 106(supplement_2):19644-19650.

Team, R. C. 2021. R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. 2012.

Valavi, R., G. Guillera‐Arroita, J. J. Lahoz‐Monfort, and J. Elith. 2021. Predictive performance of presence‐only species distribution models: a benchmark study with reproducible code. Ecological Monographs 92(0):e01486.

Valdes, P. J., E. Armstrong, M. P. Badger, C. D. Bradshaw, F. Bragg, M. Crucifix, T. Davies-Barnard, J. J. Day, A. Farnsworth, and C. Gordon. 2017. The BRIDGE HadCM3 family of climate models: HadCM3@ Bristol v1. 0. Geoscientific Model Development 10(10):3715-3743.

Valdes, P. J., C. R. Scotese, and D. J. Lunt. 2021. Deep ocean temperatures through time. Climate of the Past 17(4):1483-1506.

van Proosdij, A. S., M. S. Sosef, J. J. Wieringa, and N. Raes. 2016. Minimum required number of specimen records to develop accurate species distribution models. Ecography 39(6):542-552.

Warren, D. L., and S. N. Seifert. 2011. Ecological niche modeling in Maxent: the importance of model complexity and the performance of model selection criteria. Ecological Applications 21(2):335-342.

Yan, H., L. Feng, Y. Zhao, L. Feng, D. Wu, and C. Zhu. 2020. Prediction of the spatial distribution of Alternanthera philoxeroides in China based on ArcGIS and MaxEnt. Global Ecology and Conservation 21:e00856.