Molecular Phylogenetics and Evolution 36 (2005) 536–545 www.elsevier.com/locate/ympev
Resolving tylenchid evolutionary relationships through multiple gene analysis derived from EST data Elizabeth H. Scholl, David McK. Bird ¤ Center for the Biology of Nematode Parasitism, North Carolina State University, Raleigh, NC 27695, USA Bioinformatics Research Center, North Carolina State University, Raleigh, NC 27695, USA Received 20 July 2004; revised 11 March 2005 Available online 3 May 2005
Abstract Sequence-based phylogenetic analyses typically are based on a small number of character sets and report gene trees which may not reXect the true species tree. We employed an EST mining strategy to suppress such incongruencies, and recovered the most robust phylogeny for Wve species of plant-parasitic nematode (Meloidogyne arenaria, M. chitwoodi, M. hapla, M. incognita, and M. javanica), three closely related tylenchid taxa (Heterodera glycines, Globodera pallida, and G. rostochiensis) and a distant taxon, Caenorhabditis elegans. Our multiple-gene approach is based on sampling more than 80,000 publicly available tylenchid EST sequences to identify phylum-wide orthologues. Bayesian inference, minimum evolution, maximum likelihood and protein distance methods were employed for phylogenetic reconstruction and hypothesis tests were constructed to elucidate diVerential selective pressures across the phylogeny for each gene. Our results place M. incognita and M. javanica as sister taxa, with M. arenaria as the next closely related nematode. SigniWcant diVerences in selective pressure were revealed for some genes under some hypotheses, though all but one gene are exclusively under purifying selection, indicating conservation across the orthologous groups. This EST-based multigene analysis is a Wrst step towards accomplishing genome-wide coverage for tylenchid evolutionary analyses. 2005 Elsevier Inc. All rights reserved. Keywords: Bayesian; Caenorhabditis elegans; COG; Cyst nematode; Orthologues; Root-knot nematode
1. Introduction Meloidogyne is a large genus of plant-parasitic, rootknot nematodes (RKN) containing more than 60 species in the order Tylenchida (Eisenbeck and Triantaphyllou, 1991). Of these, four (M. incognita, M. arenaria, M. javanica, and M. hapla) are major agricultural pests worldwide. Other RKN, such as M. chitwoodi, are important pests regionally. As a genus, Meloidogyne has a very broad host-range (Sasser, 1980) but individual species may have a more restricted host preference. Because this diversity of speciWcity has obvious practical consequences in agriculture, there has been a long history of
*
Corresponding author. Fax: +1 919 515 9500. E-mail address:
[email protected] (D. McK. Bird).
1055-7903/$ - see front matter 2005 Elsevier Inc. All rights reserved. doi:10.1016/j.ympev.2005.03.016
research aimed at identifying RKN species and establishing evolutionary relationships within the genus and between closely related genera such as the cyst nematodes (Heterodera and Globodera). The traditional approach to identifying individual Meloidogyne species (reviewed by Eisenbeck and Triantaphyllou, 1991) is based on diVerentiation of morphological characteristics, which also formed the initial basis for evolutionary studies. A feature of this approach is that morphology potentially reXects the integrated outcome of the expression of many genes. Unfortunately, it is not known what those genes are or what contribution each gene makes to the organism as a whole. Subsequent attempts to establish Meloidogyne phylogenies have been based on cytogenetic criteria (Triantaphyllou, 1966, 1985), DNA hybridization techniques (Castagnone-Sereno et al., 1993; Xue et al., 1992) and
E.H. Scholl, D.McK. Bird / Molecular Phylogenetics and Evolution 36 (2005) 536–545
isozyme analysis (Dalmasso and Bergé, 1978; Dickson et al., 1971; Esbenshade and Triantaphyllou, 1987). Studies based on mitochondrial genes (Powers et al., 1986; Powers and Sandall, 1988) showed that mitochondria seem to be hypervariable, with the resultant tree reXecting only the relationship of the genes sampled for the analysis rather than the underlying species relationships. Recently, the high variability found in nuclear ribosomal internal transcribed spacers has been exploited for phylogenetic inference between closely related Meloidogyne species (Hugall et al., 1999), and the portion of this locus encoding the highly conserved 18S rRNA subunit (SSU) was exploited to examine more diverged relationships within Meloidogyne and between Meloidogyne and other tylenchid genera (DeLey et al., 2002). A limitation of using single, individual loci is that locus history within a collection of species will not always accurately reXect the relationship of the species as a whole (Dorris et al., 1999; Doyle, 1992; Pamilo and Nei, 1988; Takahata, 1989; Wu, 1991). Issues such as gene duplication and lineage sorting may contribute to varying degrees of discordance between a gene tree and a species tree. Furthermore, analysis of sequences with high numbers of tandem repeats, such as ribosomal genes, creates a risk of choosing genes that have undergone gene conversion and are therefore not actually true orthologues (Slowinski and Page, 1999). One solution to increase the reliability of capturing the correct underlying phylogenetic signal from gene sequences is to sample a larger gene set; any biases presented by a single gene with a history not reXecting that of the species will hopefully be oVset with the larger selection of genes that better reXect that proper relationship (Takahata, 1989). Ideally, entire genomes would be compared, but complete genomes remain unavailable for any tylenchids. However, large numbers of tylenchid sequences have been generated as expressed sequence tags (ESTs) (Bird et al., 2002; McCarter et al., 2000, 2003), and we have exploited these sequences for a multi-gene approach to resolve evolutionary relationships between this important group of plant-parasitic nematodes. The extensively annotated and completed C. elegans genome sequence (The C. elegans Sequencing Consortium, 1998) serves as an essential reference. Thus, we are able to initiate a program of comparative genomics before the actual completion of any tylenchid genome. Rather than select genes based on some preconceived notion of what might be the “best” Tylenchida candidates, we designed a data mining approach that was blind to gene function, based on three minimal requirements. First, the gene must be expressed and so potentially present in a cDNA library. Second, the gene must have an identiWable orthologue in C. elegans. The presence of a gene across diverse nematode orders presumably reXects some evolutionary constraint (Blaxter et al., 1998); such slowly evolving genes are prime can-
537
didates for studying evolutionary relationships (Hillis and Dixon, 1991; Kumazawa and Nishida, 1993; Mindell and Honeycutt, 1990; Smith, 1989). Hence, this requirement will tend to eliminate genes which may have undergone accelerated evolution during some speciWc aspect of the ecology of the particular nematode species as well as eliminate genes recently acquired by horizontal transfer (Scholl et al., 2003), which clearly would bias the species trees. Genes meeting these two criteria were compared to identify orthologues of each in each species. Such groups of orthologues present in multiple species are termed “clusters of orthologous groups of proteins” or COGs (Tatusov et al., 1997), and are identiWed based on an iterative series of reciprocal, pairwise comparisons of the inferred product of each gene in each species, leading to a network deWning the COG (Fig. 1). COGs thus constructed were subjected to a third minimal criterion, that each selected gene be present in the available dataset either of three cyst nematodes plus C. elegans, or three Meloidogyne plus C. elegans, or Wve species total. The analysis we present recovers the most robust phylogeny for nine nematode species (Meloidogyne arenaria, M. chitwoodi, M. hapla, M. incognita, M. javanica, Heterodera glycines, Globodera pallida, G. rostochiensis, and C. elegans). We used a set of 47 orthologous genes initially identiWed from genome-wide EST analysis. Hypothesis tests have also been established for each of the 47 genes to investigate diVerential selective pressure amongst the species represented by each gene.
2. Materials and methods 2.1. Available data All sequences were obtained from resources publicly available on April 30, 2003. The full Caenorhabditis elegans genome and the full set of predicted proteins (WormPep) were downloaded as FASTA Wles from the National Center for Biological Information (NCBI) GenBank database (Wheeler et al., 2002). G. pallida ESTs (1878 entries) were similarly obtained. Meloidogyne hapla (13,869) and M. chitwoodi (3456) ESTs were obtained from the nematode.net ftp server (McCarter et al., 2003). These ESTs had been subjected to a rigorous, automated quality Wltering pipeline prior to being made publicly available, but may retain some errors. The G. pallida ESTs had been processed through an equivalent quality pipeline. M. arenaria (3321), M. incognita (5662), M. javanica (5574), H. glycines (4307) and G. rostochiensis (5039) data were kindly provided as cluster consensus sequences by the Plant Parasitic Nematode sequencing project (Bird et al., 2002). C. elegans genome sequence is eVectively error-free (The C. elegans Sequencing Consortium, 1998).
538
E.H. Scholl, D.McK. Bird / Molecular Phylogenetics and Evolution 36 (2005) 536–545
Fig. 1. Schemata conceptualizing design and construction of clusters of orthologous groups (COGs). The tail of each arrow represents a query sequence for a BLAST analysis, and the head deWnes the best match within the target species. Double-headed arrows summarize reciprocal best matches. (A) A hypothetical, four-species COG depicting the relationship of gene A in each species to its orthologue in each of the other three species. The network is closed and complete because all potential pairwise matches are present and concordant. This represents one orthologous group. (B) A COG network exhibiting reciprocal match disagreement. The orthologues of gene A from Species 1 are identiWable in Species 3 and 4, but the relationship of this gene between Species 1 and 2 is ambiguous, owing to a best match of gene A in Species 2 to gene B in Species 1. Genes A and B in Species 1 may be paralogous. (C) Linking COGs. Solid arrows depict relationships within an initial COG, and dotted arrows reXect links between pre-assembled COGs. Dot-dash arrows reXect the incorporation of C. elegans orthologues. Because all genes have the same reciprocal best match to C. elegans, the cluster is closed and complete. (D) Some analyses as in (C) result in a gene in one species having a best match to a gene in C. elegans (dashed line) diVerent from the rest of the COG; such genes (RKN 1, Gene A in this example) were stricken from the COG. RKN, Meloidogyne spp.; Cyst, Heterodera and Globodera spp.
2.2. EST clustering Meloidogyne chitwoodi and M. hapla EST sequences were clustered into consensus sequences using the incremental clustering algorithm (ICA) as implemented by the n2tool program from the ICA-tools suite (Parsons, 1995), with the subsequence similarity score threshold set to 20; both clone orientations were compared. Each cluster was conWrmed by manual examination of the alignments. The smaller G. pallida EST dataset was clustered by pairwise BLASTN analysis (Altschul et al., 1990). The sequences for the remaining Wve species were processed through the semi-automatic NemaGene pipeline (McCarter et al., 2003), which yields cluster consensus sequences with a very low error rate. Because our cluster consensus sequences (along with those from nematode.net) can be considered to deWne genes, we refer to them as such for the subsequent analyses. 2.3. Identifying orthologues Deliberately ignoring gene ontology, genes were gathered into clusters of orthologous groups (COGs) with a procedure that is essentially the same as described (Tatusov et al., 1997). Gene-for-gene, pairwise comparisons were made between diVerent species using the Hardware Accelerated, Tera-BLAST algorithm (Time-Logic, Crystal Bay, NV) executed on three TimeLogic DeCypher nodes. As the initial Wlter, single FASTA Wles containing the genes for each parasitic species were used as queries for a six-phase translated BLASTX comparison to the
C. elegans WormPep database. Parasite genes that failed to produce a match to C. elegans at an e value of less than 1.0e¡20 were eliminated. Remaining sequences were used in additional six-phase to six-phase TeraTBLASTX searches with each of the genes from the chosen parasitic species, and the top match for each gene in each species for each search was recorded. Matches with e values of less than 1.0e¡15 were considered signiWcant. We began COG assembly within the genus Meloidogyne by comparing the reciprocal best hits between all pairs of Meloidogyne species. An initial pairwise grouping was arbitrarily selected (M. incognita and M. arenaria), and orthologues between these two species were established. The nascent COG was extended to a third species (M. javanica) by examining reciprocal best hits for the COG members established in the initial two species. The procedure was continued for each Meloidogyne species, such that multi-dimensional networks of orthologues were developed that incorporated each of the BLAST results (Fig. 1A). Once all species were exhausted based on the initial pairing of M. incognita and M. arenaria, a new “initial” pair was chosen (in this case M. incognita and M. javanica) and the procedure repeated. Similar procedures were continued until all pairs of orthologues between Meloidogyne species were incorporated into a networked grouping. Once completed for Meloidogyne, the same procedure was implemented for the three cyst nematode species. Meloidogyne COGs and cyst nematode COGs were linked by examining pairwise reciprocal best matches between the Meloidogyne genes within one COG and
E.H. Scholl, D.McK. Bird / Molecular Phylogenetics and Evolution 36 (2005) 536–545
their cyst nematode orthologues to form a new COG spanning the family Heteroderidae. To add C. elegans genes to this COG, each gene in a Heteroderidae COG was used in a six-phase to six-phase Tera-BLASTX search against the C. elegans nucleotide database (Fig. 1C). Any gene in a group of parasitic nematodes that did not have the same best match to C. elegans as the majority of the other genes within that same group was removed from the COG (Fig. 1D). If multiple sequences from a single species were represented in a COG, those sequences were aligned to the C. elegans nucleotide sequence using CLUSTAL-W (Thompson et al., 1994). When the alignment revealed distinct parasite EST clusters that clearly deWned a single gene (i.e., had not correctly been assigned to a single cluster), those EST clusters were concatenated to reXect the true relationship of the ESTs to a gene and to facilitate restoration of the COG network (Fig. 1A). In cases where the alignment did not reveal EST clustering anomalies, multiple genes from a single species were deemed to be putative paralogues. Such COGs were removed from subsequent analyses. 2.4. Phylogenetic analysis Based on deduction of the correct reading frames for the EST clusters from Tera-BLASTX alignments to C. elegans WormPep, individual genes (one COG at a time) were conceptually translated and aligned using MegAlign (DNASTAR, Madison, WI) to implement the CLUSTAL-W algorithm (Thompson et al., 1994) with default parameter values. Thus, reXecting a codon alignment, original nucleotide sequences were recovered. Adjustments to each alignment were made by hand and individual nucleotide alignments were concatenated to create one multi-gene alignment (Brown et al., 2001; Kiontke et al., 2004). Bayesian inference was employed for phylogenetic reconstruction using version 3.0b4 of the MrBayes software (Huelsenbeck and Ronquist, 2001). In accordance with the likelihood ratio test (LRT) and Akaike Information Criteria (AIC) results from ModelTest3.60 (Posada and Crandall, 1998), a general time-reversible model was selected allowing for -distributed, amongsite rate variation with four categories. Four chains were run for 100,000 iterations, sampling every 20 iterations for a total of 5000 samples. The Wrst 1000 iterations were suYcient to reach stationarity and were thus discarded as burn-in. Three additional runs to conWrm convergence were performed for 50,000 iterations each, sampling every 50 steps with the Wrst 1000 iterations discarded as burn-in. Additional phylogenetic analyses were performed using methods of minimum evolution, maximum likelihood of combined data, and maximum likelihood of partitioned data. Minimum evolution trees were
539
constructed based on nucleotide distances, using the Kimura 2-parameter and F84 models as calculated by the dnadist program in PHYLIP (Felsenstein, 1993). Nucleotide sequence alignments were further examined as combined data as well as partitioned for each of the 47 individual genes, allowing for sequence heterogeneity, in a maximum likelihood analysis using the HKY85 + model with eight categories and stepwise addition for tree searching, as implemented in the PAML v3.13a package of programs (Yang, 1997). In addition, all sequences in the alignment were conceptually translated and PHYLIP (Felsenstein, 1993) was employed to calculate protein distances with the Jones–Taylor–Thornton substitution matrix (Jones et al., 1992). Tree topologies were inferred from the protein distances using neighbor joining. Non-parametric bootstrapping with 10,000 replicates provided an estimate of overall support for the topology. Gaps were treated as missing data in all analyses. 2.5. Testing for selective pressure Three hypotheses tests were established to determine diVerences in selective pressure in pre-deWned partitions of the phylogeny for each of the 47 genes in the dataset. The null hypothesis for each test allows for only one value for the ratio of non-synonymous to synonymous substitution sites, , (Messier and Stewart, 1997), for the whole phylogeny. The Wrst alternative hypothesis (a tworatio test) allowed for two ratios, one along the branch leading to C. elegans and one for the remainder of the phylogeny. The second alternative hypothesis (a threeratio test) allowed for three ratios, one along the branch leading to C. elegans, one for the branches contained in the set of cyst nematodes (G. pallida, G. rostochiensis, and H. glycines) and one for the branches contained in the set of root-knot nematodes (Meloidogyne spp.). The Wnal alternative hypothesis (a free-ratio test) allowed each branch in the phylogeny to have its own ratio. These hypothesis establish a series of branch-model tests (Yang, 1998) which partition the phylogeny into biologically relevant categories. The ratio was estimated for each gene for each of the four models (one-ratio, two-ratio, three-ratio, and free-ratio) using the CODONML program from PAML 3.13 (Yang, 1997). Topology was set for each gene in congruence with the topology recovered from the multigene alignment and Bayesian analysis. Codon frequencies were calculated from the average nucleotide frequency at each of the three codon positions. The transition/transversion ratio and branch lengths were estimated from the data. The maximum likelihood of the topology, given the alignments for the individual genes and the diVerent ratio models, was calculated and likelihood ratio tests were used to determine signiWcant diVerence between the models.
540
E.H. Scholl, D.McK. Bird / Molecular Phylogenetics and Evolution 36 (2005) 536–545
3. Results 3.1. EST clustering Some nematode EST data sets have previously been subjected to clustering (see Section 2.1), but this was not the case for the M. hapla, M. chitwoodi, or G. pallida EST sets. Consequently, we undertook clustering, generating 1283 consensus sequences from the M. hapla EST set and 839 cluster consensuses for M. chitwoodi. As is typically observed (McCarter et al., 2003), the distribution of the cluster sizes approximately followed that expected for a Poisson process. However, an attempt to cluster the smaller, G. pallida EST set using ICA-tools revealed deviation from the anticipated Poisson distribution. Rather than use these clusters, we performed a manual clustering based on BLAST queries of each sequence against the database of all G. pallida sequences (i.e., a self BLAST), yielding 989 clusters. 3.2. COG assembly Seventy-four clusters of orthologous groups containing genes from at least three plant-parasitic nematodes plus C. elegans were deduced (Table 1). In four instances (in three COGs), the best BLAST match of an individual gene to C. elegans contradicted the remainder of the COG, so those genes were disqualiWed from their COG (Fig. 1D). Based on the best match to C. elegans, it was further apparent that some of the COGs actually deWned the same gene; manual resolution of these redundancies reduced the number of COGs to 58 (Table 1). Three additional COGs were eliminated because we considered the global alignment to be too short to be phylogenetically informative. Table 1 COG enrichment during the clustering process Step 1c Step 2d Step 3e Step 4f Step 5g
COG numbera
Putative paraloguesb
74 58 56 55 47
31 15 15 14 0
a The number of clusters of orthologous groups (COGs) remaining following sequential steps to eliminate spurious genes. b Number of COGs containing putative paralogues at each step of the COG enrichment. c Step 1. Initial COG number as determined by pair-wise BLAST results. d Step 2. COG number following merging diVerent COGs based on identical best matches to C. elegans. e Step 3. COG number following correction of EST clustering anomalies. f Step 4. COG number following removal of members with poor global alignment. g Step 5. Final number of COGs after removal of sets harboring putative paralogues.
Fig. 2. Summary of ontology classiWcations of the 47 C. elegans genes deWning the COGs established for this study, as curated by WormBase release WS100, May 10, 2003; http://www.wormbase.org (Harris et al., 2003).
One consequence of our COGing strategy is that nonoverlapping, parasite EST clusters which actually deWne a single gene will be Xagged as putative paralogues (Fig. 1B). Indeed, alignment with C. elegans revealed six such instances. Correct concatenation of those EST clusters reXected the true relationship to their C. elegans orthologue and reduced the total number of COGs with putative paralogues from 14 to 8. Because it is not apparent which, if either, of the paralogues in such COGs truly is an orthologue, these eight COGs were disregarded. The remaining 47 COGs (Table 1) were deconstructed by species, and the sequences concatenated into a multi-gene alignment of 22,920 bp. An average of six species was represented for each gene, with an average aligned sequence length of 511 bp per gene. The G + C content, which averaged 45%, and ranged from 40% in M. chitwoodi to 52% in G. rostochiensis, showed no strong species bias. Ontology was determined for each of the 47 genes post hoc (Fig. 2), and the WormPep identiWcation numbers, which are searchable at WormBase (Harris et al., 2003), serve as a convenient way to identify each COG. 3.3. Phylogenetic analysis Phylogenetic reconstruction of the multi-gene alignment using Bayesian inference recovered only four topologies with a non-zero approximated posterior probability (Fig. 3). The topology which places M. incognita and M. javanica as sister taxa (Fig. 3A) has a posterior probability of approximately 97.4%. Further, the topologies provide 99.3% bipartition support for grouping M. incognita and M. javanica as sister taxa. The next best topology (Fig. 3B), with a posterior probability of approximately 1.89%, retains the relationships between M. incognita and M. javanica, and M. arenaria, diVering from the Wrst topology only in the inferred
E.H. Scholl, D.McK. Bird / Molecular Phylogenetics and Evolution 36 (2005) 536–545
541
Fig. 3. The four topologies with non-zero posterior probability (A–D) reported by the Bayesian analysis of 47 concatenated genes assembled from individual gene alignments of COG members. The posterior probability for each topology is indicated in percentage.
relationship between M. hapla and M. chitwoodi. The third topology (Fig. 3C) groups M. incognita as most closely related to M. arenaria. The Wnal topology (Fig. 3D), which places M. javanica and M. arenaria together, as did the phylogeny reported from an analysis of 18S rRNA data (DeLey et al., 2002), has a posterior probability of only approximately 0.32%. Convergence was conWrmed by the fact that each independent run of the Bayesian analysis resulted in topologies and parameters with similar posterior estimates, regardless of the randomly chosen initial state for the Markov chain. Importantly, similar analysis using a pared-down set of genes present in M. arenaria and M. incognita and M. javanica datasets resolved the same relationship between these three species. The congruency between the two phylogenies and the use of independently evolving genes in the analysis argues against the presence of artifacts caused by lineage-sorting (Rokas et al., 2003). Based on the posterior probabilities from the Bayesian analyses, we Wnd strong support for the Wrst topology (Fig. 3A). In addition, we undertook minimum evolution and maximum likelihood analyses of the nucleotide alignment as well as protein distance analyses
of a conceptual translation of our data. All three approaches produced phylogenies congruent with our Bayesian Wndings. Bootstrapping, conducted as part of the minimum evolution and protein distance analyses, indicates high support for the relationships. Fig. 4 summarizes the phylogeny we recovered, and indicates branch lengths, posterior probabilities and bootstrap support for the protein distance analysis for each bifurcation in the tree. 3.4. Testing for selective pressure The phylogeny recovered from the multi-gene alignment was partitioned into three diVerent sub-trees reXecting parasitic and non-parasitic nematodes; nonparasitic, root-knot and cyst nematodes; and individual species. These were used to establish two-ratio, threeratio, and free ratio hypothesis tests, respectively, to investigate diVerential selective pressure across the subpartitions of the phylogeny. Likelihood ratio tests revealed signiWcantly diVerent ratios at the 95% level in 23.40% of the two-ratio tests, 45.24% of the three-ratio tests and 48.94% of the free-ratio tests. All but one of the
542
E.H. Scholl, D.McK. Bird / Molecular Phylogenetics and Evolution 36 (2005) 536–545
Fig. 4. Phylogenetic reconstruction of tyenchid species, using C. elegans as an outgroup. Branch lengths, determined by Bayesian analysis, indicate number of nucleotide substitutions per site along the length of the branch. Scale bar reports estimated number of nucleotide substitutions per site. Bootstrap support for 10,000 replicates in protein distance analysis is indicated in italics and posterior probability of clade support from Bayesian analysis of nucleotide data is indicated in bold.
ratios (across all models) was less than one, an indication of purifying, or negative, selection. The sole divergence from this was on a branch leading to M. incognita in the gene identiWed as hel-1, based on the C. elegans annotation. A value of D 1.0885 was determined in the free-ratio model. Overall, though statistically signiWcant diVerences in selective pressure are indicated under each model for some of the genes, the genes are generally under purifying selection. As the process used to establish the COGs is likely to select more conserved genes, these results were as expected. Likelihoods for each model for each gene and p values for each hypothesis test are summarized in Table 2.
4. Discussion The goal of many phylogenetic analyses is to establish the evolutionary relationship between a chosen set of species. However, species trees can be obscured by reporting relationships speciWc to the chosen trait(s). We were able to exploit EST data to employ a strategy that simultaneously examines many genes that are unlikely to have experienced diVerential evolutionary pressures. Through this strategy, we have potentially suppressed sampling bias that may result in the inferred history of a chosen gene improperly reXecting the evolutionary history of the species (Rokas et al., 2003). Traits widely conserved across the phylum are less likely to have experienced accelerated evolution, and thus are more likely to successfully recover the true underlying species relationships. A convenient way to sample such traits is through direct sampling of gene sequences. We required that each tylenchid gene we sampled also be present in C. elegans. Tylenchida and Rhabditidae are in diVerent major clades of Nematoda
(Blaxter et al., 1998). Therefore, these genes are more likely to be evolutionarily conserved, as conWrmed by our selective pressure analysis. Further, by selecting groups of genes based solely on sequence availability rather than some preconceived notion of which genes are more likely to resolve the phylogeny, we minimize biasing the results towards a speciWc relationship between genes. That is not to say that any individual gene widely distributed across the phylum might not demonstrate a history which does not reXect the true relationship between the species, but the additional use of multiple genes will tend to counteract any spurious phylogenetic signal from such genes (Rokas et al., 2003). An additional requirement of gene-based phylogenetic analysis is that truly orthologous genes be compared between species. By removing any gene in a COG that did not have the same best BLAST match to C. elegans as the other members of that COG, we were able to distinguish parasitic genes that might be suYciently related to genes in a particular COG but which in reality are not true orthologues (Fig. 1D). Examples include members of gene families. We similarly were able to Xag errors in the original EST clustering and also further identify paralogous genes which might spuriously join the COG by identifying situations in which more than one gene from a given species was represented in a COG (Fig. 1B). Thus, stringent criteria were employed to identify networks of orthologues across species (COGs) using the available EST data. Although we constructed our COGs with no a priori knowledge of gene function, it is instructive to examine the identity of the genes deWning each COG. With only one exception (F40E10.6, specifying WormPep predicted protein CE31508), all the genes encode proteins with identiWable biochemical functions. Collectively, the COGs deWned here reXect a cross-section of what are
E.H. Scholl, D.McK. Bird / Molecular Phylogenetics and Evolution 36 (2005) 536–545
543
Table 2 Summary of log-likelihood and p values for tests of diVerential selective pressure Gene IdentiWcation
Likelihoods under each hypothesis
p values
WormPep IDa
Gene nameb
Nullc
Two-ratioc
Three-ratioc
Free-ratio c
H1d
H2d
H3d
CE17986 CE03482 CE01236 CE09650 CE02618 CE09156 CE01543 CE06577 CE00664 CE31867 CE04821 CE06090 CE03025 CE05860 CE07033 CE12898 CE31508 CE06652 CE03278 CE16954 CE17154 CE02255 CE01253 CE02883 CE02249 CE20218 CE05598 CE09901 CE14956 CE09654 CE05441 CE01030 CE04561 CE06107 CE07370 CE10884 CE32364 CE11052 CE12918 CE01225 CE27706 CE09682 CE07537 CE16260 CE30997 CE20547 CE03567
rpl-12 let-70 mlc-3 F25H2.5 ifb-1 egl-21 rpl-10 rps-7 rps-1 ida-1 kin-2 K04D7.1 hel-1 rps-11 rpl-11.2 T01B11.4 F40E10.6 ZK829.4 rpl-26 pas-6 K0F2.2 rpl-5 rpt-3 T24H7.1 iV-2 Y37D8A.14 rpl-3 F33D11.10 rps-18 pas-5 daf-21 rpl-16 rps-8 K06A4.3 gpd-3 rps-22 F09E5.2 pat-10 rps-16 F01F1.12 dim-1 hsp-1 T25F10.6 tra-3 F48E8.5 unc-60 asp-4
¡1441.6987 ¡1190.1912 ¡1619.2371 ¡1676.9238 ¡2096.9080 ¡1860.1260 ¡1811.7487 ¡1907.7945 ¡1694.4378 ¡1114.2571 ¡1444.9008 ¡1381.2516 ¡1377.0339 ¡1524.6079 ¡1680.2398 ¡1774.1306 ¡935.8562 ¡1926.5655 ¡1408.1271 ¡2061.9630 ¡2344.6213 ¡1288.9829 ¡1613.8837 ¡1647.5781 ¡1318.7377 ¡1735.9355 ¡1637.8118 ¡1725.6755 ¡1387.3232 ¡2181.7693 ¡1677.4047 ¡1868.9563 ¡2211.1566 ¡1858.1318 ¡1522.4629 ¡1088.8771 ¡1917.4083 ¡1412.1055 ¡1631.6038 ¡2318.7086 ¡980.7353 ¡1760.5109 ¡1684.6186 ¡1502.3538 ¡1414.3384 ¡1482.0056 ¡1416.6127
¡1437.9624 ¡1189.1719 ¡1619.1764 ¡1674.1831 ¡2095.7031 ¡1858.2977 ¡1810.0967 ¡1906.3177 ¡1693.8161 ¡1111.7537 ¡1443.6502 ¡1372.2543 ¡1376.4903 ¡1524.1984 ¡1675.7578 ¡1773.9540 ¡935.2987 ¡1918.6859 ¡1408.1242 ¡2061.1022 ¡2343.5623 ¡1283.4461 ¡1611.9508 ¡1643.2708 ¡1313.0705 ¡1730.5813 ¡1637.4795 ¡1724.0283 ¡1387.2731 ¡2179.6595 ¡1676.2391 ¡1868.4888 ¡2210.8180 ¡1855.9615 ¡1521.3238 ¡1086.3794 ¡1916.1375 ¡1410.6432 ¡1630.8937 ¡2314.7270 ¡975.5845 ¡1758.3270 ¡1684.5532 ¡1502.0081 ¡1412.9982 ¡1481.8806 ¡1411.2025
¡1436.7090 ¡1189.0728 ¡1599.0035 ¡1671.9306 ¡2094.3284 ¡1856.4090 ¡1831.2059 ¡1901.1285 ¡1693.4103 — ¡1443.8847 — ¡1375.4805 ¡1524.0162 ¡1675.7506 ¡1765.7778 ¡929.0803 ¡1895.7788 ¡1400.6398 ¡2059.7730 ¡2343.8069 ¡1273.9159 ¡1611.5208 ¡1638.2868 — ¡1724.8014 ¡1634.1055 ¡1724.0230 ¡1386.6595 ¡2179.5159 ¡1676.2278 ¡1850.8068 ¡2208.8867 ¡1843.5355 ¡1518.3289 ¡1082.9085 ¡1908.9712 ¡1410.5506 ¡1628.1119 ¡2298.9189 — ¡1758.0897 ¡1678.8639 — ¡1412.9082 ¡1481.6447 ¡1411.0704
¡1431.5299 ¡1181.9702 ¡1583.5526 ¡1666.8454 ¡2093.4172 ¡1845.2612 ¡1805.2866 ¡1899.3056 ¡1688.6247 ¡1102.5697 ¡1434.4974 ¡1367.8524 ¡1363.8048 ¡1518.0034 ¡1668.1298 ¡1762.7760 ¡925.5441 ¡1886.7618 ¡1398.2828 ¡2058.0235 ¡2339.0120 ¡1268.5006 ¡1610.0229 ¡1628.2144 ¡1312.4221 ¡1711.2847 ¡1623.9537 ¡1721.0340 ¡1383.6139 ¡2178.9332 ¡1674.6638 ¡1842.2206 ¡2200.0424 ¡1840.5869 ¡1515.7450 ¡1081.4587 ¡1890.1827 ¡1407.6630 ¡1624.2792 ¡2291.1690 ¡971.5712 ¡1754.3719 ¡1673.2876 ¡1499.2226 ¡1407.4397 ¡1478.4533 ¡1406.2421
0.0063 0.1534 0.4254 0.0192 0.1206 0.0558 0.0691 0.0857 0.2648 0.0252 0.1389 0.0000 0.2971 0.3655 0.0028 0. 5523 0.2910 0.0001 0.9393 0.1895 0.1456 0.0009 0.0493 0.0033 0.0008 0.0011 0.4149 0.0695 0.7515 0.0400 0.1268 0.3336 0.4106 0.0372 0.1312 0.0254 0.1109 0.0872 0.2334 0.0048 0.0013 0.0366 0.7355 0.4057 0.1016 0.6171 0.0010
0.0068 0.3268 0.0000 0.0068 0.0758 0.0243 0.0000 0.0013 0.3579 — 0.4693 — 0.2115 0.5534 0.0112 0.0002 0.0011 0.0000 0.0006 0.1119 0. 4429 0.0000 0.0941 0.0001 — 0.0000 0.0246 0. 1916 0.5149 0.1050 0.3082 0.0000 0.1033 0.0000 0.0160 0.0026 0.0002 0.2112 0.0304 0.0000 — 0.0888 0. 0786 — 0. 2392 0.6970 0.0039
0.0024 0.0877 0.0000 0.0026 0.3226 0.0000 0.3746 0.0094 0.1687 0.0001 0.0876 0.0002 0.0017 0.1049 0.0070 0.0038 0.0021 0.0000 0.0031 0.4454 0.0819 0.0000 0.2592 0.0000 0.0132 0.0000 0.0005 0.3190 0.1913 0.6839 0.7051 0.0000 0.0350 0.0000 0.0977 0.0216 0.0000 0.5430 0.1454 0.0000 0.0004 0.0560 0.1143 0.3944 0.0320 0.3113 0.0230
0.2340
0.4524
0.4894
% with signiWcant diVerences under each hypothesis a
As identiWed from http://www.wormbase.org. IdentiWcation of each gene used in the analysis and gene name, according to annotation of the C. elegans member of the COG. c Log-likelihood calculated for each branch model by Wxing the topology according to the results from multiple-gene alignment, but using the alignments for the individual genes. d p value for the log likelihood test for each of the three hypotheses. p values in bold are signiWcant. b
mainly core metabolism genes (Fig. 2). As such, they are likely to be slowly evolving genes. Taken together, our method of gene selection in conjunction with the use of multiple genes for phylogenetic reconstruction is a
powerful strategy to yield a tree reXecting the actual species relationships. Resolution of the correct phylogeny for the Wve species of Meloidogyne allows for further dissection of
544
E.H. Scholl, D.McK. Bird / Molecular Phylogenetics and Evolution 36 (2005) 536–545
biological traits found within the individual species. For example, there are at least two hypotheses which attempt to explain the triploid chromosome number observed in M. arenaria. The Wrst argues that the increase in chromosome number is a recently acquired trait occurring as a result of a hybridization of two other Meloidogyne species, leading to the genesis of M. arenaria (Hugall et al., 1999; Triantaphyllou, 1985). Another postulates that the triploid chromosome number is an ancestral trait and chromosome loss and/or fusion has occurred in more recently evolved species (Triantaphyllou, 1985). If M. arenaria is a hybrid of two other Meloidogyne species, our phylogeny suggests that this hybridization involves neither M. incognita nor M. javanica, as neither are ancestral to M. arenaria. Alternatively, the idea of triploid chromosome number being an ancestral trait and chromosome fusion creating the hypotriploid species, such as M. incognita and M. javanica, is consistent with our phylogenetic framework. A fully resolved phylogeny provides us with the framework needed to develop sound hypotheses for various characteristics of the diVerent Meloidogyne species in a biologically relevant context. Issues as diverse as pathogen-binding to Meloidogyne (Davies et al., 2001) and host choice (Bird and Kaloshian, 2003) are most fruitfully examined in light of a robust phylogeny.
5. Conclusion Despite the lack of a fully sequenced tylenchid genome, we were able to establish a set of 47 COGs for use in phylogenetic data, thus illustrating that such genomic analyses can be initiated using EST data mining methods. Ongoing sequencing projects will enhance COG construction through completion of partial COGs as well as establishment of new COGs. A goal of the morphological approach to understanding evolution was ultimately to sample the entire genome using physical attributes. This EST-based multi-gene analysis is a Wrst step towards accomplishing genomewide coverage, at a Wner resolution than morphological characters are able to provide. We expect that in the near future complete genomes will be available for direct sampling, permitting the morphologists’ goals Wnally to be realized.
Acknowledgments The authors thank Drs. B. Adams, S. Aris-Brosou, K.G. Davies, J.P. McCarter, and J.L. Thorne and two anonymous reviewers for their constructive criticism and advice. This research was supported by NSF Plant Genome award DBI0077503.
References Altschul, S.F., Gish, W., Miller, W., Myers, E.W., Lipman, D.J., 1990. Basic local alignment search tool. J. Mol. Biol. 215, 403–410. Bird, D.McK., Clifton, S.W., Kepler, T., Kieber, J.J., Thorne, J., Opperman, C.H., 2002. Genomic dissection of a nematode-plant interaction: a tool to study plant biology. Plant Physiol. 129, 394–395. Bird, D.McK., Kaloshian, I., 2003. Are roots special? Nematodes have their say. Physiol. Mol. Plant P. 62, 115–123. Blaxter, M.L., DeLey, P., Garey, J., Liu, L.X., Scheldeman, P., Vierstraete, A., VanXetern, J., Mackey, L.Y., Dorris, M., Frisse, L.M., Vida, J.T., Thomas, W.K., 1998. A molecular evolutionary framework for the phylum Nematoda. Nature 392, 71–75. Brown, J.R., Douady, C.J., Italia, M.J., Marshall, W.E., Stanhope, M.J., 2001. Universal trees based on large combined protein sequence data sets. Nat. Genet. 28, 281–285. Castagnone-Sereno, P., Piotte, C., Uijthof, J., Abad, P., Wajnberg, E., Vanlerberghe-Masutti, F., Bongiovanni, M., Dalmasso, A., 1993. Phylogenetic relationships between amphimictic and parthenogenetic nematodes of the genus Meloidogyne as inferred from repetitive DNA analysis. Heredity 70, 195–204. The C. elegans Sequencing Consortium, 1998. Genome sequence of the nematode C. elegans: a platform for investigating biology. Science 282, 2012–2018. Dalmasso, A., Bergé, J.G., 1978. Molecular polymorphism and phylogenetic relationship in some Meloidogyne spp.: application to the taxonomy of Meloidogyne. J. Nematol. 10, 323–332. Davies, K.G., Fargette, M., Balla, G., Daud, A.I., Duponnois, R., Gowen, S.R., Mateille, T., Phillips, M.S., Sawadogo, A., Trivino, C., Vouyoukalou, E., Trudgill, D.L., 2001. Cuticle heterogeneity as exhibited by Pasteuria spore attachment is not linked to the phylogeny of parthenogenetic root-knot nematodes (Meloidogyne spp.). Parasitology 122, 111–120. DeLey, I.T., DeLey, P., Vierstraete, A., Karssen, G., Moens, M., VanXeteren, J., 2002. Phylogenetic analyses of Meloidogyne small subunit rDNA. J. Nematol. 34, 319–327. Dickson, D.W., Huisingh, D., Sasser, J.N., 1971. Dehydrogenases, acid and alkaline phosphatases, and esterases for chemotaxonomy of selected Meloidogyne, Ditylenchus, Heterodera, and Aphelenchus spp. J. Nematol. 3, 1–16. Dorris, M., DeLey, P., Blaxter, M.L., 1999. Molecular analysis of nematode diversity and the evolution of parasitism. Parasitol. Today 15, 188–193. Doyle, J.J., 1992. Gene trees and species tree: molecular systematics as one-character taxonomy. Syst. Bot. 17, 144–163. Eisenbeck, J.D., Triantaphyllou, H.H., 1991. Root-knot nematodes: Meloidogyne species and races. In: Nickle, W.R. (Ed.), Manual of Agricultural Nematology. Marcel Dekker, New York, pp. 191–274. Esbenshade, P.R., Triantaphyllou, A.C., 1987. Enzymatic relationships and evolution in the genus Meloidogyne (Nematoda: Tylenchida). J. Nematol. 19, 8–18. Felsenstein, J., 1993. PHYLIP (Phylogeny Inference Package) version 3.6a2. Department of Genetics, University of Washington, Seattle, WA. Harris, T.W., Lee, R., Schwarz, E., Bradnam, K., Lawson, D., Chen, W., Blasier, D., Kenny, E., Cunningham, F., Kishore, R., Chan, J., Muller, H.-M., Petcherski, A., Thorisson, G., Day, A., Bieri, T., Rogers, A., Chen, C.-K., Spieth, J., Sternberg, P., Durbin, R., Stein, L.D., 2003. WormBase: a cross-species database for comparative genomics. Nucleic Acids Res. 31, 133–137. Hillis, D.M., Dixon, M.T., 1991. Ribosomal DNA: molecular evolution and phylogenetic inference. Q. Rev. Biol. 66, 411–453. Huelsenbeck, J.P., Ronquist, F., 2001. MRBAYES: Bayesian inference of phylogeny. Bioinformatics 17, 754–755. Hugall, A., Stanton, J., Mortiz, C., 1999. Reticulate evolution and the origins of ribosomal internal transcribed spacer diversity in apomictic Meloidogyne. Mol. Biol. Evol. 16, 157–164.
E.H. Scholl, D.McK. Bird / Molecular Phylogenetics and Evolution 36 (2005) 536–545 Jones, D.T., Taylor, W.R., Thornton, J.M., 1992. The rapid generation of mutation data matrices from protein sequences. Comp. Appl. Biosci. 8, 275–282. Kiontke, K., Gavin, N.P., Raynes, Y., Roehrig, C., Piano, F., Fitch, D.H.A., 2004. Caenorhabditis phylogeny predicts convergence of hermaphroditism and extensive intron loss. PNAS 101, 9003–9008. Kumazawa, Y., Nishida, M., 1993. Sequence evolution of mitochondrial tRNA genes and deep-branch animal phylogenetics. J. Mol. Evol. 37, 380–398. McCarter, J., Abad, P., Jones, J.T., Bird, D.McK., 2000. Rapid gene discovery in plant parasitic nematodes via expressed sequence tags. Nematology 2, 719–731. McCarter, J.P., Mitreva, M.D., Martin, J., Dante, M., Wylie, T., Rao, U., Pape, D., Bowers, Y., Theising, B., Murphy, C., Kloek, A.P., Chiapelli, B., Clifton, S.W., Bird, D.McK., Waterston, R., 2003. Analysis and functional classiWcation of transcripts from the root-knot nematode Meloidogyne incognita. Genome Biol. 4, R26.1–R26.19. Messier, W., Stewart, C.-B., 1997. Episodic adaptive evolution of primate lysozymes. Nature 385, 151–154. Mindell, D.P., Honeycutt, R.L., 1990. Ribosomal RNA in vertebrates: evolution and phylogenetic applications. Annu. Rev. Ecol. Syst. 21, 541–566. Pamilo, P., Nei, M., 1988. Relationships between gene trees and species trees. Mol. Biol. Evol. 5, 568–583. Parsons, J.D., 1995. Improved tools for DNA comparison and clustering. Bioinformatics 11, 603–613. Posada, D., Crandall, K., 1998. MODELTEST: testing the model of DNA substitution. Bioinformatics 14, 817–818. Powers, T.O., Platzer, E.G., Hyman, B.C., 1986. Species-speciWc restriction site polymorphism in root-knot nematode mitochondrial DNA. J. Nematol. 18, 288–293. Powers, T.O., Sandall, L.J., 1988. Estimation of genetic divergence in Meloidogyne mitochondrial DNA. J. Nematol. 20, 505–511. Rokas, A., Williams, B.L., King, N., Carroll, S.B., 2003. Genome-scale approaches to resolving incongruence in molecular phylogenies. Nature 425, 798–804. Sasser, J.N., 1980. Root-knot nematodes: a global menace to crop production. Plant Dis. 64, 36–41.
545
Scholl, E.H., Thorne, J.L., McCarter, J.P., Bird, D.McK., 2003. Horizontally transferred genes in plant-parasitic nematodes: a highthroughput genomic approach. Genome Biol. 4, R39. Slowinski, J.B., Page, R.D.M., 1999. How should species phylogenies be inferred from sequence data? Syst. Biol. 48, 814–825. Smith, A.B., 1989. RNA sequence data in phylogenetic reconstruction: testing the limits of its resolution. Cladistics 5, 321–344. Takahata, N., 1989. Gene genealogy in three related populations: consistency probability between gene and population trees. Genetics 122, 957–966. Tatusov, R.L., Natale, D.A., Garkavtsev, I.V., Tatusova, T.A., Shankavaram, U.T., Rao, B.S., Kiryutin, B., Galperin, M.Y., Fedorova, N.D., Koonin, E.V., 1997. The COG database: new developments in phylogenetic classiWcation of proteins from complete genomes. Nucleic Acids Res. 29, 22–28. Thompson, J.D., Higgins, D.G., Gibson, T.J., 1994. CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-speciWc gap penalties and weight matrix choice. Nucleic Acids Res. 22, 4673–4680. Triantaphyllou, A.C., 1966. Polyploidy and reproductive patterns in the root-knot nematode Meloidogyne hapla. J. Morphol. 118, 403–414. Triantaphyllou, A.C., 1985. Cytogenetics, cytotaxonomy, and phylogeny of root-knot nematodes. In: Sasser, J.N., Carter, C.C. (Eds.), An Advanced Treatise on Meloidogyne, Vol. I. Biology and Control. North Carolina State University Graphics, Raleigh, NC, pp. 113–126. Wheeler, D.L., Church, D.M., Lash, A.E., Leipe, D.D., Madden, T.L., Pontius, J.U., Schuler, G.D., Schriml, L.M., Tatusova, T.A., Wagner, L., Rapp, B.A., 2002. Database resources of the National Center for Biotechnology Information: 2002 update. Nucleic Acids Res. 30, 13–16. Wu, C., 1991. Inferences of species phylogeny in relation to segregation of ancient polymorphisms. Genetics 127, 429–435. Xue, B., Baillie, D., Beckenbach, K., Webster, J., 1992. DNA hybridization probes for studying the aYnities of three Meloidogyne populations. Fundam. Appl. Nematol. 15, 35–41. Yang, Z., 1997. PAML: a program package for phylogenetic analysis by maximum likelihood. Comp. Appl. Biosci. 13, 555–556. Yang, Z., 1998. Likelihood ratio tests for detecting positive selection and application to primate lysozyme evolution. Mol. Biol. Evol. 15, 568–573.