OUP user menu

A Dated Phylogeny of Marsupials Using a Molecular Supermatrix and Multiple Fossil Constraints

Robin M. D. Beck
DOI: http://dx.doi.org/10.1644/06-MAMM-A-437.1 175-189 First published online: 19 February 2008


Phylogenetic relationships within marsupials were investigated based on a 20.1-kilobase molecular supermatrix comprising 7 nuclear and 15 mitochondrial genes analyzed using both maximum likelihood and Bayesian approaches and 3 different partitioning strategies. The study revealed that base composition bias in the 3rd codon positions of mitochondrial genes misled even the partitioned maximum-likelihood analyses, whereas Bayesian analyses were less affected. After correcting for base composition bias, monophyly of the currently recognized marsupial orders, of Australidelphia, and of a clade comprising Dasyuromorphia, Notoryctes, and Peramelemorphia, were supported strongly by both Bayesian posterior probabilities and maximum-likelihood bootstrap values. Monophyly of the Australasian marsupials, of Notoryctes + Dasyuromorphia, and of Caenolestes + Australidelphia were less well supported. Within Diprotodontia, Burramyidae + Phalangeridae received relatively strong support. Divergence dates calculated using a Bayesian relaxed molecular clock and multiple age constraints suggested at least 3 independent dispersals of marsupials from North to South America during the Late Cretaceous or early Paleocene. Within the Australasian clade, the macropodine radiation, the divergence of phascogaline and dasyurine dasyurids, and the divergence of perameline and peroryctine peramelemorphians all coincided with periods of significant environmental change during the Miocene. An analysis of “unrepresented basal branch lengths” suggests that the fossil record is particularly poor for didelphids and most groups within the Australasian radiation.

Key words
  • Ameridelphia
  • Australidelphia
  • Bayesian analysis
  • fossil record
  • marsupials
  • phylogenetic fuse
  • phylogeny
  • relaxed molecular clock
  • supermatrix

Molecular analyses of the higher-level phylogeny of marsupials have a long history (Nuttall 1904; Weymss 1953) and have been more common than equivalent morphological studies. Several of these molecular analyses have added considerably to our current understanding of marsupial phylogeny and classification (e.g., Kirsch 1968, 1977; Sarich et al. 1982). However, even with recent advances in sequencing and phylogenetic methods, different molecular data sets continue to support incon-gruent topologies (e.g., Amrine-Madsen et al. 2003; Asher et al. 2004; Baker et al. 2004; Nilsson et al. 2004); hence, the higher-level phylogeny of marsupials remains uncertain.

The affinities of several marsupial groups have proven particularly difficult to resolve. It is unclear whether the South American didelphid opossums (Didelphimorphia) and caenolestid shrew opossums (Paucituberculata) comprise a monophyletic group (Ameridelphia) or are a paraphyletic assemblage at the base of modern marsupials. Current evidence supports the monophyly of Australidelphia (Amrine-Madsen et al. 2003; Cardillo et al. 2004; Horovitz and Sanchez-Villagra 2003; Szalay 1982, 1994), a clade comprising the modern Australasian marsupial orders Peramelemorphia, Notoryctemorphia, Dasyuromorphia, and Diprotodontia and the South American monito del monte Dromiciops (Microbiotheria). The position of Dromiciops within Australidelphia remains unresolved and yet is crucial to our understanding of marsupial biogeography: if Dromiciops is the sister group to all Australasian marsupials (Amrine-Madsen et al. 2003; Phillips et al. 2006), then only a single dispersal of marsupials from South America to Australasia is implied; if Dromiciops is nested within the Australasian radiation (Asher et al. 2004; Cardillo et al. 2004; Nilsson et al. 2004), there must have been multiple dispersals of marsupials to Australasia or back-migration of microbiotheres from Australasia to South America. Within Australasian marsupials, the affinities of the highly autapomorphic, fossorial marsupial mole Notoryctes (Notoryctemorphia) are uncertain, as are those of the bandicoots (Peramelemorphia), with some studies (e.g., Kirsch et al. 1997) suggesting that the latter may not be members of Australidelphia.

Together with a robust phylogeny, divergence dates are essential for understanding the patterns, processes, and tempo of marsupial evolution and biogeography. For example, they enable assessment of whether divergences coincided with periods of environmental change, such as climatic variations (e.g., Delsuc et al. 2004; Douady and Douzery 2003) and fluctuations in sea level (Johnson et al. 2006; Mercer and Roth 2003). However, robust dates for divergences within Marsupialia have been unavailable because of incompleteness of the fossil record, notably a lack of sites from the Late Cretaceous and early Tertiary when most of the deep divergences are thought to have occurred.

Late Cretaceous marsupials from North America were relatively diverse and include taxa that have been tentatively assigned to extant lineages by some authors (Case et al. 2005; Kielan-Jaworowska et al. 2004), but they are known only from fragmentary material (largely isolated teeth) and their relationships to the marsupial crown-group remain controversial. In South America, there are no Late Cretaceous fossil sites containing marsupials; the oldest known South American marsupial is a single tooth from the earliest Paleocene (Goin et al. 2006). The Tiupampa fauna in Bolivia (Marshall and Muizon 1988; Muizon 1991) contains exceptionally well-preserved, articulated skeletons of several plesiomorphic marsupials, as well as more fragmentary remains of several other marsupial taxa. However, the Tiupampa fauna may be somewhat younger than usually assumed (Benton and Donoghue 2007; Gradstein et al. 2004), being middle (59.2–60.4 million years ago [mya]) rather than early (63–64.5 mya) Paleocene. In Australia, there is only a single mammalian fauna from the early Tertiary (the early Eocene [54.6 mya] Tingamarra fauna—Godthelp et al. 1992), and marsupial fossils from this site are highly fragmentary. The next oldest Australian sites contain diverse and often well-preserved marsupial fossils, but are some 30 million years younger (late Oligocene—Archer et al. 1999; Long et al. 2002). Attempts to determine divergence dates from fossil evidence also may be confounded by the existence of “phylogenetic fuses” (Cooper and Fortey 1998): morphological apomorphies that distinguish 2 lineages may not have evolved until a considerable time after those lineages diverged.

These factors suggest that molecular methods, which are not strictly tied to the known fossil record (although they do usually require fossils to act as calibration points), seem more likely to provide accurate divergence dates. However, previous molecular dates for marsupials have been based on dubious assumptions, such as rate homogeneity or the use of a single (often questionable) fossil calibration. Indeed, in a number of instances, molecular dates for particular marsupial clades have been considerably younger than their oldest fossil members (e.g., Nilsson et al. 2004; Sarich et al. 1982), assuming, of course, that these fossils have been correctly identified.

In an attempt to resolve these problems, I have combined molecular data from 3 recent phylogenetic analyses of marsupials (Amrine-Madsen et al. 2003; Asher et al. 2004; Nilsson et al. 2004) to produce a supermatrix of 7 nuclear and 15 mitochondrial genes. This is, to my knowledge, the largest molecular data set used so far to investigate higher-level marsupial relationships, and it is the 1st to use a Bayesian relaxed molecular clock approach to calculate divergence dates between the major groups of marsupials.

Materials and Methods

Initial data set.—Sequence data were taken from the following studies: Amrine-Madsen et al. (2003; matrix supplied by M. Springer), who sampled the nuclear genes apolipoprotein B (APOB), breast cancer 1 (BRCA1), interphotoreceptor retinoid binding-protein (IRBP), recombination activating gene 1 (RAG1), and Von Willebrand factor (VWF); Asher et al. (2004; matrix downloaded from doi:10.1016/j.ympev.2004. 05.004), who sampled the nuclear genes IRBP, phosphoglycerate kinase 1 (PGK1), and protamine 1 (P1) and the mitochondrial genes cytochrome b (CYTB), 12S rRNA, 16S rRNA, tRNA valine, and reduced nicotinamide adenine dinucleotide dehydrogenase 2 (NADH2); and Nilsson et al. (2004; matrix supplied by M. Nilsson), who sampled the 12 mitochondrial H-strand protein-coding genes. The original alignments of all 3 matrices were maintained; for the data set of Asher et al. (2004), the original authors’ static alignment produced by the program CLUSTAL was preferred over any of their direct-optimization alignments produced by the program POY (De Laet and Wheeler 2003). Portions of the matrix of Asher et al. (2004) that overlapped with those of the other 2 studies (i.e., IRBP, CYTB, and NADH2) were deleted, because taxon sampling was sparser in the former study. In total, sequences from 7 nuclear genes (APOB, BRCA1, IRBP, PGK1, P1, RAG1, and VWF) and 15 mitochondrial loci (12S rRNA, 16S rRNA, tRNA valine, and the 12 H-strand protein-coding genes) were included in this study (Appendix I). The “Macropodinae” terminal of Amrine-Madsen et al. (2003) was decomposed into separate Macropus and Dendrolagus sequences and combined with sequences for these 2 genera from the other 2 studies. Sequences of 3 genes, P1, CYTB, and 12S rRNA, for the thylacine Thylacinus were downloaded from GenBank and added manually to the alignment. In total, 22 marsupial ingroup taxa were included, with 2 monotreme (Ornithorhynchus and Tachyglossus) and 2 placental (Eulipotyphla and Xenarthra) outgroups. Five of the marsupial taxa comprised sequences from more than 1 genus: Didelphini (Didelphis and Lutreolina), Peramelinae (Perameles and Isoodon), Petauridae (Petaurus and Dactylopsila), Phascogalinae (Phascogale and Antechinus), and Pseudocheiridae (Pseudocheirus and Pseudochirops).

Choice of taxa for this study essentially replicates that of Asher et al. (2004). Some taxa suffer from relatively large amounts of missing data (notably Dasyuroides and Thylacinus, which are represented by only 2 and 3 genes, respectively), but all share at least 1 nuclear (P1) and 1 mitochondrial (12S rRNA) gene in common, thus avoiding the problem of non-overlapping sequences identified by Springer et al. (2004). The supermatrix also can be easily combined with the morphological data set of Horovitz and Sanchez-Villagra (2003) for total evidence analyses. Because recent work by Phillips et al. (2006) suggests that base composition bias might not be accounted for by current likelihood models, 2 versions of the supermatrix were analyzed, 1 including 3rd codon positions of mitochondrial protein-coding genes (“full”; 20,121 base pairs [bp]) and 1 that excludes them (“no mt3”; 16,662 bp).

Data partitioning.—Three different partitioning strategies were applied to both the full and no mt3 data sets: a single partition for the entire supermatrix (“unpartitioned”); each gene assigned its own partition (“gene-partitioned”); and each codon position within each protein-coding gene and stems and loops within each ribosomal gene assigned their own partitions (“genes, codons, stems, loops [GCSL]-partitioned”). In the GCSL-partitioned analyses, the 250-bp noncoding region of P1 identified by Queralt et al. (1995; see also Asher et al. 2004) also was assigned its own partition, and stem and loop regions of the ribosomal genes were defined using secondary structure masks from the Organellar Genome Retrieval (OGRe) database (Jameson et al. 2003). The 12 mitochondrial protein-coding genes were treated as a single supergene (which was further partitioned by codon position in the GCSL-partitioned analyses), rather than individual genes to avoid a very large number of partitions. Too many partitions risk overparameterization and greatly increase computation time and difficulty in reaching convergence in Bayesian analyses. Furthermore, base composition does not vary significantly among the 12 mitochondrial protein-coding genes (Nilsson 2006), so it seems reasonable to treat them as a single supergene.

Phylogenetic analysis.—Maximum-likelihood analyses were carried out using RAxML VI (Stamatakis 2006), which can implement partitioned analysis by applying 1 of 2 models to each partition: GTRGAMMA (general time reversible model with rate heterogeneity accommodated by a gamma distribution) or GTRCAT (as for GTRGAMMA, but rate heterogeneity is accommodated by a number of discrete rate categories). To maximize computational efficiency, the GTRMIX option of RAxML was used; this assumes the GTRCAT model (which is faster and less memory intensive than GTRGAMMA) when searching for tree topologies, but assumes the GTRGAMMA model when calculating the likelihood score of each topology. The RAxML analyses each comprised 100 tree search replicates (assuming default parameters) and were implemented using the perl script batchRAxML.

Bayesian analyses used MrBayes 3.1.2 (Ronquist and Huelsenbeck 2003), which, unlike RAxML, allows specification of 1 of a wide range of models for each partition. Models for the Bayesian analyses were identified using the model selection program MrModeltest 2.2 (Nylander 2004), and models preferred by the Akaike information criterion were implemented (following Posada and Buckley [2004]; a list of models used for each partition is available from the author). For the MrBayes analyses, all parameters except topology were unlinked across partitions, and 2 independent runs (each comprising 1 “cold” and 3 “heated” chains) were run simultaneously, with trees sampled every 100th generation. The MrBayes analyses were run for either 1 × 106 (unpartitioned), 2.5 × 106 (gene-partitioned), or 12 × 106 (GCSL-partitioned) generations for both the full and no mt3 data sets. In all cases, stationarity had been reached by the end of the analysis (average standard deviation of split frequencies <0.01). Majority-rule consensus trees were constructed, with a burn-in period of either 1 × 105 (unpartitioned), 2.5 × 105 (gene-partitioned), or 10 × 106 (GCSL-partitioned) generations excluded. Bayes factors (Nylander et al. 2004), calculated as 2 times the difference in the harmonic means of the log likelihoods of the post–burn-in trees (obtained using the “sump” command in MrBayes), were examined to investigate whether particular partitioning schemes show a better fit to the data than do others.

Support values.—Support for different clades was calculated by 1,000 bootstrap replicates for the maximum-likelihood analyses (again using RAxML and batchRAxML together with the perl script bootStrip), and by posterior probabilities for the Bayesian analyses. Studies indicate that posterior probabilities may be misleadingly high, whereas bootstrap values are more conservative (Suzuki et al. 2002; Taylor and Piel 2004). However, Hall and Salipante (2007) have argued recently that neither measure reflects the probability of a particular clade. Statistical tests of alternative topologies may therefore be more informative than either of these support measures (Shimodaira and Hasegawa 1999).

Testing alternative topologies.—The different topologies found in the various analyses, as well as a number of others taken from other studies, were tested using the approximately unbiased (AU—Shimodaira 2002) and Kishino–Hasegawa (KH—Kishino and Hasegawa 1989) tests as implemented in CONSEL (Shimodaira and Hasegawa 2001). Both the full and no mt3 data sets were used for these tests. PAUP* version 4.0b10 (Swofford 2002) was used to calculate the site likelihoods for each of the test topologies, with the GCSL-partitioning scheme assumed and the appropriate model for each partition specified using the output from MrModeltest 2.2. The CONSEL analyses employed 10 batches of 1 × 105 bootstrap replicates and 1 of 1 × 106 bootstrap replicates.

Molecular dating analysis.—A Bayesian relaxed molecular clock method was used, as implemented by the program Multidivtime (Kishino et al. 2001; Thorne et al. 1998). General methodology follows Rutschmann (2005; see also Inoue et al. 2005), with a single F84 + gamma model applied to the entire supermatrix (incomplete overlap of sequences prevented a separate model being applied to each partition) and maximum-likelihood parameters estimated under PAML version 3.15 (Yang 1997). Because a partitioned approach could not be implemented in Multidivtime, only the no mt3 data set (which excludes a partition with a known base composition bias) was used for the calculation of dates. Following Benton and Donoghue (2007), a high number of constraints was specified (23 minima and 3 maxima) and a conservative approach was used when assigning dates: if based on absolute (radiometric) dating, dates were assumed to be the lower end of the 95% confidence interval (95% CI), where given; if based on relative (biostratigraphic) dating, dates were assumed to represent the lower bound, for example, late Oligocene was assumed to be the Oligocene–Miocene boundary. Ages of boundaries were taken from Gradstein et al. (2004). The root prior rttm (the mean of the prior distribution for the time from the ingroup root to the tips) was set at 162.5 mya (the age of the monotreme–therian split) following the arguments of Benton and Donoghue (2007). A full list of calibrations is given in Table 1. Other Multidivtime parameters were calculated as follows (following Rutschmann 2005): rtrate (mean of prior distribution for the rate at the root node) = X/rttm, where X is the median amount of evolution from the root to tips; rtratesd (standard deviation of rtrtate) = 0.5 × rtrate; brownmean (mean of the prior distribution for the autocorrelation parameter, v) = 1/rttm; and brownsd (standard deviation of brownmean) = 1/rttm. Three independent Multidivtime analyses were run for 1 × 106 cycles with samples taken every 100 cycles after a burn-in period of 1 × 105 cycles. The estimated dates for each node varied <1% across the 3 runs in all analyses indicating that stationarity had been reached; the dates presented here are mean values for the 3 runs.

View this table:
Table 1

Minimum and maximum time constraints used in the Multidivtime dating analysis. Node numbers correspond to nodes in Figs. 3 and 4.

Node numberAge (mya)Justification
1162.5 (root prior)Amphilestes broderipii and Phascolotherium bucklandi—oldest known theriiomorphs (Benton and Donoghue 2007)
2124.6 (minimum)Sinodelphys szalayi and Eomaia scansoria—oldest known metatherian (Luo et al. 2003) and eutherian (Ji et al. 2002), respectively
3124.6 (maximum)All splits within placental crown-group assumed to be younger than oldest known eutherian, Eomaia
358.5 (minimum)Riostegotherium yanei—oldest known xenarthran (Rose et al. 2005)
4124.6 (maximum)All splits within marsupial crown-group assumed to be younger than oldest metatherian, Sinodelphys
4a, 6a, 7a59 (minimum)Khasia cordillerensis—oldest known microbiothere (Marshall and Muizon 1988)
4b, 6b58.5 (minimum)Oldest known tarsals with the apomorphic didelphimorphian proximal calcaneocuboid facet (Szalay 1994) and the oldest known paucituberculate (Oliveira et al. 1996)
511.8 (minimum)Micoureus laventicus (Goin 1997)—Micoureus is more closely related to Monodelphis than Didelphis (Jansa and Voss 2005; Steiner et al. 2005)
8c71.2 (maximum)Apparent absence of marsupials from South America, and no evidence of australidelphians in North America
7, 8, 923.03 (minimum)Oldest known Australasian marsupials referable to modern orders (Archer et al. 1999; Long et al. 2002)
103.53 (minimum)Cf. Peroryctes sp.—oldest known peroryctids (Turnbull et al. 2003)
11, 1223.03 (minimum)Badjcinus turnbulli—oldest known thylacinid (Muirhead and Wroe 1998)
134.45 (minimum)Antechinus sp.—oldest known phascogaline (Archer 1982; Turnbull et al. 2003)
144.45 (minimum)Cf. Dasyurus sp.—oldest known Dasyurus (Turnbull et al. 2003)
1523.03 (minimum)Oldest known vombatiforms and phalangeridans (Archer et al. 1999; Long et al. 2002)
1623.03 (minimum)Oldest known koalas and vombatoids (Long et al. 2002)
1723.03 (minimum)Oldest known petauroids, burramyids, phalangerids, and macropodoids (Archer et al. 1999; Long et al. 2002)
1823.03 (minimum)Oldest known petaurids and pseudocheirids (Archer et al. 1999; Long et al. 2002)
1923.03 (minimum)Oldest known burramyids, phalangerids, and macropodoids (Archer et al. 1999; Long et al. 2002)
204.45 (minimum)Dorcopsis wintercookorum—oldest known Dorcopsis (Flannery et al. 1992)
214.45 (minimum)Thylogale ignis—oldest known Thylogale (Flannery et al. 1992)
224.45 (minimum)Thylogale ignis—oldest known Thylogale (Flannery et al. 1992)
2323.03 (minimum)Oldest known burramyids and phalangerids (Archer et al. 1999; Long et al. 2002)
2423.03 (minimum)Oldest known trichosurin (Crosby 2004)
  • a Used in the “all constraints,” but not in the “no Khasia constraint” MultiDivTime analyses.

  • b Used in the “no Khasia constraint,” but not in the “all constraints” MultiDivTime analyses.

  • c Not used in the “no 71.2 mya constraint” MultiDivTime analyses.

An additional analysis was carried out with the lower constraint based on Khasia cordillerensis (Table 1) excluded, as doubts have been raised (e.g., Szalay 1994) as to whether Khasia is indeed a microbiothere. When Khasia was excluded, the lower constraint on the Didelphimoiphia–(Caenolesto + Australidelphia) and the Caenolestes–Australidelphia splits was 58.5 mya, based on the presence of didelphimorphian-type tarsals (Szalay 1994) and an undescribed paucituberculate (Oliveira et al. 1996) at Itaborai. Likewise, in this analysis the lower constraint on the Microbiotheria–Australasian radiation split was 23.03 mya, based on the oldest known members of the Australasian crown-group (Table 1). A further analysis was carried out with the 71.2 upper constraint on divergences within the Australasian radiation excluded, because the apparent absence of marsupials from South America before this date may be an artifact of the incomplete fossil record.

All 3 analyses (i.e., all constraints, no Khasia constraint, and no 71.2 mya constraint) were repeated using a reduced data set excluding Burramyidae, Dasyuroides, Dorcopsis, Thylacinus, and Thylogale; these taxa were represented by fewer than 6 of the 22 genes (Appendix I), which may lead to inappropriate branch lengths if these genes are evolving particularly slowly or rapidly relative to the average rate across all 22 genes.

Analysis of unrepresented basal branch lengths (UBBLs).—Following Teeling et al. (2005), I calculated unrepresented basal branch lengths (UBBLs) for several lineages to quantify the incompleteness of the fossil record implied by the molecular divergence dates (Appendix II). UBBL was calculated as 1 – (age of oldest fossil assignable to a particular branch divided by age of divergence of that branch) and converted into a percentage by multiplying by 100. Unlike the dating analysis, in which it is preferable to specify conservative minimum divergence dates (Benton and Donoghue 2007), the point estimate age (rather than the lower bound of the 95% CI) of fossils from radiometrically dated sites was used in the UBBL analysis, and the ages of sites estimated by biocorrelation were assumed to be the midpoint of the assigned age range (following the dates given by Gradstein et al. [2004]), for example, late Oligocene = 25.715 mya ( of 23.03–28.4 mya).


The unpartitioned Bayesian analysis (Fig. 1A) and all maximum-likelihood analyses (Fig. 1B) of the full data set recover a broadly similar topology in which Peramelemorphia is sister to all other australidelphians (= Eometatheria) and Dromiciops forms a clade with either Dasyuromorphia (unpartitioned Bayesian and unpartitioned maximum-likelihood analyses; Fig. 1A) or Dasyuromorphia + Notoryctes (gene- and GCSL-partitioned maximum-likelihood analyses; Fig. 1B). By contrast, gene-partitioned and GCSL-partitioned Bayesian analyses of the full data set (Fig. 2A) and all Bayesian and maximum-likelihood analyses of the no mt3 data set (Figs. 2B and 3) support Peramelemorphia as the sister to Dasyuromorphia + Notoryctes and place Dromiciops as the sister to all other australidelphians. The only differences between the topologies presented in Fig. 2A (no mt3 data set, GCSL-partitioning, Bayesian), Fig. 2B (full data set, GCSL-partitioning, maximum-likelihood), and Fig. 3 (no mt3 data set, GCSL-partitioning, Bayesian) concern the position of Caenolestes (sister to Didelphimorphia in Fig. 2B, but sister to Australidelphia in Figs. 2A and 3) and relationships within macropodines (Thylogale is sister to Dendrolagus in Figs. 2B and 3, but sister to Macropus in Fig. 2A). These alternatives are statistically indistinguishable (Table 2), and I focus my discussion on the tree in Fig. 3 (no mt3 data set, GCSL-partitioning, Bayesian), because this represents the topology and data set used in the molecular dating analysis.

Fig. 1

Topologies that result from A) Bayesian analysis of the full data set without partitioning, and B) maximum-likelihood analysis of the full data set using the GCSL-partitioning scheme. The Bayesian tree is a 50% majority rule consensus tree. Numbers in A are Bayesian posterior probabilities, and in B are maximum-likelihood bootstrap values.

Fig. 2

Topologies that result from A) Bayesian analysis of the full data set using the GCSL-partitioning scheme, and B) maximum-likelihood analysis of the no mt3 data set using the GCSL-partitioning scheme. Other details are as in Fig. 1.

Fig. 3

Topology that results from Bayesian analysis of the no mt3 data set using the GCSL-partitioning scheme. Numbers to the left of the nodes correspond to constraints used in the Multidivtime analyses (Table 1). All Bayesian posterior probability values are ≥0.98, unless otherwise indicated. Other details are as in Fig. 1.

View this table:
Table 2

Results of approximately unbiased (AU—Shimodaira 2002) and Kishino–Hasegawa (KH—Kishino and Hasegawa 1989) tests of alternative topologies as calculated by CONSEL (Shimodaira and Hasegawa 2001) assuming both the full and no mt3 (excluding 3rd codon positions of mitochondrial protein-coding genes) data sets. P-values less than 0.05 (indicating that the topology is statistically rejected by the data set) are indicated with asterisks. The first 5 topologies are from Figs. 13, with the remainder from the studies cited.

Full data setNo mt3 data set
Full, unpartitioned, BayesianFig. 1A0.0840.0760.031*0.026*
Full, GCSL-partitioned, MLFig. IB0.1180.0710.047*0.036*
No mt3, GCSL-partitioned, MLFig. 2A0.480.2990.5290.371
Full, GCSL-partitioned, BayesianFig. 2B0.5770.380.6080.405
No mt3, Bayesian, GCSL-partitionedFig. 30.7080.4140.7741
Didelphimorphia + Caenolestes (= Ameridelphia)0.4730.2980.5260.372
Caenolestes sister to rest of MarsupialiaBaker et al. 2004; Kirsch et al. 1997; Szalay and Sargis 2001, 20060.7320.5860.6540.481
Peramelemorphia sister to rest of Australidelphia (= Eometatheria)Asher et al. 2004: figure 1; Kirsch et al. 19970.4080.2670.2610.154
Diprotodontia + Peramelemorphia (= Syndactyla excluding Notoryctes)Szalay 1994; Szalay and Sargis 2001, 20060.012*0.024*0.016*0.011*
Diprotodontia + (Peramelemorphia + Notoryctes) (= Syndactyla including Notoryctes)Szalay 1994; Szalay and Sargis 2001, 20060.038*0.036*0.037*0.023*
Peramelemorphia + NotoryctesHorovitz and Sanchez-Villagra 2003; Warburton 20030.30.1860.30.143
Dromiciops + DiprotodontiaHorovitz and Sanchez-Villagra 2003, Cardillo et al. 20040.3390.2360.270.159
Dromiciops + Dasyuromorphia (= Gondwanadelphia)Szalay 1994; Szalay and Sargis 20010.022*0.021*0.009*0.008*
Dromiciops + (Dasyuromorphia + Notoryctes)0.004*0.014*0.029*0.004*
Dromiciops + PeramelemorphiaAsher et al. 2004: figure 2—right0.0750.0630.0740.048*
Dromiciops + (Peramelemorphia + (Dasyuromorphia + Notoryctes))Nilsson et al. 20040.025*0.04*0.023*0.017*
Macropodoidea + PetauroideaAmrine-Madsen et al. 20030.2080.1750.2270.159
Petauroidea + (Burramyidae + Phalangeridae) (= possum monophyly)0.4980.360.5150.376
Burramyidae sister to rest of DiprotodontiaAsher et al. 2004: figure 10.036*0.0630.041*0.02*

When the probable microbiothere Khasia is not used as a fossil constraint, molecular dates become 1–6% younger, whereas removal of the 71.2 mya upper constraint on the diversification of Australasian marsupials results in an 8–12% increase in age (Table 3). Exclusion of the 5 taxa with the most missing data (Burramyidae, Dasyuroides, Dorcopsis, Thylacinus, and Thylogale) generally results in somewhat (∼2–8%) younger dates, although some nodes instead show a slight increase in age (Table 3). In the absence of compelling evidence that Khasia is not a microbiothere or that marsupials were present in South America before 71.2 mya, I will assume for the purposes of discussion the divergence dates that result when all constraints (Table 1; Fig. 4) are enforced.

View this table:
Table 3

Divergence dates estimated using Multidivtime (node numbers are given in Figs. 3 and 4). Two different data sets were used, 1 including all sequences and the other excluding the 5 taxa with most incomplete data (Burramyidae, Dasyuroides, Dorcopsis, Thylacinus, and Thylogale). Multidivtime analyses were carried out using all constraints listed in Table 1 (“All constraints”), excluding Khasia as the oldest known microbiothere (“No Khasia constraint”), or excluding the 71.2 mya upper constraint on divergences within the Australasian radiation (“No 71.2 mya constraint”). Divergence dates from the analysis using all sequences and all constraints are illustrated in Fig. 4.

No mt3 data set, all sequencesNo mt3 data set, 5 most incomplete taxa deleted
NodeAll constraintsNo Khasia constraintNo 71.2 mya constraintAll constraintsNo Khasia constraintNo 71.2 mya constraint
2197.2 ± 23.5 (154.9–246.9)193.8 ± 24.5 (150.1–245.5)215.4 ± 33.1 (159.2–285.5)198.9 ± 22.6 (159.6–247.2)190.6 ± 23.6 (151.5–242.6)216.2 ± 32.8 (163.4–287.8)
378.2 ± 11.5 (60.3–104.0)77.0 ±11.5 (59.8–103.4)85.1 ± 14.7 (61.3–116.7)77.9 ± 11.3 (60.3–103.3)75.0 ± 11.1 (59.5–101.1)84.2 ± 14.7 (61.1–116.5)
480.6 ± 5.2 (70.7–90.0)79.0 ± 6.3 (65.8–89.7)90.3 ± 12.7 (71.5–118.8)79.3 ± 5.5 (70.2–90.1)75.0 ± 7.5 (62.6–89.1)88.4 ± 13.0(70.9–118.0)
544.4 ± 4.5 (35.9–53.4)43.4 ± 4.9 (34.1–53.1)50.4 ± 8.7 (37.1–70.3)42.4 ± 4.4 (34.4–51.5)39.8 ± 5.3 (30.6–50.5)47.9 ± 8.5 (35.4–67.5)
676.5 ± 4.7 (67.5–84.7)74.9 ± 5.8 (62.7–84.4)85.9 ± 12.1 (68.1–113.0)75.1 ± 5.0 (67.0–84.5)70.9 ± 7.0 (59.4–83.8)83.9 ± 12.4(67.6–112.3)
767.4 ± 3.9 (59.8–73.6)66.0 ± 4.9 (55.5–73.5)75.9 ± 10.8 (60.3–100.2)65.8 ± 4.2 (59.3–73.3)62.0 ± 6.2 (51.7–73.0)73.6 ± 11.0(59.5–99.1)
865.2 ± 3.7 (57.9–70.9)63.8 ± 4.7 (53.7–70.8)73.3 ± 10.4 (58.3–96.8)63.5 ± 4.0(57.1–70.7)59.9 ± 6.0 (49.9–70.5)71.2 ± 10.7 (57.5–95.8)
963.4 ± 3.7 (56.3–69.2)62.1 ± 4.6 (52.3–69.1)71.4 ± 10.2 (56.8–94.3)61.6 ± 4.0 (55.2–68.8)58.1 ± 5.8 (48.3–68.5)69.1 ± 10.4 (55.7–93.1)
1011.7 ± 1.5 (9.1–15.1)11.5 ± 1.6(8.7–14.9)13.2 ± 2.5 (9.4–19.0)11.3 ± 1.5 (8.8–14.5)10.6 ± 1.6 (7.9–14.0)12.8 ± 2.5 (9.1–18.7)
1160.1 ± 3.6 (53.2–66.1)58.9 ± 4.4 (49.6–65.9)67.6 ± 9.7 (53.8–89.8)58.3 ± 3.9 (51.8–65.5)54.9 ± 5.6 (45.4–65.1)65.4 ± 10.0 (52.3–88.3)
1226.0 ± 2.2 (23.2–31.4)25.7 ± 2.2 (23.1–31.1)28.6 ± 4.5 (23.2–39.7)
1311.6 ± 1.6(8.8–15.2)11.4 ± 1.6(8.7–15.1)12.8 ± 2.5 (9.1–18.9)10.6 ± 1.3 (8.3–13.5)10.0 ± 1.5 (7.5–13.2)12.0 ± 2.3 (8.6–17.5)
147.8 ± 1.5 (5.3–11.1)7.7 ± 1.5 (5.2–11.0)8.6 ± 2.1 (5.5–13.5)
1557.1 ± 3.4 (50.5–62.9)55.9 ± 4.3 (47.1–62.8)64.3 ± 9.3 (50.9–85.2)55.4 ± 3.7 (49.2–62.3)52.2 ± 5.3 (43.3–61.8)62.1 ± 9.5 (49.7–84.0)
1640.7 ± 3.3 (34.4–47.0)39.8 ± 3.7 (32.5–46.9)45.8 ± 7.1 (35.0–62.1)39.2 ± 3.4 (33.0–45.8)36.9 ± 4.2 (29.5–45.2)44.0 ± 7.2 (33.7–60.6)
1751.1 ± 3.2 (44.9–57.0)50.0 ± 4.0 (41.8–56.8)57.4 ± 8.4 (45.2–76.4)49.0 ± 3.5 (43.0–55.8)46.2 ± 4.8 (38.2–55.2)55.0 ± 8.5 (43.4–74.7)
1832.4 ± 3.1 (26.6–38.7)31.7 ± 3.4 (25.3–38.5)36.4 ± 6.0 (27.0–50.1)30.8 ± 3.1 (25.3–37.2)29.1 ± 3.5 (23.6–36.4)34.6 ± 6.0 (25.7–48.5)
1947.1 ± 3.2 (41.1–53.0)46.0 ± 3.8 (38.3–52.8)52.9 ± 7.9 (41.3–70.7)45.0 ± 3.4 (39.1–51.6)42.4 ± 4.5 (34.8–51.0)50.5 ± 8.0 (39.5–68.9)
2016.7 ± 2.3 (12.6–21.7)16.3 ± 2.4 (12.0–21.4)18.7 ± 3.7 (12.9–27.2)
2114.7 ± 2.1 (10.9–19.2)14.3 ± 2.2(10.3–19.0)16.4 ± 3.3 (11.1–24.0)17.1 ± 2.5 (12.6–22.4)16.1 ± 2.6(11.5–21.8)19.2 ± 3.9 (13.1–28.2)
2213.1 ± 2.0 (9.4–17.4)12.7 ± 2.1 (9.0–17.2)14.6 ± 3.1 (9.7–21.6)
2340.0 ± 3.7 (33.0–47.4)38.8 ± 4.2 (30.6–46.9)44.5 ± 7.4 (33.0–61.3)
2428.1 ± 3.2 (23.3–35.3)26.5 ± 4.1 (18.7–34.9)30.4 ± 6.2 (20.2–44.4)27.1 ± 4.1 (19.5–35.4)25.5 ± 4.3 (17.8–34.5)30.4 ± 6.2 (20.4–44.7)
Fig. 4

Modified version of Fig. 3 with dates calculated by Multidivtime assuming full constraints and including all taxa (estimated dates from all the Multidivtime analyses are given in Table 3). Monotreme and placental outgroups have been pruned. Dark gray bars represent 1 SD either side of the point estimate, and light gray bars represent the 95% CI. Circles on branches and their associated letters represent the earliest known fossils assignable to that lineage and used to calculate unrepresented basal branch lengths (UBBLs; Appendix II). Black circle = fossil confidently assigned to that lineage and used to calibrate the Multidivtime analysis (Table 1); gray circle = fossil confidently assigned to that lineage but not used to calibrate the Multidivtime analysis; white circle = fossil only tentatively assigned to that lineage.


Although relatively strongly supported by both bootstrapping and posterior probabilities, the positions of Peramelemorphia and Dromiciops in Figs. 1A and 1B are most likely anomalous and result from misleading signals in the 3rd codon positions of the mitochondrial protein-coding genes. Increased partitioning of the full data set, which is strongly supported by Bayes factors ≫ 10 (Nylander et al. 2004), results in a different topology (Fig. 2A) that is largely identical to that from the no mt3 data set, albeit only when Bayesian analysis is used. The topologies in Fig. 1 also are rejected by both the AU and KH tests when the no mt3 data set is assumed (Table 2). Finally, Phillips et al. (2006) found clear support for a close relationship between Peramelemorphia and Dasyuromorphia (and for Dromiciops as the sister to the Australasian radiation) when the analysis was corrected for base composition bias. The inability of maximum-likelihood analysis (as implemented by RAxML) to recover the inferred correct topology, even under the GCSL-partitioning scheme (Fig. 1B), was unexpected and may relate to the fact that maximum-likelihood model parameters must be specified before the analysis, whereas they are estimated during the analysis in the Bayesian approach. Also, RAxML only allows use of the GTR + gamma model for each partition, rather than other models that may be more appropriate.

Figure 3 is fully resolved, and all but 1 of the interordinal divergences are strongly supported by Bayesian posterior probabilities (BPP ≥ 0.98), although equivalent bootstrap values from the maximum-likelihood analysis using the same data set and partitioning scheme (Fig. 2A) are low for some nodes. The deepest split within marsupials (BPP = 1.00) is between Didelphimorphia and a clade comprising Caenolestes and Australidelphia. This arrangement is congruent with recent molecular phylogenies (Amrine-Madsen et al. 2003; Nilsson et al. 2004), the morphological analysis of Horovitz and Sanchez-Villagra (2003), and the combined analysis of Asher et al. (2004). However, the Caenolestes + Australidelphia clade is relatively weakly supported (BPP = 0.72), and neither the AU nor KH test rejects a monophyletic Ameridelphia (recovered in the equivalent maximum-likelihood analysis with 60% bootstrap support; Fig. 2A) or Didelphimorphia + Australidelphia. The latter clade was favored by Szalay and Sargis (2001, 2006) on morphological grounds and has been found in some molecular studies, for example, Baker et al. (2004) and Kirsch et al. (1997). Thus, the relationships of caenolestids remain uncertain, and inclusion of sequence data from 1 or both of the other extant caenolestid genera (Lestoros and Rhyncholestes) may help increase resolution by breaking up the long branch leading to Caenolestes. At least part of this uncertainty also may be due to rooting difficulties, because studies suggest that the position of the root is often the most difficult region of a molecular phylogeny to resolve robustly (Sanderson and Shaffer 2002). Additional sources of molecular data, such as retrotransposon insertions, insertions–deletions, or other rare genomic changes (Rokas and Holland 2000), may be required to root the marsupial tree with confidence.

Within Australidelphia, the South American Dromiciops is sister to a monophyletic Australasian clade (BPP = 1.00; bootstrap = 100). The alternative positions for Dromiciops recovered by Asher et al. (2004—sister to all Australasian marsupials except Peramelemorphia) and Nilsson et al. (2004—sister to Notoryctes + Dasyuromorphia + Peramelemorphia) are most likely erroneous, because both studies used mitochondrial genes without correcting for the biases identified by Phillips et al. (2006). The topology of Nilsson et al. (2004) is rejected by both the AU and KH tests, whereas the topology of Asher et al. (2004) is not (Table 2).

Within the Australasian radiation, the 1st divergence is between the diverse order Diprotodontia and a clade comprising the 3 Australasian “polyprotodont” orders, Peramelemorphia, Notoryctemorphia, and Dasyuromorphia (BPP = 1.00; bootstrap = 40). This topology agrees with that of Amrine-Madsen et al. (2003) and Phillips et al. (2006), and it indicates that syndactyly (fusion of digits II and III and their integration into a single functional digit—Weisbecker and Nilsson 2006) evolved independently in these 2 orders (contra Szalay 1994). Although syndactyly is a derived feature present in all diprotodontians and peramelemorphians, a monophyletic Syndactyla (Szalay 1982, 1994), whether including or excluding Notoryctes (which may or may not be syndactylous), is rejected by both the KH and AU tests (Table 2). Notoryctes is the sister-group of Dasyuromorphia (BPP = 0.98; bootstrap = 70). Recent morphological studies (Horovitz and Sanchez-Villagra 2003; Warburton 2003) have linked Notoryctes with Peramelemorphia, and a Notoryctes + Peramelemorphia clade is not rejected by KH and AU tests (Table 2). Tarsal material of a fossil notoryctid from early Miocene Faunal Zone B (previously “System B”—Arena 2004; Travouillon et al. 2006) sites at Riversleigh also shares distinctive apomorphies (not present in the more derived Notoryctes) with peramele-morphians (R. M. D. Beck, in litt.). Notoryctes is the sole extant representative of the order Notoryctemorphia and is extremely autapomorphic; it thus represents both a molecular and morphological long branch that may be difficult to place robustly with either morphological or sequence data. Other kinds of molecular data may be required to resolve the affinities of Notoryctes conclusively. Within Dasyuromorphia, Dasyuridae is monophyletic relative to the thylacine Thylacinus (BPP = 1.00; bootstrap = 100).

In agreement with morphological data (e.g., Aplin and Archer 1987), Diprotodontia comprises 2 major clades: Vombatiformes (the koala and wombats; BPP = 1.00; bootstrap = 100) and Phalangerida (the various families of extant “possums” and kangaroos; BPP = 1.00; bootstrap = 90). An alternative topology, with Burramyidae sister to other diprotodontians (Asher et al. 2004), is rejected by both the AU and KH tests with the no mt3 data set, but only by the AU test for the full data set (Table 2). Petauroidea, a clade comprising pseudocheirids and petaurids, is monophyletic (BPP = 1.00; bootstrap = 100) and also is well supported by morphology (Aplin and Archer 1987). A 2nd possum clade consisting of Burramyidae and Phalangeridae (BPP = 0.98; bootstrap = 70) also was seen in the DNA-hybridization trees of Springer and Kirsch (1991) and Kirsch et al. (1997), a maximum-likelihood analysis of RAG1 sequences by Baker et al. (2004), and morphological analyses of phalangeridans by Crosby (2004). “Possums,” as a whole, appear to be a paraphyletic assemblage (although possum monophyly cannot be rejected; Table 2) because kangaroos are the sister-group of the Burramyidae + Phalangeridae clade (BPP = 0.98; bootstrap = 50). Plesio-morphic macropodoids, burramyids, and phalangerids all possess a hypertrophied 3rd premolar that may represent a morphological synapomorphy of this clade.

Divergence dates.—On initial examination, the divergence dates in Fig. 4 and Table 3 seem congruent with the known fossil record; the dates are not implausibly ancient (i.e., older than the Late Cretaceous), nor are any lineages younger than the oldest known fossils referred to them. In contrast, at least some of the molecular dates of Nilsson et al (2004), who used a penalized log-likelihood method and a relatively young estimate of 135 mya for the divergence between marsupials and placentals (compared to 186–193 mya calculated by van Rheede et al. [2006]), appear incongruent with the fossil record. For example, their estimated age for the divergence of Microbiotheria is 46 mya, but the oldest known fossil microbiothere, Khasia cordillerensis from Tiupampa (Marshall and Muizon 1988), is approximately 60 million years old (however, there are some doubts as to whether Khasia is indeed a microbiothere—Szalay 1994). Below, I undertake a more detailed examination of the divergence dates to identify the degree of congruence with the fossil record and to determine whether particular divergences coincided with major environmental changes.

The 2 deepest splits within marsupials are here estimated as occurring 80.6 mya (Didelphimorphia–other marsupials) and 76.5 mya (Caenolestes–Australidelphia), that is, in the Campanian of the Late Cretaceous. If the apparent absence of marsupials from the Alamitian (∼71.2 mya) faunas of South America is accurate, then these 2 divergences must have occurred in Laurasia (probably North America, given that marsupials appear to have been rare in Eurasia during the Late Cretaceous), and members of 3 separate marsupial lineages (didelphimorphians, paucituberculates, and australidelphians) dispersed from North America to South America during the latest Cretaceous–earliest Paleocene. This total does not include fossil groups, such as borhyaenoids and polydolopi-morphians, that may lie outside the marsupial crown-group. Fossil evidence indicates that as many as 4 separate placental lineages also dispersed from North to South America during the latest Cretaceous–earliest Paleocene (Muizon 1991; Muizon and Cifelli 2000), as did other terrestrial vertebrates such as hadrosaurs and snakes (Case et al. 2000, 2005), suggesting the existence of a relatively robust connection between the 2 continents during this period.

Evidence of crown marsupials in the Late Cretaceous of North America is equivocal given the limited fossil material currently known, but the dates presented here are congruent with the assignment of Glasbius to Paucituberculata (Kielan-Jaworowska et al. 2004) and Nortedelphys to Didelphimorphia (Case et al. 2005), in which case UBBL is 17.4% and 13% for Didelphimorphia and Paucituberculata, respectively (Appendix II). Unequivocal didelphimorphians (identified by the presence of the distinctive proximal calcaneocuboid facet of the tarsus seen in all extant didelphids—Szalay 1994) and paucituberculatans (Oliveira et al. 1996) are 1st known from the late Paleocene (58.7–59.2 mya) Itaborai fauna in Brazil, but have not been identified in the slightly older (59.2–60.4 mya) Tiupampa fauna. If the Itaboraian fossils are preferred as the earliest known members of these orders, then UBBL increases to 26.9% and 22.9% for Didelphimorphia and Paucituberculata, respectively (Appendix II). Efforts to discover earlier records of these orders are hampered by a lack of South American marsupials older than those from Tiupampa, with the exception of the highly fragmentary Punta Peligro fauna (possibly early Paleocene) and a single tooth from the early Paleocene Lefipan Formation (Goin et al. 2006). UBBLs for the crown didelphimorphians Monodelphis (68.4%) and Didelphini (87.8%) are high, probably because didelphids are poorly defined on dental grounds (the family is characterized by a highly plesiomorphic dentition), making referral of isolated teeth to the family or to specific genera difficult.

The split between microbiotheres and the Australasian marsupials is estimated as 67.4 mya. There is no evidence of microbiotheres or any member of the Australasian radiation in the Late Cretaceous North American deposits, which may be indirect evidence that this split occurred in Gohdwana. However, stem australidelphians most likely had a relatively plesiomorphic tribosphenic dentition (perhaps similar to that of the early Eocene Australian marsupial Djarthia murgonensis—R. M. D. Beck, in litt.), and unequivocal australidelphian synapomorphies are restricted to the tarsus, suggesting that it may be impossible to identify them from isolated teeth alone. The oldest known australidelphian is the probable micro-biothere Khasia cordillerensis from Tiupampa; hence, UBBL for Microbiotheria is only 11.3% (Appendix II).

Although the estimated age of the 1st split within the Australasian marsupials (65.2 mya) is almost immediately after the Cretaceous–Tertiary boundary, this is not necessarily evidence of a causal link between these 2 events. Whether this split occurred in South America, Australia, or Antarctica is uncertain. However, given that the modern Australasian radiation appears to be monophyletic, and that there is no fossil evidence of members of any Australasian order in South America or the middle Eocene of Seymour Island in western Antarctica, an Australasian or eastern Gondwanan center of origin seems likely. If so, the date of this split represents the minimum age for the arrival of marsupials in Australasia–eastern Gondwana, which therefore must have occurred between 65.2 (± 3.7; 95% CI = 57.9–70.9) and 71.2 mya. All 4 extant Australasian orders—Peramelemorphia, Dasyuromorphia, Notoryctemorphia, and Diprotodontia—are estimated to have diverged by the age of the early Eocene Tingamarra fauna (54.6 mya). However, no undoubted members of any of these orders have been identified among the marsupial specimens recovered from this site; hence UBBL for these 4 orders is 59.4%, 57.2%, 67.5%, and 60.6%, respectively (Appendix II). A single upper molar from Tingamarra was tentatively identified as a possible bandicoot (Archer et al. 1999), but it lacks unequivocal peramelemorphian dental synapomorphies. A 2nd upper molar from Tingamarra reported as an undoubted peramelemorphian (Beck et al. 2006) is in fact a misattributed specimen from a younger site. Woodburne and Case (1996) suggested that the bunodont Tingamarran marsupial Thylacotinga is a peramelemorphian, but it also lacks unequivocal apomorphies that would distinguish it from other bunodont forms such as polydolopimorphians. Given that the most distinctive feature of peramelemorphian molars—a highly invasive, incomplete centrocrista that breaches the ectoloph—is not present in the late Oligocene Yarala burchfieldi (Muirhead and Filan 1995), identifiable dental synapomorphies of Peramelemorphia may not have evolved until long after the early Eocene. If so, fossil peramelemorphians may be present at Tingamarra but cannot be recognized as such. This explanation also may apply to Notoryctemorphia; although notoryctemorphians are estimated as originating at 60.0 mya, the only known fossil notoryctid (from the middle Miocene of Riversleigh—Long et al. 2002) is considerably more plesiomorphic than Notoryctes, suggesting that many characteristic notoryctid specializations may have evolved much later than the early Eocene. The absence of notoryctids from the late Oligocene Faunal Zone A of Riversleigh is more surprising, given that a relatively derived notoryctid is present in the middle Miocene Faunal Zone B. However, Faunal Zone A sites appear biased (whether taphonomically or ecologically) toward large-bodied taxa such as thylacinids and larger diprotodontians (H. Godthelp, University of New South Wales, pers. comm.).

A detailed analysis of marsupial dental characters by Godthelp et al. (1999) failed to identify dental apomorphies that distinguish dasyuromorphians from other marsupials with relatively unspecialized tribosphenic dentitions. Thus, dasyuromorphians may be represented by some of the hundreds of isolated teeth now recovered from Tingamarra but cannot be distinguished unequivocally from other dentally plesiomorphic groups. Within Dasyuromorphia, the dasyurid–thylacinid divergence date of 26.0 mya is congruent with presence of the plesiomorphic thylacinids Nimbacinus and Badjcinus in late Oligocene Riversleigh Faunal Zone A deposits and the dasyurid Barinya (which shows relatively derived dasyurid apomorphies of the petrosal—Wroe 1999) in middle Miocene Faunal Zone B deposits. UBBL for Thylacinidae is only 1.1%, yet Muirhead and Wroe (1998) identified 6 thylacinid apomorphies in Badjcinus. This suggests either that morphological evolution in the earliest thylacinids was quite rapid, that Riversleigh Faunal Zone A sites may be somewhat younger than the 25.7 mya assumed here, or that the molecular divergence date calculated here is an underestimate resulting from the high proportion of missing data for Thylacinus.

The estimated divergence dates for the split between phascogaline and dasyurine dasyurids (11.6 mya) and between peramelid and peroryctid bandicoots (11.7 mya) are almost identical and coincide with the middle–late Miocene boundary, which saw a major drop in sea level (Haq et al. 1987). These dates explain the absence of these groups from Faunal Zones A–C (late Oligocene–middle Miocene) sites at Riversleigh; they are 1st seen in the fossil record only in the Pliocene, probably because of a lack of Australian sites from the late Miocene.

Within Diprotodontia, the point estimate for the Vombati-formes–Phalangerida split is 57.1 mya, implying that these lineages had diverged by the estimated age of the Tingamarra fauna. Both vombatiforms and phalangeridans possess apomorphically enlarged “diprotodont” 1st lower incisors that were presumably present in their last common ancestor, but no such teeth (or any other specimens that can be referred unequivocally to Diprotodontia) have been found at Tingamarra. However, the standard deviation for the Vombati-formes-Phalangerida split is 3.4 million years and the 95% CI is 50.5–62.9 mya, and so it may in fact postdate Tingamarra. Drummond et al. (2006) used a different relaxed molecular clock method that (unlike Multidivtime) does not employ rate autocorrelation and estimated a slightly younger age for this split, at about 48 mya.

Both koalas and vombatoids (but not vombatids), estimated here as having diverged 40.5 mya, are known from the late Oligocene deposits of Riversleigh, as are petaurids and pseudocheirids (divergence date = 32.2 mya), and burramyids and phalangerids (divergence date = 39.4 mya). Crosby (2004) identified trichosurins among the late Oligocene Riversleigh phalangerids, congruent with the 26.9 mya divergence date estimated for the trichosurin–phalangerin split presented here, hence UBBL for trichosurins is only 8.5%. However, undoubted phalangerins have not been found at Riversleigh or indeed any fossil site older than the Pleistocene, resulting in a very large UBBL of 96.6%. Phalangerins may have evolved in New Guinea or Sulawesi and invaded Australia during the Pleistocene (Raterman et al. 2006). The earliest undoubted macropodoids are from the late Oligocene (Archer et al. 1999; Long et al. 2002), and an undescribed macropodine is known from a middle Miocene Faunal Zone C site at Riversleigh (B. N. Cooke, Queensland University of Technology, pers. comm.). The radiation of the macropodines, which appears to have started in the middle Miocene (16.7 mya), coincided with the mid-Miocene climatic optimum (Gradstein et al. 2004; Zachos et al. 2001) and the 1st appearance of grasslands in Australia (Martin 2006). Late Miocene macropodines have yet to be found, but, as already mentioned, few Australian sites of this age are known.


I thank M. Springer and M. Nilsson for kindly providing me with molecular matrices. My particular gratitude goes to O. Bininda-Emonds, who carried out all of the RAxML analyses on my behalf. The perl scripts batchRAxML and bootStrip for automating RAxML analyses can be downloaded from his Web site (http://www2.uni-jena.de/∼b6biol2/). H. Godthelp, M. Archer, B. Cooke, and R. Voss provided extremely helpful information and advice, and Y. Gurovich assisted me in the preparation of the figures. V. Weisbecker, R. Tomlins, M. Hafner, and 2 anonymous reviewers gave insightful comments that have greatly improved this manuscript. Funding for this work is provided by the Leverhulme Trust through Study Abroad Studentship SAS/30110.

Appendix I

View this table:

GenBank accession numbers for all gene sequences used in this analysis. All sequences and alignments except those for Thylacinus(which were added manually) are taken from Amrine-Madsen et al. (2003), Asher et al. (2004), and Nilsson et al. (2004). Abbreviations for gene names are given in text.

APOBBRCA1IRBPP1PGK1RAG1VWF12S rRNA16S rRNAtRNA valine12 mt protein-coding genes
ThylacinusU87140U87405M99452 (CYTBonly)

Appendix II

View this table:

Unrepresented basal branch lengths (UBBLs—Teeling et al. 2005) of selected lineages, with fossils used to calculate UBBLs indicated. The UBBL value for a particular lineage is the proportion of the age of that lineage for which fossils are unknown. Age ranges are taken from Gradstein et al. (2004). Fossils are indicated on Fig. 4 by the letters in parentheses after the fossil names.

LineageEarliest fossilAge (mya)ReferenceUBBL (%)
DidelphimorphiaNortedelphysspp. (A)66.55 (late Maastrichtian = 65.5–67.6)Case et al. 200517.4
IMG V and IMG XII (B)58.95 (Itaboraian = 58.7–59.2)Szalay 199426.9
MonodelphisMicoureus laventicus(C)14.05 (Laventan–Friasian = 11.8–16.3)Marshall 1976; Goin 199768.4
DidelphisDidelphis pattersoni(D)5.4 (Montehermosan = 4–6.8)Simpson 197487.8
PaucituberculataGlasbiusspp. (E)66.55 (Late Maastrichtian = 65.5–67.6)Kielan-Jaworowska et al. 200413.0
Unnamed paucituberculate (F)58.95 (Itaboraian = 58.7–59.2)Oliveira et al. 199622.9
MicrobiotheriaKhasia cordillerensis(G)59.8 (Tiupampan = 59.2–60.4)Marshall and Muizon 198811.3
PeramelemorphiaYarala burchfieldi(H)25.715 (late Oligocene = 23.03–28.4)Muirhead and Filan 199559.4
PeramelidaePerameles allinghamensis(I)3.615 (Allingham Formation = 3.6–3.63)Archer and Wade 197669.1
PeroryctidaeCf. Peroryctesspp. (J)4.46 (Hamilton Fauna = 4.46)Turnbull et al. 200361.9
NotoryctemorphiaUnnamed notoryctid (K)19.5 (early Miocene = 15.97–23.03)Archer et al. 1999; Long et al. 200267.6
ThylacinidaeBadjcinus turnbulli(L)25.715 (late Oligocene = 23.03–28.4)Muirhead and Wroe 19981.1
DasyuridaeBarinya wangala(M)19.5 (early Miocene = 15.97–23.03)Wroe 199925
PhascogalinaeAntechinussp. (N)4.46 (Hamilton Fauna = 4.46)Archer 1982; Turnbull et al. 200361.6
DasyurusCf. Dasyurussp. (O)4.46 (Hamilton Fauna = 4.46)Turnbull et al. 200342.8
DasyuroidesDasyuroides achilpatna(P)1.81 (late Pliocene–early Pleistocene = 1.81)Archer 198276.9
VombatidaeVombatoid spp. (Q)25.715 (late Oligocene = 23.03–28.4)Archer et al. 1999; Long et al. 200236.8
PhascolarctidaePhascolarctid spp. (R)25.715 (late Oligocene = 23.03–28.4)Archer et al. 1999; Long et al. 200236.8
PseudocheiridaePseudocheirid spp. (S)25.715 (late Oligocene = 23.03–28.4)Archer et al. 1999; Long et al. 200220.6
PetauridaePetaurid spp. (T)25.715 (late Oligocene = 23.03–28.4)Archer et al. 1999; Long et al. 200220.6
MacropodidaeMacropodoid spp. (U)25.715 (late Oligocene = 23.03–28.4)Archer et al. 1999; Long et al. 200245.4
DorcopsisDorcopsis wintercookorum(V)4.46 (Hamilton Fauna = 4.46)Flannery et al. 199273.3
MacropusMacropus(Osphranter) pavana(W)3.615 (Allingham Formation = 3.6–3.63)Bartholomai 197875.4
ThylogaleThylogale ignis(X)4.46 (Hamilton Fauna = 4.46)Flannery et al. 199266.0
DendrolagusCf. Dendrolagussp. (Y)4.46 (Hamilton Fauna = 4.46)Flannery et al. 199266.0
BurramyidaeBurramys brutyi(Z)25.715 (late Oligocene = 23.03–28.4)Brammall and Archer 199735.7
TrichosuriniUnnamed trichosurins (AA)25.715 (late Oligocene = 23.03–28.4)Crosby 20048.5
PhalangeriniPleistocene phalangerins (AB)0.955 (Pleistocene = 0.01–1.81)Crosby et al. 200496.6

Literature Cited

View Abstract