Phylogenomic perspectives on speciation and reproductive isolation in a North American biodiversity hotspot: an example using California sages (Salvia subgenus Audibertia: Lamiaceae)

BACKGROUND AND AIMS

The California Floristic Province (CA-FP) is the most species-rich region of North America north of Mexico. One of several proposed hypotheses explaining the exceptional diversity of the region is that the CA-FP harbours myriad recently diverged lineages with nascent reproductive barriers. Salvia subgenus Audibertia is a conspicuous element of the CA-FP, with multiple sympatric and compatible species.

METHODS

Using 305 nuclear loci and both organellar genomes, we reconstruct species trees, examine genomic discordance, conduct divergence-time estimation, and analyse contemporaneous patterns of gene flow and mechanical reproductive isolation.

KEY RESULTS

Despite strong genomic discordance, an underlying bifurcating tree is supported. Organellar genomes capture additional introgression events not detected in the nuclear genome. Most interfertility is found within clades, indicating that reproductive barriers arise with increasing genetic divergence. Species are generally not mechanically isolated, suggesting that it is unlikely to be the primary factor leading to reproductive isolation.

CONCLUSIONS

Rapid, recent speciation with some interspecific gene flow in conjunction with the onset of a Mediterranean-like climate is the underlying cause of extant diversity in Salvia subgenus Audibertia. Speciation has largely not been facilitated by gene flow. Its signal in the nuclear genome seems to mostly be erased by backcrossing, but organellar genomes each capture different instances of historical gene flow, probably characteristic of many CA-FP lineages. Mechanical reproductive isolation appears to be only part of a mosaic of factors limiting gene flow.

Summary Tests of Introgression Are Highly Sensitive to Rate Variation Across Lineages

The evolutionary implications and frequency of hybridization and introgression are increasingly being recognized across the tree of life. To detect hybridization from multi-locus and genome-wide sequence data, a popular class of methods are based on summary statistics from subsets of 3 or 4 taxa. However, these methods often carry the assumption of a constant substitution rate across lineages and genes, which is commonly violated in many groups. In this work, we quantify the effects of rate variation on the D test (also known as ABBA-BABA test), the D3 test, and HyDe. All 3 tests are used widely across a range of taxonomic groups, in part because they are very fast to compute. We consider rate variation across species lineages, across genes, their lineage-by-gene interaction, and rate variation across gene-tree edges. We simulated species networks according to a birth-death-hybridization process, so as to capture a range of realistic species phylogenies. For all 3 methods tested, we found a marked increase in the false discovery of reticulation (type-1 error rate) when there is rate variation across species lineages. The D3 test was the most sensitive, with around 80% type-1 error, such that D3 appears to more sensitive to a departure from the clock than to the presence of reticulation. For all 3 tests, the power to detect hybridization events decreased as the number of hybridization events increased, indicating that multiple hybridization events can obscure one another if they occur within a small subset of taxa. Our study highlights the need to consider rate variation when using site-based summary statistics, and points to the advantages of methods that do not require assumptions on evolutionary rates across lineages or across genes.

Ancient Rapid Radiation Explains Most Conflicts Among Gene Trees and Well-Supported Phylogenomic Trees of Nostocalean Cyanobacteria

Prokaryotic genomes are often considered to be mosaics of genes that do not necessarily share the same evolutionary history due to widespread horizontal gene transfers (HGTs). Consequently, representing evolutionary relationships of prokaryotes as bifurcating trees has long been controversial. However, studies reporting conflicts among gene trees derived from phylogenomic data sets have shown that these conflicts can be the result of artifacts or evolutionary processes other than HGT, such as incomplete lineage sorting, low phylogenetic signal, and systematic errors due to substitution model misspecification. Here, we present the results of an extensive exploration of phylogenetic conflicts in the cyanobacterial order Nostocales, for which previous studies have inferred strongly supported conflicting relationships when using different concatenated phylogenomic data sets. We found that most of these conflicts are concentrated in deep clusters of short internodes of the Nostocales phylogeny, where the great majority of individual genes have low resolving power. We then inferred phylogenetic networks to detect HGT events while also accounting for incomplete lineage sorting. Our results indicate that most conflicts among gene trees are likely due to incomplete lineage sorting linked to an ancient rapid radiation, rather than to HGTs. Moreover, the short internodes of this radiation fit the expectations of the anomaly zone, i.e., a region of the tree parameter space where a species tree is discordant with its most likely gene tree. We demonstrated that concatenation of different sets of loci can recover up to 17 distinct and well-supported relationships within the putative anomaly zone of Nostocales, corresponding to the observed conflicts among well-supported trees based on concatenated data sets from previous studies. Our findings highlight the important role of rapid radiations as a potential cause of strongly conflicting phylogenetic relationships when using phylogenomic data sets of bacteria. We propose that polytomies may be the most appropriate phylogenetic representation of these rapid radiations that are part of anomaly zones, especially when all possible genomic markers have been considered to infer these phylogenies. [Anomaly zone; bacteria; horizontal gene transfer; incomplete lineage sorting; Nostocales; phylogenomic conflict; rapid radiation; Rhizonema.].

Integrated phylogenomic analyses unveil reticulate evolution in Parthenocissus (Vitaceae), highlighting speciation dynamics in the Himalayan-Hengduan Mountains

Hybridization caused by frequent environmental changes can lead both to species diversification (speciation) and to speciation reversal (despeciation), but the latter has rarely been demonstrated. Parthenocissus, a genus with its trifoliolate lineage in the Himalayan-Hengduan Mountains (HHM) region showing perplexing phylogenetic relationships, provides an opportunity for investigating speciation dynamics based on integrated evidence. We investigated phylogenetic discordance and reticulate evolution in Parthenocissus based on rigorous analyses of plastome and transcriptome data. We focused on reticulations in the trifoliolate lineage in the HHM region using a population-level genome resequencing dataset, incorporating evidence from morphology, distribution, and elevation. Comprehensive analyses confirmed multiple introgressions within Parthenocissus in a robust temporal-spatial framework. Around the HHM region, at least three hybridization hot spots were identified, one of which showed evidence of ongoing speciation reversal. We present a solid case study using an integrative methodological approach to investigate reticulate evolutionary history and its underlying mechanisms in plants. It demonstrates an example of speciation reversal through frequent hybridizations in the HHM region, which provides new perspectives on speciation dynamics in mountainous areas with strong topographic and environmental heterogeneity.

A Phylogenomic Assessment of Processes Underpinning Convergent Evolution in Open-Habitat Chats

Insights into the processes underpinning convergent evolution advance our understanding of the contributions of ancestral, introgressed, and novel genetic variation to phenotypic evolution. Phylogenomic analyses characterizing genome-wide gene tree heterogeneity can provide first clues about the extent of ILS and of introgression and thereby into the potential of these processes or (in their absence) the need to invoke novel mutations to underpin convergent evolution. Here, we were interested in understanding the processes involved in convergent evolution in open-habitat chats (wheatears of the genus Oenanthe and their relatives). To this end, based on whole-genome resequencing data from 50 taxa of 44 species, we established the species tree, characterized gene tree heterogeneity, and investigated the footprints of ILS and introgression within the latter. The species tree corroborates the pattern of abundant convergent evolution, especially in wheatears. The high levels of gene tree heterogeneity in wheatears are explained by ILS alone only for 30% of internal branches. For multiple branches with high gene tree heterogeneity, D-statistics and phylogenetic networks identified footprints of introgression. Finally, long branches without extensive ILS between clades sporting similar phenotypes provide suggestive evidence for the role of novel mutations in the evolution of these phenotypes. Together, our results suggest that convergent evolution in open-habitat chats involved diverse processes and highlight that phenotypic diversification is often complex and best depicted as a network of interacting lineages.

Origin and diversification of the saguaro cactus (Carnegiea gigantea): a within-species phylogenomic analysis

Reconstructing accurate historical relationships within a species poses numerous challenges, not least in many plant groups in which gene flow is high enough to extend well beyond species boundaries. Nonetheless, the extent of tree-like history within a species is an empirical question on which it is now possible to bring large amounts of genome sequence to bear. We assess phylogenetic structure across the geographic range of the saguaro cactus, an emblematic member of Cactaceae, a clade known for extensive hybridization and porous species boundaries. Using 200 Gb of whole genome resequencing data from 20 individuals sampled from 10 localities, we assembled two data sets comprising 150,000 biallelic single nucleotide polymorphisms (SNPs) from protein coding sequences. From these we inferred within-species trees and evaluated their significance and robustness using five qualitatively different inference methods. Despite the low sequence diversity, large census population sizes, and presence of wide-ranging pollen and seed dispersal agents, phylogenetic trees were well resolved and highly consistent across both data sets and all methods. We inferred that the most likely root, based on marginal likelihood comparisons, is to the east and south of the region of highest genetic diversity, which lies along the coast of the Gulf of California in Sonora, Mexico. Together with striking decreases in marginal likelihood found to the north, this supports hypotheses that saguaro's current range reflects post-glacial expansion from the refugia in the south of its range. We conclude with observations about practical and theoretical issues raised by phylogenomic data sets within species, in which SNP-based methods must be used rather than gene tree methods that are widely used when sequence divergence is higher. These include computational scalability, inference of gene flow, and proper assessment of statistical support in the presence of linkage effects.

Comparative plastomes of Carya species provide new insights into the plastomes evolution and maternal phylogeny of the genus

Carya, in the Juglandiodeae subfamily, is to a typical temperate-subtropical forest-tree genus for studying the phylogenetic evolution and intercontinental disjunction between eastern Asia (EA) and North America (NA). Species of the genus have high economic values worldwide for their high-quality wood and the rich healthy factors of their nuts. Although previous efforts based on multiple molecular markers or genome-wide SNPs supported the monophyly of Carya and its two EA and NA major subclades, the maternal phylogeny of Carya still need to be comprehensively evaluated. The variation of Carya plastome has never been thoroughly characterized. Here, we novelly present 19 newly generated plastomes of congeneric Carya species, including the recently rediscovered critically endangered C. poilanei. The overall assessment of plastomes revealed highly conservative in the general structures. Our results indicated that remarkable differences in several plastome features are highly consistent with the EA-NA disjunction and showed the relatively diverse matrilineal sources among EA Carya compared to NA Carya. The maternal phylogenies were conducted with different plastome regions and full-length plastome datasets from 30 plastomes, representing 26 species in six genera of Juglandoideae and Myrica rubra (as root). Six out of seven phylogenetic topologies strongly supported the previously reported relationships among genera of Juglandoideae and the two subclades of EA and NA Carya, but displayed significant incongruencies between species within the EA and NA subclades. The phylogenetic tree generated from full-length plastomes demonstrated the optimal topology and revealed significant geographical maternal relationships among Carya species, especially for EA Carya within overlapping distribution areas. The full-length plastome-based phylogenetic topology also strongly supported the taxonomic status of five controversial species as separate species of Carya. Historical and recent introgressive hybridization and plastid captures might contribute to plastome geographic patterns and inconsistencies between topologies built from different datasets, while incomplete lineage sorting could account for the discordance between maternal topology and the previous nuclear genome data-based phylogeny. Our findings highlight full-length plastomes as an ideal tool for exploring maternal relationships among the subclades of Carya, and potentially in other outcrossing perennial woody plants, for resolving plastome phylogenetic relationships.

Andean uplift, drainage basin formation, and the evolution of plants living in fast-flowing aquatic ecosystems in northern South America

Northern South America is a geologically dynamic and species-rich region. Fossil and stratigraphic data show that mountain uplift in the tropical Andes reconfigured river drainages. These landscape changes shaped the evolution of the flora in the region, yet the impacts on aquatic taxa have been overlooked. We explore the role of landscape change on the evolution of plants living strictly in rivers across drainage basins in northern South America by conducting population structure, phylogenetic inference, and divergence-dating analyses for two species in the genus Marathrum (Podostemaceae). Mountain uplift and drainage basin formation isolated populations of M. utile and M. foeniculaceum in northern South America and created barriers to gene flow across river drainages. Sympatric species hybridize and the hybrids show the phenotype of one parental line. We propose that the pattern of divergence of populations reflects the formation of river drainages, which was not complete until < 4.1 million yr ago (Ma). Our study provides a clear picture of the role of landscape change on the evolution of plants living strictly in rivers in northern South America. By shifting the focus to aquatic taxa, we provide a novel perspective on the processes shaping the evolution of the Neotropical flora.

Leaf form diversification in an ornamental heirloom tomato results from alterations in two different HOMEOBOX genes

Domesticated plants display diverse phenotypic traits. However, the influence of breeding effort on this phenotypic diversity remains unknown. Here, we demonstrate that a single nucleotide deletion in the homeobox motif of BIPINNATA, a BEL-LIKE HOMEODOMAIN gene, led to a highly complex leaf phenotype in an heirloom tomato (Solanum lycopersicum), Silvery Fir Tree (SiFT), which is used as a landscaping and ornamental plant. A comparative gene network analysis revealed that repression of SOLANIFOLIA, the ortholog of WUSCHEL RELATED HOMEOBOX 1, caused the narrow leaflet phenotype seen in SiFT. Comparative genomics indicated that the bip mutation in SiFT likely arose de novo and is unique to SiFT and not introgressed from other tomato genomes. These results provide new insights into the natural variation in phenotypic traits introduced into crops during improvement processes after domestication and establish homeobox genes as evolutionary hotspots.

Leaf form diversification in an ornamental heirloom tomato results from alterations in two different HOMEOBOX genes

Domesticated plants display diverse phenotypic traits. However, the influence of breeding effort on this phenotypic diversity remains unknown. Here, we demonstrate that a single nucleotide deletion in the homeobox motif of BIPINNATA, a BEL-LIKE HOMEODOMAIN gene, led to a highly complex leaf phenotype in an heirloom tomato (Solanum lycopersicum), Silvery Fir Tree (SiFT), which is used as a landscaping and ornamental plant. A comparative gene network analysis revealed that repression of SOLANIFOLIA, the ortholog of WUSCHEL RELATED HOMEOBOX 1, caused the narrow leaflet phenotype seen in SiFT. Comparative genomics indicated that the bip mutation in SiFT likely arose de novo and is unique to SiFT and not introgressed from other tomato genomes. These results provide new insights into the natural variation in phenotypic traits introduced into crops during improvement processes after domestication and establish homeobox genes as evolutionary hotspots.

Plotting for change: an analytical framework to aid decisions on which lineages are candidate species in phylogenomic species discovery

A recent study argued that coalescent-based models of species delimitation mostly delineate population structure, not species, and called for the validation of candidate species using biological information additional to the genetic information, such as phenotypic or ecological data. Here, we introduce a framework to interrogate genomic datasets and coalescent-based species trees for the presence of candidate species in situations where additional biological data are unavailable, unobtainable or uninformative. For de novo genomic studies of species boundaries, we propose six steps: (1) visualize genetic affinities among individuals to identify both discrete and admixed genetic groups from first principles and to hold aside individuals involved in contemporary admixture for independent consideration; (2) apply phylogenetic techniques to identify lineages; (3) assess diagnosability of those lineages as potential candidate species; (4) interpret the diagnosable lineages in a geographical context (sympatry, parapatry, allopatry); (5) assess significance of difference or trends in the context of sampling intensity; and (6) adopt a holistic approach to available evidence to inform decisions on species status in the difficult cases of allopatry. We adopt this approach to distinguish candidate species from within-species lineages for a widespread species complex of Australian freshwater fishes (Retropinna spp.). Our framework addresses two cornerstone issues in systematics that are often not discussed explicitly in genomic species discovery: diagnosability and how to determine it, and what criteria should be used to decide whether diagnosable lineages are conspecific or represent different species.

Taxonomic Uncertainty and the Anomaly Zone: Phylogenomics Disentangle a Rapid Radiation to Resolve Contentious Species (Gila robusta Complex) in the Colorado River

Species are indisputable units for biodiversity conservation, yet their delimitation is fraught with both conceptual and methodological difficulties. A classic example is the taxonomic controversy surrounding the Gila robusta complex in the lower Colorado River of southwestern North America. Nominal species designations were originally defined according to weakly diagnostic morphological differences, but these conflicted with subsequent genetic analyses. Given this ambiguity, the complex was re-defined as a single polytypic unit, with the proposed "threatened" status under the U.S. Endangered Species Act of two elements being withdrawn. Here we re-evaluated the status of the complex by utilizing dense spatial and genomic sampling (n = 387 and >22 k loci), coupled with SNP-based coalescent and polymorphism-aware phylogenetic models. In doing so, we found that all three species were indeed supported as evolutionarily independent lineages, despite widespread phylogenetic discordance. To juxtapose this discrepancy with previous studies, we first categorized those evolutionary mechanisms driving discordance, then tested (and subsequently rejected) prior hypotheses which argued phylogenetic discord in the complex was driven by the hybrid origin of Gila nigra. The inconsistent patterns of diversity we found within G. robusta were instead associated with rapid Plio-Pleistocene drainage evolution, with subsequent divergence within the "anomaly zone" of tree space producing ambiguities that served to confound prior studies. Our results not only support the resurrection of the three species as distinct entities but also offer an empirical example of how phylogenetic discordance can be categorized within other recalcitrant taxa, particularly when variation is primarily partitioned at the species level.

Integrating Sequence Capture and Restriction Site-Associated DNA Sequencing to Resolve Recent Radiations of Pelagic Seabirds

The diversification of modern birds has been shaped by a number of radiations. Rapid diversification events make reconstructing the evolutionary relationships among taxa challenging due to the convoluted effects of incomplete lineage sorting (ILS) and introgression. Phylogenomic data sets have the potential to detect patterns of phylogenetic incongruence, and to address their causes. However, the footprints of ILS and introgression on sequence data can vary between different phylogenomic markers at different phylogenetic scales depending on factors such as their evolutionary rates or their selection pressures. We show that combining phylogenomic markers that evolve at different rates, such as paired-end double-digest restriction site-associated DNA (PE-ddRAD) and ultraconserved elements (UCEs), allows a comprehensive exploration of the causes of phylogenetic discordance associated with short internodes at different timescales. We used thousands of UCE and PE-ddRAD markers to produce the first well-resolved phylogeny of shearwaters, a group of medium-sized pelagic seabirds that are among the most phylogenetically controversial and endangered bird groups. We found that phylogenomic conflict was mainly derived from high levels of ILS due to rapid speciation events. We also documented a case of introgression, despite the high philopatry of shearwaters to their breeding sites, which typically limits gene flow. We integrated state-of-the-art concatenated and coalescent-based approaches to expand on previous comparisons of UCE and RAD-Seq data sets for phylogenetics, divergence time estimation, and inference of introgression, and we propose a strategy to optimize RAD-Seq data for phylogenetic analyses. Our results highlight the usefulness of combining phylogenomic markers evolving at different rates to understand the causes of phylogenetic discordance at different timescales. [Aves; incomplete lineage sorting; introgression; PE-ddRAD-Seq; phylogenomics; radiations; shearwaters; UCEs.].

Assessing the fit of the multi-species network coalescent to multi-locus data

MOTIVATION

With growing genome-wide molecular datasets from next-generation sequencing, phylogenetic networks can be estimated using a variety of approaches. These phylogenetic networks include events like hybridization, gene flow or horizontal gene transfer explicitly. However, the most accurate network inference methods are computationally heavy. Methods that scale to larger datasets do not calculate a full likelihood, such that traditional likelihood-based tools for model selection are not applicable to decide how many past hybridization events best fit the data. We propose here a goodness-of-fit test to quantify the fit between data observed from genome-wide multi-locus data, and patterns expected under the multi-species coalescent model on a candidate phylogenetic network.

RESULTS

We identified weaknesses in the previously proposed TICR test, and proposed corrections. The performance of our new test was validated by simulations on real-world phylogenetic networks. Our test provides one of the first rigorous tools for model selection, to select the adequate network complexity for the data at hand. The test can also work for identifying poorly inferred areas on a network.

AVAILABILITY AND IMPLEMENTATION

Software for the goodness-of-fit test is available as a Julia package at https://github.com/cecileane/QuartetNetworkGoodnessFit.jl.

SUPPLEMENTARY INFORMATION

Supplementary data are available at Bioinformatics online.

Gene tree discord, simplex plots, and statistical tests under the coalescent

A simple graphical device, the simplex plot of quartet concordance factors, is introduced to aid in the exploration of a collection of gene trees on a common set of taxa. A single plot summarizes all gene tree discord, and allows for visual comparison to the expected discord from the multispecies coalescent model (MSC) of incomplete lineage sorting on a species tree. A formal statistical procedure is described that can quantify the deviation from expectation for each subset of four taxa, suggesting when the data is not in accord with the MSC, and thus that either gene tree inference error is substantial or a more complex model such as that on a network may be required. If the collection of gene trees is in accord with the MSC, the plots reveal when substantial incomplete lineage sorting is present. Applications to both simulated and empirical multilocus data sets illustrate the insights provided.

Introgression across evolutionary scales suggests reticulation contributes to Amazonian tree diversity

Hybridization has the potential to generate or homogenize biodiversity and is a particularly common phenomenon in plants, with an estimated 25% of plant species undergoing interspecific gene flow. However, hybridization in Amazonia's megadiverse tree flora was assumed to be extremely rare despite extensive sympatry between closely related species, and its role in diversification remains enigmatic because it has not yet been examined empirically. Using members of a dominant Amazonian tree family (Brownea, Fabaceae) as a model to address this knowledge gap, our study recovered extensive evidence of hybridization among multiple lineages across phylogenetic scales. More specifically, using targeted sequence capture our results uncovered several historical introgression events between Brownea lineages and indicated that gene tree incongruence in Brownea is best explained by reticulation, rather than solely by incomplete lineage sorting. Furthermore, investigation of recent hybridization using ~19,000 ddRAD loci recovered a high degree of shared variation between two Brownea species that co-occur in the Ecuadorian Amazon. Our analyses also showed that these sympatric lineages exhibit homogeneous rates of introgression among loci relative to the genome-wide average, implying a lack of selection against hybrid genotypes and persistent hybridization. Our results demonstrate that gene flow between multiple Amazonian tree species has occurred across temporal scales, and contrasts with the prevailing view of hybridization's rarity in Amazonia. Overall, our results provide novel evidence that reticulate evolution influenced diversification in part of the Amazonian tree flora, which is the most diverse on Earth.

Biased Gene Retention in the Face of Introgression Obscures Species Relationships

Phylogenomic analyses are recovering previously hidden histories of hybridization, revealing the genomic consequences of these events on the architecture of extant genomes. We applied phylogenomic techniques and several complementary statistical tests to show that introgressive hybridization appears to have occurred between close relatives of Arabidopsis, resulting in cytonuclear discordance and impacting our understanding of species relationships in the group. The composition of introgressed and retained genes indicates that selection against incompatible cytonuclear and nuclear-nuclear interactions likely acted during introgression, whereas linkage also contributed to genome composition through the retention of ancient haplotype blocks. We also applied divergence-based tests to determine the species branching order and distinguish donor from recipient lineages. Surprisingly, these analyses suggest that cytonuclear discordance arose via extensive nuclear, rather than cytoplasmic, introgression. If true, this would mean that most of the nuclear genome was displaced during introgression whereas only a small proportion of native alleles were retained.

Out of Sight, Out of Mind: Widespread Nuclear and Plastid-Nuclear Discordance in the Flowering Plant Genus Polemonium (Polemoniaceae) Suggests Widespread Historical Gene Flow Despite Limited Nuclear Signal

Phylogenomic data from a rapidly increasing number of studies provide new evidence for resolving relationships in recently radiated clades, but they also pose new challenges for inferring evolutionary histories. Most existing methods for reconstructing phylogenetic hypotheses rely solely on algorithms that only consider incomplete lineage sorting (ILS) as a cause of intra- or intergenomic discordance. Here, we utilize a variety of methods, including those to infer phylogenetic networks, to account for both ILS and introgression as a cause for nuclear and cytoplasmic-nuclear discordance using phylogenomic data from the recently radiated flowering plant genus Polemonium (Polemoniaceae), an ecologically diverse genus in Western North America with known and suspected gene flow between species. We find evidence for widespread discordance among nuclear loci that can be explained by both ILS and reticulate evolution in the evolutionary history of Polemonium. Furthermore, the histories of organellar genomes show strong discordance with the inferred species tree from the nuclear genome. Discordance between the nuclear and plastid genome is not completely explained by ILS, and only one case of discordance is explained by detected introgression events. Our results suggest that multiple processes have been involved in the evolutionary history of Polemonium and that the plastid genome does not accurately reflect species relationships. We discuss several potential causes for this cytoplasmic-nuclear discordance, which emerging evidence suggests is more widespread across the Tree of Life than previously thought. [Cyto-nuclear discordance, genomic discordance, phylogenetic networks, plastid capture, Polemoniaceae, Polemonium, reticulations.].

Leaf shape is a predictor of fruit quality and cultivar performance in tomato

Commercial tomato (Solanum lycopersicum) is one of the most widely grown vegetable crops worldwide. Heirloom tomatoes retain extensive genetic diversity and a considerable range of fruit quality and leaf morphological traits. Here the role of leaf morphology was investigated for its impact on fruit quality. Heirloom cultivars were grown in field conditions, and BRIX by yield (BY) and other traits were measured over a 14-wk period. The complex relationships among these morphological and physiological traits were evaluated using partial least-squares path modeling, and a consensus model was developed. Photosynthesis contributed strongly to vegetative biomass and sugar content of fruits but had a negative impact on yield. Conversely leaf shape, specifically rounder leaves, had a strong positive impact on both fruit sugar content and yield. Cultivars such as Stupice and Glacier, with very round leaves, had the highest performance in both fruit sugar and yield. Our model accurately predicted BY for two commercial cultivars using leaf shape data as input. This study revealed the importance of leaf shape to fruit quality in tomato, with rounder leaves having significantly improved fruit quality. This correlation was maintained across a range of diverse genetic backgrounds and shows the importance of leaf morphology in tomato crop improvement.

Divergence-Based Introgression Polarization

Introgressive hybridization results in the transfer of genetic material between species, often with fitness implications for the recipient species. The development of statistical methods for detecting the signatures of historical introgression in whole-genome data has been a major area of focus. Although existing techniques are able to identify the taxa that exchanged genes during introgression using a four-taxon system, most methods do not explicitly distinguish which taxon served as donor and which as recipient during introgression (i.e., polarization of introgression directionality). Existing methods that do polarize introgression are often only able to do so when there is a fifth taxon available and that taxon is sister to one of the taxa involved in introgression. Here, we present divergence-based introgression polarization (DIP), a method for polarizing introgression using patterns of sequence divergence across whole genomes, which operates in a four-taxon context. Thus, DIP can be applied to infer the directionality of introgression when additional taxa are not available. We use simulations to show that DIP can polarize introgression and identify potential sources of bias in the assignment of directionality, and we apply DIP to a well-described hominin introgression event.

Reticulate Evolution Helps Explain Apparent Homoplasy in Floral Biology and Pollination in Baobabs (Adansonia; Bombacoideae; Malvaceae)

Baobabs (Adansonia) are a cohesive group of tropical trees with a disjunct distribution in Australia, Madagascar, and continental Africa, and diverse flowers associated with two pollination modes. We used custom-targeted sequence capture in conjunction with new and existing phylogenetic comparative methods to explore the evolution of floral traits and pollination systems while allowing for reticulate evolution. Our analyses suggest that relationships in Adansonia are confounded by reticulation, with network inference methods supporting at least one reticulation event. The best supported hypothesis involves introgression between Adansonia rubrostipa and core Longitubae, both of which are hawkmoth pollinated with yellow/red flowers, but there is also some support for introgression between the African lineage and Malagasy Brevitubae, which are both mammal-pollinated with white flowers. New comparative methods for phylogenetic networks were developed that allow maximum-likelihood inference of ancestral states and were applied to study the apparent homoplasy in floral biology and pollination mode seen in Adansonia. This analysis supports a role for introgressive hybridization in morphological evolution even in a clade with highly divergent and geographically widespread species. Our new comparative methods for discrete traits on species networks are implemented in the software PhyloNetworks. [Comparative methods; Hyb-Seq; introgression; network inference; population trees; reticulate evolution; species tree inference; targeted sequence capture.]

Phylotranscriptomic Insights into the Diversification of Endothermic Thunnus Tunas

Birds, mammals, and certain fishes, including tunas, opahs and lamnid sharks, are endothermic, conserving internally generated, metabolic heat to maintain body or tissue temperatures above that of the environment. Bluefin tunas are commercially important fishes worldwide, and some populations are threatened. They are renowned for their endothermy, maintaining elevated temperatures of the oxidative locomotor muscle, viscera, brain and eyes, and occupying cold, productive high-latitude waters. Less cold-tolerant tunas, such as yellowfin tuna, by contrast, remain in warm-temperate to tropical waters year-round, reproducing more rapidly than most temperate bluefin tuna populations, providing resiliency in the face of large-scale industrial fisheries. Despite the importance of these traits to not only fisheries but also habitat utilization and responses to climate change, little is known of the genetic processes underlying the diversification of tunas. In collecting and analyzing sequence data across 29,556 genes, we found that parallel selection on standing genetic variation is associated with the evolution of endothermy in bluefin tunas. This includes two shared substitutions in genes encoding glycerol-3 phosphate dehydrogenase, an enzyme that contributes to thermogenesis in bumblebees and mammals, as well as four genes involved in the Krebs cycle, oxidative phosphorylation, β-oxidation, and superoxide removal. Using phylogenetic techniques, we further illustrate that the eight Thunnus species are genetically distinct, but found evidence of mitochondrial genome introgression across two species. Phylogeny-based metrics highlight conservation needs for some of these species.

Differences in flowering time maintain species boundaries in a continental radiation of Viburnum

PREMISE

We take an integrative approach in assessing how introgression and Pleistocene climate fluctuations have shaped the diversification of the core Lentago clade of Viburnum, a group of five interfertile species with broad areas of sympatry. We specifically tested whether flowering time plays a role in maintaining species isolation.

METHODS

RAD-seq data for 103 individuals were used to infer the species relationships and the genetic structure within each species. Flowering times were compared among species on the basis of historical flowering dates documented by herbarium specimens.

RESULTS

Within each species, we found a strong relationship between flowering date and latitude, such that southern populations flower earlier than northern ones. In areas of sympatry, the species flower in sequence rather than simultaneously, with flowering dates offset by ≥9 d for all species pairs. In two cases it appears that the offset in flowering times is an incidental consequence of adaptation to differing climates, but in the recently diverged sister species V. prunifolium and V. rufidulum, we find evidence that reinforcement led to reproductive character displacement. Long-term trends suggest that the two northern-most species are flowering earlier in response to recent climate change.

CONCLUSIONS

We argue that speciation in the Lentago clade has primarily occurred through ecological divergence of allopatric populations, but differences in flowering time were essential to maintain separation of incipient species when they came into secondary contact. This combination of factors may underlie diversification in many other plant clades.

Repeated evolution of a morphological novelty: a phylogenetic analysis of the inflated fruiting calyx in the Physalideae tribe (Solanaceae)

PREMISE OF THE STUDY

The evolution of novel fruit morphologies has been integral to the success of angiosperms. The inflated fruiting calyx, in which the balloon-like calyx swells to completely surround the fruit, has evolved repeatedly across angiosperms and is postulated to aid in protection and dispersal. We investigated the evolution of this trait in the tomatillos and their allies (Physalideae, Solanaceae).

METHODS

The Physalideae phylogeny was estimated using four regions (ITS, LEAFY, trnL-F, waxy) with maximum likelihood (ML) and Bayesian inference. Under the best-fitting ML model of trait evolution, we estimated ancestral states along with the numbers of gains and losses of fruiting calyx accrescence and inflation with Bayesian stochastic mapping. Also, phylogenetic signal in calyx morphology was examined with two metrics (parsimony score and Fritz and Purvis's D).

KEY RESULTS

Based on our well-resolved and densely sampled phylogeny, we infer that calyx evolution has proceeded in a stepwise and directional fashion, from non-accrescent to accrescent to inflated. In total, we inferred 24 gains of accrescence, 24 subsequent transitions to a fully inflated calyx, and only two reversals. Despite this lability, fruiting calyx accrescence and inflation showed strong phylogenetic signal.

CONCLUSIONS

Our phylogeny greatly improves the resolution of Physalideae and highlights the need for taxonomic work. The comparative analyses reveal that the inflated fruiting calyx has evolved many times and that the trajectory toward this phenotype is generally stepwise and irreversible. These results provide a strong foundation for studying the genetic and developmental mechanisms responsible for the repeated origins of this charismatic fruit trait.

A multilocus phylogeny of the fish genus Poeciliopsis: Solving taxonomic uncertainties and preliminary evidence of reticulation

The fish genus Poeciliopsis constitutes a valuable research system for evolutionary ecology, whose phylogenetic relationships have not been fully elucidated. We conducted a multilocus phylogenetic study of the genus based on seven nuclear and two mitochondrial loci with a thorough set of analytical approaches, that is, concatenated (also known as super-matrix), species trees, and phylogenetic networks. Although several relationships remain unresolved, the overall results uncovered phylogenetic affinities among several members of this genus. A population previously considered of undetermined taxonomic status could be unequivocally assigned to P. scarlli; revealing a relatively recent dispersal event across the Trans-Mexican Volcanic Belt (TMVB) or Pacific Ocean, which constitute a strong barrier to north-south dispersal of many terrestrial and freshwater taxa. The closest relatives of P. balsas, a species distributed south of the TMVB, are distributed in the north; representing an additional north-south split in the genus. An undescribed species of Poeciliopsis, with a highly restricted distribution (i.e., a short stretch of the Rio Concepcion; just south of the US-Mexico border), falls within the Leptorhaphis species complex. Our results are inconsistent with the hypothesis that this species originated by "breakdown" of an asexual hybrid lineage. On the other hand, network analyses suggest one or more possible cases of reticulation within the genus that require further evaluation with genome-wide marker representation and additional analytical tools. The most strongly supported case of reticulation occurred within the subgenus Aulophallus (restricted to Central America), and implies a hybrid origin for P. retropinna (i.e., between P. paucimaculata and P. elongata). We consider that P. balsas and P. new species are of conservation concern.

Modern Phylogenomics: Building Phylogenetic Trees Using the Multispecies Coalescent Model

The multispecies coalescent (MSC) model provides a compelling framework for building phylogenetic trees from multilocus DNA sequence data. The pure MSC is best thought of as a special case of so-called "multispecies network coalescent" models, in which gene flow is allowed among branches of the tree, whereas MSC methods assume there is no gene flow between diverging species. Early implementations of the MSC, such as "parsimony" or "democratic vote" approaches to combining information from multiple gene trees, as well as concatenation, in which DNA sequences from multiple gene trees are combined into a single "supergene," were quickly shown to be inconsistent in some regions of tree space, in so far as they converged on the incorrect species tree as more gene trees and sequence data were accumulated. The anomaly zone, a region of tree space in which the most frequent gene tree is different from the species tree, is one such region where many so-called "coalescent" methods are inconsistent. Second-generation implementations of the MSC employed Bayesian or likelihood models; these are consistent in all regions of gene tree space, but Bayesian methods in particular are incapable of handling the large phylogenomic data sets currently available. Two-step methods, such as MP-EST and ASTRAL, in which gene trees are first estimated and then combined to estimate an overarching species tree, are currently popular in part because they can handle large phylogenomic data sets. These methods are consistent in the anomaly zone but can sometimes provide inappropriate measures of tree support or apportion error and signal in the data inappropriately. MP-EST in particular employs a likelihood model which can be conveniently manipulated to perform statistical tests of competing species trees, incorporating the likelihood of the collected gene trees on each species tree in a likelihood ratio test. Such tests provide a useful alternative to the multilocus bootstrap, which only indirectly tests the appropriateness of competing species trees. We illustrate these tests and implementations of the MSC with examples and suggest that MSC methods are a useful class of models effectively using information from multiple loci to build phylogenetic trees.

A New Quartet-Based Statistical Method for Comparing Sets of Gene Trees Is Developed Using a Generalized Hoeffding Inequality

Extracting the strength of the tree signal that is encompassed by a collection of gene trees is an exceptionally challenging problem in phylogenomics. Often, this problem not only involves the construction of individual phylogenies based on different genes, which may be a difficult endeavor on its own, but is also exacerbated by many factors that create conflicts between the evolutionary histories of different gene families, such as duplications or losses of genes; hybridization events; incomplete lineage sorting; and horizontal gene transfer, the latter two play central roles in the evolution of eukaryotes and prokaryotes, respectively. In this work, we tackle the aforementioned problem by focusing on quartet trees, which are the most basic unit of information in the context of unrooted phylogenies. In the first part, we show how a theorem of Janson that generalizes the classical Hoeffding inequality can be used to develop a statistical test involving quartets. In the second part, we study real and simulated data using this theoretical advancement, thus demonstrating how the significance of the differences between sets of quartets can be assessed. Our results are particularly intriguing since they nonstandardly require the analysis of dependent random variables.

The Biogeography of Deep Time Phylogenetic Reticulation

Most phylogenies are typically represented as purely bifurcating. However, as genomic data have become more common in phylogenetic studies, it is not unusual to find reticulation among terminal lineages or among internal nodes (deep time reticulation; DTR). In these situations, gene flow must have happened in the same or adjacent geographic areas for these DTRs to have occurred and therefore biogeographic reconstruction should provide similar area estimates for parental nodes, provided extinction or dispersal has not eroded these patterns. We examine the phylogeny of the widely distributed New World kingsnakes (Lampropeltis), determine if DTR is present in this group, and estimate the ancestral area for reticulation. Importantly, we develop a new method that uses coalescent simulations in a machine learning framework to show conclusively that this phylogeny is best represented as reticulating at deeper time. Using joint probabilities of ancestral area reconstructions on the bifurcating parental lineages from the reticulating node, we show that this reticulation likely occurred in northwestern Mexico/southwestern US, and subsequently, led to the diversification of the Mexican kingsnakes. This region has been previously identified as an area important for understanding speciation and secondary contact with gene flow in snakes and other squamates. This research shows that phylogenetic reticulation is common, even in well-studied groups, and that the geographic scope of ancient hybridization is recoverable.

Plastomes of nine hornbeams and phylogenetic implications

Poor phylogenetic resolution and inconsistency of gene trees are major complications when attempting to construct trees of life for various groups of organisms. In this study, we addressed these issues in analyses of the genus Carpinus (hornbeams) of the Betulaceae. We assembled and annotated the chloroplast (cp) genomes (plastomes) of nine hornbeams representing main clades previously distinguished in this genus. All nine plastomes are highly conserved, with four regions, and about 158-160 kb long, including 121-123 genes. Phylogenetic analyses of whole plastome sequences, noncoding sequences, and the well-aligned coding genes resulted in high resolution of the sampled species in contrast to the failure based on a few cpDNA markers. Phylogenetic relationships in a few clades based only on the coding genes are slightly inconsistent with those based on the noncoding and total plastome datasets. Moreover, these plastome trees are highly incongruent with those based on bi-parentally inherited internal transcribed spacer (ITS) sequence variations. Such high inconsistencies suggest widespread occurrence of incomplete lineage sorting and hybrid introgression during diversification of these hornbeams.

Pleistocene glacial cycles drive isolation, gene flow and speciation in the high-elevation Andes

Mountain ranges are amongst the most species-rich habitats, with many large and rapid evolutionary radiations. The tempo and mode of diversification in these systems are key unanswered questions in evolutionary biology. Here we study the Andean Lupinus radiation to understand the processes driving very rapid diversification in montane systems. We use genomic and transcriptomic data of multiple species and populations, and apply phylogenomic and demographic analyses to test whether diversification proceeded without interspecific gene flow - as expected if Andean orogeny and geographic isolation were the main drivers of diversification - or if diversification was accompanied by gene flow, in which case other processes were probably involved. We uncover several episodes of gene flow between species, including very recent events likely to have been prompted by changes in habitat connectivity during Pleistocene glacial cycles. Furthermore, we find that gene flow between species was heterogeneously distributed across the genome. We argue that exceptionally fast diversification of Andean Lupinus was partly a result of Late Pleistocene glacial cycles, with associated cycles of expansion and contraction driving geographic isolation or secondary contact of species. Furthermore, heterogeneous gene flow across the genome suggests a role for selection and ecological speciation in rapid diversification in this system.

Comparative transcriptomics and genomic patterns of discordance in Capsiceae (Solanaceae)

The integration of genomics and phylogenetics allows new insight into the structure of gene tree discordance, the relationships among gene position, gene history, and rate of evolution, as well as the correspondence of gene function, positive selection, and gene ontology enrichment across lineages. We explore these issues using the tribe Capsiceae (Solanaceae), which is comprised of the genera Lycianthes and Capsicum (peppers). In combining the annotated genomes of Capsicum with newly sequenced transcriptomes of four species of Lycianthes and Capsicum, we develop phylogenies for 6747 genes, and construct a backbone species tree using both concordance and explicit phylogenetic network approaches. We quantify phylogenetic discordance among individual gene trees, measure their rates of synonymous and nonsynonymous substitution, and test whether they were positively selected along any branch of the phylogeny. We then map these genes onto the annotated Capsicum genome and test whether rates of evolution, gene history, and gene ontology vary significantly with gene position. We observed substantial discordance among gene trees. A bifurcating species tree placing Capsicum within a paraphyletic Lycianthes was supported over all phylogenetic networks. Rates of synonymous and nonsynonymous substitution varied 41-fold and 130-fold among genes, respectively, and were significantly lower in pericentromeric regions. We found that results of concordance tree analyses vary depending on the subset of genes used, and that genes within the pericentromeric regions only capture a portion of the observed discordance. We identified 787 genes that have been positively selected throughout the diversification history of Capsiceae, and discuss the importance of these genes as targets for investigation of economically important traits in the domesticated peppers.

Multi-locus phylogenetics, lineage sorting, and reticulation in Pinus subsection Australes

PREMISE OF THE STUDY

Both incomplete lineage sorting and reticulation have been proposed as causes of phylogenetic incongruence. Disentangling these factors may be most difficult in long-lived, wind-pollinated plants with large population sizes and weak reproductive barriers.

METHODS

We used solution hybridization for targeted enrichment and massive parallel sequencing to characterize low-copy-number nuclear genes and high-copy-number plastomes (Hyb-Seq) in 74 individuals of Pinus subsection Australes, a group of ~30 New World pine species of exceptional ecological and economic importance. We inferred relationships using methods that account for both incomplete lineage sorting and reticulation.

KEY RESULTS

Concatenation- and coalescent-based trees inferred from nuclear genes mainly agreed with one another, but they contradicted the plastid DNA tree in recovering the Attenuatae (the California closed-cone pines) and Oocarpae (the egg-cone pines of Mexico and Central America) as monophyletic and the Australes sensu stricto (the southern yellow pines) as paraphyletic to the Oocarpae. The plastid tree featured some relationships that were discordant with morphological and geographic evidence and species limits. Incorporating gene flow into the coalescent analyses better fit the data, but evidence supporting the hypothesis that hybridization explains the non-monophyly of the Attenuatae in the plastid tree was equivocal.

CONCLUSIONS

Our analyses document cytonuclear discordance in Pinus subsection Australes. We attribute this discordance to ancient and recent introgression and present a phylogenetic hypothesis in which mostly hierarchical relationships are overlain by gene flow.

Quartet Sampling distinguishes lack of support from conflicting support in the green plant tree of life

PREMISE OF THE STUDY

Phylogenetic support has been difficult to evaluate within the green plant tree of life partly due to a lack of specificity between conflicted versus poorly informed branches. As data sets continue to expand in both breadth and depth, new support measures are needed that are more efficient and informative.

METHODS

We describe the Quartet Sampling (QS) method, a quartet-based evaluation system that synthesizes several phylogenetic and genomic analytical approaches. QS characterizes discordance in large-sparse and genome-wide data sets, overcoming issues of alignment sparsity and distinguishing strong conflict from weak support. We tested QS with simulations and recent plant phylogenies inferred from variously sized data sets.

KEY RESULTS

QS scores demonstrated convergence with increasing replicates and were not strongly affected by branch depth. Patterns of QS support from different phylogenies led to a coherent understanding of ancestral branches defining key disagreements, including the relationships of Ginkgo to cycads, magnoliids to monocots and eudicots, and mosses to liverworts. The relationships of ANA-grade angiosperms (Amborella, Nymphaeales, Austrobaileyales), major monocot groups, bryophytes, and fern families are likely highly discordant in their evolutionary histories, rather than poorly informed. QS can also detect discordance due to introgression in phylogenomic data.

CONCLUSIONS

Quartet Sampling is an efficient synthesis of phylogenetic tests that offers more comprehensive and specific information on branch support than conventional measures. The QS method corroborates growing evidence that phylogenomic investigations that incorporate discordance testing are warranted when reconstructing complex evolutionary histories, in particular those surrounding ANA-grade, monocots, and nonvascular plants.

Reticulation, divergence, and the phylogeography-phylogenetics continuum

Phylogeography, and its extensions into comparative phylogeography, have their roots in the layering of gene trees across geography, a paradigm that was greatly facilitated by the nonrecombining, fast evolution provided by animal mtDNA. As phylogeography moves into the era of next-generation sequencing, the specter of reticulation at several levels-within loci and genomes in the form of recombination and across populations and species in the form of introgression-has raised its head with a prominence even greater than glimpsed during the nuclear gene PCR era. Here we explore the theme of reticulation in comparative phylogeography, speciation analysis, and phylogenomics, and ask how the centrality of gene trees has fared in the next-generation era. To frame these issues, we first provide a snapshot of multilocus phylogeographic studies across the Carpentarian Barrier, a prominent biogeographic barrier dividing faunas spanning the monsoon tropics in northern Australia. We find that divergence across this barrier is evident in most species, but is heterogeneous in time and demographic history, often reflecting the taxonomic distinctness of lineages spanning it. We then discuss a variety of forces generating reticulate patterns in phylogeography, including introgression, contact zones, and the potential selection-driven outliers on next-generation molecular markers. We emphasize the continued need for demographic models incorporating reticulation at the level of genomes and populations, and conclude that gene trees, whether explicit or implicit, should continue to play a role in the future of phylogeography.

Anchoring quartet-based phylogenetic distances and applications to species tree reconstruction

Background

Inferring species trees from gene trees using the coalescent-based summary methods has been the subject of much attention, yet new scalable and accurate methods are needed.

Results

We introduce DISTIQUE, a new statistically consistent summary method for inferring species trees from gene trees under the coalescent model. We generalize our results to arbitrary phylogenetic inference problems; we show that two arbitrarily chosen leaves, called anchors, can be used to estimate relative distances between all other pairs of leaves by inferring relevant quartet trees. This results in a family of distance-based tree inference methods, with running times ranging between quadratic to quartic in the number of leaves.

Conclusions

We show in simulated studies that DISTIQUE has comparable accuracy to leading coalescent-based summary methods and reduced running times.

Electronic supplementary material

The online version of this article (doi:10.1186/s12864-016-3098-z) contains supplementary material, which is available to authorized users.

Whole genome duplication in coast redwood (Sequoia sempervirens) and its implications for explaining the rarity of polyploidy in conifers

Polyploidy is common and an important evolutionary factor in most land plant lineages, but it is rare in gymnosperms. Coast redwood (Sequoia sempervirens) is one of just two polyploid conifer species and the only hexaploid. Evidence from fossil guard cell size suggests that polyploidy in Sequoia dates to the Eocene. Numerous hypotheses about the mechanism of polyploidy and parental genome donors have been proposed, based primarily on morphological and cytological data, but it remains unclear how Sequoia became polyploid and why this lineage overcame an apparent gymnosperm barrier to whole-genome duplication (WGD). We sequenced transcriptomes and used phylogenetic inference, Bayesian concordance analysis and paralog age distributions to resolve relationships among gene copies in hexaploid coast redwood and close relatives. Our data show that hexaploidy in coast redwood is best explained by autopolyploidy or, if there was allopolyploidy, it happened within the Californian redwood clade. We found that duplicate genes have more similar sequences than expected, given the age of the inferred polyploidization. Conflict between molecular and fossil estimates of WGD can be explained if diploidization occurred very slowly following polyploidization. We extrapolate from this to suggest that the rarity of polyploidy in gymnosperms may be due to slow diploidization in this clade.

Inferring Phylogenetic Networks with Maximum Pseudolikelihood under Incomplete Lineage Sorting

Phylogenetic networks are necessary to represent the tree of life expanded by edges to represent events such as horizontal gene transfers, hybridizations or gene flow. Not all species follow the paradigm of vertical inheritance of their genetic material. While a great deal of research has flourished into the inference of phylogenetic trees, statistical methods to infer phylogenetic networks are still limited and under development. The main disadvantage of existing methods is a lack of scalability. Here, we present a statistical method to infer phylogenetic networks from multi-locus genetic data in a pseudolikelihood framework. Our model accounts for incomplete lineage sorting through the coalescent model, and for horizontal inheritance of genes through reticulation nodes in the network. Computation of the pseudolikelihood is fast and simple, and it avoids the burdensome calculation of the full likelihood which can be intractable with many species. Moreover, estimation at the quartet-level has the added computational benefit that it is easily parallelizable. Simulation studies comparing our method to a full likelihood approach show that our pseudolikelihood approach is much faster without compromising accuracy. We applied our method to reconstruct the evolutionary relationships among swordtails and platyfishes (Xiphophorus: Poeciliidae), which is characterized by widespread hybridizations.

Phylogenetic networks display the evolutionary history of groups of individuals (species or populations) including reticulation events such as hybridization, horizontal gene transfer or migration. Here, we present a likelihood method to learn networks from molecular sequences at multiple genes. Our model accounts for several biological processes: mutations, incomplete lineage sorting of alleles in ancestral populations, and reticulations in the network. The likelihood is decomposed into 4-taxon subsets to make the analyses scale to many species and many genes. Our work makes it possible to learn large phylogenetic networks from large data sets, with a statistical approach and a biologically relevant model.