| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
© 2005 Society of Systematic Biologists
Intraspecific Variability and Timing in Ancestral Ecology Reconstruction: A Test Case from the Cape Flora
Edited by Todd Oakley: Associate Editor
1 Institute of Systematic Botany Zollikerstrasse 107, CH-8008 Zurich, Switzerland E-mail: chardy{at}systbot.unizh.ch (C.R.H.)
| Abstract |
|---|
|
|
|---|
Thamnochortus (ca. 32 species) is an ecologically diverse genus of Restionaceae. Restionaceae comprise a major component of the southern African Cape flora, wherein eco-diversification might have been important in the generation of high levels of species richness. In an attempt to reconstruct the macroecological history of Thamnochortus, it was found that standard procedures for character state optimization make two inappropriate assumptions. The first is that ancestors are monomorphic (i.e., ecologically uniform) and the second is that eco-diversification follows, or is slower than, lineage diversification. We demonstrate a variety of coding schemes with which the assumption of monomorphy can be avoided. For unordered discrete ecological characters, presence coding and generalized frequency coding (GFC) are suboptimal because they occasionally yield illogical assignments of no state to ancestors. Polymorphism coding or use of the program DIVA are preferable in this respect but are applicable only with parsimony. For continuous eco-characters (e.g., a rainfall gradient, where individual species occur in ranges), GFC and MaxMin coding provide equally valid solutions to optimizing ranges with parsimony. However, MaxMin can be extended to likelihood approaches and is therefore preferable. With respect to rates and timing, all algorithms currently employed for ancestral ecology reconstruction bias toward slow rates of eco-diversification relative to lineage diversification. An alternative to this bias is provided by DIVA, which biases toward accelerated rates of eco-diversification and thus inferences of ecology-driven speciation. We see no way of choosing between these biases; however, phylogeneticists should be aware of them. Applying these methods to Thamnochortus, we find there to be important differences in details, yet general congruence, regarding the historical ecology of this clade. We infer the most recent common ancestor of Thamnochortus to have been a post-fire resprouting species distributed on rocky, well-drained, sandstone-derived soils at lower-middle elevations, in regions of moderate levels of yearly (primarily winter) rainfall. This species would have been distributed in habitats much like those of the southwestern Cape mountains today. Major ecological trends include shifts to lower rainfall regimes and shifts from sandstone to limestone-derived alkaline soils at lower altitudes.
Keywords: Ancestor reconstruction; Cape flora; evolutionary ecology; generalized frequency coding; macroecology; optimization; Restionaceae; Thamnochortus
Received February 6, 2004; Revised June 20, 2004; Accepted August 24, 2004
In addition to genetic and morphological diversification, a clade's evolutionary history includes ecological diversification. Our ability to reconstruct a clade's ecological history and to test related hypotheses relies on our ability to infer the habitats and other ecological attributes of ancestral species. However, commonly employed methods such as Fitch and Wagner optimizations were originally developed as step-counting algorithms in order to calculate and compare tree lengths during tree searches (Farris, 1970; Fitch, 1971; Swofford and Maddison, 1987). Consequently, these may not be appropriate for ecologically meaningful reconstructions. In fact, all coding and optimization procedures most commonly used for ancestor reconstruction share two fundamental implicit assumptions that are often violated when applied to ecological characters.
The most prominent assumption is that of ancestral monomorphy. That is, optimization algorithms allow just one state per character per internal node. This assumption may not be valid for many ecological characters (e.g., soil type), because individual species are often ecologically variable (Sterelny, 1999). Assuming such species are phylogenetic species (sensu Nixon and Wheeler, 1990), subdivision into less inclusive "monomorphic" units (as per Nixon and Davis, 1991) is not appropriate. The existence of ecologically variable (polymorphic) extant species demands the possibility of reconstructing polymorphic ancestors. Related to this problem of polymorphism is the fact that most species are distributed along portions of habitat gradients. A species' range along some habitat gradient is a special form of intraspecific polymorphism, and when this occurs for a potentially ordered character such as altitude, optimization methods that account for not only this special polymorphism but also the potentially additive nature of possible change may be preferable.
The second fundamental assumption concerns the timing of character transformations. For example, Fitch optimization for discretely coded characters parsimoniously assigns character state transformations to the internodes following the cladogenic events that lead to the lineages displaying the character differences. We demonstrate here that all standard methods for the optimization of both discrete and continuous characters, whether parsimony or likelihood based, show this same general bias with respect to timing. This precludes the possibility that the ecological variations found in descendant lineages arose prior to speciation, and that speciation was associated with a subsequent split along those ecological boundaries (i.e., ecological vicariance). Optimization methods that allow for this latter alternative are desired, and the sensitivity of one's conclusions to such issues of timing might be investigated.
Here we draw attention to these normally unspoken assumptions and explore a number of possible solutions in an empirical framework provided by the Restionaceae (Poales) genus Thamnochortus (32 species; Linder, 1991, 2002). Thamnochortus (Fig. 1A) is an ecologically diverse clade of wind-pollinated graminoids largely restricted to the Cape Floristic Region (Goldblatt, 1978) of southern Africa, a region where ecological factors may have played an important role in the origin and maintenance of unusually high levels of plant species richness (Linder, 1985, 2003; Cowling, 1987; Cowling et al., 1996). Thus, Thamnochortus represents a suitable framework with which to investigate these important issues in the inference of ancestral ecology. Other critical issues for ancestor reconstruction in general, such as accounting for uncertainty in phylogenetic hypotheses or in the estimates of particular states at particular nodes, have already received considerable attention (e.g., Frumhoff and Reeve, 1994; Maddison, 1995; Donoghue and Ackerly, 1996; Schluter et al., 1997; Cunningham et al., 1998; Sharkey, 1999; Huelsenbeck et al., 2003). Although there is still much work to be done in these areas, we do not pursue them here. Likewise, similar attention has been given to the fact that different optimization procedures (e.g., parsimony versus likelihood, etc.) may yield incongruent results (e.g., Swofford and Maddison, 1992; Butler and Losos, 1997; Martins and Hansen, 1997; Cunningham et al., 1998; Omland, 1999; Losos, 1999). Although we apply both likelihood and parsimony methods here, the differences between such methods are well understood and a reiteration of these differences is beyond the scope of this article. Instead, we focus here on the more neglected challenges of dealing with polymorphism, continuous range data, and issues of relative timing in the reconstruction of ancestral habitats and ecologies.
|
| MATERIALS AND METHODS |
|---|
|
|
|---|
Thamnochortus Phylogeny
The Thamnochortus phylogeny utilized for this study includes 30 of 32 recognized species, each represented by a single sample (Table 1). Of the two unsampled species, Th. amoena H.P. Linder is known only from a single population killed by a recent fire and Th. ellipticus Pillans is known only from the type collection, is possibly not distinct from Th. guthrieae, and attempts at polymerase chain reaction (PCR) amplification from herbarium specimens have failed. The family-level analyses of Linder (1984), Linder et al. (2000), and Eldenäs and Linder (2000) have resolved Rhodocoma as sister to Thamnochortus, a relationship supported by their shared possession of pendulous male inflorescences (Fig. 1A), scattered cavities in the culm central ground tissue, and a distinctively thin-walled and only slightly raised circular pollen aperture (Linder, 1984). Accordingly, outgroup sampling included five of the eight species of Rhodocoma (Table 1). However, preliminary results from an ongoing phylogenetic study of the African Restionaceae as a whole indicate that while Thamnochortus and Rhodocoma are closely related, they may not be sister taxa (Hardy and Linder, unpublished), and additional outgroup taxa (Table 1) were chosen in accordance with these results.
|
DNA sequencing
DNA sequences were generated from the plastid regions spanning the trnL intron and the trnL-trnF intergenic spacer (Taberlet et al., 1991), the complete gene encoding rbcL (Chase and Albert, 1998), the complete atpB-rbcL intergenic spacer (Manen et al., 1994), and matK plus the flanking trnK intron (Hilu and Liang, 1997). Total DNA was isolated from silica gel-dried culms using the Dneasy Plant Mini Kit (Qiagen, Inc., Valencia, California, USA). Sequences for trnL-F were obtained as described in Eldenäs and Linder (2000). The two regions spanning the contiguous atpB-rbcL spacer plus rbcL, as well as matK and the flanking trnK intron, were each amplified from a single polymerase chain reaction using the primers designated in Table 2. Sequences were generated using standard methods for automated sequencing using the primers designated in Table 2.
|
Phylogenetic analysis
Raw sequence data files were analyzed with the ABI Prism 377 Software Collection 2.1. Contigs were constructed in Sequencher and alignments were performed using the default alignment parameters in Clustal X (Thompson et al., 1997), followed by adjustment by eye. These sequences were assembled into a single matrix in WinClada (Nixon, 2002). Indels were coded at the end of the matrix using Simple Indel Coding (Simmons and Ochoterena, 2000) as implemented in the program GapCoder (Young and Healy, 2001). The data matrix used in the analysis is available from the authors and is deposited at TreeBASE (http://www.treebase.org; accession number SN1765). Sequences were deposited in GenBank (Table 1). Heuristic parsimony searches were conducted with NoNa (Goloboff, 1993), run as a daughter process from WinClada. Eleven thousand tree searches were conducted, with each search initiated with the generation of a Wagner tree, using a random taxon entry sequence, and followed by tree bisection and reconnection (TBR) swapping on the Wagner tree, with one shortest tree retained and subjected to branch swapping. All most parsimonious trees accumulated during these searches then were subjected to TBR swapping, including swapping on all trees propagated during this phase of the search. Nonparametric bootstrap support values (Felsenstein, 1985) were obtained using NoNa spawned as a daughter process in WinClada using 1000 replicates with 100 TBR searches per replicate, holding one tree per replicate, followed by "max*" to swap to completion. Percentages were then based on the strict consensus tree of each of the 1000 replicates.
A tree search was also conducted using maximum likelihood (ML) in PAUP* 4.0 (Swofford, 2002). We used ModelTest 3.06 (Posada and Crandall, 1998) in tandem with PAUP* to select a statistically adequate model of sequence evolution from a possible 56 models for a modified matrix (i.e., excluding indel characters). The ModelTest default of a neighbor-joining tree was used. A general time reversible model was selected in which the proportion of invariant sites is estimated and among site rate variation is modeled using a discrete approximation of a gamma distribution with four rate categories. Using this model (GTR+I+G), the tree with the highest likelihood was estimated.
The parsimony analysis yielded two fully resolved most parsimonious trees under ACCTRAN optimization (i.e., some nodes had no character support under DELTRAN). One of these two trees (Fig. 1B) placed Th. pluristachyus (Fig. 1A) sister to the clade (Th. paniculatus (Th. fraternus, Th. spicigerus)) and the other placed Th. paniculatus sister to the clade (Th. pluristachyus (Th. fraternus, Th. spicigerus)). From the ML analysis, the tree with the highest likelihood (ln likelihood –13,032.26) was also congruent with the second parsimony-derived tree (i.e., that in which Th. paniculatus is sister to the clade comprising Th. fraternus and Th. spicigerus). The only differences were that four of the five Thamnochortus nodes circled in Fig. 1B were resolved by indel characters with parsimony, yet were unresolved (of zero length) by the ML analysis (which excluded indel characters). On the basis of this agreement between ML and the second most parsimonious tree, and because the parsimony tree was fully resolved through the addition of indel characters, the parsimony tree (Fig. 1B) was chosen for subsequent analysis.
Dating of nodes
A likelihood ratio test (Felsenstein, 1981) was conducted to assess interlineage rate heterogeneity. Using the model determined above and PAUP*, likelihoods with and without assuming a clock were calculated for the selected most parsimonious topology (Fig. 1B). As a clock was rejected (
2 = 79.8, df = 37, P < 0.05), ML branch lengths were estimated without a clock assumption. These branches were then made ultrametric using penalized likelihood as implemented in r8s (Sanderson, 2002a, 2002b). The optimal rate-smoothing parameter (10–2.50) was estimated using the cross-validation procedure described by Sanderson (2000a). Node ages were estimated using the Truncated Newton algorithm with no internal age constraints. As we were only interested in relative node ages, the root node was arbitrarily set at 10.0. The ultrametric tree (Fig. 1B) was viewed and printed from TreeView 1.6.6 (Page, 1996).
Macroecological Parameters
The set of ecological parameters described by Linder and Hardy (in press) were scored for each species (Table 3) and are summarized in Appendix 1 (available at the Society of Systematic Biologists Web site, http://systematicbiology.org). Data were obtained from label data in the Bolus Herbarium, in addition to field observations of numerous populations over a period of 20 years.
|
Coding and Optimizing Macroecological Parameters
AsIs coding
What we have termed "AsIs coding" treats each ecological parameter as one character. It creates a data matrix that is a straightforward, untransformed version of Table 3. For continuously variable parameters, this is equivalent to continuous coding. For discrete parameters, it is equivalent to multistate (including binary) coding, which is the normal procedure for morphological and molecular characters for cladistic analysis and ancestor reconstruction.
For unordered discrete parameters, the applicable methods are Fitch (Fitch, 1971; here implemented in WinClada) and ML optimizations (Lewis, 2001; as implemented in Mesquite 1.0 by Maddison and Maddison, 2003). We call these combinations AsIs+Fitch and AsIs+ML (Table 4). For continuous parameters, the applicable algorithm is linear parsimony (Farris, 1970; as generalized by Swofford and Maddison, 1987). We refer to this latter combination as AsIs+LP and implement it using Mesquite. These algorithms allow for polymorphic terminals. The polymorphism parsimony method implemented in PHYLIP (Felsenstein, 1989) optimizes polymorphisms to internal nodes to avoid a second, independent acquisition of the derived state (Felsenstein, 1979), an assumption that is akin to Dollo parsimony and will not generally be appropriate for ecological characters. Squared-change parsimony (Maddison, 1991), ML, and the generalized least-squares (GLS) approaches implemented in programs such as Continuous (Pagel, 1999) and COMPARE (Martins, 2003), although designed for continuous data, are not applicable to AsIs coding of our continuously variable parameters (Altitude and Average Rainfall). This is because the algorithms used accept only single values for terminals, which means that intraspecific variation must be converted into single-valued descriptive statistics (e.g., a mean) for each terminal species (e.g., a species' range in altitude cannot be coded). Reducing intraspecific variation to means is a commonly followed procedure for continuous characters in general (e.g., Archie, 1985; Thiele, 1993; Rae, 1998; Rice et al., 2003), but doing so for ancestor reconstruction assumes that the mean encompasses more information about the historical ecology of the species than any other value within the variation range of the species. Furthermore, although means may be easily applied to morphometric data, such a statistical treatment of our altitude and rainfall data (taken primarily from herbarium labels, distributional data, and collectors' notes) may not be appropriate. COMPARE does enable an accounting for intraspecific variability in the form of standard error, but this is yet an additional descriptive statistic.
|
Presence coding
Presence coding was first developed by Mickevich and coworkers for the cladistic analysis of allozyme data (as "independent-allele" coding; Mickevich and Johnson, 1976; Mickevich and Mitter, 1981). It is applicable to coding polymorphisms in unordered, discrete parameters and is accomplished by coding each state as a separate presence-absence character. Consequently, each parameter is represented by as many characters as it has states. Although presence coding is superficially similar to binary coding, it differs in that an ordinary binary character (e.g., FLOWER COLOR: blue/red) is transformed into two presence characters (BLUE FLOWER COLOR: absent/present; RED FLOWER COLOR:absent/present). By doing so, each state can be optimized to a node regardless of the number of other states of the same character also assigned to that node.
Fitch, ML, and DIVA (Ronquist, 1996) optimizations are applicable to presence-coded parameters, and each was employed here for demonstrative purposes (Presence+Fitch, Presence+ML, and Presence+DIVA, respectively). Although DIVA was designed for historical biogeographical analysis, it is readily adapted for the optimization of discrete, unordered ecological parameters by substituting ecological parameter states for "areas" (Linder and Hardy, in press).
Polymorphism coding
Maddison and Maddison (1987, as cited in Maddison and Maddison, 1992) proposed a stepmatrix approach to allow for both extant and ancestral polymorphisms and we applied it here to our unordered discrete parameters (Table 4). The Maddisons' method is to code the possible polymorphisms as separate states along with monomorphic states in a single multistate character. A stepmatrix is then employed to assign the desired transition costs (steps) between the states. We assigned equal costs for gains and losses of states (Fig. 2) and then determined the most parsimonious solutions using Sankoff optimization (Sankoff and Rousseau, 1975), as implemented in Mesquite (Polymorphism+Sankoff). However, because the maximum allowable number of states is 56 (Mesquite 1.0), the parameter Bedrock with its eight states, could not be optimized algorithmically because the number of possible combinations of these eight component states exceeds 56. Thus, the Polymorphism+Sankoff procedure for Bedrock was carried out manually (as per Felsenstein, 2004:13–15, 67–69).
|
Generalized frequency coding
Developed by Smith and Gutberlet (2001) for phylogeny reconstruction, Generalized frequency coding, or GFC, allows for the coding and optimization of polymorphic binary or multistate characters. It is also applicable to continuously variable characters, where the additive axis is subdivided into discrete increment states. The precision with which this is done is limited only by the precision with which character data were measured (here, Average Rainfall into seven 200-mm/year increments). The Altitude data presented in Table 1, although measured approximately to the nearest 50-m increment, was recoded into ten 200-m increments to enable comparison with the next coding method described (range coding), which is limited to just 10 increment states.
GFC is similar to presence coding in that these states are then converted into the same number of subcharacters that report the occurrence of a particular state in the species. GFC is different than all of the other methods investigated here in that a state's occurrence is weighted by its frequency. In this respect, GFC captures the essential component of the tree-searching program FreqPars (Swofford and Berlocher, 1987) and the "frequency bins"–coding approach of Wiens (1995). However, we treat all states within a species here at equal frequencies because of the probable collector bias inherent with systematics specimen collections. Smith and Gutberlet thoroughly described GFC and its use of "cumulative" and "noncumulative" frequencies for ordered (and continuous) and unordered characters, respectively. GFC codes were obtained automatically using the program FastMorphologyGFC ver. 1.0 (Chang and Smith, 2001). Optimizations of each subcharacter were carried out under Wagner optimization in Mesquite (GFC+Wagner).
Range coding
We propose a modification of polymorphism coding to preserve both range and order information in continuously variable ecological parameters. As for GFC, gradient parameters such as altitude and rainfall are subdivided into increment states. Thus, the continuous characters are transformed into ordered multistate characters. However, instead of coding all possible polymorphisms (i.e., all possible combinations of monomorphic states) as states themselves, we only consider contiguous combinations of increment states (i.e., continuous ranges) (Fig. 3). A stepmatrix is then employed to assign a cost to each possible transformation and Sankoff optimization is applied (Range+Sankoff). For comparative purposes, two types of stepmatrices were employed, one in which no additivity or order is assumed (e.g., no cost distinctions are made between transitions of 1000 or 100 m in altitude), and the other in which additivity is assumed (Appendix 2, available at the Society of Systematic Biologists Web site, http://systematicbiology.org).
|
While this article was in press, what we've termed "MaxMin coding" was independently applied in the inference of ancestral temperature and rainfall regimes for a clade of Ecuadorian dendrobadid frogs (Graham et al. 2004). MaxMin Coding codes ranges for a terminal by two characters: one giving the maximum value and the second the minimum value for that parameter in the species. The range for the terminal is then the set of values between the maximum and the minimum. These two characters are then optimized independently to the internal nodes using linear parsimony (MaxMin+LP) or squared-change parsimony (MaxMin+SqCP), the latter may be weighted by branch lengths (Maddison, 1991; as implemented in Mesquite). Weighted SqCP produces results identical to the simplest models implemented in Schluter et al.'s (1997) ML method, and the GLS-linear model of COMPARE (Martins, 2003). The results of these independent optimizations of maxima and minima are then combined, and the ancestral ranges inferred to lie between them. In order to allow a direct comparison with AsIs+LP, GFC+Wagner, and Range+Sankoff results, MaxMin Coding was conducted with the same coarse 200 m/200 mm scale employed with these other methods.
| RESULTS AND DISCUSSION |
|---|
|
|
|---|
Thamnochortus Phylogeny and Age Estimation
The parsimony and ML analyses of the combined data matrix converged toward a single topology (Fig. 1B). Although most nodes are strongly supported (as inferred from the bootstrap values and parsimony branch lengths), there are others for which support is very low; thus, some of the relationships depicted in Fig. 1B may not be robust to additional data. The branches circled in Fig. 1B received very short or zero length in the penalized likelihood analysis.
The monophyly of Thamnochortus received strong bootstrap support: the genus was recovered in 99% of the replicates. The basal split in the clade is that between a lineage represented by Th. levynsiae plus Th. pulcher and the clade including all other Thamnochortus species, each strongly supported. Within this larger clade, the lineage comprising Th. gracilis and Th. nutans is sister to the remaining Thamnochortus species.
Within this remaining group, the clade including both Th. muirii and Th. spicigerus received strong bootstrap support (94%) and is the most distinctive of all clades in Thamnochortus in terms of habitat, occurring at low elevations and primarily on limestone soils or the alkaline coastal sands of the western Cape (Fig. 6E). The one exception to this is Th. karooica, which is found on sandstone in the dry fynbos in the low to middle elevations of the Little Karoo Mountains and along the inland side of the Langeberg. The sister-group relationship of this alkaline-soil clade to the remaining Thamnochortus received no bootstrap support and was resolved only by a single (substitution) character.
Of the remaining four major clades, the one including both Th. cinereus and Th. fruticosus is sister, albeit with no bootstrap support, to a clade that includes the other three. Within the latter, the clade comprising Th. guthrieae, Th. erectus, and Th. insignis is sister to the lineage comprising the remaining two major clades. This latter clade comprises the sister groups (Th. arenarius(Th. stokoei(Th. pellucidus(Th. dumosus, Th. lucens)))) and (Th. schlechteri(Th. platypteris((Th. bachmannii, Th. punctatus)(Th. obtusus, Th. sporadicus)))). The basal position of Th. arenarius in the former, although recovered in all most parsimonious trees, received only 27% bootstrap support.
Polymorphism in Ancestral Habitat and Ecology Reconstruction
Discrete parameters
The combination of multistate coding and Fitch or ML optimizations (AsIs+Fitch and AsIs+ML) precluded the reconstruction of polymorphic ancestors (Fig. 4). This is because these procedures treat states of a multistate character as mutually exclusive. The occasional ambiguity at nodes (e.g., differences between ACCTRAN and DELTRAN with Fitch, or equally likely state assignments with ML) cannot be interpreted as polymorphism. When applied to ancestral ecology reconstruction, the unavoidable monomorphy restriction of AsIs coding is analogous to an inappropriate assumption of ecological uniformity within species. Of the species included in this study, 30.8% possess more than one state for one or more unordered ecological parameters, leaving 11.5% of the AsIs-coded data matrix cells polymorphic for these parameters (Table 3). Such intraspecific ecological variability may be due to ecological plasticity, differentially adapted conspecifics, or different community-level species interactions in different parts of a species' range. Regardless of the biological basis for such intraspecific ecological variability, optimization procedures that account for this are preferable.
|
In contrast to AsIs coding, Presence Coding and GFC (as applied to unordered characters) allow ancestral polymorphism (uncircled reconstructions below the branches in Fig. 4) because the possession of each ecological state is treated independently of all other states. In fact, these two methods (Presence+Fitch and GFC+Wagner) produce identical results when state occurrences are not weighted by their frequencies. Ironically, however, the strength of these methods for allowing polymorphic reconstructions (i.e., the independent optimization of each state) is also their weakness. Although usually providing decisive optimizations (e.g., all of Fig. 4B, Fig. 5A, B, D, Fig. 6B; most of Fig. 4A, Fig. 6A, E, F), Presence+Fitch, Presence+ML, and GFC+Wagner occasionally find nodes with no state (i.e., the "absence" state for each presence character; nodes with a "?" in Fig. 4A, Fig. 6A, E, F), an artifact in clades where shared occurrences of particular ecological states are few or lacking. Using presence coding, "no-state" assignments were first recognized by Mickevich and Mitter (1981, for allozyme data) and are artifacts of the assumption of character independence, an assumption that quite correctly lies at the core of normal cladistic theory and practice. However, presence coding and GFC for nonadditive parameters requires a special conditional relationship, which is lacking, between the separate subcharacters that comprise an ecological parameter. This conditional dependence is that at least one "presence" must be optimized for each parameter.
|
|
This conditional dependence among presence characters is automatically implemented in the program DIVA (Ronquist, 1996); thus, "no-state" assignments are avoided (Fig. 5C, D, Fig. 6C, G). Although also utilizing Presence coding and the parsimony criterion, DIVA optimization has certain explicit assumptions concerning potential character transformations that are different than those implicit with Presence+Fitch. If applied to ecological parameters, Fitch unconditionally counts both gains and losses of parameter states and attempts to minimize these, with the consequent implicit assumption that new states in just one descendant lineage are acquired after the speciation event leading to it (Fig. 7A). In contrast, DIVA utilizes a three-dimensional stepmatrix to assign no cost to losses associated with the partitioning of ancestral states into sister lineages following a speciation event, an explicit assumption that biases toward the acquisition of descendent traits to before cladogenic events (Fig. 7B). These different underlying assumptions may translate into qualitatively different ancestral reconstructions at certain nodes (Fig. 5, Fig. 6, Fig. 7), indicating that DIVA must not be viewed merely as an alternative to the use of Presence+Fitch, Presence+ML, or GFC+Wagner in order to avoid "no-state" assignments.
|
Polymorphism coding also avoids the "no-state" assignments that arise with presence coding and GFC (Fig. 5, Fig. 6) while also allowing polymorphic ancestors. Through the application of Sankoff optimization and the assignment of equal and unconditional costs for gains and losses of parameter states (Fig. 2), Polymorphism+Sankoff is similar in its cost scheme to Presence+Fitch and therefore yields similar, if not identical, optimizations (Fig. 5A, D, Fig. 6D, H). The differences between Polymorphism+Sankoff and Presence+Fitch arise not only where the latter returns "no-state" assignments, but also at nodes for which Presence+Fitch yields ambiguity with respect to two or more states (compare Fig. 6A and E with D and H, respectively). Given that the cost scheme is, in principal, identical for both Presence+Fitch and Polymorphism+Sankoff, all differences in optimization that arise between the two methods should be attributable to the latter's more restricted model of possible change (i.e., parsimonious histories that proceed via nodes with no state are not allowed). Unfortunately, the utility of Polymorphism+Sankoff is limited by the number of states allowed by programs such as Mesquite 1.0 (56), PAUP*4.0 (32), and MacClade 4.0 (26). Because polymorphism coding entails the coding of every possible combination of monomorphic states as separate polymorphic states, parameters with more than five monomorphic states exceed the maximum of 56 states allowed for a single character in available phylogenetics packages and must be optimized by the tedious manual implementation of Sankoff optimization (e.g., Bedrock; Fig. 6H). However, we suspect that more than five states will not be needed for most discrete ecological parameters.
Continuous parameters
The direct (AsIs) coding and optimization of Altitude and Average Rainfall was possible only with linear parsimony, as the algorithms implemented in SqCP, ML, and GLS approaches do not accept ranges (e.g., 200 to 600 m) to be coded for species. LP's treatment of ranges in terminals (e.g., AsIs+LP; Fig. 8A, D), however, is of limited value because the algorithm treats these ranges as ambiguous polymorphisms; i.e., sets of possible values, any of which may be chosen by the algorithm to achieve a most parsimonious solution of character evolution (Maddison and Slatkin, 1990). Furthermore, LP is designed to yield single values at ancestral nodes (e.g., Fig. 8A, D) and these provide no hypothesis as to how widespread ancestral species may have been along a particular habitat gradient. What look like "ranges" recovered at some nodes by LP (Fig. 8A, D) do not represent estimates of ancestral ranges, but ambiguity (Maddison and Slatkin, 1990; Losos, 1999). Accordingly, we consider the use of AsIs coding with any of the algorithms mentioned above to be inappropriate.
|
Our modification of polymorphism coding, which we term range coding, provides one means of optimizing ranges to internal nodes (Fig. 3). The stepmatrix utilized here for Sankoff optimization of range-coded parameters (Appendix 2) employs a uniform cost scheme (e.g., all incremental increases or decreases of 200 m receive the same cost, no matter where along the gradient they occur). It is also possible to differentially weight transitions along different parts of the gradient (e.g., transitions from 1800 to 2000 m may be assigned a higher cost than transitions at lower altitudes). The implementation of stepmatrix costs that confer order and, consequently, additivity along the parameter axis is significant not only because it captures a key attribute of the algorithms designed for continuous data (e.g., LP, SqCP), but also because not doing so may impact the reconstructions inferred. The differences between unordered and ordered range coding are illustrated by the optimizations of Average Rainfall and Altitude (Fig. 9). Although similar, ordered range coding generally reconstructs ancestors with larger ranges than with unordered coding, with the unordered results often being subsets of the ordered reconstructions. However, the unordered and ordered reconstructions for Altitude in the clade comprising Th. papyraceus,Th. acuminatus, and Th. fruticosus (compare Fig. 9D, E) reveal that the reconstructions obtained under both methods are not always so similar. These examples demonstrate that the assumption of additivity in relevant habitat parameters may have important consequences in the inference of ancestral habitats.
|
Although the Range+Sankoff procedure provides a viable and, in terms of possible weighting schemes, flexible alternative to the assumption of "monomorphy" that is implemented through the application of standard algorithms such as SqCP and LP to continuous parameters, there remains the potential for loss of information when transforming continuously variable parameters into discrete characters. For example, our application of Range Coding to the ca. 2000-m altitudinal gradient observed for Thamnochortus was restricted to the recognition of just ten increment states of 200 m each. The increment state "a" (0–200 m) in Fig. 9, for example, did not distinguish between observed ranges of 100–200 m, 50–200 m, or 0–200 m (Table 3). It is theoretically possible to improve resolution by subdividing continuously variable parameters more finely. However, Mesquite 1.0's 56 allowable states sets the upper limit to just 10 monomorphic increment states in range coding, which are combined to produce an additional 45 polymorphic (range) states, for a total of 55 states (Appendix 2).
Generalized frequency coding is superior to range coding in at least one respect: GFC has no limit on how finely parameters such as Altitude and Rainfall are subdivided for coding. GFC's prerequisite subdivision of continuous axes into multiple discrete increment states is limited only by the precision with which data are measured. Furthermore, GFC is readily implemented in the program FastMorphologyGFC (Chang and Smith, 2001), thereby greatly facilitating the implementation of an otherwise laborious coding procedure. Reconstructions implied by GFC+Wagner are, however, tedious to interpret, because reconstructions must be inferred by summing across multiple subcharacters on a node-by-node basis. This is because GFC was developed as a way only to extract information from polymorphic (and continuously variable) characters for cladistic analysis: we have merely co-opted it here for the inference of ancestral states. The labor involved with ancestor inference using GFC+Wagner increases drastically with increasing subdivision of continuously variable character axes. Although we have not done so here on account of a perceived bias in the habitat distributions implied by systematics specimen collections, the capacity of GFC to accommodate species' frequency distributions may be seen by some investigators as an additional desirable attribute of the method.
MaxMin coding has several advantages over both range coding and GFC. First, it is not necessary to transform continuously variable data into discrete increment states. As demonstrated above, the number of increment states (and therefore precision) possible using range coding is limited. Although the precision of axis subdivision using GFC is not limiting, inference of ancestral ranges becomes increasingly tedious with increasing subdivision for the reason described above. MaxMin coding requires that inferences be made simply by summing across just two subcharacters (the minimum and the maximum). Our results suggest that the added tedium associated with GFC comes without added benefit: GFC and MaxMin produced identical results under parsimony when ranges, rather than frequency distributions, were coded (Fig. 9C, F).
A second advantage that MaxMin has over both range coding and GFC is that it allows for the application of SqCP or more powerful GLS or ML approaches (Pagel, 1999; Martins and Hansen, 1997; Martins, 1999) designed to account for branch lengths during optimization. These algorithms/models implement variations upon an implicit (SqCP) or explicit (ML, GLS) assumption that continuous character evolution can be modeled as a stochastic process (Maddison, 1991; Martins, 1999). The probability of change along any branch is proportional to its length. One possible objection to modeling ecological change as a stochastic process is the possibility that many ecological attributes are under selective pressures and evolve under the influence of other forces such as divergent or stabilizing selection (Losos, 1999; Martins et al., 2002 and references cited therein). Indeed, LP's tendency of reconstructing most branches with no change and relatively few branches with large amounts of change has been likened to an implicit, and perhaps appropriate, model of stabilizing selection with occasional adaptive shifts (Losos, 1999). As appealing as this assumption implicit with LP might be, recent theoretical and algorithmic progress has introduced the potential for modeling the influence of selection (e.g., GLS-exponential and the Phylogenetic Mixed Model of Martins, 2003, and Housworth et al., 2004). The efficacy with which such forces can be modeled and used for accurate reconstructions is not yet clear (compare Martins, 1999 with Martins et al., 2002). However, the increasing availability of a variety of parameter-driven, statistically based evolutionary models points toward the powerful prospect of objectively testing and selecting a model appropriate for one's study group.
Timing of Ecological Diversification
A neglected aspect of ancestor reconstruction is the impact that alternative algorithms have on the rate and timing of character evolution inferred (Linder and Hardy, in press). The rate and timing of eco-diversification relative to lineage (species) diversification are of particular interest in the pursuit of identifying the forces driving speciation. If eco-diversification can be shown to precede lineage diversification in particular clades, then there is support for the hypothesis that ecological factors may be driving speciation. The converse would be that eco-diversification follows and, by implication, does not drive lineage diversification. Our results reveal that, for the same data and phylogeny, either scenario can be demonstrated, depending on the optimization method employed, and we see no a priori means of choosing between these biases. Using Presence+DIVA on unordered parameters, the accumulation of ecological diversity precedes that of lineage diversity (Fig. 10A). In contrast, all other combinations produce ecological diversity curves that lag behind lineage diversity (Fig. 10A, B). These biases are especially evident when trend lines are fitted to the data plotted in Fig. 10 (not shown). Among the parsimony-based optimizations, this dichotomy can be understood as a consequence of different transition cost schemes. Schemes that assign equal and unconditional costs to gains and losses of states (Fitch and, as implemented here, Sankoff) will tend to place the acquisition of new ecological states on branches following speciation events. In contrast, DIVA employs a three-dimensional stepmatrix to assign no cost to "losses" associated with the partitioning of ancestral states into sister lineages following a speciation event. Thus, DIVA biases toward inferences of ecological vicariance or ecological speciation (Ehrlich and Raven, 1969; Andersson, 1990; Schluter, 2000; Levin, 2000, 2003), where the acquisition of states exclusive to each of two descendant sister lineages is accelerated to the ancestral lineage, followed by the splitting of that lineage along these ecological boundaries reflected in the descendants (Fig. 7). Our results also demonstrate that ML behaves like other non-DIVA parsimony algorithms with respect to the relative rates of ecological and lineage diversification. In fact, within coding types (e.g., Presence, MaxMin), optimizations that account for branch lengths (ML and, for continuous parameters, weighted SqCP) produce the slowest ecological diversity curves, all of which lay below that of lineage diversity (Fig. 10). Furthermore, the delay in eco-diversification inferred by procedures that use branch lengths may be exacerbated by branch length estimates based on saturated genetic loci when, despite employing reasonable models of sequence evolution, ML may fail to fully correct the bias toward the artifactual inference of shorter internal nodes and longer terminal nodes (Graybeal, 1994; Cunningham et al., 1998; Cunningham, 1999). As ML and weighted SqCP procedures assign higher probabilities of change to longer branches, optimizations on trees with artifactually long terminal branches will be biased toward the inference of even slower rates of eco-diversification relative to lineage diversification (cf. Cunningham et al., 1998; Cunningham, 1999). Thus, it is clear from an analysis of each of these algorithms' optimality criteria and transition cost schemes that the ecological-lineage diversity curves observed here probably depict general properties of the algorithms themselves, and that DIVA seems uniquely capable of resolving putative instances of ecological speciation.
|
Unfortunately, it is not yet possible to implement DIVA-like optimization with ordered and continuous characters. The optimization procedures employed here (LP, Range+Sankoff, GFC+Wagner, MaxMin+LP, and MaxMin+SqCP) all produce ecological diversity curves that lag behind that of lineage diversity (Fig. 10B). Although parameters such as Altitude and Average Rainfall could be presence-coded and optimized using DIVA, the additive and continuous nature of possible change along these gradient axes would then be lost. It may be possible, however, to re-code two of the three dimensions of the stepmatrix employed by DIVA, to reflect the additivity modeled in the stepmatrix utilized here for range coding (Appendix 2). DIVA assumptions might then be implemented. Then, as for unordered and discrete parameters, it would be possible to explore for additive and continuous parameters the impact that biases toward accelerated and delayed character evolution have on ancestor reconstructions.
| CONCLUSIONS |
|---|
|
|
|---|
Methodological Problems and Solutions
Methods that account for ancestral polymorphism and continuously variable or additive range characters are prerequisites for meaningful reconstructions of the macroecological history of clades. Using a molecular phylogenetic hypothesis for Thamnochortus, the advantages and disadvantages of available methods were demonstrated. For discrete, unordered ecological parameters, polymorphism coding is superior to AsIs coding, GFC, and non-DIVA implementations of presence coding because it avoids illogical "no-state" assignments. Presence+DIVA also avoids these no-state assignments, yet DIVA is biased towards accelerated rates of eco-diversification and therefore inferences of ecology-driven speciation. All other parsimony and ML methods share the opposite bias. We see no objective means of choosing between these biases; however, investigators must be aware of them. Although both DIVA and Polymophism coding methods are optimal solutions for use with parsimony, they cannot be implemented with ML options available in existing software. For continuously variable parameters, the ranges along a given habitat gradient that individual species occupy cannot be optimized to internal nodes using standard procedures. Range coding provides a solution to this problem and is highly flexible in terms of possible transition-weighting schemes along a given habitat axis such as rainfall or altitude. Range coding does, however, require the discrete coding of continuously variable parameters and is therefore limited by the number of states allowed by phylogenetic software. GFC and MaxMin coding have no such limit and are therefore superior to range coding when a strict linear additivity is assumed along continuously variable parameter axes. When ranges were coded for Thamnochortus, GFC and MaxMin gave identical results; however, MaxMin is preferable to GFC in such cases because it is easier to implement and is applicable to both parsimony and ML or GLS approaches. For continuous parameters, there is nothing analogous to DIVA's accelerated eco-diversification bias with discrete, unordered parameters.
Ecological Diversification in Thamnochortus
Although inferences differed in details, the methods converged upon a single basal optimization and are congruent regarding the major trends of macroecological diversification. Accordingly, the ancestral lineage may have been a post-fire resprouting species distributed on rocky and well-drained, sandstone-derived soils at lower-middle elevations in regions characterized by moderate levels of yearly (particularly winter) rainfall and negligible contributions of moisture from clouds during the summer months (Fig. 5, Fig. 6, Fig. 7, Fig. 8 and Fig. 9). This species would have been distributed in habitats much like those existing throughout much of the southwestern Cape mountains today.
With respect to fire survival, post-fire reseeders are derived within Thamnochortus. DIVA infers the same polarity to this transformation, but the transformations to reseeding proceeded via a polymorphic intermediate lineage on at least eight separate occasions, each invoking a hypothesis of vicariant speciation with respect to these two life-history traits. Presence+ML inferences are contrary to parsimony inferences in suggesting that all ancestral species were polymorphic for Fire Survival Mode and that loss of either of the two traits occurred recently and independently in several lineages. However, non-DIVA implementations of presence coding, such as Presence+ML, are known to be prone to illogical, artifactual reconstructions are therefore suboptimal.
With regards to moisture regime, all methods point to a shift from regions of moderate yearly rainfall (ca. 600–1000 mm/year) to regions of lower rainfall (200–600 mm/year) in the Th. muirii—Th. spicigerus lineage. Independently, the Th. platypteris–Th. bachmannii lineage experienced a shift to a lower and narrower rainfall regime (ca. 400–700 mm/year). The occurrences of species in winter rainfall regions that receive substantial amounts of summer moisture from clouds off the Indian Ocean are recently and separately derived in several lineages. Expectedly, DIVA suggests that these summer (southeast) cloud habitat shifts took place via intermediate lineages that existed in both habitats. All parsimony methods suggest that lineages found on deep, rockless, and summer-dry soils (Appendix 1) are derived from lineages occurring on well-drained, rocky soils.
Finally, all methods also point to occurrences at lower-middle elevations (< 1000 m) on sandstone-derived soils as plesiomorphic. Subsequent shifts to lower, including coastal, elevations have occurred independently in both the Th. guthrieae–Th. erectus and Th. muirii–Th. spicigerus lineages. Both were associated with a shift from acidic sandstone-derived soils to limestone or otherwise alkaline soils.
| ACKNOWLEDGMENTS |
|---|
|
|
|---|
Support for this research was generously provided by the Swiss National Fund, the Swiss Academy of Natural Sciences, Georges and Antoine Claraz-Schenkung, and The National Geographic Society. We thank The Western Cape Nature Conservation Board for the permission to collect the plants, Frank Rutschmann for assistance with the penalized likelihood analysis, Philip Moline for providing unpublished primer sequences, and Philip Moline, Terry Trinder-Smith, and the staff at the Bolus Herbarium for facilitating fieldwork in the Cape. Chloé Galley, Timo van der Niet, Roderic Page, Todd Oakley, and two anonymous reviewers provided helpful comments on an earlier draft of this article. We also thank Emilia Martins for helpful correspondence regarding COMPARE.
|
| Notes |
|---|
|
|
|---|
2 Current Address: Biology Department and Herbarium, Millersville University, Millersville PA, 17551, USA
| REFERENCES |
|---|
|
|
|---|
-
Andersson L. The driving force: Species concepts and ecology. Taxon (1990) 39:375–382.[CrossRef][Web of Science]
Archie J. W. Methods for coding variable morphological features for numerical taxonomic analysis. Syst. Zool. (1985) 34:326–345.
Asmussen C. B., Chase M. W. Coding and noncoding plastid DNA in palm systematics. Am. J. Bot. (2001) 88:1103–1117.
Butler M. A., Losos J. B. Testing for unequal amounts of evolution in a continuous character on different branches of a phylogenetic tree using linear and squared-change parsimony: An example using Lesser Antillean Anolis lizards. Evolution (1997) 51:1623–1635.[CrossRef][Web of Science]
Campbell B. M. Montane plant environments in the Fynbos Biome. Bothalia (1983) 14:283–298.
Change V., Smith E. N. FastMorphologyGFC (2001) version 1.0. http://www3.uta.edu/faculty/ensmith.
Chase M. W., Albert V. A. A perspective on the contribution of plastid rbcL DNA sequences to angiosperm phylogenetics. In: Molecular systematics of plants II: DNA sequencing—Soltis D. E., Soltis P. S., Doyle J. J., eds. (1998) Boston, Massachusetts: Kluwer Academic Publishers. Pages 488–507.
Cowling R. M. Fire and its role in coexistence and speciation in Gondwanan shrublands. S. African J. Sci. (1987) 83:106–112.
Cowling R. M., Rundel P. W., Lamont B. B., Arroyo M. K., Arianoutsou M. Plant diversity in Mediterranean-climate regions. Trends Ecol. Evol. (1996) 11:362–366.[CrossRef]
Cunningham C. W. Some limitations of ancestral character-state reconstruction when testing evolutionary hypotheses. Syst. Biol. (1999) 48:665–674.
Cunningham C. W., Omland K. E., Oakley T. H. Reconstructing ancestral character states: A critical reappraisal. Trends Ecol. Evol. (1998) 13:361–366.[CrossRef]
Donoghue M. J., Ackerly D. D. Phylogenetic uncertainties and sensitivity analysis in comparative biology. Philos. Trans. R. Soc. Lond. B (1996) 351:1241–1249.
Ehrlich P. R., Raven P. H. Differentiation in populations: Gene flow seems to be less important in speciation than the neo-Darwinians thought. Science (1969) 165:1228–1232.
Eldenäs P. K., Linder H. P. Congruence and complementarity of morphological and trnL-trnF sequence data and the phylogeny of the African Restionaceae. Syst. Bot. (2000) 25:692–707.[CrossRef]
Farris J. S. Methods for computing Wagner trees. Syst. Zool. (1970) 19:83–92.[Abstract]
Felsenstein J. Alternative methods of phylogenetic inference and their interrelationship. Syst. Zool. (1979) 28:49–62.
Felsenstein J. Evolutionary trees from DNA sequences: A maximum likelihood approach. J. Mol. Evol. (1981) 17:368–376.[CrossRef][Web of Science][Medline]
Felsenstein J. Confidence limits on phylogenies: An approach using the bootstrap. Evolution (1985) 39:783–791.[CrossRef][Web of Science]
Felsenstein J. PHYLIP phylogeny inference package. Cladistics (1989) 5:164–166.
Felsenstein J. Inferring phylogenies (2004) Sunderland, Massachusetts: Sinauer Associates.
Fitch W. Toward defining the course of evolution: Minimum change for a specific tree topology. Syst. Zool. (1971) 20:406–416.[Abstract]
Frumhoff P. C., Reeve H. K. Using phylogenies to test hypotheses of adaptation: A critique of some current proposals. Evolution (1994) 48:172–180.[CrossRef][Web of Science]
Goldblatt P. An analysis of the flora of southern Africa: Its characteristics, relationships, and origins. Ann. MO. Bot. Gard. (1978) 65:369–436.[CrossRef]
Goldblatt P., Manning J. C. Cape plants: A conspectus of the Cape flora of South Africa (2000) St. Louis, Missouri, USA: National Botanical Institute of South Africa, Pretoria, South Africa, and MBG Press.
Goloboff P. NONA (1993) Buenos Aires, Argentina. version 1.6. Distributed by the author.
Graham C. H., Ron S. R., Santos J. C., Schneider C. J., Moritz C. Integrating phylogenetics and environmental niche models to explore speciation mechanisms in dendrobatid frogs. Evolution (2004) 58:1781–1793.[Web of Science][Medline]
Graybeal A. Evaluating the phylogenetic utility of genes: A search for genes informative about deep divergences among vertebrates. Syst. Biol. (1994) 43:174–193.
Hilu K. W., Liang H. The matK gene: Sequence variation and application in plant systematics. Am. J. Bot. (1997) 84:830–839.[Abstract]
Housworth E. A., Martins E. P., Lynch M. The phylogenetic mixed model. Am. Nat. (2004) 163:84–96.[CrossRef][Web of Science][Medline]
Huelsenbeck J. P., Nielsen R., Bollback J. P. Stochastic mapping of morphological characters. Syst. Biol. (2003) 53:131–158.
Kocyan A., Y.-L. Qiu, Endress P. K., Conti E. A phylogenetic analysis of Apostasioideae (Orchidaceae) based on ITS, trnL-F and matK sequences. Plant Syst. Evol. (2004) 247:203–213.
Lambrechts J. J. N. Geology, geomorphology and soils. In: Fynbos ecology: A preliminary synthesis—Day J., Siegfried W. R., Louw G. N., Jarman M. L., eds. (1979) Pretoria, South Africa: Council for Scientific and Industrial Research. Pages 16–26.
Levin D. A. The origin, expansion, and demise of plant species (2000) New York: Oxford University Press.
Levin D. A. Ecological speciation: Lessons from invasive species. Syst. Bot. (2003) 28:643–650.
Lewis P. O. A likelihood approach to estimating phylogeny from discrete morphological character data. Syst. Biol. (2001) 50:913–925.
Linder H. P. A phylogenetic classification of the genera of the African Restionaceae. Bothalia (1984) 15:11–76.
Linder H. P. Gene flow, speciation, and species diversity patterns in a species-rich area: The Cape flora. In: Species and speciation—Vrba E. S., ed. (1985) South Africa. Pages 53–57Transvaal. Museum Monograph 4.
Linder H. P. A review of the southern African Restionaceae. Contr. Bolus Herb. (1991) 13:209–264.
Linder H. P. The African Restionaceae: An IntKey identification and description system. (2002) version 2. Contr. Bol. Herb. 20. Available at http://www.systbot.unizh.ch/datenbanken/restionaceae/index.htm.
Linder H. P. The radiation of the Cape flora, southern Africa. Biol. Rev. (2003) 78:597–638.[Medline]
Linder H. P., Briggs B. G., Johnson L. A. S. Restionaceae: A morphological phylogeny. In: Monocots: Systematics and evolution—Wilson K. L., Morrison D. A., eds. (2000) Collingwood, Victoria, Australia: CSIRO Publishing. Pages 653–660.
Linder H. P., Hardy C. R. Speciation in the Cape flora: A macroevolutionary and macroecological perspective. In: Plant species-level systematics: New perspectives on pattern and process—Bakker F. T., Chatrou L. W., Gravendeel B., Pelser P. B., eds. In press. Koeltz, Königstein.
Losos J. B. Uncertainty in the reconstruction of ancestral character states and limitations on the use of phylogenetic comparative methods. Anim. Behav. (1999) 58:1319–1324.[CrossRef][Web of Science][Medline]
Maddison W. P. Squared-change parsimony reconstructions of ancestral states for continuous-valued characters on a phylogenetic tree. Syst. Zool. (1991) 40:304–314.[Abstract]
Maddison W. P. Calculating the probability distributions of ancestral states reconstructed by parsimony on phylogenetic trees. Syst. Biol. (1995) 44:474–481.
Maddison W. P., Maddison D. R. MacClade (1987) Cambridge, Massachusetts. version 2.1. Distributed by the authors.
Maddison W. P., Maddison D. R. Analysis of phylogeny and character evolution (1992) Sunderland, Massachusetts: Sinauer Associates. MacClade version 3.
Maddison W. P., Maddison D. R. Mesquite: A modular system for evolutionary analysis. (2003) version 1.0. Available at: http://mesquiteproject.org.
Maddison W. P., Slatkin M. Parsimony reconstructions of ancestral states do not depend on the relative distances between linearly-ordered character states. Syst. Zool. (1990) 39:175–178.
Manen J., Natali A., Ehrendorfer F. Phylogeny of Rubiaceae-Rubieae inferred from the sequence of a cpDNA intergene region. Pl. Syst. Evol. (1994) 190:195–211.[CrossRef]
Marloth R. Results of experiments on Table Mountain for ascertaining the amount of moisture deposited from the southeast clouds. Trans. S. Afr. Phil. Soc. (1904) 14:403–408.
Martins E. P. Estimation of ancestral states of continuous characters: A computer simulation study. Syst. Biol. (1999) 48:642–650.
Martins E. P. Computer programs for the statistical analysis of comparative data (2003) Bloomington, Indiana: Department of Biology, Indiana University. COMPARE, version 4.5. Distributed by the author at http://compare.bio.indiana.edu/.
Martins E. P., Diniz-Filho J. A., Housworth E. A. Adaptation and the phylogenetic comparative method. Evolution (2002) 21:317–340.
Martins E. P., Hansen T. F. Phylogenies and the comparative method: A general approach to incorporating phylogenetic information into the analysis of interspecific data. Am. Nat. (1997) 149:646–667.[CrossRef][Web of Science]
Mickevich M. F., Johnson M. F. Congruence between morphological and allozyme data in evolutionary inference and character evolution. Syst. Zool. (1976) 25:260–270.
Mickevich M. F., Mitter C. Treating polymorphic characters in systematics: A phylogenetic treatment of electrophoretic data. Funk V. A., Brooks D. R., eds. Advances in cladistics: Proceedings of the first meeting of the Willi Hennig Society: New York. New York Botanical Garden Press. Pages 45–58.
Nixon K. C. WinClada (2002) Ithaca, New York. ver. 1.00.08. Published by the author.
Nixon K. C., Davis J. I. Polymorphic taxa, missing values and cladistic analysis. Cladistics (1991) 7:233–241.[CrossRef][Web of Science]
Nixon K. C., Wheeler Q. D. An amplification of the phylogenetic species concept. Cladistics (1990) 6:211–223.[CrossRef][Web of Science]
Omland K. E. The assumptions and challenges of ancestral state reconstructions. Syst. Biol. (1999) 48:604–611.
Page R. D. M. TREEVIEW: An application to display phylogenetic trees on personal computers. Comp. Appl. Biosci. (1996) 12:357–358.[Medline]
Pagel M. Inferring the historical patterns of biological evolution. Nature (1999) 401:877–884.[CrossRef]
Posada D., Crandall K. A. MODELTEST: Testing the model of DNA substitution. Bioinformatics (1998) 14:917–818.
Rae T. The logical basis for the use of continuous characters in phylogenetic systematics. Cladistics (1998) 14:221–228.[CrossRef][Web of Science]
Rice N. H., Martínez-Meyer E., Peterson A. T. Ecological niche differentiation in the Aphelocoma jays: A phylogenetic perspective. Biol. J. Linn. Soc. (2003) 80:369–383.[CrossRef][Web of Science]
Ronquist F. DIVA (1996) version 1.1. Computer program and manual available by anonymous FTP from Uppsala University ftp.uu.se or ftp.systbot.uu.se.
Sanderson M. J. Estimating absolute rates of molecular evolution and divergence times: A penalized likelihood approach. Mol. Biol. Evol. (2002a) 19:101–109.
Sanderson M. J. Computer program distributed by the author. (2002b) r8s, version 1.5http://ginger.ucdavis.edu/r8s.
Sankoff D., Rousseau P. Locating the vertices of a Steiner tree in arbitrary space. Math. Program. (1975) 9:240–246.[CrossRef]
Schluter D. T. The ecology of adaptive radiation (2000) New York: Oxford University Press.
Schluter D., Price T., Mooers A. Ø, Ludwig D. Likelihood of ancestor states in adaptive radiation. Evolution (1997) 51:1699–1711.[CrossRef][Web of Science]
Schutte A. L., Vlok J. H. J., Van Wyk B. W. Fire survival strategy: A character of taxonomic, ecological and evolutionary importance in fynbos legumes. Pl. Syst. Evol. (1995) 195:243–259.[CrossRef]
Sharkey M. J. Transition confidence and modified mean values: Confidence measures for hypotheses of character state transition between nodes and ancestral state optimizations. Cladistics (1999) 15:113–120.[CrossRef][Web of Science]
Simmons M. P., Ochoterena H. Gaps as characters in sequence-based phylogenetic analyses. Syst. Biol. (2000) 49:369–381.
Smith E. N., Gutberlet R. L. Generalized frequency coding: A method of preparing polymorphic multistate characters for phylogenetic analysis. Syst. Biol. (2001) 50:156–169.
Sterelny K. Species as ecological mosaics. In: Species: New interdisciplinary essays—Wilson R. A., ed. (1999) Cambridge, Massachusetts: Massachusetts Institute of Technology. Pages 119–138.
Swofford D. L. PAUP*. (2002) Sunderland, MA, USA: Sinauer Associates. version 4.0.
Swofford D. L., Berlocher S. H. Inferring evolutionary trees from gene frequency data under the principle of maximum parsimony. Syst. Zool. (1987) 36:293–325.
Swofford D. L., Maddison W. P. Reconstructing ancestral character states under Wagner parsimony. Math. Biosci. (1987) 87:199–229.[CrossRef][Web of Science]
Swofford D. L., Maddison W. P. Parsimony, character-state reconstructions, and evolutionary inferences. In: Systematics, historical ecology, and North American freshwater fishes—Mayden R. L., ed. (1992) Stanford, California: Stanford University Press. Pages 187–223.
Taberlet P., Gielly L., Pautou G., Bouvet J. Universal primers for amplification of three non-coding regions of chloroplast DNA. Pl. Mol. Bio. (1991) 17:1105–1109.[CrossRef]
Thiele K. The holy grail of the perfect character: The cladistic treatment of morphometric data. Cladistics (1993) 9:275–304.[CrossRef][Web of Science]
Thompson J. D., Gibson T. J., Plewniak F., Jeanmougin F., Higgins D. G. The CLUSTAL_X windows interface: Flexible strategies for multiple sequence alignment aided by quality analysis tools. Nucleic Acids Res. (1997) 25:4876–4882.
van Wilgen B. W. Fire regimes in the Fynbos biome. In: Disturbance and dynamics of fynbos biome communities—Cowling R. M., Le Maitre C. D., McKenzie B., Prys-Jones R. P., Van Wilgen B. W., eds. (1987) Pretoria, South Africa: South African National Scientific Programmes. Pages 6–14. report number 135.
Wiens J. J. Polymorphic characters in phylogenetic systematics. Syst. Biol. (1995) 47:482–500.
Young N. D., Healy J. GapCoder: A computer program for including indels in phylogenetic analysis. (2001) Available from the authors at http://www.home.cr.duq.edu/
youngnd/GapCoder.
This article has been cited by other articles:
![]() |
J. R. Clark, R. H. Ree, M. E. Alfaro, M. G. King, W. L. Wagner, and E. H. Roalson A Comparative Study in Ancestral Range Reconstruction Methods: Retracing the Uncertain Histories of Insular Lineages Syst Biol, October 1, 2008; 57(5): 693 - 707. [Abstract] [Full Text] [PDF] |
||||
![]() |
C. Galley, B. Bytebier, D. U Bellstedt, and H Peter Linder The Cape element in the Afrotemperate flora: from Cape to Cairo? Proc R Soc B, February 22, 2007; 274(1609): 535 - 543. [Abstract] [Full Text] [PDF] |
||||
![]() |
C. Yesson and A. Culham Phyloclimatic Modeling: Combining Phylogenetics and Bioclimatic Modeling Syst Biol, October 1, 2006; 55(5): 785 - 802. [Abstract] [Full Text] [PDF] |
||||
![]() |
H. P. Linder, T. Dlamini, J. Henning, and G. A. Verboom The evolutionary history of Melianthus (Melianthaceae) Am. J. Botany, July 1, 2006; 93(7): 1052 - 1064. [Abstract] [Full Text] [PDF] |
||||
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||

50% below. The five circled nodes were subtended by branches of effectively zero length in the penalized likelihood tree; they are resolved here with short branches for illustrative purposes, in accordance with the topology of this most parsimonious tree.











