<title>Introduction</title> <p>
<italic>Aspergillus oryzae</italic> is a filamentous ascomycete used in traditional Japanese soy sauce production, and has important industrial applications as a producer of hydrolytic enzymes in solid-state fermentations ##UREF##0##[1]##, ##REF##17341818##[2]##. When the <italic>A. oryzae</italic> genome was sequenced it was found to be significantly larger and to contain more genes than the genomes of other species in the genus <italic>Aspergillus</italic>
##REF##16372010##[3]##–##REF##16372009##[5]##. At 37 Mb, the <italic>A. oryzae</italic> genome assembly is 23% larger than the <italic>A. nidulans</italic> assembly and 32% larger than the <italic>A. fumigatus</italic> assembly. It is predicted to code for 12,074 proteins of &gt;100 amino acids, which is 1,412 more proteins than <italic>A. nidulans</italic> and 2,444 more than <italic>A. fumigatus</italic>
##REF##16372000##[4]##. Although some of the differences in the predicted gene count between species may be due to differences in the bioinformatics methods used by the three sequencing laboratories to annotate the genomes, the larger DNA content of <italic>A. oryzae</italic> is incontrovertible. Moreover, even if one excludes singleton genes (those with no homologs in completed fungal genome databases), which are the ones most likely to be annotation artifacts, <italic>A. oryzae</italic> still contains 16–26% more genes than the other species: the numbers of non-singletons are 12,044, 10,425 and 9,574 for <italic>A. oryzae</italic>, <italic>A. nidulans</italic> and <italic>A. fumigatus</italic> respectively. There is also indirect evidence that most of the predicted genes in <italic>A. oryzae</italic> are functional: sequencing the genome of <italic>A. flavus</italic>, which is a very close relative of <italic>A. oryzae</italic> and has a similarly large genome, shows that 8,953 genes (85% of all orthologs between <italic>A. oryzae and A. flavus</italic>) have a conservative pattern of nucleotide substitution (<italic>K<sub>A</sub></italic>&lt;<italic>K<sub>S</sub></italic>) consistent with purifying selection to retain protein-coding capacity (N. Khaldi and G. Payne, unpublished results).</p> <p>Comparisons of gene order revealed the presence of syntenic blocks common to <italic>A. oryzae</italic>, <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>, as well as genomic blocks specific to <italic>A. oryzae</italic> that lack synteny with the other two species ##REF##16372010##[3]##. These non-syntenic blocks share the characteristics of, first, containing genes that seem to exist only in <italic>A. oryzae</italic> (that is, they lack orthologs in <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>) and, second, appearing in a mosaic manner throughout the genome of <italic>A. oryzae</italic>. Machida <italic>et al.</italic>
##REF##16372010##[3]## reported that, surprisingly, phylogenetic analysis of some genes from these non-syntenic regions indicated that they are distantly related paralogs of other genes that have a syntenic and apparently orthologous relationships among <italic>A. oryzae</italic>, <italic>A. nidulans</italic> and <italic>A. fumigatus.</italic> Thus, the extra genes in <italic>A. oryzae</italic> did not seem to be the result of recent gene duplications that happened specifically in <italic>A. oryzae</italic> after it diverged from the other two species; rather, they seemed to be the products of much older divergence event(s) ##REF##16372010##[3]##. This is particularly surprising because comprehensive genome-based phylogenetic analysis shows that <italic>A. nidulans</italic> is an outgroup to <italic>A. fumigatus</italic> and <italic>A. oryzae</italic>
##REF##16372000##[4]##, so on that basis the increased gene number in <italic>A. oryzae</italic> would most parsimoniously be explained by a more recent species-specific expansion.</p> <p>In this study we aim to find out the origin of the <italic>A. oryzae</italic> genes giving rise to this unusual phylogenetic tree topology. This topology was first noticed, but not named, by Machida <italic>et al.</italic>
##REF##16372010##[3]##. We refer to it as the “S topology” . We show that the discovery of these trees is not just a subjective observation: there is a genuine statistical excess of trees with the S topology in <italic>A. oryzae</italic> as compared to other <italic>Aspergillus</italic> species. We then explore several possible explanations for these trees. First, we test whether the pattern could have been caused by an ancient whole-genome duplication (WGD) in the common ancestor of the three <italic>Aspergillus</italic> species, similar to the WGD that occurred in an ancestor of <italic>Saccharomyces cerevisiae</italic>
##REF##9192896##[6]##, ##REF##15004568##[7]##, followed by a lower rate of gene loss in <italic>A. oryzae</italic> than in <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>. We find no evidence for WGD or similar events such as whole-chromosome duplication (aneuploidy). Second, we consider the possibility that the Topology S trees might be artifacts. If the paralogous pairs in <italic>A. oryzae</italic> were in fact the products of recent (species-specific) gene duplications, but one gene copy subsequently underwent rapid sequence divergence, the ensuing asymmetry of evolutionary rates could cause phylogenetic methods to infer an incorrect tree due to the phenomenon of long-branch attraction between the accelerated branch and the outgroup ##UREF##1##[8]##, ##REF##16207937##[9]##. Third, we consider the possibility that the extra genes were added to the <italic>A. oryzae</italic> genome by horizontal gene transfer (HGT).</p>
<title>Methods</title> <title>Data</title> <p>Genome sequences and annotations of <italic>A. nidulans</italic>, <italic>F. graminearum</italic> and <italic>U. maydis</italic> were downloaded from the Broad Institute (<ext-link ext-link-type="uri" xlink:href="http://www.broad.mit.edu/annotation/fgi/">http://www.broad.mit.edu/annotation/fgi/</ext-link>); <italic>A. fumigatus</italic> from TIGR (<ext-link ext-link-type="uri" xlink:href="http://www.tigr.org/tdb/e2k1/afu1/">http://www.tigr.org/tdb/e2k1/afu1/</ext-link>); and <italic>A. oryzae</italic> from the DDBJ database (accession numbers AP007150-AP007177). BLAST searches against the <italic>A. flavus</italic> genome were performed on <ext-link ext-link-type="uri" xlink:href="http://gaplabg5.cfr.ncsu.edu/blast/blast.html">http://gaplabg5.cfr.ncsu.edu/blast/blast.html</ext-link>.</p> <title>Identification of Sets S<sub>AO</sub>, S<sub>AN</sub> and S<sub>AFU</sub>
</title> <p>We accepted as triplet sets of orthologous genes in <italic>A. oryzae</italic>, <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>, those genes that were best mutual BLASTP hits (<italic>i.e</italic>., the protein pairs AO1-AN, AO1-AFU, and AFU-AN were all bidirectional best BLASTP hits between their respective genomes, using NCBI BLASTP with default settings). We then tested whether a paralog of the <italic>A. oryzae</italic> gene was present. To identify pairs of paralogs in <italic>A. oryzae</italic> that satisfy topology S, we required that (<italic>i</italic>) the <italic>A. nidulans</italic> and <italic>A. fumigatus</italic> genes in the triplet do not have paralogs produced by species-specific duplications, and (<italic>ii</italic>) the amino acid sequence divergences between the three orthologs in the triplet were all lower than the distance between each of the orthologs and the <italic>A. oryzae</italic> paralog. If these conditions were met, the genes were retained for phylogenetic analysis and the <italic>A. oryzae</italic> orthologous and paralogous copies were termed <italic>AO1</italic> and <italic>AO2</italic> respectively. An automated search using these criteria resulted in the Set S<sub>AO</sub> of 456 loci that are duplicated only in <italic>A. oryzae</italic> and have topology S. We then used phylogenetic methods to verify the topology of the trees in Set S<sub>AO</sub>. Using the likelihood ratio test, we tested whether the unrooted topology ((AO2,AN)(AO1,AFU)) has a significantly higher likelihood than both of the two possible alternative topologies for the four sequences. We found that at a cutoff of α = 0.05, the expected tree was significantly more likely at 443 of the 456 loci, and at the remaining 13 loci no other topology was significantly better than the expected. We retained all 456 loci for further analysis. Analogous methods were used to identify Sets S<sub>AN</sub> and S<sub>AFU</sub>. Synonymous and nonsynonymous sequence divergence was calculated using the yn00 program in the PAML package ##REF##17483113##[35]##. Fop values were calculated using CodonW (<ext-link ext-link-type="uri" xlink:href="http://mobyle.pasteur.fr">http://mobyle.pasteur.fr</ext-link>) with its predefined set of optimal codons for <italic>A. nidulans</italic>.</p> <title>Testing for WGD</title> <p>To evaluate the hypothesis that the pairs of genes in Set S<sub>AO</sub> could be the products of a whole-genome duplication, we tested whether these pairs tend to be located on sister genomic regions. If a WGD occurred and produced two pairs of duplicates, A′ and A″, B′ and B″, we would expect that if A′ and B′ are physically close together in the genome, their paralogs A″ and B″ should also be close together somewhere else in the genome. We examined whether the Set S<sub>AO</sub> paralogs in the <italic>A. oryzae</italic> genome tend to be arranged in this pattern, relative to random expectations. In this analysis we did not take into account the distinction between the <italic>AO1</italic> and <italic>AO2</italic> members of the <italic>A. oryzae</italic> pairs (which have and do not have, respectively, orthologs in <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>). We used a sliding window approach to find regions of relatively high paralog density in <italic>A. oryzae</italic>. We considered windows of size <italic>j</italic> genes (<italic>j</italic> = 200 genes in the examples shown in ##FIG##3##Figure 4##), and counted how many genes in Set S<sub>AO</sub> they contained. We then retained only the denser windows for further analysis: those with some number <italic>i</italic> genes in Set S<sub>AO</sub> (<italic>i</italic> = 6 in ##FIG##3##Figure 4C##). For each window with a particular value of <italic>i</italic>, we then tried to measure the physical distance (in units of genes) spanned by their paralogs elsewhere in the <italic>A. oryzae</italic> genome. In many cases the paralogs were not all on the same chromosome. For those windows whose paralogs were all on the same chromosome, we computed the total chromosomal distance <italic>d</italic> occupied by the paralogous windows. We then evaluated the statistical significance of the observed value of <italic>d</italic> by comparing it to the empirical distribution of values obtained in 10,000 computer simulations where the locations of the paralogous genes were randomized. We carried out a parallel analysis on <italic>S. cerevisiae</italic> as a positive control, using the 551 gene pairs from Byrne and Wolfe ##REF##16169922##[36]## in place of Set S<sub>AO</sub>. In a second analysis, we also included data from Set S<sub>AO</sub> genes whose paralogs were located on different <italic>A. oryzae</italic> chromosomes, by adding a fixed penalty to the <italic>d</italic> value in these cases (##FIG##3##Figure 4D##). The penalty was set at 3000, which is approximately the average number of genes on an <italic>A. oryzae</italic> chromosome. Experimenting with different values of <italic>j</italic> and <italic>i</italic> (data not shown) did not reveal any statistically significant evidence of WGD for <italic>A. oryzae</italic>, whereas the results for <italic>S. cerevisiae</italic> were consistently significant.</p> <title>Testing for HGT</title> <p>We used the basidiomycete <italic>U. maydis</italic> as an outgroup for all phylogenetic analyses. Among the 456 loci in Set S<sub>AO</sub>, we retained a subset (Q<sub>AO</sub>) of 122 loci for which the same <italic>U. maydis</italic> protein was identified as the reciprocal best BLASTP hit of the S<sub>AO</sub> set members from all three <italic>Aspergillus</italic> species. We then used BLASTP to identify additional fungal homologs of the Set Q<sub>AO</sub> genes pairs in the NCBI nonredundant protein sequence database. Proteins were aligned using ClustalW ##REF##7984417##[37]##, and Gblocks ##REF##10742046##[38]## was used to remove poorly aligned positions and regions. Maximum likelihood trees were constructed using PHYML ##REF##14530136##[39]## with the JTT amino acid substitution matrix and four categories of substitution rates. Bootstrapping was done using PHYML's default options with 100 replicates per run. The trees for each gene were inspected manually to detect common themes and classified into the three possible mutually exclusive topologies (##FIG##4##Figure 5##). Analogous methods were used to make and study Sets Q<sub>AN</sub> and Q<sub>AFU</sub>.</p>
<title>Results</title> <title>An excess of topology S trees</title> <p>Machida <italic>et al.</italic>
##REF##16372010##[3]## noted that many gene families are expanded in <italic>A. oryzae</italic> and that phylogenetic trees of these families often show the topology that we have called Topology S. This topology is remarkable because it implies that a pair of paralogous genes (<italic>AO1</italic> and <italic>AO2</italic>; ##FIG##0##Figure 1a##) that are present in <italic>A. oryzae</italic> originated by a gene duplication that occurred before the speciation events that separated <italic>A. oryzae</italic>, <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>. If that interpretation is correct, then orthologs of the <italic>AO2</italic> gene must subsequently have been deleted from the <italic>A. nidulans</italic> and <italic>A. fumigatus</italic> genomes, and the shape of the species phylogeny necessitates that two independent deletions of <italic>AO2</italic> must have occurred in <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>.</p> <p>We quantified this phenomenon by searching systematically for loci that form a paralogous pair in <italic>A. oryzae</italic> but are single-copy in <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>, and where one of the <italic>A. oryzae</italic> copies is divergent from the other three sequences as expected under Topology S (see <xref ref-type="sec" rid="s4">Methods</xref>). We found 456 such pairs, which accounts for 19% of all pairs of paralogous genes in <italic>A. oryzae</italic>. We refer to these loci as Set S<sub>AO</sub>. For comparison, we likewise defined two analogous sets of loci that contain extra copies in other species. Set S<sub>AFU</sub> consists of genes that are double-copy only in <italic>A. fumigatus</italic> and have the equivalent of Topology S for that species (##FIG##0##Figure 1b##). This set contains 202 pairs, which is 11% of all paralog pairs in <italic>A. fumigatus</italic>. Similarly, Set S<sub>AN</sub> consists of 219 pairs (12% of all paralog pairs) that are duplicated only in <italic>A. nidulans</italic> and have the topology shown in ##FIG##0##Figure 1c##.</p> <p>The number of loci in Set S<sub>AO</sub> is significantly greater than the numbers in Set S<sub>AN</sub> (<italic>P</italic> = 2×10<sup>−10</sup> by Fisher test) and Set S<sub>AFU</sub> (<italic>P = </italic>7×10<sup>−12</sup>), whereas there is no significant difference between Sets S<sub>AN</sub> and S<sub>AFU</sub> (<italic>P</italic> = 0.43). This result confirms Machida <italic>et al.</italic>'s report that there is an excess of divergent paralogous gene copies in the <italic>A. oryzae</italic> genome.</p> <title>Topology S trees are unlikely to be artifacts of rate acceleration after gene duplication</title> <p>Gene duplication is often followed by a period of accelerated sequence evolution in one or both gene copies ##REF##11073452##[10]##–##REF##12952876##[12]##. If the rates are sufficiently unequal, the resulting phylogenetic trees can suffer from long branch attraction ##UREF##1##[8]##, ##REF##15861208##[13]##, an artifact of tree reconstruction methods that causes the longer branches to clump together in the tree regardless of their true phylogenetic relationship. Where one member of a duplicated gene pair has accelerated, attraction between this gene's branch and the outgroup used to root the tree (which is usually also a long branch) can result in an incorrect topology that makes the gene duplication appear older than it actually was. In the case of <italic>Saccharomyces cerevisiae</italic>, this artifact affected more than half of all the duplicate gene pairs that were formed by its WGD ##REF##16207937##[9]##, ##REF##12594514##[14]##.</p> <p>To test whether the <italic>A. oryzae</italic> Topology S trees are the result of a similar artifact, we compared the levels of synonymous nucleotide divergence (<italic>K<sub>S</sub></italic>) between each of the <italic>A. oryzae</italic> genes (<italic>AO1</italic> and <italic>AO2</italic>) and their single homolog (<italic>AN</italic>) in <italic>A. nidulans</italic>. Rate acceleration after gene duplication is expected to affect synonymous divergence to a much lesser extent than nonsynonymous divergence, if at all. If the duplication that produced <italic>AO1</italic> and <italic>AO2</italic> was a species-specific event in <italic>A. oryzae</italic> (<italic>i.e.</italic>, the Topology S trees are artifacts) then the synonymous distances |<italic>AO1</italic>-<italic>AN</italic>| and |<italic>AO2</italic>-<italic>AN</italic>| should be equal, and the excess amino acid sequence divergence of <italic>AO2</italic> could be attributed to accelerated protein evolution. Alternatively, if <italic>AO2</italic> truly branched off before the <italic>A. oryzae</italic>-<italic>A. nidulans</italic> divergence (<italic>i.e.</italic>, the Topology S trees are correct) then the synonymous distance |<italic>AO1</italic>-<italic>AN</italic>| should be less than the synonymous distance |<italic>AO2</italic>-<italic>AN</italic>| (Note that <italic>AO2</italic> is always defined as the <italic>A. oryzae</italic> gene that lacks an apparent <italic>A. nidulans</italic> ortholog; ##FIG##0##Figure 1a##). We find that the synonymous divergence is indeed higher for the <italic>AO2</italic> genes: the median <italic>K<sub>S</sub></italic> values are 3.35 for |<italic>AO2</italic>-<italic>AN</italic>| and 2.87 for |<italic>AO1</italic>-<italic>AN</italic>| (<italic>P</italic>&lt;10<sup>−4</sup> by Wilcoxon test; ##FIG##1##Figure 2##). This result indicates that the excess of genes in Set S<sub>AO</sub> is not an artifact of long branch attraction.</p> <p>One possible caveat with the above analysis is that, if the <italic>AO2</italic> genes are expressed at lower levels than their <italic>AO1</italic> homologs, then codon bias in the <italic>AO1</italic> genes might decrease the apparent level of synonymous divergence from <italic>A. nidulans</italic>
##REF##3328816##[15]##. It is known that <italic>A. oryzae</italic> genes located in genomic regions that are not syntenic to other species have lower average expression levels than those in syntenic blocks ##REF##16372010##[3]##. However, this expression difference is due primarily to the low expression of genes in non-syntenic blocks that lack homologs elsewhere in the genome, and not to <italic>AO2</italic>-type genes ##REF##17967552##[16]##. Using codon bias as a proxy for gene expression levels, we found that the frequency of optimal codons (Fop) is slightly higher in <italic>AO1</italic> genes than in their <italic>AO2</italic> homologs (median Fop values 0.4290 and 0.3945 respectively), which is suggestive of higher expression but the difference is not statistically significant (Wilcoxon test of the hypothesis that the ratio Fop(AO1)/Fop(AO2) in each locus is greater than 1; P = 0.25). We doubt that the difference in Fop values is sufficient to cause the difference in <italic>K<sub>S</sub></italic> values observed in ##FIG##1##Figure 2##, but are unable to test this directly.</p> <title>No evidence for ancient WGD in <italic>Aspergillus</italic>
</title> <p>Another possible scenario that might explain the extra genes in <italic>A. oryzae</italic> is an ancient whole genome duplication (WGD) prior to the speciation of <italic>A. oryzae</italic>, <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>, followed by parallel independent losses of most of the duplicated genes in the latter two species (##FIG##2##Figure 3##). Although this scenario might seem unlikely – and indeed Machida <italic>et al</italic>. ##REF##16372010##[3]## did not observe any long duplicated segments in the <italic>A. oryzae</italic> genome – it is a formal possibility so we searched rigorously for evidence of it.</p> <p>The WGD events described so far in other genomes were detected using a number of different methods ##REF##15510166##[reviewed in 17]##. Patterns in the chromosomal locations of paralogs give by far the clearest signals in situations where a WGD is present in a sequenced genome but no genome sequence from a closely-related unduplicated outgroup species is available ##REF##9192896##[6]##, ##REF##12374856##[18]##–##REF##17721507##[21]##. Immediately after a WGD, the genome consists of pairs of identical chromosomes containing identical genes in the same order. This pattern will gradually become eroded by gene losses and chromosomal rearrangements, but we would expect that even an old WGD should leave a signal: the order of genes along a given chromosomal region should be correlated with the order of their paralogs along some other region in the genome. In other words, paralogs should not just be randomly distributed with regard to one another.</p> <p>To compare the distribution of locations of paralogous genes in the <italic>A. oryzae</italic> genome to random expectations we defined a distance measure, <italic>d</italic>, of the extent to which genes that are close together in a genome also have paralogs that are close together somewhere else in the genome (see <xref ref-type="sec" rid="s4">Methods</xref>). We estimated the statistical significance of the observed distances in <italic>A. oryzae</italic> by comparing it to the distribution of distances obtained in computer simulations where gene locations were shuffled randomly (##FIG##3##Figure 4##).</p> <p>We applied this method to both <italic>S. cerevisiae</italic> and <italic>A. oryzae</italic>. The method successfully detects a strong signal of non-random locations of paralogous gene pairs in <italic>S. cerevisiae</italic>, a species known to have undergone WGD (empirical <italic>P</italic>&lt;10<sup>−3</sup>; ##FIG##3##Figure 4A,B##), but there is no significant signal of WGD in the <italic>A. oryzae</italic> genome (<italic>P</italic> = 0.79; ##FIG##3##Figure 4C##). This result not only argues against WGD but also against segmental duplication which would be expected to leave a similar, though weaker, signal. We also used a modified measure to permit consideration of situations where two genes on one chromosome have paralogs on two different chromosomes (see <xref ref-type="sec" rid="s4">Methods</xref>) but again found no evidence of ancient WGD in <italic>A. oryzae</italic> (##FIG##3##Figure 4D##).</p> <p>In an independent analysis we carried out all-against-all BLASTP searches within the protein sets of each of <italic>A. oryzae</italic>, <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>. We then made a dot matrix plot for each species showing the chromosomal locations of reciprocal-best-hit pairs within each genome. In previous studies on <italic>S. cerevisiae</italic> and plant species, approaches similar to this revealed strong diagonals corresponding to paralogous chromosomal segments that were formed by WGD ##REF##9192896##[6]##, ##REF##17721507##[21]##, ##REF##11006329##[22]##. In contrast, inspection of our plots for <italic>Aspergillus</italic> species did not reveal any indication of WGD or segmental duplications, in agreement with the previous report for <italic>A. oryzae</italic>
##REF##16372010##[3]##.</p> <title>Can the excess of Topology S genes be explained by Horizontal Gene Transfer?</title> <p>Taken together, the above results show that the excess of <italic>A. oryzae</italic> genes in Set S<sub>AO</sub> cannot be explained satisfactorily by recent gene duplications in <italic>A. oryzae</italic>, nor by an ancient WGD or segmental duplications in an <italic>Aspergillus</italic> ancestor. We therefore examined the third possibility of HGT into <italic>A. oryzae</italic>.</p> <p>To detect a possible donor taxon, we first analyzed the phylogenetic trees of a random sample of genes from Set S<sub>AO</sub>. This showed that the extra (<italic>AO2</italic>) copies of most genes are nested in the Ascomycota kingdom, so if HGT occurred the donor should also be in the Ascomycota. To provide a common reference point for phylogenetic reconstructions, we imposed the basidiomycete <italic>Ustilago maydis</italic> [whose genome is completely sequenced; 23], as the outgroup for all trees. This reduced our Set S<sub>AO</sub> to a subset Q<sub>AO</sub> of 122 genes for which an ortholog in <italic>U. maydis</italic> could be identified unambiguously by a reciprocal-best-BLAST-hits (RBH) approach.</p> <p>We constructed rooted phylogenetic trees for the 122 genes in Set Q<sub>AO</sub> with their homologs in other fungi (primarily Sordariomycetes), and sorted them manually into three mutually exclusive topology types depending on the position of the <italic>AO2</italic> gene (Types A, B and C; ##FIG##4##Figure 5##; Supplementary ##SUPPL##1##Figure S1##). In all three topologies the <italic>AO1</italic> gene of <italic>A. oryzae</italic> is clustered in a monophyletic group with its orthologs from <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>, and in all three the species in the subphylum Sordariomycetes (<italic>e.g.</italic>, <italic>Neurospora crassa</italic>, <italic>Fusarium graminearum</italic>, <italic>Magnaporthe grisea</italic>) form a monophyletic group. Type A trees are those where the second <italic>A. oryzae</italic> copy <italic>AO2</italic> clusters with Sordariomycete species to the exclusion of the genus <italic>Aspergillus</italic> (##FIG##4##Figure 5##, upper). Type B trees have a topology in which <italic>AO2</italic> lies outside both the <italic>Aspergillus</italic> genus and the Sordariomycetes (##FIG##4##Figure 5##, center). Type C trees are those where <italic>AO2</italic> forms a sister group to the clade of <italic>Aspergillus</italic> orthologs that includes <italic>AO1</italic>, with Sordariomycete sequences outside this pair (##FIG##4##Figure 5##, lower). Of these three topologies, only the Type A trees are directly suggestive of HGT because they implicate an identifiable donor lineage.</p> <p>There are only three possible topologies for a phylogenetic tree that consists of an outgroup and three ingroup taxa. The A, B, and C trees correspond to these three topologies; the ingroup clades are <italic>AO2</italic>, the <italic>Aspergillus</italic> clade that includes <italic>AO1</italic>, and the Sordariomycete clade (##FIG##4##Figure 5##, <italic>left</italic>). No other topologies can exist for our data unless the Sordariomycete or <italic>Aspergillus</italic> clades are broken up.</p> <p>The sorting process for trees in Set Q<sub>AO</sub> was done manually and with reference to the bootstrap support values for critical branches on the tree (we required a bootstrap value of ≥70% for topology-defining branches). To keep the process simple we rejected trees from complex gene families or with poor bootstrap support, and we permitted minor deviations from the expected species phylogeny for the <italic>Aspergillus</italic> (<italic>AO1</italic>) and Sordariomycete clades. We also classified trees where the <italic>AO2</italic> gene was placed within the Sordariomycetes, rather than as sister to Sordariomycetes, as Type A. These additional restrictions left only 35 classifiable trees from Set Q<sub>AO</sub>, of which 12 were Type A, 9 were Type B, and 14 were Type C (##TAB##0##Table 1## and Supplementary ##SUPPL##0##Table S1##). Likelihood Ratio Tests on these trees confirmed that monophyly of the <italic>AO1</italic> and <italic>AO2</italic> genes could be rejected at α&lt;0.0001. The fact that we find approximately equal numbers of trees of the three possible types is not consistent with the hypothesis that all <italic>AO2</italic> genes originated by HGT from the same donor. Instead, the result suggests that the <italic>AO2</italic> genes are heterogeneous in terms of their origins: if they all arose by HGT then the HGT events involved multiple donors, or if they all arose by gene duplication then those duplications occurred at multiple times.</p> <p>For comparison, we made similar analyses in <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>. We compiled two other sets of genes that are analogous to Set Q<sub>AO</sub> but have duplications only in <italic>A. nidulans</italic> (Set Q<sub>AN</sub>, 40 genes), or only in <italic>A. fumigatus</italic> (Set Q<sub>AFU</sub>, 28 genes), and that give phylogenetic trees that place the second gene (<italic>AN2</italic> or <italic>AFU2</italic>) outside a monophyletic group containing one gene from each of the three <italic>Aspergillus</italic> species. As expected, these gene sets are smaller than the corresponding set with duplication in <italic>A. oryzae</italic>. After phylogenetic tree construction we find that sets Q<sub>AN</sub> and Q<sub>AFU</sub> both present an absence of Topology A trees (##TAB##0##Table 1##). Although the numbers of classified trees are small (11 in Set Q<sub>AN</sub> and 6 in Set Q<sub>AFU</sub>) this deficit contrasts with the number of Topology A trees seen in the <italic>A. oryzae</italic> set. Topology A trees are the only ones where a putative HGT donor lineage can be identified (##FIG##4##Figure 5##), so it is striking that trees with this topology are found only in <italic>A. oryzae</italic>. If HGT from Sordariomycetes to <italic>A. oryzae</italic> accounts for the type A topologies, we can estimate that it maximally accounts for one-third of the extra genes in <italic>A. oryzae</italic> (that is, by extrapolation from the 12/35 proportion found in the loci that were classifiable).</p> <p>It is also notable that among the 12 <italic>A. oryzae</italic> gene pairs that give Topology A trees, the trees for 10 of these pairs place the <italic>AO2</italic> gene on a branch <italic>within</italic> the Sordariomycetes as opposed to lying sister to the Sordariomycetes. An example is the Topology A gene in ##FIG##4##Figure 5##, where the <italic>AO2</italic> sequence (NCBI identifier 83775924) clusters specifically with <italic>F. graminearum</italic>. Of the 10 genes that showed this property, five clustered with <italic>F. graminearum</italic>, three with <italic>M. grisea</italic> and two with <italic>N. crassa</italic>. A single consistent candidate donor lineage is again lacking, even within the Sordariomycetes.</p> <title>Functions and physical clustering of <italic>AO2</italic> genes with Topology A within the <italic>A. oryzae</italic> genome</title> <p>The set of 12 <italic>A. oryzae</italic> genes with Topology A is enriched in particular functions, most notably the hydrolytic traits that are characteristic of <italic>A. oryzae</italic> (Supplementary ##SUPPL##0##Table S1##). This species breaks down starch into sugars, and proteins into amino acids, by producing vast amounts of hydrolytic enzymes. The enrichment of some families in these classes among duplicated <italic>A. oryzae</italic> genes is striking, as previously also noted by Machida <italic>et al.</italic>
##REF##16372010##[3]##. For example, gene AO090103000299 (gi 83775625) is a carboxypeptidase that our method identified as an <italic>AO2</italic> gene giving Topology A. Further exploration of this gene family (##FIG##5##Figure 6##) shows that it is a member of a six-gene family in <italic>A. oryzae</italic>, only two of which have orthologs in <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>. Indeed this family consists of only two genes in the latter two species. Of the four other <italic>A. oryzae</italic> genes in the family, three cluster with Sordariomycete genes in the phylogenetic tree (two with different <italic>F. graminearum</italic> genes, and one (the gene AO090103000299 originally identified as an <italic>AO2</italic> candidate) with a <italic>M. grisea-N. crassa</italic> pair. In this <italic>A. oryzae</italic> gene family, only two genes were inherited vertically and are shared with other <italic>Aspergillus</italic> species, two (possibly three) originated by HGT from Sordariomycetes, and the sixth gene has an undetermined origin. Thus, HGT appears to have played a major role in the expansion of this carboxypeptidase family in <italic>A. oryzae</italic>.</p> <p>Comparison to the genome sequence of <italic>A. flavus</italic> (G. Payne <italic>et al</italic>., unpublished), which has a recent and clonal relationship to <italic>A. oryzae</italic>
##REF##9419385##[24]##, shows that many of the <italic>AO2</italic> genes have orthologs in that species. Therefore, the <italic>AO2</italic> genes were gained by <italic>A. oryzae</italic> after it diverged from <italic>A. fumigatus</italic> but before the split between <italic>A. flavus</italic> and <italic>A. oryzae</italic>.</p> <p>If foreign DNA enters a recipient genome in a single event, the foreign genes might be expected to exhibit some clustering in the recipient genome. Conversely if they arrived in independent events, or if the hypothesis of HGT is incorrect, they should exhibit random genomic locations. If an HGT event is very recent, the transferred genes may also still be clustered in the donor genome. We applied the method of Lee and Sonnhammer ##REF##12695325##[25]## to test whether the locations of the <italic>AO2</italic> genes that may have undergone HGT from Sordariomycetes show significant physical clustering in the <italic>A. oryzae</italic> genome. For the 12 <italic>AO2</italic> genes with tree topology A (Supplementary ##SUPPL##0##Table S1##), the clustering score is 56, as compared to a mean score of 21.7 for 100,000 replicates with randomized data (empirical <italic>P</italic> = 0.02). In contrast, the <italic>AO1</italic> partners of these genes do not show physical clustering (score = 6.6; <italic>P</italic> = 0.95), nor did we find significant clustering in comparisons of the <italic>A. oryzae AO2</italic> genes that gave B, C, or both B and C topologies (<italic>P</italic> = 0.19; <italic>P</italic> = 0.06 and <italic>P</italic> = 0.19 respectively). This result again highlights the distinctive nature of the <italic>AO2</italic> genes with tree topology A, in terms of both their phylogenetic and genomic placement. The physical clustering of the Topology A <italic>AO2</italic> genes is indicative of either a one-time transfer event, or a preference for integrating new genes into the genome at particular sites. We did not detect significant clustering in the <italic>F. graminearum</italic> genome of homologs of the <italic>AO2</italic> genes that form branches specifically with <italic>F. graminearum</italic> (though the dataset is small; only 5 genes).</p>
<title>Discussion</title> <p>A frustrating feature of our study was the extent of data loss that we encountered as the analysis progressed. The <italic>A. oryzae</italic> genome contains more than one thousand extra genes as compared to <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>, but when we attempted to use phylogenetic methods to determine the origin of these genes the dataset collapsed to only 35 informative loci (##TAB##0##Table 1##). One major cause of the reduction of the dataset is that we were only able to work with genes that are duplicated in <italic>A. oryzae</italic>. The reason for this is pragmatic: we cannot actually identify which, out of the 12,074 genes in <italic>A. oryzae</italic>, are the ‘extra’ ones, except in the 456 cases where the <italic>AO1</italic> gene was retained too (Set S<sub>AO</sub>). Moreover, our method ignored recent species-specific duplications in all species. A second problem was the lack of identifiable basidiomycete outgroup sequences for many of the genes that are duplicated in <italic>A. oryzae</italic>. It is difficult to speculate why these sequences do not exist, but one possibility is that the set of genes that has been added to the <italic>A. oryzae</italic> genome tends to be derived from the fastest-evolving, and hence ascomycete-specific, portion of the fungal genome ##REF##16755356##[26]##, ##REF##17090923##[27]##.</p> <p>We were able to rule out annotation error and WGD as possible sources of the extra genes in <italic>A. oryzae</italic> relative to other <italic>Aspergillus</italic> species. Our phylogenetic analysis also rules out species-specific gene duplication events in <italic>A. oryzae</italic> as the source of many (456) of these genes. These results leave only two viable possibilities for the origin of the extra genes: HGT, or older single-gene duplications (in the common ancestor of the three <italic>Aspergillus</italic> species) followed by independent parallel losses in <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>. If multiple single-gene duplications occurred, then the level of retention of ancient duplicated genes must have been much higher in <italic>A. oryzae</italic> than in the other two <italic>Aspergillus</italic> species, and parallel losses of the same genes must have occurred in <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>. This scenario seems unparsimonious, but it is not impossible and it may have occurred in a systematic manner if – for example – there were large differences in the effective population sizes of the species ##REF##14631042##[28]##.</p> <p>Parallel losses in <italic>A. nidulans</italic> and <italic>A. fumigatus</italic> may therefore account for some of the extra content in <italic>A. oryzae</italic>. However, barring phylogenetic tree error, parallel loss is specifically ruled out in the case of the 12 loci with Topology A (##FIG##4##Figure 5## and Supplementary ##SUPPL##0##Table S1##). This result seems to leave HGT as the most viable option for the source of much of the increased gene content of <italic>A. oryzae</italic>, which makes it all the more puzzling that we were unable to identify a consistent donor lineage. Instead, if we accept that all the tree topologies are correct, we must postulate that <italic>A. oryzae</italic> has received genes from at least five different donor lineages: one lineage that diverged prior to the Eurotiomycetes/Sordariomycetes split (giving Topology B trees), one lineage that is in Eurotiomycetes but outside <italic>Aspergillus</italic> (Topology C), and three separate Sordariomycete lineages that have specific affinities to <italic>F. graminearum</italic>, <italic>N. crassa</italic> and <italic>M. grisea</italic> respectively (accounting for the heterogeneity within the Topology A trees). We can speculate that some aspect of the lifestyle or physiology of <italic>A. oryzae</italic> (and <italic>A. flavus</italic>) could make it more able to take up foreign DNA by HGT than other <italic>Aspergillus</italic> species.</p> <p>Many previous published examples of HGT have involved genes that were transferred between distantly related kingdoms, such as from bacteria to eukaryotes ##REF##15947202##[29]##–##REF##16482157##[34]##. In many of these cases the gained genes were inferred to have changed the lifestyle of the recipient species markedly. In our case we are studying a more local exchange between reasonably closely related genomes that already share quite similar gene contents and functions. Our method of analysis necessitates that all the genes we identify as putative gains by HGT result in gene family expansion rather than in the acquisition of new biochemical functions. It is tempting to speculate that HGT helped <italic>A. oryzae</italic> adapt to some new evolutionary niche, but based on our analysis we can only assert that HGT enhanced the hydrolytic and fermentation abilities of this species. We can infer that, because these foreign genes were taken up and maintained in the genome, they contributed in some way to increasing <italic>A. oryzae</italic>'s evolutionary fitness, but we do not know exactly what that contribution was.</p>
<p>Conceived and designed the experiments: NK KHW. Performed the experiments: NK. Analyzed the data: NK. Wrote the paper: NK KHW.</p> <p>The genome sequence of <italic>Aspergillus oryzae</italic> revealed unexpectedly that this species has approximately 20% more genes than its congeneric species <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>. Where did these extra genes come from? Here, we evaluate several possible causes of the elevated gene number. Many gene families are expanded in <italic>A. oryzae</italic> relative to <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>, but we find no evidence of ancient whole-genome duplication or other segmental duplications, either in <italic>A. oryzae</italic> or in the common ancestor of the genus <italic>Aspergillus</italic>. We show that the presence of divergent pairs of paralogs is a feature peculiar to <italic>A. oryzae</italic> and is not shared with <italic>A. nidulans</italic> or <italic>A. fumigatus.</italic> In phylogenetic trees that include paralog pairs from <italic>A. oryzae</italic>, we frequently find that one of the genes in a pair from <italic>A. oryzae</italic> has the expected orthologous relationship with <italic>A. nidulans</italic>, <italic>A. fumigatus</italic> and other species in the subphylum Eurotiomycetes, whereas the other <italic>A. oryzae</italic> gene falls outside this clade but still within the Ascomycota. We identified 456 such gene pairs in <italic>A. oryzae</italic>. Further phylogenetic analysis did not however indicate a single consistent evolutionary origin for the divergent members of these pairs. Approximately one-third of them showed phylogenies that are suggestive of horizontal gene transfer (HGT) from Sordariomycete species, and these genes are closer together in the <italic>A. oryzae</italic> genome than expected by chance, but no unique Sordariomycete donor species was identifiable. The postulated HGTs from Sordariomycetes still leave the majority of extra <italic>A. oryzae</italic> genes unaccounted for. One possible explanation for our observations is that <italic>A. oryzae</italic> might have been the recipient of many separate HGT events from diverse donors.</p>
<title>Supporting Information</title>
<p>We thank Gary Payne and members of his lab for sharing <italic>A. flavus</italic> data, and are grateful to Gavin Conant, Kevin Byrne and Meg Woolfit for discussion.</p>
<fig id="pone-0003036-g001" position="float"><object-id pub-id-type="doi">10.1371/journal.pone.0003036.g001</object-id><label>Figure 1</label><caption><title>Definition of Topology S loci in (a) <italic>A. oryzae</italic>, (b) <italic>A. fumigatus</italic>, and (c) <italic>A. nidulans</italic>.</title><p>Yellow boxes indicate the genes that are orthologous in the three species. The red box contains the extra paralog (labeled <italic>AO2</italic>, <italic>AFU2</italic> or <italic>AN2</italic>). Loci having Topology S were identified in each species by an automated process without the use of an outgroup, but the root of the tree is likely to lie on the upper (longest) branch.</p></caption></fig> <fig id="pone-0003036-g002" position="float"><object-id pub-id-type="doi">10.1371/journal.pone.0003036.g002</object-id><label>Figure 2</label><caption><title>Boxplots of the distributions of synonymous nucleotide divergence (<italic>K</italic>
<sub>S</sub>) values between <italic>A. nidulans</italic> genes (<italic>AN</italic>) and <italic>A. oryzae</italic> gene copies <italic>AO1</italic> and <italic>AO2</italic>, for all genes in Set S<sub>AO</sub>.</title><p>The boxes show the median and interquartile values, with 95% of the data falling within the whiskers. <italic>K</italic>
<sub>S</sub> values were calculated using yn00 from the PAML package ##REF##9367129##[40]##. <italic>K</italic>
<sub>S</sub> values for the <italic>AO1</italic> genes are significantly lower than for the corresponding <italic>AO2</italic> genes (<italic>P</italic>&lt;10<sup>−4</sup>; Wilcoxon test).</p></caption></fig> <fig id="pone-0003036-g003" position="float"><object-id pub-id-type="doi">10.1371/journal.pone.0003036.g003</object-id><label>Figure 3</label><caption><title>Model of a WGD hypothesis for the origin of extra genes in <italic>A. oryzae</italic>.</title><p>A WGD in the common ancestor of the three <italic>Aspergillus</italic> species, followed by extensive gene losses in both <italic>A. nidulans</italic> and <italic>A. fumigatus</italic> could potentially account for the large number of genes in Set S<sub>AO</sub>. In this schematic diagram black squares represent genes, double vertical lines represent chromosome breaks, and dashed horizontal lines represent gene deletions or intrachromosomal rearrangements. Our method of testing for WGD involves sliding a window (blue box) through the <italic>A. oryzae</italic> genome to find clusters of nearby genes that are enriched in genes from Set S<sub>AO</sub> (red box). If a WGD occurred, the paralogs of the genes in the blue window are expected to be closer together on the sister <italic>A. oryzae</italic> genomic region than expected by chance.</p></caption></fig> <fig id="pone-0003036-g004" position="float"><object-id pub-id-type="doi">10.1371/journal.pone.0003036.g004</object-id><label>Figure 4</label><caption><title>Testing for whole genome duplication in <italic>A. oryzae</italic> and <italic>S. cerevisiae</italic>.</title><p>Black histograms show the distribution of distance scores (<italic>d</italic>) obtained in 10,000 simulations where paralog locations were randomized. Blue and green lines show the cutoff points for the most-clustered 1% and 5% of these data, respectively. Red lines show the observed <italic>d</italic> values for <italic>S. cerevisiae</italic> or <italic>A. oryzae</italic>. All the graphs are for windows of size <italic>j</italic> = 200. Panels A-C consider only windows where the paralogous genes are all on the same chromosome. (A) Results for <italic>S. cerevisiae</italic> for windows with <italic>i</italic> = 10 paralogs. (B) Results for <italic>S. cerevisiae</italic> for windows with <italic>i</italic> = 12 paralogs. (C) Results for <italic>A. oryzae</italic> for windows with <italic>i</italic> = 6 paralogs. (D) Results for <italic>A. oryzae</italic> for windows with <italic>i</italic> = 3 paralogs, including a penalty of 3000 for each interchromosomal rearrangement event.</p></caption></fig> <fig id="pone-0003036-g005" position="float"><object-id pub-id-type="doi">10.1371/journal.pone.0003036.g005</object-id><label>Figure 5</label><caption><title>Topologies A, B and C for rooted trees.</title><p>
<italic>Left</italic>: The three topologies correspond to the three possible points for placement of a root (arrows labeled A, B, C) onto an unrooted tree consisting of a gene pair from <italic>A. oryzae</italic> and their homolog from a Sordariomycete such as <italic>Fusarium graminearum</italic>. <italic>Center</italic>: Rooted representations of the topologies. <italic>Right</italic>: real examples of <italic>A. oryzae</italic> genes classified in each of the three topologies. Numbers in parentheses are NCBI protein identifiers. Trees were constructed using PHYML with bootstrap percentages (form 500 replicates) shown. The main groups differentiating the topologies from each other are colored.</p></caption></fig> <fig id="pone-0003036-g006" position="float"><object-id pub-id-type="doi">10.1371/journal.pone.0003036.g006</object-id><label>Figure 6</label><caption><title>A case of family expansion due to apparent HGT.</title><p>The gene family is a member of the metal dependent hydrolytic enzyme superfamily (amidases, aminoacyclases, and carboxypeptidases). The <italic>A. oryzae</italic> genes <italic>AO1</italic> and <italic>AO2</italic> originally identified in our automated search are labeled. Four of the six <italic>A. oryzae</italic> genes (83775625, 83765251, 83766521, 83770136) in this family lack orthologs in <italic>A. nidulans</italic> and <italic>A. fumigatus</italic>, and for the first two of these there is strong bootstrap support for HGT from Sordariomycete. Of these four genes, three show Topology A and one (83766521) shows Topology C relative to their nearest Sordariomycete homologs.</p></caption></fig>
<table-wrap id="pone-0003036-t001" position="float"><object-id pub-id-type="doi">10.1371/journal.pone.0003036.t001</object-id><label>Table 1</label><caption><title>Numbers of phylogenetic trees with Topologies A, B and C obtained from genes in Sets Q<sub>AO</sub>, Q<sub>AN</sub> and Q<sub>AFU</sub>.</title></caption><alternatives><table frame="hsides" rules="groups"><colgroup span="1"><col align="left" span="1"/><col align="center" span="1"/><col align="center" span="1"/><col align="center" span="1"/></colgroup><thead><tr><td align="left" rowspan="1" colspan="1">Topology</td><td colspan="3" align="left" rowspan="1">Number of trees</td></tr><tr><td align="left" rowspan="1" colspan="1"/><td align="left" rowspan="1" colspan="1">
<italic>A. oryzae</italic> (Set Q<sub>AO</sub>)</td><td align="left" rowspan="1" colspan="1">
<italic>A. nidulans</italic> (Set Q<sub>AN</sub>)</td><td align="left" rowspan="1" colspan="1">
<italic>A. fumigatus</italic> (Set Q<sub>AFU</sub>)</td></tr></thead><tbody><tr><td align="left" rowspan="1" colspan="1">A</td><td align="left" rowspan="1" colspan="1">12</td><td align="left" rowspan="1" colspan="1">0</td><td align="left" rowspan="1" colspan="1">0</td></tr><tr><td align="left" rowspan="1" colspan="1">B</td><td align="left" rowspan="1" colspan="1">9</td><td align="left" rowspan="1" colspan="1">9</td><td align="left" rowspan="1" colspan="1">2</td></tr><tr><td align="left" rowspan="1" colspan="1">C</td><td align="left" rowspan="1" colspan="1">14</td><td align="left" rowspan="1" colspan="1">2</td><td align="left" rowspan="1" colspan="1">4</td></tr></tbody></table></alternatives></table-wrap>
<supplementary-material content-type="local-data" id="pone.0003036.s001"><label>Table S1</label><caption><p>Names and putative functions of A. oryzae gene pairs showing Topologies A, B or C.</p><p>(0.03 MB DOC)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pone.0003036.s002"><label>Figure S1</label><caption><p>Trees classified as Types A, B and C in each Aspergillus species. Trees were constructed using PHYML as described in <xref ref-type="sec" rid="s4">Methods</xref>. In each tree, the sequences identified as AO1 and AO2 (for duplications A. oryzae), AN1 and AN2 (for duplications in A. nidulans), or AFU1 and AFU2 (for duplications in A. fumigatus) are labeled. NCBI identifier (GI) numbers for each sequence are shown.</p><p>(0.93 MB PDF)</p></caption></supplementary-material>
<fn-group><fn fn-type="COI-statement"><p><bold>Competing Interests: </bold>The authors have declared that no competing interests exist.</p></fn><fn fn-type="financial-disclosure"><p><bold>Funding: </bold>This study was supported by Science Foundation Ireland. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.</p></fn></fn-group>
<graphic xlink:href="pone.0003036.g001"/> <graphic xlink:href="pone.0003036.g002"/> <graphic xlink:href="pone.0003036.g003"/> <graphic xlink:href="pone.0003036.g004"/> <graphic xlink:href="pone.0003036.g005"/> <graphic id="pone-0003036-t001-1" xlink:href="pone.0003036.t001"/> <graphic xlink:href="pone.0003036.g006"/>
<media xlink:href="pone.0003036.s001.doc"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pone.0003036.s002.pdf"><caption><p>Click here for additional data file.</p></caption></media>
[{"label": ["1"], "element-citation": ["\n"], "surname": ["Hamada", "Fukushima", "Motai"], "given-names": ["T", "Y", "H"], "year": ["1992"], "article-title": ["Continuous production of soy sauce in a bioreactor."], "italic": ["Applications of Biotechnology to Traditional Fermented Foods"], "publisher-loc": ["Washington DC"], "publisher-name": ["National Academy Press"], "fpage": ["114"], "lpage": ["118"]}, {"label": ["8"], "element-citation": ["\n"], "surname": ["Felsenstein"], "given-names": ["J"], "year": ["1978"], "article-title": ["Cases in which parsimony and compatibility methods will be positively misleading."], "source": ["Syst Zool"], "volume": ["27"], "fpage": ["401"], "lpage": ["410"]}]
PLoS One. 2008 Aug 22; 3(8):e3036