Abstract

Madagascar's flora and fauna are remarkable both for their diversity and supraspecific endemism. Moreover, many taxa contain large numbers of species with limited distributions. Several hypotheses have been proposed to explain this high level of microendemism, including 1) riverine barrier, 2) mountain refuge, and 3) watershed contraction hypotheses, the latter 2 of which center on fragmentation due to climatic shifts associated with Pliocene/Pleistocene glaciations. The Malagasy leaf chameleon genus Brookesia is a speciose group with a high proportion of microendemic taxa, thus making it an excellent candidate to test these vicariance scenarios. We used mitochondrial and nuclear sequence data to construct a Brookesia phylogeny, and temporal concordance with Pliocene/Pleistocene speciation scenarios was tested by estimating divergence dates using a relaxed-clock Bayesian method. We strongly reject a role for Pliocene/Pleistocene climatic fluctuations in species-level diversification of Brookesia. We also used simulations to test the spatial predictions of the watershed contraction model in a phylogenetic context, independent of its temporal component, and found no statistical support for this model. The riverine barrier model is likewise a qualitatively poor fit to our data, but some relationships support a more ancient mountain refuge effect. We assessed support for the 3 hypotheses in a nonphylogenetic context by examining altitude and species richness and found a significant positive correlation between these variables. This is consistent with a mountain refuge effect but does not support the watershed contraction or riverine barrier models. Finally, we find repeated higher level east-west divergence patterns 1) between the 2 sister clades comprising the Brookesia minima group and 2) within the clade of larger leaf chameleons, which shows a basal divergence between western and eastern/northern sister clades. Our results highlight the central role of phylogeny in any meaningful tests of species-level diversification theories.

Once part of the southern supercontinent of Gondwana, Madagascar, is currently situated some 400 km off the southeastern shore of Africa and is famous for a remarkable flora and fauna that are increasingly threatened by loss of habitat and other human-induced changes. Madagascar's last subaerial connections to Africa and India date to approximately 160 and 90 million years ago, respectively (Storey et al. 1995), although limited dispersal to and from other continents may have been possible over land bridges or small sea barriers in the Late Cretaceous or Paleocene (Briggs 2003; Noonan and Chippindale 2006; Van Bocxlaer et al. 2007; van der Meijden et al. 2007). This long history of isolation has contributed greatly to the remarkable degree of floral and faunal endemism of Madagascar (Goodman and Benstead 2003), which amounts to 100% in native amphibians and 92% in nonmarine nonavian reptiles (Glaw and Vences 2007).

In terrestrial vertebrates, the majority of species diversity corresponds to a limited number of endemic clades that colonized most of the available bioclimatic regions, including the eastern rainforests, central highlands, western dry deciduous forests, and southern arid spiny forests. There are several sister-taxon pairs that follow an east-west pattern and possibly arose by specialization to these different climatic conditions (the “ecographic constraint hypothesis”; see Yoder and Heckman 2006), from groups such as skinks, snakes, boophine treefrogs, and spiders (Nussbaum and Raxworthy 1998; Nussbaum et al. 1998; Vences and Glaw 2002, 2003; Wood et al. 2007). In other groups, the major evolutionary splits consistently separate a clade occurring roughly in the northern fifth of Madagascar from a more southern clade, with examples found in chameleons, geckos, skinks, and mouse lemurs (Yoder and Heckman 2006; Boumans et al. 2007). As a general pattern, the spatial distribution of species richness in some higher taxonomic groups may have been shaped by a latitudinal and altitudinal middomain effect (Lees et al. 1999). Nonetheless, even within the major bioclimatic regions, species-level local endemism is high, and a major goal of researchers has been to uncover the causes of this pattern (Goodman and Benstead 2003 and references therein; Vences et al. 2009).

Potential Climatic and Geographic Barriers

In their study of biogeographic patterns in the leaf chameleons (genus Brookesia) of northern Madagascar, Raxworthy and Nussbaum (1995) recognized 5 rainforest regions delimited largely by altitudinal differences and intervening dry forests and characterized by a high degree of endemism in these lizards. Raxworthy and Nussbaum (1995) hypothesized that Pleistocene climatic fluctuations had caused fragmentation of the rainforest, resulting in multiple allopatrically distributed sister-taxon pairs (hereafter referred to as the mountain refuge [MR] hypothesis). These conclusions were not based on explicit phylogenetic hypotheses; sister-taxon pairs were assumed based on general morphological similarities.

Wilmé et al. (2006) analyzed Madagascar's striking microendemicity in the context of watersheds. Using museum data for more than 35,000 georeferenced land vertebrate specimens, they found that watersheds with low-elevation headwaters tended to define centers of endemism (COEs), whereas those with connections to the 3 highest summits in Madagascar (each at > 2000 m) tended to contain more widespread species. Wilmé et al. (2006) reasoned that during periods of Late Tertiary and Pleistocene glacial maxima, aridification caused contraction of previously widespread mesic environments. Lower elevations were more severely affected than higher elevations, leading to fragmentation and isolation of low-elevation watersheds. In contrast, areas with riverine connections to high-elevation source waters were buffered from this effect and thus served as “retreat–dispersion watersheds” (RDWs), providing a means of refuge and recolonization during dry and wet cycles, respectively. The authors presented this model of Pliocene–Pleistocene fragmentation and refugia (hereafter referred to as the watershed contraction [WC] hypothesis) as the main generator of Madagascar's famously microendemic biota. The Wilmé et al. (2006) study was based entirely on extant species distributions coupled with climatological, hydrological, and topographical variables and incorporated no phylogenetic data. Pearson and Raxworthy (2009) found some support for the WC hypothesis in lemurs and day geckos of the genus Phelsuma (a climatological gradient effect was also cited for these lizards), but once again relied on distribution patterns alone.

Finally, the presence of large rivers flowing from the highlands either toward the west or the east has also been discussed as a factor influencing species formation and microendemism (hereafter referred to as the riverine barrier [RB] hypothesis). In western Madagascar, these rivers coincide with the limits of distribution areas of species or phylogeographic lineages within species of lemurs and are therefore likely to constitute significant barriers to gene flow for these animals (Pastorini et al. 2003). A similar situation seems to exist in eastern Madagascar and may especially influence species occurring at low elevations where rivers are widest (Goodman and Ganzhorn 2004; Louis et al. 2005). For a more detailed overview of the processes inherent to the MR, WC, and RB hypotheses, see Vences et al. (2009).

Brookesia Leaf Chameleons

The chameleon genus Brookesia is an excellent group to test hypotheses of species-level diversification and microendemism on Madagascar. These lizards, which appear to form the sister taxon of all other chameleons (Rieppel 1987; Townsend and Larson 2002), constitute one of the largest Malagasy reptile groups (Raxworthy and Nussbaum 1995; Glaw and Vences 2007). Most Brookesia species have very small ranges (Raxworthy and Nussbaum 1995; Raselimanana and Rakotomamalala 2003), with almost half known essentially from single localities (Carpenter and Robson 2005). Brookesia can be divided into 3 main groups based on morphology. One is represented by just 2 species (B. nasus and B. lolontany) that are highly divergent from other Brookesia by both molecular (Raxworthy et al. 2002; Townsend and Larson 2002) and morphological measures (long snouts, laterally compressed bodies; Raxworthy and Nussbaum 1995). A second larger group is composed of approximately 18 species with more robust bodies and blunted snouts. The remaining 6 described species have taken the already diminutive body form of Brookesia to an extreme, with total lengths of about 28–45 mm, making them some of the world's smallest vertebrate species (Glaw and Vences 2007). For ease of reference, these 3 groups will hereafter be referred to as the B. nasus group, the “typical” Brookesia, and the B. minima group, respectively.

Hypothesis Testing

In this study, we use Brookesia to test the temporal and spatial predictions of 3 species-level diversification hypotheses for Madagascar (MR, WC, and RB; Vences et al. 2009). The first step in this process is to infer a phylogeny of Brookesia to identify statistically supported sister-species pairs using DNA sequence data from multiple mitochondrial and nuclear protein-coding genes. Next, we use divergence dating to statistically test the temporal prediction of both the MR and the WC hypotheses that recent (Pleistocene or possibly Pliocene) climatic cycles are a major force promoting Brookesia species diversification.

Major climatic cycles have of course occurred repeatedly throughout the earth's history, and each of these 3 general hypotheses make spatial predictions that can be tested independently of any temporal predictions (e.g., Miocene WCs could generate the same general fragmentation patterns as the proposed Pleistocene contractions). Specifically, the Wilmé et al. (2006) WC model states that contraction of mesic habitats within adjacent lowland watersheds during periods of glaciation created gaps in ancestral species’ distributions, leading to allopatric differentiation/speciation. Assuming stability in the relative geographic positions of populations, this model therefore predicts that sister species should generally occupy adjacent COEs or possibly a COE and an adjacent RDW (although the species in the RDW would be expected to have a wider geographic range). Because most of the major watersheds are composed of several smaller drainages, this model could also be construed to predict that sister-species pairs would occupy the same drainage. Likewise, the MR model of Raxworthy and Nussbaum (1995) predicts that sister taxa will tend to occupy montane forested areas separated by lower altitude dry forests. Finally, the RB model (Pastorini et al. 2003) predicts species’ ranges to be delimited by major lowland rivers. To evaluate these hypotheses, we first reconstruct species ranges using a georeferenced database of Brookesia distribution records and adjust these ranges (when sample sizes permit) through environmental niche modeling. We then use simulations to test the fit of our data to the spatial predictions of the WC hypothesis, and we assess the qualitative fit of our data to the MR and RB hypotheses (all in a phylogenetic context).

Finally, each of the 3 hypotheses suggests predictions about the altitudinal distribution of species richness. The WC and RB models both predict a greater number of lowland species; low-headwater watersheds are more likely than high-headwater watersheds to become fragmented from their neighbors during periods of increased aridity, and rivers are widest at lower elevations and therefore more likely to constitute significant barriers to gene flow as altitude decreases. Conversely, the MR model predicts higher species diversity in montane areas because it is here that populations tended to become isolated during periods of forest contraction; lowlands were subsequently recolonized by some species. We use these predictions to evaluate relative support for these hypotheses in a nonphylogenetic context by correlation of species richness with altitudinal range.

MATERIALS AND METHODS

Taxa and Gene Regions Sampled

A previous phylogenetic study (Townsend and Larson 2002) identified 9 clades with uncertain interrelationships that diverged from each other early in the history of Chamaeleonidae. Because current taxonomy places 2 of these major clades within the genus Brookesia, outgroup samples from all 7 non-Brookesia chamaeleonid clades were included in this study. Sampling within Brookesia was broad, including more than three-quarters of the named species as well as several undescribed forms, and species were sampled from multiple localities whenever possible. To facilitate divergence-dating calibrations, outgroup sampling included representatives of all major squamate lineages (Table 1).

TABLE 1.

Nonchamaeleonid outgroup sampling for all phylogenetic analysesa

Higher taxon  Representative taxonb 
Rhyncocephalia  Sphenodon punctatus 
Squamata   
 Dibamidae Dibamus 
 Gekkonidae Gekko gecko 
 Xantusiidae Xantusia vigilis 
 Scincidae Mabuya 
 Cordylidae Cordylus 
 Teiidae Teiinae 
 Amphisbaenia Trogonophidae 
 Lacertidae Lacertidae 
 Serpentes Dinodon 
 Shinisauridae Shinisaurus crocodilurus 
 Iguanidae Liolaemus 
 Uromastycinae Uromastyx 
 Leiolepidinae Leiolepis 
Higher taxon  Representative taxonb 
Rhyncocephalia  Sphenodon punctatus 
Squamata   
 Dibamidae Dibamus 
 Gekkonidae Gekko gecko 
 Xantusiidae Xantusia vigilis 
 Scincidae Mabuya 
 Cordylidae Cordylus 
 Teiidae Teiinae 
 Amphisbaenia Trogonophidae 
 Lacertidae Lacertidae 
 Serpentes Dinodon 
 Shinisauridae Shinisaurus crocodilurus 
 Iguanidae Liolaemus 
 Uromastycinae Uromastyx 
 Leiolepidinae Leiolepis 
a

All phylogenetic analyses included these outgroup taxa as well as all non-Brookesia chamaeleonid taxa from Figure 2.

b

Composite taxa with more than 1 species represented among the different gene partitions are called by the name of the most exclusive higher taxon possible.

All specimens were sequenced for 2 nuclear (RAG1, ∼1500 bp and CMOS, ∼800 bp) and 3 mitochondrial (ND1, ∼ 70 bp; ND2, ∼ 1035 bp; and ND4, ∼ 700 bp) protein-coding genes; we included several ND4 sequences from Raxworthy et al. (2002) that expanded our taxonomic (3 species) or geographic (6 species) coverage. Our final Brookesia sampling covered approximately 28 species, including all but 2 of the named species. Museum/collection and GenBank numbers of specimens can be found in the online Supplementary Material (available from http://sysbio.oxfordjournals.org).

Molecular Data and Phylogenetic Analyses

Genomic DNA was extracted and amplified using standard protocols (see Townsend et al. 2004). Polymerase chain reaction products were sequenced with ABI 3100 PRISM™ automated sequencers, and contigs were assembled using Sequencher 4.1.2 (Gene Codes Corporation, Ann Arbor, MI). Primer sequences are given in Table 2. We used the Clustal algorithm (Thompson et al. 1994) implemented in the program DAMBE (Xia and Xie 2001) to align all sequences by their amino acid translations and adjusted alignments by eye using MacClade 4.03 (Maddison D.R. and Maddison W.P. 2000). Ambiguously aligned regions were excluded from all analyses (see Results section). Gaps were treated as a separate (binary) character partition in the MrBayes analyses and as missing data in the BEAST and maximum-likelihood (ML) analyses.

TABLE 2.

Amplification and sequencing primers used in this study

Directiona Location Sequence Reference 
Forward ND1 5′-CAACTAATACACCTACTATGAAA-3′ Macey et al. (1997) 
Forward ND1 5′-CGATTCCGATATGACCARCT-3′ Kumazawa and Nishida (1993) 
Reverse tRNAAla 5′-AAAATRTCTGRGTTGCATTCAG-3′ Macey et al. (1997) 
Forward ND4 5′-TGACTACCAAAAGCTCATGTAGAAGC-3′ Raxworthy et al. (2002) 
Reverse tRNALeu 5′-CATTACTTTTACTTGGATTTGCACCA-3′ Raxworthy et al. (2002) 
Forward RAG1 5′-TCTGAATGGAAATTCAAGCTGTT-3′ Groth and Barrowclough (1999) 
Forward RAG1 5′-CCACTTGGAAAAATACTCCCTGA-3′ This study 
Reverse RAG1 5′-GTCATCAACCAAATGTTGTATGCCTG-3′ This study 
Reverse RAG1 5′-GTGTCYACTGGGTARTCATC-3′ Townsend et al. (2004) 
Forward CMOS 5′-ATTATGCGATCMCCTMTTCC-3′ This study 
Forward CMOS 5′-TCTGGAATTTTCTCCWTCTGT-3′ This study 
Reverse CMOS 5′-GCTACCACAGARTASAGTACA-3′ This study 
Directiona Location Sequence Reference 
Forward ND1 5′-CAACTAATACACCTACTATGAAA-3′ Macey et al. (1997) 
Forward ND1 5′-CGATTCCGATATGACCARCT-3′ Kumazawa and Nishida (1993) 
Reverse tRNAAla 5′-AAAATRTCTGRGTTGCATTCAG-3′ Macey et al. (1997) 
Forward ND4 5′-TGACTACCAAAAGCTCATGTAGAAGC-3′ Raxworthy et al. (2002) 
Reverse tRNALeu 5′-CATTACTTTTACTTGGATTTGCACCA-3′ Raxworthy et al. (2002) 
Forward RAG1 5′-TCTGAATGGAAATTCAAGCTGTT-3′ Groth and Barrowclough (1999) 
Forward RAG1 5′-CCACTTGGAAAAATACTCCCTGA-3′ This study 
Reverse RAG1 5′-GTCATCAACCAAATGTTGTATGCCTG-3′ This study 
Reverse RAG1 5′-GTGTCYACTGGGTARTCATC-3′ Townsend et al. (2004) 
Forward CMOS 5′-ATTATGCGATCMCCTMTTCC-3′ This study 
Forward CMOS 5′-TCTGGAATTTTCTCCWTCTGT-3′ This study 
Reverse CMOS 5′-GCTACCACAGARTASAGTACA-3′ This study 

Note: tRNA = transfer RNA.

a

Mitochondrial forward and reverse primers extend the light and heavy strands, respectively.

Partitioned Bayesian and likelihood methods were used to infer phylogenies. BEAST (Drummond et al. 2006) implements a Bayesian relaxed-clock method that allows the simultaneous estimation of topology and divergence times. The incorporation of a relaxed-clock constraint into the phylogenetic inference procedure has intuitive appeal, as it seems likely to model biological reality better than the 2 extreme alternatives of either a strict clock or no clock at all (i.e., the unrooted method of Felsenstein 1981). The relaxed clock is also expected to have greater statistical power due to the smaller number of estimated parameters relative to the unrooted method (Drummond et al. 2006). MrBayes v3.1.2 (Huelsenbeck and Ronquist 2001) and BEAST v1.4.6 (Drummond and Rambaut 2007) were used to conduct Bayesian analyses under unrooted and relaxed-clock models, respectively. The data were divided a priori into partitions by gene and codon position, with models chosen using the Akaike information criterion as implemented in MrModelTest (Nylander 2004). Analyses were performed under several a priori partitioning schemes, and the harmonic means of likelihood scores from the posterior distribution were compared using Bayes factors to choose a final partitioning scheme (see Brandley et al. 2005).

We used the ML program RAxML (Stamatakis 2006) via its Web server (Stamatakis et al. 2008) to estimate a phylogeny and conduct nonparametric bootstrap analyses under the same partitioning scheme favored in the Bayesian analyses, except that all partitions were assigned a general time reversible model (currently the only option) with a proportion of invariable sites estimated and branch lengths optimized separately for each data partition. For simplicity, further references to the unrooted Bayesian, relaxed-clock Bayesian, and maximum-likelihood analyses will be referred to by the initials UB, RCB, and ML, respectively.

Divergence Time Estimates

In most cases, fossil constraints should place the maximum probability near the estimated age of the fossil, with probabilities dropping off rapidly for younger ages and more gradually for older ages (Benton and Donoghue 2007), effectively placing a hard bound on the minimum age and a soft bound (Yang and Rannala 2006) on the maximum age. A translated lognormal (TL) distribution models this situation well (Hedges and Kumar 2004; Drummond et al. 2006; Ho 2007) and was used for most of the calibrations in this study. Because overestimation of divergence times could lead to false rejection of Pliocene/Pleistocene speciation scenarios, we also ran analyses under a more conservative set of normally distributed fossil constraints in place of the more geologically appropriate lognormal constraints.

Fossils or geologic data suitable for constraining nodes within Brookesia are not available, and we therefore used several external calibration points for our analyses. For fossil-based calibrations, we constructed TL distributions in BEAST with the lower bound of the 95% highest probability density (HPD) at a point 1% more recent than the estimated fossil age, and we left (somewhat arbitrarily) long probability tails for the soft maxima. The following fossil-based calibration points were used (Table 3): 1) the rhynchocephalian–squamate split within lepidosaurs (Sues and Olsen 1990; Evans et al. 2001), 2) the oldest known stem scincomorphs (here defined as the clade containing the crown taxa Scincidae, Cordylidae, and Xantusiidae) (Evans 1993, 1998), 3) the oldest known anguimorph (Evans 1994, 1998), 4) the oldest known amphisbaenian (Gao and Nessov 1998), and 5) the oldest teiids (Winkler et al. 1990; Nydam and Cifelli 2002). Finally, the Comoros archipelago is home to the chameleon species Furcifer polleni and Furcifer cephalolepis, which are endemic to Mayotte and Grand Comoro, respectively. The oldest of the 2 islands (Mayotte) formed 10–15 million years ago (Nougier et al. 1986), and the older of these dates served as 6) a gamma-distributed maximum age constraint for the most recent common ancestor (mrca) of these taxa.

TABLE 3.

Divergence date calibration priors

Node Relevant fossila Median (95% CI) TL zero offset 
  (million years ago) (mean, SD) 
Rhyncocephalian–squamate Indeterminate taxon 231 (225–301) 224 (2.0, 1.2) 
Scincomorph–iguanian Balnealacerta 167 (162–204) 161 (1.8, 1.0) 
Anguimorph–iguanian Parviraptor 167 (162–204) 161 (1.8, 1.0) 
Lacertid–amphisbaenian Hodzhhakulia 102 (97–182) 97 (1.7, 1.4) 
Teiid–(lacertid, amphisbaenian) Ptilotodon 116 (110–187) 110 (1.8, 1.3) 
Furcifer polleni–Furcifer cephalolepis N/A 6.3 (1.7–16) 0.0 (3.5, 2.0)b 
Node Relevant fossila Median (95% CI) TL zero offset 
  (million years ago) (mean, SD) 
Rhyncocephalian–squamate Indeterminate taxon 231 (225–301) 224 (2.0, 1.2) 
Scincomorph–iguanian Balnealacerta 167 (162–204) 161 (1.8, 1.0) 
Anguimorph–iguanian Parviraptor 167 (162–204) 161 (1.8, 1.0) 
Lacertid–amphisbaenian Hodzhhakulia 102 (97–182) 97 (1.7, 1.4) 
Teiid–(lacertid, amphisbaenian) Ptilotodon 116 (110–187) 110 (1.8, 1.3) 
Furcifer polleni–Furcifer cephalolepis N/A 6.3 (1.7–16) 0.0 (3.5, 2.0)b 

Notes: SD, standard deviation; N/A, not available.

a

See text for references.

b

Nonfossil calibration based on geologic data and modeled with a gamma distribution: zero offset (γ-shape, γ-scale) (see text).

We were concerned about the effects that substitutional saturation at deeper levels might have on divergence time estimates, especially with regard to mitochondrial DNA (mtDNA) data (Hugall et al. 2007). However, the mtDNA data were needed to confidently resolve more recent branching points (see below), suggesting also that branch lengths within Brookesia would be better estimated by including mtDNA data. We therefore conducted analyses using 1) all data, 2) nuclear data plus first and second positions from the mtDNA data, and 3) nuclear data alone to test the robustness of our divergence time estimates. We also tested our fossil calibrations using the cross-validation method of Near et al. (2005), which identifies potentially problematic fossils that cause incongruent age estimates of other dated nodes in the tree.

All BEAST analyses were run for a sufficient number of generations to achieve an effective sample size of at least 200 for all estimated parameters, and 3 replicate runs were conducted for each analysis. Initial analyses were run without data to check the influence of the priors on the results. BEAST output was examined for evidence of proper mixing and convergence using Tracer (Rambaut and Drummond 2004), runs were combined using LogCombiner (part of BEAST package), and maximum credibility trees with divergence time means and 95% HPDs were produced using Tree Annotator (part of BEAST package). All BEAST files used are included in the online Supplementary Material.

Distribution Mapping and Diversity Analyses

Locality data for all Brookesia species of Madagascar were gathered from museum data, our own global positioning system readings, and pertinent literature. Single-species maps that almost fully agree with the data used herein are reproduced in Glaw and Vences (2007). Small sample sizes of locality records are not appropriate to obtain reliable estimates of potential distribution by modeling. For 9 species, more than 6 (range 7–39) locality records were available. For these species, we defined distributions by potential distribution models (pruned for overprediction), and for the rest of the species, we used point locality data.

For the models, we used 23 variables as predictors: potential evapotranspiration, yearly water balance (annual rainfall minus annual evapotranspiration), number of months with a positive water balance (a measure of drought), percentage of forest cover in 2000, and 19 climatic variables from the WorldClim database version 1.4 (Hijmans et al. 2005): annual mean temperature; mean diurnal temperature range; isothermality (monthly/annual temperature range); temperature seasonality (standard deviation across months); maximum temperature of warmest month; minimum temperature of coldest month; annual temperature range; mean temperature of wettest, driest, warmest, and coldest quarters; annual precipitation; precipitation of wettest and driest months; precipitation seasonality (coefficient of variation); and precipitation of wettest, driest, warmest, and coldest quarters.

We used Maxent v2.3 (Phillips et al. 2006) for distribution modeling, as it performs well in predicting species distributions (Elith et al. 2006), even when sample sizes are small (Hernandez et al. 2006). Analyses followed Vieites et al. (2008), using 1000 randomly selected data points across Madagascar as background pseudoabsence data. All models had areas under the receiver operating characteristic curve (Hanley and McNeil 1982) higher than 0.7, suggesting that the models were good at discriminating between presence and absence sites (Fielding and Bell 1997). We applied a pruning algorithm to remove areas of overprediction from the mean model. This algorithm thresholds the Maxent output models using a user-defined probability of occurrence (t), defining a convex hull around all occupied regions, and buffering this hull by a given number of grid cells (b) and width (f), within which the model values are reduced until they reach zero. We used the following parameters: t = 40, b = 40, and f = 80, as they provide a conservative scenario to remove biogeographic overprediction areas (see Kremen et al. 2008).

To calculate values of spatial species richness and endemism, the distributional data (potential distribution models and point distribution data) were transformed into a grid cell data set and plotted on a one-quarter degree square grid covering the entire island. When more than 2 records were available per species, we drew a minimum convex polygon between localities, and grid cells within the polygon were considered to contain the species. We used Worldmap v. 4.20.18. (Williams 2002) to calculate species richness as the total number of species per grid cell. Endemism was calculated as a measure of range-size rarity, expressed as the percentage aggregated reciprocal range size for all species per grid cell.

Hypothesis Testing

Temporal concordance with the Wilmé et al. (2006) and the Raxworthy and Nussbaum (1995) Pliocene/ Pleistocene speciation scenarios (WC and MR hypotheses, respectively) was statistically tested by comparing the most recent tail of the 95% HPD for each sister-species pair in our phylogeny to the oldest limits of these geological epochs. The RB hypothesis was not formulated with an explicit time frame and thus could not be tested in this manner.

The Wilmé et al. (2006) WC model predicts that sister species should generally occupy adjacent COEs or possibly a COE and an adjacent RDW. To test the fit of our Brookesia data to the spatial aspect of this model, we first plotted the distributions of all Brookesia sister-species pairs identified in our phylogenetic analyses onto the Wilmé et al. (2006) watershed map (Fig. 1) and counted the number of pairs matching the model's predictions. Next, we randomly assigned each species from all sister-species pairs to one of the watersheds collectively occupied by them and counted the number of pairs fitting the model's predictions. We repeated this 10,000 times to generate a null distribution of expected number of fits to the model if sister species are actually distributed randomly with respect to relative watershed positions. If less than 5% of the simulation rounds resulted in a number of matches equal to or greater than the number of matches from the real data, we rejected the null in favor of the Wilmé et al. (2006) WC hypothesis. We repeated this analysis under the interpretation that allopatric sister species occupying the same watershed also fit the model.

FIGURE 1.

Map showing proposed centers of endemism (COEs; dark gray), retreat-dispersion watersheds (RDWs; white), and the 3 highest peaks in Madagascar. Light gray areas are designated as RDWs, with possibly some function as COEs. Modified from Wilmé et al. (2006).

FIGURE 1.

Map showing proposed centers of endemism (COEs; dark gray), retreat-dispersion watersheds (RDWs; white), and the 3 highest peaks in Madagascar. Light gray areas are designated as RDWs, with possibly some function as COEs. Modified from Wilmé et al. (2006).

The MR hypothesis of Raxworthy and Nussbaum (1995) states that rainforest contraction around high-elevation massifs during cooler periods isolated populations on either side of unsuitable intervening habitat, leading to allopatric speciation. Referencing Figure 1, the 5 proposed areas of endemism (AOEs) of Raxworthy and Nussbaum (1995, see their fig. 7) correspond to 1) the Montagne d'Ambre mountain range in the far north (straddling the border between COEs 1 and 12); 2) a northwestern region west from Tsaratanana to Nosy Be and Presqu'Ile d'Ampasindava; 3) the Tsaratanana massif and surrounding high-altitude environs; 4) a northeast region extending east from Tsaratanana and encompassing parts of the RDWs A and a2, the southern tip of COE 1, and the northern section of COE 2, with the Antongil Bay watershed as its southwestern border; and finally 5) an eastern region corresponding to at least the southern section of COE 2.

Several factors make it impractical to test this MR hypothesis with our data using the statistical framework outlined above for the WC hypothesis. First, these AOEs comprise only a fraction of the land area of Madagascar (see Fig. 1), which would necessitate the exclusion of several sister-species pairs that are not confined to them (see Results section) and thus reduce statistical power. Also, there are multiple potential vicariant relationships among most of these areas (e.g., the northeastern region could have once been connected by suitable habitat to all the other AOEs), so the number of sister-species pairs consistent with the MR hypothesis under the null model (i.e., no real effect of forest contraction on speciation) would be high. Furthermore, as stated in Raxworthy and Nussbaum (1995), there is no clear biogeographic boundary between the northeastern and the eastern AOEs, making interpretation of relationships among taxa from these areas unclear. However, there is intuitive appeal to arguments that forest-fragment contraction promoted vicariant speciation in geologically complex northern Madagascar. We thus examine congruence of several potential examples with our new phylogenetic hypothesis (some of the proposed examples of Raxworthy and Nussbaum were flawed because they compared nonsister taxa; see Results section).

The RB hypothesis (e.g., Pastorini et al. 2003) predicts sister taxa to be separated by major rivers. A few Malagasy rivers can be unambiguously defined as “major” (e.g., the Betsiboka, Mangoro, and Mananara rivers) because they are quite long, with headwaters at relatively high elevations (Wilmé et al. 2006), and because previous works (e.g., Pastorini et al. 2003) have identified them as potential barriers for lemurs. However, for most rivers, an objective assessment of their potential as barriers would require combining parameters such as length, width, headwater elevation, and permanence. These factors, combined with the difficulty of determining the size and shape of the test areas, preclude a formalized test of the RB hypothesis. However, we examine the spatial relationship of Brookesia sister taxa with several of the largest rivers to look for any evidence of a major role in species-level diversification. Our geographic sampling within some species also allows some evaluation of the role that rivers might play in intraspecific genetic structuring within Brookesia.

The WC, RB, and MR hypotheses all make general predictions about relative species abundance at different altitudes. Following Wollenberg et al. (2008), we divided Madagascar into a coarser grid with cell sizes of 82 × 63 km=5166 km2 and used the ArcView extension “Endemicity Tools” (provided by N. Danho) to calculate Brookesia species richness and corrected weighted endemism (Crisp et al. 2001) for each cell. We then tested for nonparametric correlation of these values with altitudinal range, defined as maximum–minimum elevations per grid cell (excluding grids with no occurrence of any Brookesia species).

RESULTS

Data Characteristics and Phylogenetic Analyses

Alignments were unambiguous across chameleons for all partitions. However, alignment of the 3’ ends of ND1 (∼ 70 bp) and ND2 (∼ 90 bp), as well as part of ND4 (∼60 bp), was problematic across more distant outgroups due to high levels of amino acid replacement and multiple indels. We therefore coded these regions as missing data for nonchamaeleonid taxa. Maximum-likelihood analyses of near-complete ND1 and ND2 data across several chameleon species found very similar evolutionary model parameter values for these 2 genes (results not shown). To avoid the statistical problem of low data/parameter ratios, the ND1 sequence (∼ 70 bp) was combined with the ND2 for partitioning purposes. For the same reason, data from the 3 transfer RNA genes separating ND1 and ND2 (and from which the loop regions were impossible to align confidently) were not used in this study. Bayes factor analysis on various partitioning schemes (2, 3, 6, 9, and 12 partitions) strongly favored the 12-partition scheme (nuclear and mtDNA genes each partitioned by codon position; Bayes factor ∼ 178 in all comparisons). All results reported here are from this partitioning scheme.

Within each of the ML, UB, and RCB analysis sets, phylogenetic results were broadly congruent across the different genomes and genes, with topological conflicts involving only poorly supported nodes (i.e., ML bootstraps < 70%, Bayesian posterior probabilities [PPs] < 95%). One exception to this pattern involves strongly conflicting nuclear and mitochondrial results for the placement of 1 specimen of Brookesia brygooi. This discrepancy is probably due to ancient mtDNA introgression, and the nuclear topology almost certainly reflects the true organismal history (see online Supplementary Material). Mitochondrial data for this specimen were therefore excluded in all further phylogenetic analyses. In general, the mtDNA data alone gave strong support for relationships within more terminal clades, but poor support for basal nodes, and the nuclear data showed the opposite pattern. All results presented here are from analyses of the combined data set (Fig. 2). This data set and trees from this study can be found at www.treebase.org under accession number SN4455.

FIGURE 2.

Combined mitochondrial (ND1/ND2/ND4) and nuclear (RAG1/CMOS) data, 12-partition relaxed-clock Bayesian (RCB) cladogram with maximum-likelihood (ML) phylogram insert. Analyses included all ougroup taxa, although for clarity, most nonchameleon taxa are not shown. Major clades are color coded, and major clades of “typical” Brookesia are numbered (see text). RCB/unrooted Bayesian (UB) PP > 95% are denoted by asterisks (PP between 90% and 95% are shown as actual values) above branches, and ML bootstrap values > 50% are shown below branches. Black diamond indicates alternate attachment point for Brookesia stumpffi clade in the UB and ML analyses. Dashed terminal branches mark individuals represented by ND4 data only.

FIGURE 2.

Combined mitochondrial (ND1/ND2/ND4) and nuclear (RAG1/CMOS) data, 12-partition relaxed-clock Bayesian (RCB) cladogram with maximum-likelihood (ML) phylogram insert. Analyses included all ougroup taxa, although for clarity, most nonchameleon taxa are not shown. Major clades are color coded, and major clades of “typical” Brookesia are numbered (see text). RCB/unrooted Bayesian (UB) PP > 95% are denoted by asterisks (PP between 90% and 95% are shown as actual values) above branches, and ML bootstrap values > 50% are shown below branches. Black diamond indicates alternate attachment point for Brookesia stumpffi clade in the UB and ML analyses. Dashed terminal branches mark individuals represented by ND4 data only.

The mrca of the B. nasusB. lolontany clade and all other Brookesia appears early in the history of the group (although monophyly of neither Brookesia nor this subclade is significantly supported; Fig. 2). The remainder of the genus is clearly monophyletic and is divided into 2 well-supported clades. One of these clades corresponds to the extremely miniaturized B. minima group. The other major clade, representing the typical Brookesia morphotype, consists of a series of 5 maximally supported subclades (one of which is a single species) whose interrelationships are strongly supported by the RCB analysis but less so by the UB and ML analyses (Fig. 2).

Referencing the numbered clades in Figure 2, Clade 1 contains Brookesia perarmata, Brookesia decaryi, Brookesia bonsi, and B. brygooi. These species are restricted to the dry deciduous forests of the west and/or southwest and together form the sister taxon of the remaining typical Brookesia. The eastern rainforest sister species Brookesia therezieni and Brookesia superciliaris (Clade 2) and the northern rainforest species Brookesia ebenaui (Clade 3) are sister taxa to progressively less inclusive clades. Clade 4 contains several populations that show substantial sequence divergence from each other (11.1% average uncorrected ND4 divergence among specimens from distinct localities) but whose interrelationships are mostly poorly supported. The eastern rainforest species Brookesia thieli is paraphyletic, although basal divergences are poorly supported. One subclade of B. thieli appears to be the sister taxon of the morphologically distinctive Brookesia vadoni, suggesting the possibility of cryptic species within B. thieli. Finally, Clade 5 comprises several species (Brookesia antakarana, Brookesia ambreensis, Brookesia valerieae, Brookesia griveaudi, and Brookesia stumpffi) largely restricted to evergreen rainforest of the north, northeast, and northwest (B. stumpffi is exceptional among Brookesia in its adaptation to both dry deciduous and evergreen forest). Brookesia valerieae is strongly supported by all analyses as the sister taxon of B. griveaudi. Brookesia stumpffi is weakly supported as the sister taxon of a strongly monophyletic B. antakarana–B. ambreensis clade in the RCB analysis (Fig. 2). However, the other 2 methods find weak (ML, 54% bootstrap) to strong (UB, 96% PP) support for B. stumpffi as the sister taxon of the B. griveaudi–B. valerieae clade, with this larger clade as the sister group to the partially sympatric species pair B. antakarana and B. ambreensis (Fig. 2). These latter 2 species exhibit very low levels of divergence from one another (0.6–1.3% uncorrected ND4 distance; see Fig. 2 inset) and may not be reciprocally monophyletic. There do appear to be 2 distinct morphologies in this lineage (see photos in Glaw and Vences 2007), but mtDNA data from several additional specimens suggest that either incomplete lineage sorting or mitochondrial introgression is a major factor at this site (results not shown).

Generally, the RCB method gave a more precise phylogenetic estimate than did the UB method (83% vs. 72% of nodes significantly supported in Fig. 2). The increase in support may be a consequence of the expected increase in statistical power inherent to the relaxed-clock approach. Alternatively, model misspecification may play a role in the difference.

Geographical Centers of Diversity and Endemism in Brookesia

For the 9 modelable species, all distribution models predicted known localities, although comparison with museum records indicated some overprediction. The pruning algorithm removed some (but not all) areas appearing to have suitable environmental niches for the species, but which were located far from the known distribution ranges. This suggests that nonautecological (e.g., biogeographical, community ecological) factors are needed to explain the restricted microendemic distribution of some species.

Brookesia distributions span the island, although they are absent from many large areas in the central highlands and the south. Three main centers of diversity (Fig. 3a) are found in the northeast, north, and northwest, respectively. These areas are also COEs, with the highest percent values of range-size rarity (Fig. 3b). However, there are several areas, mainly in the west and northwest regions of the island, with a high degree of endemism corresponding to several species only known from single locality records.

FIGURE 3.

Maps illustrating (a) differential Brookesia species diversity (richness) across Madagascar, quantified as numbers of species present within one-quarter degree square grid cells covering the entire island and (b) degree of endemism in different parts of Madagascar. Endemism is calculated as a measure of range-size rarity, expressed as the percentage aggregated reciprocal range size for all species per one-quarter degree square grid cell.

FIGURE 3.

Maps illustrating (a) differential Brookesia species diversity (richness) across Madagascar, quantified as numbers of species present within one-quarter degree square grid cells covering the entire island and (b) degree of endemism in different parts of Madagascar. Endemism is calculated as a measure of range-size rarity, expressed as the percentage aggregated reciprocal range size for all species per one-quarter degree square grid cell.

The 3 main clades recovered in our phylogenetic analyses (Figs. 2 and 4) show distinct biogeographic patterns. The basally diverging B. nasus clade contains 2 species, B. nasus in the southeast and B. lolontany in the north, suggesting an old north-south connection. The B. minima clade shows the highest degree of endemism; most of the 11 species are restricted to very small areas in the north and north-central regions. The third and largest clade (typical Brookesia) spans most of the island, although centers of both diversity and endemism in this clade are in the north. Within this clade, B. brygooi, B. bonsi, B. decaryi, and B. perarmata (all western species; Clade 1, Fig. 2) form the sister taxon of the remaining species, none of which are restricted to the western forests.

FIGURE 4.

Chronogram from the 12-partition (ND1/ND2, ND4, CMOS, and RAG1, each partitioned by codon position) RCB phylogenetic analysis. Time units on scale bar in millions of years ago. See Supplementary Materials for exact divergence time means and 95% CIs for all nodes. Maps above chronogram illustrate patterns of species richness (quantified as numbers of species present within one-quarter degree square grid cells) for the 3 main clades (1 = Brookesianasus group, 2 =Brookesiaminima group, and 3=“typical” Brookesia group).

FIGURE 4.

Chronogram from the 12-partition (ND1/ND2, ND4, CMOS, and RAG1, each partitioned by codon position) RCB phylogenetic analysis. Time units on scale bar in millions of years ago. See Supplementary Materials for exact divergence time means and 95% CIs for all nodes. Maps above chronogram illustrate patterns of species richness (quantified as numbers of species present within one-quarter degree square grid cells) for the 3 main clades (1 = Brookesianasus group, 2 =Brookesiaminima group, and 3=“typical” Brookesia group).

Timing of Diversification in Brookesia

The RCB 12-partition analysis places the basal split within Brookesia at approximately 72 million years ago (95% confidence interval [CI] from ∼ 63 to ∼81 million years ago; Fig. 4), which is possibly older than any divergence within non-Brookesia chameleons (see Supplementary Material). Generally, divergences between recognized species groups and even sister taxa appear to be quite deep. Within the typical Brookesia clade, aside from the problematic B. ambreensis–B. antakarana complex, the most recent sister-species split (B. bonsi–B. decaryi) has a 95% CI extending no more recently than the Middle Miocene (Fig. 4). Within the B. minima clade, mean estimates for sister-taxon divergences are even older (the most recent divergence time mean is at the Eocene–Oligocene boundary), although the 95% CIs within the 2 clades overlap. Thus, we can confidently reject the temporal component (i.e., Pliocene–Pleistocene time frame) of both the MR (Raxworthy and Nussbaum 1995) and the WC (Wilmé et al. 2006) hypotheses.

Additional analyses using alternative calibration distributions and data sets had little effect on divergence time estimates. Analyses performed with normally distributed calibrations (with standard deviations arbitrarily set at 10% of the mean) gave average intrachameleon divergences about 10% more recent than those using lognormally distributed calibrations. Likewise, analyses that excluded mtDNA third codon positions and all mtDNA data gave respective divergence estimates on average about 15% and 18% more recent than the full-data analysis. In all these analyses, Brookesia sister-taxon divergences remained solidly Miocene in age, and thus, our above conclusions are unaffected.

Using the method of Near et al. (2005), we found no significant incongruence among our fossil calibrations with one exception; the root calibration (representing the rhynchocephalian/squamate split) was a significant outlier (P =0.017). Following Near et al. (2005), this calibration should be discarded as potentially misleading. Exclusion of this calibration point draws all extrachamaeleonid divergences to markedly more recent dates, and the estimated root age varies from approximately 170–150 million years ago. However, fossil rhynchocephalians are known from multiple Laurasian and Gondwanan localities by the Late Triassic (Evans et al. 2001). Thus, even allowing for reasonable error in fossil dates, Middle to Late Jurassic estimates for the Rhynchocephalia–Squamata split are clearly untenable (see Marshall 2008 for further discussion of the potential pitfalls of excluding outliers). We therefore favor the results from the full constraints analysis. That said, we note that even in analyses performed without the root constraint, in no case did the 95% HPDs for Brookesia sister-taxon divergences reach the Pliocene.

Fit of Spatial Hypotheses

We compared the distribution of Brookesia with the 12–15 COEs recently proposed by Wilmé et al. (2006) (Fig. 1). Of approximately 30 Brookesia species, 18 are limited to single COEs (3 of these are also found in 1 RDW), 9 are found in 2 or more separate COEs (1 is also found in an RDW), and 2 additional species are limited to single RDWs. However, these numbers are difficult to interpret for several reasons. First, Brookesia distributions are heavily skewed toward the northern end of the island, and certain COEs are occupied by disproportionate numbers of species. For example, 22 of 30 species (73%) are distributed across just 6 COEs, and 7 of these species are restricted to a single COE (S. Bemarivo/N. Mangoro [2] on the eastern and northeastern coast) (Fig. 1). Furthermore, the level of microendemism within Brookesia is extreme: 18 of 30 species are essentially known from single localities (Glaw and Vences 2007). Thus, most species of Brookesia will be confined to single COEs regardless of the mechanism of their isolation. But perhaps most importantly, simple tabulations of species counts across proposed COEs ignores phylogeny completely.

To test for spatial concordance with the WC hypothesis, we first plotted the distributions of the 10 strongly supported (PP > 95%) Brookesia sister-species pairs identified in our phylogenetic analysis (i.e., 20 species, see Fig. 2, also Supplementary Table S1) onto the Wilmé et al. watershed map (Fig. 1). Although none of these sister-species pairs are strictly confined to adjacent COEs, Brookesia tuberculata is found on Montagne d'Ambre, which straddles COEs 1 and 12, and its sister species Brookesia sp. “Montagne de Francais” occupies COE 1. Two sister-species pairs (Brookesia exarmata [8]/Brookesia dentata [G; also 9] and Brookesia sp.“Betampona” [2]/Brookesia ramanantsoai [B]) follow the adjacent COE/RDW pattern (Figs. 1 and 2). However, neither of these RDW residents is a wide-ranging species, as the Wilmé et al. (2006) hypothesis would predict. Nonetheless, accepting these 3 pairs as potential fits to the model, we performed 10,000 rounds of random allocations of the 20 relevant species (i.e., those belonging to sister-species pairs) to the watersheds collectively occupied by them and found that 12% of the simulation rounds allocated the members of 3 or more sister-species pairs to adjacent watersheds (i.e., P =0.12). When allopatric taxa occupying the same watershed were included as potential matches in the simulations (there were none of these in the real data), P =0.32. Thus, we found no support for the spatial component of the WC hypothesis.

Unlike the WC hypothesis, the MR hypothesis of Raxworthy and Nussbaum (1995) was formulated with respect to 5 AOEs in the north and east of Madagascar (see Materials and Methods section and Fig. 1), as opposed to watersheds covering the entire island. Of the 10 Brookesia sister-species pairs from our phylogeny (Fig. 2), only a subset of these taxa are distributed in these areas and therefore useful for evaluation of this hypothesis. Specifically, the species pairs B. antakarana–B. ambreensis, B. griveaudi–B. valeriae, Brookesia karchei–Brookesia peyrierasi, and B. tuberculata–B. sp. “Montagne de Francais” are wholly confined to these AOEs. Of these pairs, B. griveaudi (northeast) and B. valeriae (northwest) both inhabit forests less than about 900 m altitude and are absent from the higher altitude (> 1600 m) forests of the intervening Tsaratanana massif. As proposed by Raxworthy and Nussbaum (1995), these areas may once have been connected by an east-west corridor of even lower altitude forest that was lost during a period of aridification. Tsaratanana itself is home to B. lolontany, likely sister taxon of the only other Brookesia species with populations known to inhabit this altitudinal range, B. nasus of southeastern Madagascar. Although none of the other species pair distributions are consistent with this spatial model, 2 interesting higher level patterns within the B. minima group clade are potential fits. First, B. minima is restricted to the northwest, whereas its sister clade (B. tuberculata + B. sp. “Montagne des Francais”) is restricted to the Montagne d'Ambre region. Second, the clade formed by these 3 species is found only to the west of Tsaratanana, whereas its sister clade (B. cf. karchei + B. peyrierasi) is found only to the east of Tsaratanana.

Finally, a comparison of Brookesia distribution and phylogeny with major rivers yields no clear support for the RB hypothesis. The Betsiboka River in western Madagascar, which has been shown to be a major barrier for several lemur populations (Pastorini et al. 2003), may possibly lay between the sister species B. exarmata (to the southwest) and B. dentata (to the northeast). However, our B. dentata is the only known specimen from Ankaranfantsika National Park, which actually straddles the river. At any rate, B. exarmata is confined to the Tsingy de Bemaraha approximately 300 km (and several more rivers) away to the southwest. Similarly, B. bonsi is found at a single locality (Namoroka) about 75 and 100 km to the west of the large Mahavavy and Betsiboka rivers, respectively, and its sister species B. decaryi is known only from Ankarafantsika. Once again, the highly restricted and disjunct distributions of these species make any inference of a river barrier as the isolating mechanism highly speculative. No other sister-species pair (or higher level sister-clade pair) from our phylogeny is divided by what might reasonably be considered a major river.

River barriers can of course play a role in intraspecific genetic structuring (e.g., Pastorini et al. 2003; Kozak et al. 2006; Lemmon et al. 2007), and there are a few Brookesia species with distributions that span large rivers (B. nasus, B. superciliaris, and B. thieli in the east; B. stumpffi and possibly B. ebenaui in the north; and B. brygooi in the west). To fully study the phylogeographic effects of these rivers will require much more thorough sampling. However, our sampling within some species does allow at least some preliminary comments. The western B. brygooi has an extensive distribution, including populations on both sides of the Betsiboka River. However, the deepest split within B. brygooi separates a population to the south of the river from multiple populations spanning both sides of the river. Likewise, the distribution of B. stumpffi spans the Sambirano, Mahavavy, and Mananjeba rivers (among others) in northwestern Madagascar. Although the deepest split does separate a population to the south of the Sambirano River from several more northern populations, relationships among the remaining populations do not support the Mahavavy and Mananjeba rivers as barriers to gene flow. Finally, B. nasus is found both to the north and to the south of the large Mananara River in southeastern Madagascar. However, an individual from the northernmost population is nested within a clade of individuals from populations to the south of the potential barrier. See Supplementary Material for figures illustrating these scenarios. Thus, although further sampling could reveal a more complex pattern (possibly also involving smaller rivers), our data do not suggest a major river effect on genetic structuring within any of these species.

Species richness of Brookesia was negatively correlated with the minimum altitude (P < 0.05) and positively correlated with the maximum altitude (P < 0.001) in each grid cell. Species richness and corrected weighted endemism of Brookesia furthermore showed a significant correlation with altitudinal range per grid cell (P < 0. 001 and P < 0.05, respectively). These results are not congruent with the RB and WC hypotheses, both of which predict a higher species richness and endemism at lower elevations where rivers are widest and where river basins forming COEs have their headwaters. In contrast, the results are consistent with the MR hypothesis, which predicts more species will be restricted to montane regions where there is a higher probability of isolation during dry periods.

DISCUSSION

Our ability to evaluate past processes that produced current patterns of organismal distribution and diversity is constrained and biased by the relative availability of climatological and geological data from different periods in the earth's history. In general, Pleistocene climatic changes (specifically, glaciation cycles that have predictable effects on sea level, aridity, and temperature) are much better understood than those from the Pliocene, Miocene, and earlier (Wells 2003). This may explain their prominence in many glaciation-related diversification theories. Results from phylogeographic studies have challenged the importance of the most recent Pleistocene and Holocene glaciations as major drivers of species-level diversification for several groups (Klicka and Zink 1997; e.g., Avise et al. 1998; Rull 2008). In Madagascar, the existence of Pleistocene glacial activity and dynamic vegetational shifts has been documented (Burney 1995; Vidal-Romaní et al. 2002), but well-studied examples indicating the dependence of species diversification on these processes are lacking.

Most phylogeny-based biogeographical studies of Madagascar have centered on the relative importance of vicariance and dispersal as the source of Madagascar's many endemic higher taxa (Yoder and Nowak 2006 and references therein). Given the long-standing interest in Madagascar's high level of micorendemicity, it is surprising that only a few phylogenetic studies have explicitly estimated sibling-species divergence times (e.g., Vences et al. 2002; Raxworthy et al. 2008; Wirta et al. 2008). To our knowledge, none of these studies have found species-level divergence times compatible with Pleistocene–Pliocene radiations. Other mtDNA studies of various groups including geckos, skinks, spiders, and lemurs (Yoder et al. 2000; Vences et al. 2002; Schmitz et al. 2005; Wood et al. 2007; Jackman et al. 2008) have reported genetic divergences that are likewise incompatible with widespread Pleistocene–Pliocene speciation, using standard rates of mtDNA sequence evolution (but see Wirta and Montreuil 2008 for a possible Pliocene radiation of beetles). Thus, although well-sampled phylogenies of many more groups are needed, existing evidence (including this study) suggests that recent climatological cycles have had relatively minor effects on levels of species richness and associated microendemicity on Madagascar.

Theoretical considerations (e.g., Hewitt 1996) and empirical examples from other groups have suggested the importance of earlier glaciation cycles on the North American and European flora and fauna (Veith et al. 2003; Ayoub and Riechert 2004; Hewitt 2004; Kadereit et al. 2004; Turgeon et al. 2005). Unfortunately, historical species distributions in Madagascar have been very difficult to determine due to the paucity of Tertiary fossil deposits earlier than the Quaternary (Wells 2003). This has in turn greatly limited the ability to infer paleoclimates and historical vegetation cover that might help explain existing complex microendemic distribution patterns. Nonetheless, assuming relative stability of drainage systems and topographical features, resultant phylogenetic patterns could be very similar between earlier and more recent climatic cycles. An evaluation of the competing hypotheses thus depends largely on our knowledge of current distribution, habitat, and phylogeny.

Watershed Contractions, Mountain Refuges, and Riverine Barriers

Our test of the spatial component of the Wilmé et al. WC hypothesis is potentially biased in that phylogenetically closer species might be expected to have distributions geographically closer than more distantly related species. This could result in a nonrandom association between sister taxa and adjacent watersheds, even if WC was not the cause of their separation. Thus, our failure to reject a random watershed/sister-taxon association, combined with a negative correlation between altitude and species richness, strongly suggests that this mechanism cannot explain current Brookesia distributions. On the other hand, because species-level divergences in this group are relatively deep (Late Eocene to Late Miocene), there is the possibility that extinctions, climate/habitat changes, and range expansions have obscured patterns that would support this hypothesis. These are problems inherent to the study of any older radiation (Losos and Glor 2003). Unless new relevant fossils (including those indicative of past climates) are discovered, the first 2 problems are largely unavoidable. However, Brookesia may be less prone than many groups to the last potential problem due to the probable low vagility of these leaf litter inhabitants.

The Raxworthy and Nussbaum (1995) MR model, although restricted to northern Madagascar, seems to fit to our data reasonably well. Our study reveals that some of the sister-taxon pairs originally cited in favor of this model (B. ambreensis/B. valerieae and B. stumpffi /B. griveaudi) are actually not sister taxa but only more distantly related. Nonetheless, the Tsaratanana massif and associated montane areas do separate 2 sister taxa as well as 2 multispecies clades, and another clade is split between the Northwest and Montagne d'Ambre regions. We noted earlier that for any given forest fragment, multiple potential vicariance relationships often exist. This fact may account for some of the apparent fit to this model. However, independent support for this model is provided by the significant positive correlation of altitude and altitudinal range with species richness and endemism, which is the pattern expected if most speciation occurs by fragmentation and isolation in montane forest fragments.

Rivers appear to have played a substantial role in the diversification of several lemur groups (Pastorini et al. 2003; Goodman and Ganzhorn 2004; Louis et al. 2005), some frogs (Vieites et al. 2006), and possibly even some chameleons of the genus Furcifer (Raselimanana and Rakotomamalala 2003). However, this does not seem to be true for leaf chameleons, at least at the species level and above. It seems unlikely that individual Brookesia disperse over large distances, and they certainly cannot swim across rivers. However, given their small size and leaf litter microhabitat, they could quite plausibly float across rivers on debris washed away during heavy rains.

Higher Level Phylogenetic and Biogeographical Patterns

Wells (2003) hypothesized that the northward drift of Madagascar from its pre-Paleocene position south of the subtropical arid zone precipitated the aridification of most or all the island, causing the loss or great reduction of much of the existing biodiversity. Later entry into the “trade wind” zone of southeasterly winds (most likely sometime in the Eocene) brought rains to the east coast, causing the formation of the eastern humid forests. According to Wells (2003), this increased moisture, combined with flow of the warm southern equatorial current through the widening Mozambique Channel, also led to the formation of the western deciduous forests. Correlating this scenario with our dated phylogeny allows some speculation on drivers of diversification and distribution at higher levels within Brookesia. Specifically, aridification may have fragmented the ancient B. nasus clade into disjunct high montane areas, with expansion of B. nasus into lower altitudes of the eastern rainforest during the Eocene–Oligocene. Most diversification within both the B. minima and the typical Brookesia clades occurred during this same period and was thus likely tied to the same expansion of mesic habitats. The striking size disparity between members of these 2 clades, combined with largely overlapping distributions, suggests some sort of niche partitioning (most likely prey size and/or microhabitat usage), but data supporting or refuting this hypothesis are lacking.

Our geographical analyses confirm that northern Madagascar is the center of species richness and endemism of Brookesia (Fig. 3), as postulated by Raxworthy and Nussbaum (1995). However, within the clade of typical Brookesia, the most basal splits do not divide taxa restricted to the north. Instead, the first 2 splits separate a western and eastern clade, respectively, from the remainder, which itself contains mostly northern and a few more nested eastern species. This pattern suggests ancient east-west dispersal or vicariance with subsequent diversification within these areas and a later diversification in the north. Although the B. minima group is speciose, it is split at its base into a north-restricted clade and a more southern clade. Thus, hypotheses of northern or southern origins for this group are equally parsimonious. However, within the southern clade, there is a clear basal split between 1 western (dry deciduous) and 2 eastern (rainforest) groups. Interestingly, the repeated east-west splits within Brookesia represent a type of within-landmass biome shift in one of the lineages, a phenomenon that was found to be quite rare in a recent global-wide study of plants (Crisp et al. 2009).

In summary, higher level diversification within Brookesia was quite possibly tied to an Eocene–Oligocene shift to more mesic habitats on Madagascar, and our analyses support a substantial role for rainforest fragmentation (especially in the elevationally heterogeneous north) in species-level diversification. A true test of any model purported to generally explain the complex patterns of diversity in the Malagasy fauna will require at least a combination of distributional and current climatic data such as that provided by Wilmé et al. (2006) and a series of well-supported phylogenies (with reasonably accurate estimates of divergence times) from a variety of taxonomic groups. As paleoclimatological and paleovegetational reconstructions become available, they will also be very valuable because Malagasy terrestrial fossil deposits from the Tertiary are virtually nonexistent.

SUPPLEMENTARY MATERIAL

Supplementary material can be found at http://www.sysbio.oxfordjournals.org/.

FUNDING

This work was supported by grants from the Deutsche Forschungsgemeinschaft and the Volkswagen Foundation to M.V. and F.G. and by U.S. National Science Foundation grants (postdoctoral bioinformatics fellowship DBI-0204451 to T.M.T., Tree of Life NSF 1-10135 grant to D.R.V., and Tree of Life grant EF 0334923 to Tod Reeder).

We are grateful to Daniel Scantlebury for assistance in the lab and to Edward Louis of the Henry Doorly Zoo's Center for Conservation and Research for providing collecting opportunities for T.M.T. in Madagascar. Akira Mori generously provided important additional tissue samples. Franco Andreone, Parfait Bora, Kathrin Glaw, Fabio Mattioli, Roger Daniel Randrianiaina, and Jasmin E. Randrianirina assisted M.V., F.G., and D.R.V. during fieldwork. Katharina C. Wollenberg helped with the geographical analyses. Scott Kelley gave helpful Python coding advice. This research was carried out in the framework of collaboration accords with the Département de Biologie Animale, Université d'Antananarivo, and the Association Nationale des Aires Protegées ANGAP. We are indebted to the Malagasy authorities for issuing permits for research and export of specimens and samples. Tod Reeder, Bradford Hollingsworth, Ben Lowe, Adam Leaché, Marshal Hedin, Jack Sullivan, and an anonymous reviewer gave helpful comments on earlier versions of the manuscript.

References

Avise
JC
Walker
D
Johns
GC
Speciation durations and Pleistocene effects on vertebrate phylogeography. Proc. R. Soc. Lond
Ser. B Biol. Sci.
 , 
1998
, vol. 
265
 (pg. 
1707
-
1712
)
Ayoub
NA
Riechert
SE
Molecular evidence for Pleistocene glacial cycles driving diversification of a North American desert spider, Agelenopsis aperta
Mol. Ecol
 , 
2004
, vol. 
13
 (pg. 
3453
-
3465
)
Benton
MJ
Donoghue
PCJ
Paleontological evidence to date the tree of life
Mol. Biol. Evol.
 , 
2007
, vol. 
24
 (pg. 
26
-
53
)
Boumans
L
Vieites
DR
Glaw
F
Vences
M
Geographical patterns of deep mitochondrial differentiation in widespread Malagasy reptiles
Mol. Phylogenet. Evol.
 , 
2007
, vol. 
45
 (pg. 
822
-
839
)
Brandley
MC
Schmitz
A
Reeder
TW
Partitioned Bayesian analyses, partition choice, and the phylogenetic relationships of scincid lizards
Syst. Biol.
 , 
2005
, vol. 
54
 (pg. 
373
-
390
)
Briggs
JC
The biogeographic and tectonic history of India
J. Biogeogr
 , 
2003
, vol. 
30
 (pg. 
381
-
388
)
Burney
DA
Rakotondravony
D
, et al.  . 
Theories and facts on Holocene environmental change before and after human colonization. Environmental change in Madagascar
Fiovan'ny tontolo iainana eto Madagasikara
 , 
1995
Chicago
Field Museum of Natural History
(pg. 
108
-
109
)
Carpenter
AI
Robson
O
A review of the endemic chameleon genus
Brookesia from Madagascar, and the rationale for its listing on CITES Appendix II. Oryx
 , 
2005
, vol. 
39
 (pg. 
375
-
380
)
Crisp
MD
Arroyo
MTK
Cook
LG
Gandolfo
MA
Jordan
GJ
McGlone
MS
Weston
PH
Westoby
M
Wilf
P
Linder
HP
Phylogenetic biome conservatism on a global scale
Nature
 , 
2009
, vol. 
458
 (pg. 
754
-
756
)
Crisp
MD
Laffan
S
Linder
HP
Monro
A
Endemism in the Australian flora
J. Biogeogr
 , 
2001
, vol. 
28
 (pg. 
183
-
198
)
Drummond
AJ
Ho
SYW
Phillips
MJ
Rambaut
A
Relaxed phylogenetics and dating with confidence.
2006
 
PLoS Biol. 4:e88
Drummond
AJ
Rambaut
A
BEAST: Bayesian evolutionary analysis by sampling trees
BMC Evol. Biol.
 , 
2007
, vol. 
7
 pg. 
214
 
Elith
J
Graham
CH
Anderson
RP
Dudik
M
Ferrier
S
Guisan
A
Hijmans
RJ
Huettmann
F
Leathwick
JR
Lehmann
A
Li
J
Lohmann
LG
Loiselle
BA
Manion
G
Moritz
C
Nakamura
M
Nakazawa
Y
Overton
JM
Peterson
AT
Phillips
SJ
Richardson
K
Scachetti-Pereira
R
Schapire
RE
Soberon
J
Williams
S
Wisz
MS
Zimmermann
NE
Novel methods improve prediction of species’ distributions from occurrence data
Ecography
 , 
2006
, vol. 
29
 (pg. 
129
-
151
)
Evans
SE
Buffetaut
E
Mazin
J-M
Jurassic lizard assemblages
Proceedings of the Second Georges Cuvier Symposium
 , 
1993
(pg. 
55
-
65
Revue de Paleobiologie, Special Volume 7. Montbeliard (France)
Evans
SE
A new anguimorph lizard from the Jurassic and Lower Cretaceous of England
Palaeontology
 , 
1994
, vol. 
37
 (pg. 
33
-
49
)
Evans
SE
Crown-group lizards from the Middle Jurassic of Britain
Palaeontographica A
 , 
1998
, vol. 
250
 (pg. 
1
-
32
)
Evans
SE
Prasad
GVR
Manhas
BK
Rhynchocephalians (Diapsida: Lepidosauria) from the Jurassic Kota Formation of India
Zool. J. Linn. Soc.
 , 
2001
, vol. 
133
 (pg. 
309
-
334
)
Felsenstein
J
Evolutionary trees from DNA sequences: a maximum-likelihood approach
J. Mol. Evol.
 , 
1981
, vol. 
17
 (pg. 
368
-
376
)
Fielding
AH
Bell
JF
A review of methods for the assessment of prediction errors in conservation presence/absence models. Environ
Conserv
 , 
1997
, vol. 
24
 (pg. 
38
-
49
)
Gao
KQ
Nessov
LA
Early Cretaceous squamates from the Kyzylkum Desert
Uzbekistan. Neues Jahrb. Geol. Palaontol. Abh
 , 
1998
, vol. 
207
 (pg. 
289
-
309
)
Glaw
F
Vences
M
A fieldguide to the amphibians and reptiles of Madagascar.
2007
3rd ed
Köln (Germany)
Vences & Glaw
Goodman
SM
Benstead
JP
The natural history of Madagascar
 , 
2003
Chicago (IL)
The University of Chicago Press
Goodman
SM
Ganzhorn
JU
Biogeography of lemurs in the humid forests of Madagascar: the role of elevational distribution and rivers
J. Biogeogr
 , 
2004
, vol. 
31
 (pg. 
47
-
55
)
Groth
JG
Barrowclough
GF
Basal divergences in birds and the phylogenetic utility of the nuclear RAG-1 gene
Mol. Phylogenet. Evol.
 , 
1999
, vol. 
12
 (pg. 
115
-
123
)
Hanley
JA
McNeil
BJ
The meaning and use of the area under a receiver operating characteristic (Roc) curve
Radiology
 , 
1982
, vol. 
143
 (pg. 
29
-
36
)
Hedges
SB
Kumar
S
Precision of molecular time estimates
Trends Genet.
 , 
2004
, vol. 
20
 (pg. 
242
-
247
)
Hernandez
PA
Graham
CH
Master
LL
Albert
DL
The effect of sample size and species characteristics on performance of different species distribution modeling methods
Ecography
 , 
2006
, vol. 
29
 (pg. 
773
-
785
)
Hewitt
GM
Some genetic consequences of ice ages, and their role in divergence and speciation
Biol. J. Linn. Soc.
 , 
1996
, vol. 
58
 (pg. 
247
-
276
)
Hewitt
GM
Genetic consequences of climatic oscillations in the Quaternary. Philos. Trans. R. Soc. Lond
Ser. B Biol. Sci.
 , 
2004
, vol. 
359
 (pg. 
183
-
195
)
Hijmans
RJ
Cameron
SE
Parra
JL
Jones
PG
Jarvis
A
Very high resolution interpolated climate surfaces for global land areas
Int. J. Climatol
 , 
2005
, vol. 
25
 (pg. 
1965
-
1978
)
Ho
SYW
Calibrating molecular estimates of substitution rates and divergence times in birds
J. Avian Biol.
 , 
2007
, vol. 
38
 (pg. 
409
-
414
)
Huelsenbeck
JP
Ronquist
F
MRBAYES: Bayesian inference of phylogenetic trees
Bioinformatics
 , 
2001
, vol. 
17
 (pg. 
754
-
755
)
Hugall
AF
Foster
R
Lee
MSY
Calibration choice, rate smoothing, and the pattern of tetrapod diversification according to the long nuclear gene RAG-1
Syst. Biol.
 , 
2007
, vol. 
56
 (pg. 
543
-
563
)
Jackman
TR
Bauer
AM
Greenbaum
E
Glaw
F
Vences
M
Molecular phylogenetic relationships among species of the Malagasy-Comoran gecko genus Paroedura (Squamata: Gekkonidae)
Mol. Phylogenet. Evol.
 , 
2008
, vol. 
46
 (pg. 
74
-
81
)
Kadereit
JW
Griebeler
EM
Comes
HP
Quaternary diversification in European alpine plants: pattern and process. Philos. Trans. R. Soc. Lond
Ser. B Biol. Sci.
 , 
2004
, vol. 
359
 (pg. 
265
-
274
)
Klicka
J
Zink
RM
The importance of recent ice ages in speciation: a failed paradigm
Science
 , 
1997
, vol. 
277
 (pg. 
1666
-
1669
)
Kozak
KH
Blaine
RA
Larson
A
Gene lineages and eastern North American palaeodrainage basins: phylogeography and speciation in salamanders of the Eurycea bislineata species complex
Mol. Ecol
 , 
2006
, vol. 
15
 (pg. 
191
-
207
)
Kremen
C
Cameron
A
Moilanen
A
Phillips
SJ
Thomas
CD
Beentje
H
Dransfield
J
Fisher
BL
Glaw
F
Good
TC
Harper
GJ
Hijmans
RJ
Lees
DC
Louis
E
Nussbaum
RA
Raxworthy
CJ
Razafimpahanana
A
Schatz
GE
Vences
M
Vieites
DR
Wright
PC
Zjhra
ML
Aligning conservation priorities across taxa in Madagascar with high-resolution planning tools
Science
 , 
2008
, vol. 
320
 (pg. 
222
-
226
)
Kumazawa
Y
Nishida
M
Sequence evolution of mitochondrial transfer-RNA genes and deep-branch animal phylogenetics
J. Mol. Evol.
 , 
1993
, vol. 
37
 (pg. 
380
-
398
)
Lees
DC
Kremen
C
Andriamampianina
L
A null model for species richness gradients: bounded range overlap of butterflies and other rainforest endemics in Madagascar
Biol. J. Linn. Soc.
 , 
1999
, vol. 
67
 (pg. 
529
-
584
)
Lemmon
EM
Lemmon
AR
Cannatella
DC
Geological and climatic forces driving speciation in the continentally distributed trilling chorus frogs (Pseudacris)
Evolution
 , 
2007
, vol. 
61
 (pg. 
2086
-
2103
)
Losos
JB
Glor
RE
Phylogenetic comparative methods and the geography of speciation
Trends Ecol. Evol.
 , 
2003
, vol. 
18
 (pg. 
220
-
227
)
Louis
EE
Ratsimbazafy
JH
Razakamaharauo
VR
Pierson
DJ
Barber
RC
Brenneman
RA
Conservation genetics of black and white ruffed lemurs, Varecia variegata, from southeastern Madagascar. Anim
Conserv
 , 
2005
, vol. 
8
 (pg. 
105
-
111
)
Macey
JR
Larson
A
Ananjeva
NB
Fang
ZL
Papenfuss
TJ
Two novel gene orders and the role of light-strand replication in rearrangement of the vertebrate mitochondrial genome
Mol. Biol. Evol.
 , 
1997
, vol. 
14
 (pg. 
91
-
104
)
Maddison
DR
Maddison
WP
analysis of phylogeny and character evolution. Version 4.03
Tucson
 , 
2000
Sinauer
 
MacClade 4 (AZ)
Marshall
CR
A simple method for bracketing absolute divergence times on molecular phylogenies using multiple fossil calibration points
Am. Nat
 , 
2008
, vol. 
171
 (pg. 
726
-
742
)
Near
TJ
Meylan
PA
Shaffer
HB
Assessing concordance of fossil calibration points in molecular clock studies: an example using turtles
Am. Nat
 , 
2005
, vol. 
165
 (pg. 
137
-
146
)
Noonan
BP
Chippindale
PT
Vicariant origin of Malagasy reptiles supports Late Cretaceous Antarctic land bridge
Am. Nat
 , 
2006
, vol. 
168
 (pg. 
730
-
741
)
Nougier
J
Cantagrel
JM
Karche
JP
The Comores archipelago in the western Indian Ocean: volcanology, geochronology, and geodynamic setting
J. Afr. Earth Sci.
 , 
1986
, vol. 
5
 (pg. 
135
-
145
)
Nussbaum
RA
Raxworthy
CJ
Revision of the genus Ebenavia Boettger (Reptilia: Squamata: Gekkonidae)
Herpetologica
 , 
1998
, vol. 
54
 (pg. 
18
-
34
)
Nussbaum
RA
Raxworthy
CJ
Pronk
O
The ghost geckos of Madagascar: a further revision of the Malagasy leaf-toed geckos (Reptilia, Squamata, Gekkonidae)
Misc. Publ. Mus. Zool
 , 
1998
Univ. Michigan 186
(pg. 
1
-
26
)
Nydam
RL
Cifelli
RL
Lizards from the Lower Cretaceous (Aptian-Albian) Antlers and Cloverly formations
J. Vert. Paleontol
 , 
2002
, vol. 
22
 (pg. 
286
-
298
)
Nylander
JAA
Distributed by the author. Uppsula (Sweden): Evolutionary Biology Centre
2004
Uppsala University
 
MrModelTest v.2
Pastorini
J
Thalmann
U
Martin
RD
A molecular approach to comparative phylogeography of extant Malagasy lemurs
Proc. Natl. Acad. Sci. USA
 , 
2003
, vol. 
100
 (pg. 
5879
-
5884
)
Pearson
RG
Raxworthy
CJ
The evolution of local endemism in Madagascar: watershed versus climatic gradient hypotheses evaluated by null biogeographic models
Evolution
 , 
2009
, vol. 
63
 (pg. 
959
-
967
)
Phillips
SJ
Anderson
RP
Schapire
RE
Maximum entropy modeling of species geographic distributions
Ecol. Modell
 , 
2006
, vol. 
190
 (pg. 
231
-
259
)
 
Rambaut A., Drummond A.J. 2004. Tracer v1.3. Available from http://beast.bio.ed.ac.uk/Tracer
Raselimanana
AP
Rakotomamalala
D
Goodman
SM
Benstead
JP
Chamaeleonidae Chameleons
The natural history of Madagascar
 , 
2003
Chicago (IL)
The University of Chicago Press
(pg. 
960
-
972
)
Raxworthy
CJ
Forstner
MRJ
Nussbaum
RA
Chameleon radiation by oceanic dispersal
Nature
 , 
2002
, vol. 
415
 (pg. 
784
-
787
)
Raxworthy
CJ
Nussbaum
RA
Systematics, speciation and biogeography of the dwarf chameleons (Brookesia, Reptilia, Squamata, Chamaeleontidae) of northern Madagascar
J. Zool
 , 
1995
, vol. 
235
 (pg. 
525
-
558
)
Raxworthy
CJ
Pearson
RG
Zimkus
BM
Reddy
S
Deo
AJ
Nussbaum
RA
Ingram
CM
Continental speciation in the tropics: contrasting biogeographic patterns of divergence in the Uroplatus leaf-tailed gecko radiation of Madagascar
J. Zool
 , 
2008
, vol. 
275
 (pg. 
423
-
440
)
Rieppel
O
The phylogenetic relationships within the Chamaeleonidae, with comments on some aspects of cladistic analysis. Zool
J. Linn. Soc.
 , 
1987
, vol. 
89
 (pg. 
41
-
62
)
Rull
V
Speciation timing and neotropical biodiversity: the Tertiary-Quaternary debate in the light of molecular phylogenetic evidence
Mol. Ecol
 , 
2008
, vol. 
17
 (pg. 
2722
-
2729
)
Schmitz
A
Brandley
MC
Mausfeld
P
Vences
M
Glaw
F
Nussbaum
RA
Reeder
TW
Opening the black box: phylogenetics and morphological evolution of the Malagasy fossorial lizards of the subfamily “Scincinae”
Mol. Phylogenet. Evol.
 , 
2005
, vol. 
34
 (pg. 
118
-
133
)
Stamatakis
A
RAxML-VI-HPC: maximum likelihood-based phylogenetic analyses with thousands of taxa and mixed models
Bioinformatics
 , 
2006
, vol. 
22
 (pg. 
2688
-
2690
)
Stamatakis
A
Hoover
P
Rougemont
J
A rapid bootstrap algorithm for the RAxML Web servers
Syst. Biol.
 , 
2008
, vol. 
57
 (pg. 
758
-
771
)
Storey
M
Mahoney
JJ
Saunders
AD
Duncan
RA
Kelley
SP
Coffin
MF
Timing of hot spot-related volcanism and the breakup of Madagascar and India
Science
 , 
1995
, vol. 
267
 (pg. 
852
-
855
)
Sues
HD
Olsen
PE
Triassic vertebrates of Gondwanan aspect from the Richmond basin of Virginia
Science
 , 
1990
, vol. 
249
 (pg. 
1020
-
1023
)
Thompson
JD
Higgins
DG
Gibson
TJ
Clustal-W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice
Nucleic Acids Res.
 , 
1994
, vol. 
22
 (pg. 
4673
-
4680
)
Townsend
TM
Larson
A
Molecular phylogenetics and mitochondrial genomic evolution in the Chamaeleonidae (Reptilia, Squamata)
Mol. Phylogenet. Evol.
 , 
2002
, vol. 
23
 (pg. 
22
-
36
)
Townsend
TM
Larson
A
Louis
E
Macey
JR
Molecular phylogenetics of Squamata: the position of snakes, amphisbaenians, and dibamids, and the root of the squamate tree
Syst. Biol.
 , 
2004
, vol. 
53
 (pg. 
735
-
757
)
Turgeon
J
Stoks
R
Thum
RA
Brown
JM
McPeek
MA
Simultaneous Quaternary radiations of three damselfly clades across the Holarctic
Am. Nat
 , 
2005
, vol. 
165
 (pg. 
E78
-
E107
)
 
Van Bocxlaer I., Roelants K., Biju S., Nagaraju J., Bossuyt F. 2007. Late Cretaceous vicariance in Gondwanan amphibians. PLoS One. 1:e74. doi:10.1371/journal.pone.0000074
van der Meijden
A
Vences
M
Hoegg
S
Boistel
R
Channing
A
Meyer
A
Nuclear gene phylogeny of narrow-mouthed toads (Family: Microhylidae) and a discussion of competing hypotheses concerning their biogeographical origins
Mol. Phylogenet. Evol.
 , 
2007
, vol. 
44
 (pg. 
1017
-
1030
)
Veith
M
Schmidtler
JF
Kosuch
J
Baran
I
Seitz
A
Palaeoclimatic changes explain Anatolian mountain frog evolution: a test for alternating vicariance and dispersal events
Mol. Ecol
 , 
2003
, vol. 
12
 (pg. 
185
-
199
)
Vences
M
Andreone
F
Glaw
F
Kosuch
J
Meyer
A
Schaefer
HC
Veith
M
Exploring the potential of life-history key innovation: brook breeding in the radiation of the Malagasy treefrog genus Boophis
Mol. Ecol
 , 
2002
, vol. 
11
 (pg. 
1453
-
1463
)
Vences
M
Glaw
F
Molecular phylogeography of Boophis tephraeomystax: a test case for east-west vicariance in Malagasy anurans (Amphibia, Anura, Mantellidae)
Spixiana
 , 
2002
, vol. 
25
 (pg. 
79
-
84
)
Vences
M
Glaw
F
Phylogeography, systematics and conservation status of boid snakes from Madagascar (Sanzinia and Acrantophis)
Salamandra
 , 
2003
, vol. 
39
 (pg. 
181
-
206
)
Vences
M
Wollenberg
KC
Vieites
DR
Lees
DC
Madagascar as a model region of species diversification
Trends Ecol. Evol.
 , 
2009
, vol. 
24
 (pg. 
456
-
465
)
Vidal-Romaní
JR
Mosquera
DF
Campos
ML
A 12,000 yr BP record from Andringitra Massif, (southern Madagascar): post-glacial environmental evolution from geomorphological and sedimentary evidence. Quat
Int
 , 
2002
, vol. 
93–4
 (pg. 
45
-
51
)
Vieites
DR
Chiari
Y
Vences
M
Andreone
F
Rabemananjara
F
Bora
P
Nieto-Roman
S
Meyer
A
Mitochondrial evidence for distinct phylogeographic units in the endangered Malagasy poison frog Mantella bernhardi
Mol. Ecol
 , 
2006
, vol. 
15
 (pg. 
1617
-
1625
)
Vieites
DR
Nieto-Román
S
Vences
M
Andreone
F
Towards understanding the spatial pattern of amphibian diversity in Madagascar
A conservation strategy for the amphibians of Madagascar. Turin (Italy): Monografie del Museo Regionale di Scienze Naturali di Torino
 , 
2008
(pg. 
397
-
410
)
Wells
NA
Goodman
SM
Benstead
JP
Some hypotheses on the Mesozoic and Cenozoic paleoenvironmental history of Madagascar
The natural history of Madagascar
 , 
2003
Chicago: (IL)
The University of Chicago Press
(pg. 
16
-
34
)
Williams
PH
Using WORLDMAP
Priority areas for biodiversity v4.20. Distributed by the author
 , 
2002
London
The Natural History Museum
Wilmé
L
Goodman
SM
Ganzhorn
JU
Biogeographic evolution of Madagascar's microendemic biota
Science
 , 
2006
, vol. 
312
 (pg. 
1063
-
1065
)
Winkler
DA
Murry
PA
Jacobs
LL
Early Cretaceous (Comanchean) vertebrates of central Texas
J. Vert. Paleontol
 , 
1990
, vol. 
10
 (pg. 
95
-
116
)
Wirta
H
Montreuil
O
Evolution of the Canthonini Longitarsi (Scarabaeidae) in Madagascar. Zool
Scr
 , 
2008
, vol. 
37
 (pg. 
651
-
663
)
Wirta
H
Orsini
L
Hanski
I
An old adaptive radiation of forest dung beetles in Madagascar
Mol. Phylogenet. Evol.
 , 
2008
, vol. 
47
 (pg. 
1076
-
1089
)
Wollenberg
KC
Vieites
DR
van der Meijden
A
Glaw
F
Cannatella
DC
Vences
M
Patterns of endemism and species richness in Malagasy cophyline frogs support a key role of mountainous areas for speciation
Evolution
 , 
2008
, vol. 
62
 (pg. 
1890
-
1907
)
Wood
HM
Griswold
CE
Spicer
GS
Phylogenetic relationships within an endemic group of Malagasy `assassin spiders’ (Araneae, Archaeidae): ancestral character reconstruction, convergent evolution and biogeography
Mol. Phylogenet. Evol.
 , 
2007
, vol. 
45
 (pg. 
612
-
619
)
Xia
X
Xie
Z
DAMBE: software package for data analysis in molecular biology and evolution
J. Hered
 , 
2001
, vol. 
92
 (pg. 
371
-
373
)
Yang
ZH
Rannala
B
Bayesian estimation of species divergence times under a molecular clock using multiple fossil calibrations with soft bounds
Mol. Biol. Evol.
 , 
2006
, vol. 
23
 (pg. 
212
-
226
)
Yoder
AD
Heckman
KL
Fleagle
J
Lehman
SM
Mouse lemur phylogeography revises a model of ecographic constraint in Madagascar
Primate biogeography: progress and prospects
 , 
2006
New York
Kluwer
(pg. 
255
-
268
Press
Yoder
AD
Nowak
MD
Has vicariance or dispersal been the predominant biogeographic force in Madagascar? Only time will tell
Annu. Rev. Ecol. Evol. Syst
 , 
2006
, vol. 
37
 (pg. 
405
-
431
)
Yoder
AD
Rasoloarison
RM
Goodman
SM
Irwin
JA
Atsalis
S
Ravosa
MJ
Ganzhorn
JU
Remarkable species diversity in Malagasy mouse lemurs (Primates: Microcebus)
Proc. Natl. Acad. Sci. USA
 , 
2000
, vol. 
97
 (pg. 
11325
-
11330
)

Author notes

Associate Editor: Marshal Hedin

Supplementary data