Turgidodon lillegraveni OMNH

Brannick, Alexandria L., Fulghum, Henry Z., Grossnickle, David M. & Wilson Mantilla, Gregory P., 2023, Dental ecomorphology and macroevolutionary patterns of North American Late Cretaceous metatherians, Palaeontologia Electronica (a 48) 26 (3), pp. 1-42 : 7-11

Turgidodon lillegraveni OMNH 20117 M2(?)

Turgidodon madseni OMNH 20538 M3

Turgidodon praesagus UALVP 55849 M3

Turgidodon praesagus UCMP 122168 M3

Turgidodon praesagus UCMP 131345 M2

Turgidodon rhaister UCMP 47366 M2

Turgidodon russelli UALVP 55852 M3

Turgidodon russelli UALVP 6983 M3

Varalphadon creber UALVP 29525 M3

Varalphadon creber UALVP 29527 M3

Varalphadon creber UALVP 5529 M3

Varalphadon UALVP 5544 M3 wahweapensis

Atokatheridium boreni OMNH 61623 M2

Glasbius intricatus UCMP 102111 M3

Glasbius twitchelli UCMP 153679 M3

Glasbius twitchelli UCMP 156143 M3

Glasbius twitchelli UCMP 224090 M3

Nortedelphys intermedius UCMP 134776 M3

Nortedelphys jasoni UCMP 174506 M3

Nortedelphys jasoni UCMP 177838 M3

Nortedelphys magnus UA 2846 M3

Nortedelphys minimus UCMP 52715 M2

Nortedelphys minimus UCMP 72211 M3

Anchistodelphys archibaldi OMNH 21033 M3

Apistodon exiguus UALVP 29693 M3

Dakotadens morrowi OMNH 49450 Mx

Alphadontidae Judithian cast

Alphadontidae Judithian cast

Alphadontidae Aquilan fossil

Alphadontidae Aquilan fossil

Alphadontidae Aquilan fossil

Alphadontidae Judithian fossil

Alphadontidae Lancian fossil

Alphadontidae Lancian fossil

Alphadontidae Judithian cast

Alphadontidae Lancian fossil

Alphadontidae pre-Aquilan fossil

Alphadontidae pre-Aquilan cast

Alphadontidae pre-Aquilan fossil

Alphadontidae Lancian fossil

Alphadontidae Lancian fossil

Alphadontidae Lancian fossil

Alphadontidae Judithian cast

Alphadontidae Judithian cast

Alphadontidae Judithian fossil

Alphadontidae Judithian fossil

Alphadontidae Judithian fossil

Alphadontidae Lancian fossil

Alphadontidae Judithian fossil

Alphadontidae Judithian fossil

Alphadontidae Aquilan fossil

Alphadontidae Aquilan fossil

Alphadontidae Aquilan fossil

Alphadontidae Aquilan fossil

Deltatheridiidae pre-Aquilan cast

Glasbiidae Lancian fossil

Glasbiidae Lancian fossil

Glasbiidae Lancian fossil

Glasbiidae Lancian fossil

Herpetotheriidae Lancian fossil

Herpetotheriidae Lancian fossil

Herpetotheriidae Lancian fossil

Herpetotheriidae Lancian cast

Herpetotheriidae Lancian fossil

Herpetotheriidae Lancian fossil incertae sedis Aquilan cast incertae sedis Aquilan fossil incertae sedis pre-Aquilan cast include in our sample any extant mammal species with teeth that have, to our knowledge, secondary wear-induced functionality. Specimens with cusps missing due to breakage were also excluded.

Dietary Categories

We classified each extant species in our dataset into one of six dietary categories: carnivory (carn), animal-dominated omnivory (ado), plant-dominated omnivory (pdo), frugivory (frug), invertivory (inv; i.e., ‘insectivory’), or soft-invertebrate specialist (sis) ( Table 2). We used six specific dietary categories rather than the classic three-diet classification scheme (herbivory-omnivory-carnivory) to provide more detailed dietary information and to avoid oversimplification ( Pineda-Munoz and Alroy, 2014). Our choice of diet categories follows Smith (2017), who used these six categories (along with folivory and hard-object invertivory) in DTA analyses of lower molars. Our ‘soft-invertebrate specialist’ group includes taxa that primarily consume soft invertebrates such as earthworms and slugs, whereas our ‘invertivory’ (i.e., ‘insectivory’) group includes taxa that primarily eat relatively harder-bodied insects, such as beetles and moths. Following Pineda-Munoz and Alroy (2014), we classified diets of each species, with emphasis on its primary food resource. A species was classified as a specialist (i.e., non-omnivore) if one food resource makes up 50% or more of its total diet. For dietary information, we used online archives (EltonTraits [ Wilman et al., 2014] and Mammal DIET [Kissling et al., 2014]) and a natural history compendia ( Nowak, 1999). We supplemented each classification with primary literature sources (see Table 2 for sources), which were especially important when species-level information was extrapolated from genus-level information in the online archives (see Kissling et al., 2014; Table 2).

We acknowledge that our decision to use six dietary categories rather than the classic ‘carnivore-omnivore-herbivore’ trophic classification could lead to greater overlap of categories in the morphospace and less power to predict diet. We classified the diet of some extant taxa in our sample differently than previous studies have. For example, Nasua narica (white-nosed coati) is known to eat insects, but it is strictly frugivorous when fruit is available (e.g., Nowak, 1999). Although some studies classified its diet as plant-dominated omnivory ( Smith, 2017), we followed EltonTraits, which records its diet as 70% fruit and considered this taxon a frugivore. We recognize that in this and any large-scale study of mammalian feeding behaviors, decisions that reduce the complexity of dietary data into discrete categories could have an impact on the results.

Fossil Metatherian Sampling

We sampled 71 isolated upper molars of 42 species (22 genera; six major clades) of North American Late Cretaceous (NALK) metatherians from the Western Interior region ( Table 3). Our sample includes two stagodontids, one deltatheriid, two glasbiids, eight pediomyids, six taxa classified as incertae sedis, four herpetotheriids, and 19 alphadontids. To increase our taxonomic sampling of Cretaceous metatherians, we substituted the M2 (which tends to be morphologically very similar to the M3) for some species that did not have an available M3, and we used upper molar specimens of uncertain position (i.e., “Mx”) for some species that did not have definitive M2 or M 3 specimens available (see Table 3 for details). Our sample includes 62% of the known species of NALK metatherians (42 of 68 known species; Case et al., 2005; Williamson et al., 2014; Cohen, 2018; Cohen et al., 2020). Some species were omitted from our sample because of either a lack of a well-preserved upper molar in the fossil record or an appropriate specimen was not available for loan. Our sampling of deltatheriids and stagodontids is limited, and this likely artificially reduced both morphological disparity values and morphospace occupation (e.g., Nanocuris has been interpreted as a specialized carnivore), especially in the pre-Aquilan and Lancian time bins (see below). The absence of other taxa may have had a negligible effect on the results because their morphologies are approximated by other sampled taxa (e.g., the absence of the pediomyoid genus Aquiladelphis may be accounted for by the presence of other pediomyoid genera in our sample to some degree).

We assigned each fossil species in our sample to one (or more) of four time bins depending on the known temporal range of each species ( Williamson et al., 2014), using a range-through approach. Three bins are Cretaceous NALMAs (Woodburne, 2004): Aquilan (ca. 86–79 Ma), Judithian (79–69 Ma), and Lancian (69–66 Ma). We binned the eight specimens from geologic units that pre-date the Aquilan NALMA into a “pre-Aquilan” time bin (ca.126–86 Ma). Most taxa that we assigned to the pre-Aquilan time bin are from 100– 86 Ma, but we also include Atokatheridium , which has a range of ca. 126–100 Ma. Because the “Edmontonian” NALMA is poorly characterized and not well sampled (Cifelli et al., 2004), we lumped the “Edmontonian” taxa into the Judithian bin. We recognize that these time bins are uneven in duration and that the longer duration bins could artificially inflate measures of disparity and diversity; however, we were unable to more finely and precisely bin our data due to uneven sample sizes across time bins and the lack of high-precision ages for certain geologic units.

Collection of 3D Tooth Surface Data

Three-dimensional digital models of the sampled teeth were created using micro-computed tomography (μCT) scan data. We scanned original specimens of teeth, molds of teeth, and epoxy casts of teeth ( Tables 2–3). López-Torres et al. (2017) found that OPCR values of epoxy casts tend to be higher than those from their original specimens due to potential for artificially rougher surfaces on the casts (both DNE and RFI are more robust to this effect). Thus, we interpret OPCR results for the relatively few casts in our sample (15 of 71 specimens) with caution. Specimens were scanned using either a Bruker Skyscan 1172, Skyscan 1173, or NSI X5000 scanner, all of which are housed on the University of Washington campuses. We also downloaded image stacks (TIFF format) of scan data for eight extant specimens ( Table 2) from the MorphoSource online repository (morphosource.org) to bolster our modern comparative dataset (Appendix 1). For detailed information regarding scanner types and scan settings, see Appendix 2. Molds of extant teeth were made using Coltene President Plus polyvinylsiloxane (type 2, medium consistency), and epoxy casts were collected from the UWBM, University of California Museum of Paleontology, and Sam Noble Oklahoma Museum of Natural History collections. For specimens scanned with Bruker Skyscan scanners, scan data were reconstructed using NRecon (Bruker microCT, Belgium); scans completed using the NSI X5000 were reconstructed using efX Reconstruction (North Star Imaging, Inc.). We segmented raw scan data using Avizo Lite 9.2.0 (Thermo Fisher Scientific). We then removed artifacts (“cleaning”), cropped, and oriented tooth models using GeoMagic Studio (3DS Systems). Specimen models were cropped to include the entire enamel cap (EEC cropping method; see Berthaume et al., 2019 for details) and were oriented such that the occlusal plane is perpendicular to the Z-axis. We exported the cleaned and oriented 3D tooth models from GeoMagic Studio as PLY files. These PLY files were imported back into Avizo Lite 9.2.0, and the 3D tooth models were simplified to 20,000 faces using the Simplification Editor tool. We then used the “Remesh Surface” function to downsample the tooth models to ~10,000 faces. The remesh function was used because it reduces the chance that surfaces with extremely disparate polygon mesh face-sizes are produced during simplification (Spradley, personal comm., 2018). We then used the “Smooth Surface” function with 25 iterations and lambda = 0.6 ( Spradley et al., 2017; Spradley, personal comm., 2018). Because the consistency of model creation and processing is extremely important for producing comparable DTA results ( Spradley et al., 2017; Berthaume et al., 2019), we used the same workflow for the creation of all models in this study. The resulting smoothed tooth models were saved as PLY files and used in our DTA analyses.

Dental Topographic Analyses (DTA)

We computed RFI, DNE, and OPCR for all 3D tooth models using the molaR_Batch function from the package molaR, version 4.2 ( Pampush et al., 2016), in R version 3.3.3 (R Core Team, 2017). RFI is the ratio between the 3D surface area of a tooth crown and the 2D “footprint” area of a tooth ( Ungar and M’Kirera, 2003). We use a modified version of this ratio in which the entire tooth crown is more accurately considered (Boyer, 2008). The modified RFI calculation is:

(A3D = 3D embedded surface area of the tooth crown, A2D = 2D tooth crown footprint area in occlusal view; Boyer, 2008; López-Torres et al., 2017). DNE represents the curvature of the tooth crown by calculating the sum energy values across the entire occlusal surface (Bunn et al., 2011; Winchester et al., 2014; Winchester, 2016). OPCR measures tooth crown complexity by calculating the total number of patches, or “tools,” on the crown of a tooth. A patch is a contiguous group of pixels that face the same cardinal direction on the tooth model (Evans et al., 2007; Evans and Jernvall, 2009; Wilson et al., 2012). Parameters for each metric were set as follows: RFI—alpha = 0.15; DNE—boundary discard = “Vertex”; and OPCR—step size = 8 and minimum patch size = 3 pixels (Evans et al., 2007; Pampush et al., 2016; Smith, 2017; Spradley, 2017). We ran a second DTA with the OPCR minimum patch size = 5 pixels to minimize any “noise” that might artificially inflate values for our extant and fossil samples, which include molds and casts, respectively (Winchester, 2016; López-Torres et al., 2017).

We log-transformed our DTA data to reduce skew. We generated scatter biplots of all possible combinations of the dental metrics to visualize morphospace occupation of extant dietary groups. We then plotted our fossil metatherian DTA values within the same morphospace of the extant dataset to examine both how fossil morphospace occupation compared to extant mammal morphospace occupation and how fossil morphospace occupation changed through time. We also tested for correlation between our DTA metrics by calculating Spearman’s rho and using least-squares linear regressions.

To test for differences between DTA values of the six dietary groups, we used one-way analysis of variance (ANOVA) and Tukey’s honest significant difference (HSD) post hoc test. We also performed a MANOVA using all three DTA metrics as independent variables.

Dietary Inference of Fossil Metatherians

To quantitatively infer diet in our sample of fossil metatherians, we conducted a discriminant function analysis (DFA) using the function lda() from the package MASS (Venables and Ripley, 2002). We first used the extant comparative dataset and a leave-one-out cross validation to assess the accuracy of discriminant functions in predicting diet (see MASS package documentation for more information). We then applied this DFA to the fossil metatherian DTA data (with fossils treated as having unknown diets), and we used posterior probabilities of dietary groupings to infer fossil diets. In a second permutation, we conducted a DFA on the extant comparative dataset using both the DTA data and mean body mass (compiled from the primary literature) to test whether this would significantly improve the discriminatory power of our model (Winchester et al., 2014). Because the resulting accuracy did not significantly improve discrimination, we only report the results from the first permutation. For fossil species in which different specimens were classified differently by the DFA (e.g., Didelphodon vorax ), we based our dietary inferences on additional evidence, such as which diet was most commonly reconstructed by the specimens and evidence from previous studies, or we simply report two possible diet classifications for the species.

Dental Disparity of Fossil Metatherians

We calculated morphological disparity in our sample of fossil metatherians as: i) intra-family disparity and ii) total disparity per time bin. We did not calculate the intra-family disparity per time bin because sample sizes were too small. All disparity calculations used mean species values of each standardized, log-transformed DTA metric. We measured disparity as both the variance of each DTA metric and the sum of variances (Ciampaglio et al., 2001) using the morphol.disparity function in the geomorph package in R (Adams et al., 2020), which calculates a simulation-based P -value for statistical comparison between groups (i.e., between families or between time bins). We generated 95% confidence intervals using a custom bootstrapping function in R with 1,000 replicates.

Testing for Phylogenetic Signal

We tested for phylogenetic signal in the DTA results of our extant comparative dataset using a phylogenetic tree that we generated suing TimeTree (www.timetree.org; Kumar et al., 2017). We calculated Blomberg’s K (Blomberg et al., 2003) and Pagel’s lambda ( Pagel, 1992) using the phylosig function in the package phytools in R ( Revell, 2012). We did not test for phylogenetic signal in the DTA results of our fossil taxa because the most recent species-level phylogeny that includes all of the fossil taxa in our sample ( Williamson et al., 2014) is highly unresolved with a large polytomy.