Population genomics from pool sequencing

A Molecular Survey of the Occurrence of Coffee Berry Disease Resistant Coffee Cultivars Near the Wild Gene Pool of Arabica Coffee in Its Region of Origin in Southwest Ethiopia

Cultivation of crops close to their wild relatives may jeopardise the integrity of wild genetic resources. Detecting cultivars among wild plants is necessary to characterise crop-wild gene flow, but can be challenging if cultivars and wild plants are phenotypically highly similar. Genomics tools can be used instead, but the selection of diagnostic loci for cultivar identification can be difficult if the wild and cultivated genepools are closely related. In Ethiopia, Arabica coffee cultivars resistant to coffee berry disease (CBD) occur near wild Coffea arabica plants and local landraces. However, the abundance and distribution of these cultivars across coffee sites remains unclear. Here, we present a new module of the SMAP package called SMAP relatedness pairwise to characterise pairwise genetic relationships between individuals based on haplotype calls and to identify diagnostic loci that distinguish (sets of) individuals from each other. Next, we estimate the relative abundance of CBD-resistant cultivars across 60 Ethiopian Arabica coffee sites using a genome-wide fingerprinting approach. We confirm the presence of these cultivars in around 75% of the coffee sites with a high agreement between a field survey and our DNA fingerprinting approach. At least 20 out of 60 sites with supposedly wild C. arabica individuals contain signatures of the cultivated genepool. Overall, we conclude that CBD-resistant cultivars are widespread in Ethiopian coffee sites. The development of SMAP relatedness pairwise opens opportunities to assess the distribution of coffee cultivars in other regions in Ethiopia and to apply similar screenings near wild relatives from other crops.

Dynamic Outlier Slicing Allows Broader Exploration of Adaptive Divergence: A Comparison of Individual Genome and Pool-Seq Data Linked to Humic Adaptation in Perch

How genetic variation contributes to adaptation at different environments is a central focus in evolutionary biology. However, most free-living species still lack a comprehensive understanding of the primary molecular mechanisms of adaptation. Here, we characterised the targets of selection associated with drastically different aquatic environments-humic and clear water-in the common freshwater fish, Eurasian perch (Perca fluviatilis). By using whole-genome sequencing (WGS) on a large population dataset (n = 42 populations) and analysing 873,788 SNPs, our primary aim was to uncover novel and confirm known footprints of selection. We compared individual and pooled WGS, and developed a novel approach, termed dynamic outlier slicing, to assess how the choice of outlier-calling stringency influences functional and Gene Ontology (GO) enrichment. By integrating genome-environment association (GEA) analysis with allele frequency-based approaches, we estimated composite selection signals (CSS) and identified 2679 outlier SNPs distributed across 324 genomic regions, involving 468 genes. Dynamic outlier slicing identified robust enrichment signals in five annotation categories (upstream, downstream, synonymous, 5'UTR and 3'UTR) highlighting the crucial role of regulatory elements in adaptive evolution. Furthermore, GO analyses revealed strong enrichment of molecular functions associated with gated channel activity, transmembrane transporter activity and ion channel activity, emphasising the importance of osmoregulation and ion balance maintenance. Our findings demonstrate that despite substantial random drift and divergence, WGS of high number of population pools enabled the identification of strong selection signals associated with adaptation to both humic and clear water environments, providing robust evidence of widespread adaptation. We anticipate that the dynamic outlier slicing method we developed will enable a more thorough exploration of adaptive divergence across a diverse range of species.

Genomic Diversity Illuminates the Environmental Adaptation of Drosophila suzukii

Biological invasions carry substantial practical and scientific importance and represent natural evolutionary experiments on contemporary timescales. Here, we investigated genomic diversity and environmental adaptation of the crop pest Drosophila suzukii using whole-genome sequencing data and environmental metadata for 29 population samples from its native and invasive range. Through a multifaceted analysis of this population genomic data, we increase our understanding of the D. suzukii genome, its diversity and its evolution, and we identify an appropriate genotype-environment association pipeline for our dataset. Using this approach, we detect genetic signals of local adaptation associated with nine distinct environmental factors related to altitude, wind speed, precipitation, temperature, and human land use. We uncover unique functional signatures for each environmental variable, such as the prevalence of cuticular genes associated with annual precipitation. We also infer biological commonalities in the adaptation to diverse selective pressures, particularly in terms of the apparent contribution of nervous system evolution to enriched processes (ranging from neuron development to circadian behavior) and to top genes associated with all nine environmental variables. Our findings therefore depict a finer-scale adaptive landscape underlying the rapid invasion success of this agronomically important species.

grenedalf: population genetic statistics for the next generation of pool sequencing

SUMMARY

Pool sequencing is an efficient method for capturing genome-wide allele frequencies from multiple individuals, with broad applications such as studying adaptation in Evolve-and-Resequence experiments, monitoring of genetic diversity in wild populations, and genotype-to-phenotype mapping. Here, we present grenedalf, a command line tool written in C++ that implements common population genetic statistics such as θ, Tajima's D, and FST for Pool sequencing. It is orders of magnitude faster than current tools, and is focused on providing usability and scalability, while also offering a plethora of input file formats and convenience options.

AVAILABILITY AND IMPLEMENTATION

grenedalf is published under the GPL-3, and freely available at github.com/lczech/grenedalf.

A new method for ecologists to estimate heterozygote excess and deficit for multi-locus gene families

The fixation index, F IS, has been a staple measure to detect selection, or departures from random mating in populations. However, current Next Generation Sequencing (NGS) cannot easily estimate F IS, in multi-locus gene families that contain multiple loci having similar or identical arrays of variant sequences of ≥1 kilobase (kb), which differ at multiple positions. In these families, high-quality short-read NGS data typically identify variants, but not the genomic location, which is required to calculate F IS (based on locus-specific observed and expected heterozygosity). Thus, to assess assortative mating, or selection on heterozygotes, from NGS of multi-locus gene families, we need a method that does not require knowledge of which variants are alleles at which locus in the genome. We developed such a method. Like F IS, our novel measure, 1 H IS, is based on the principle that positive assortative mating, or selection against heterozygotes, and some other processes reduce within-individual variability relative to the population. We demonstrate high accuracy of 1 H IS on a wide range of simulated scenarios and two datasets from natural populations of penguins and dolphins. 1 H IS is important because multi-locus gene families are often involved in assortative mating or selection on heterozygotes. 1 H IS is particularly useful for multi-locus gene families, such as toll-like receptors, the major histocompatibility complex in animals, homeobox genes in fungi and self-incompatibility genes in plants.

Association of a rapidly selected 4.3kb transposon-containing structural variation with a P450-based resistance to pyrethroids in the African malaria vector Anopheles funestus

Deciphering the evolutionary forces controlling insecticide resistance in malaria vectors remains a prerequisite to designing molecular tools to detect and assess resistance impact on control tools. Here, we demonstrate that a 4.3kb transposon-containing structural variation is associated with pyrethroid resistance in central/eastern African populations of the malaria vector Anopheles funestus. In this study, we analysed Pooled template sequencing data and direct sequencing to identify an insertion of 4.3kb containing a putative retro-transposon in the intergenic region of two P450s CYP6P5-CYP6P9b in mosquitoes of the malaria vector Anopheles funestus from Uganda. We then designed a PCR assay to track its spread temporally and regionally and decipher its role in insecticide resistance. The insertion originates in or near Uganda in East Africa, where it is fixed and has spread to high frequencies in the Central African nation of Cameroon but is still at low frequency in West Africa and absent in Southern Africa. A marked and rapid selection was observed with the 4.3kb-SV frequency increasing from 3% in 2014 to 98% in 2021 in Cameroon. A strong association was established between this SV and pyrethroid resistance in field populations and is reducing pyrethroid-only nets' efficacy. Genetic crosses and qRT-PCR revealed that this SV enhances the expression of CYP6P9a/b but not CYP6P5. Within this structural variant (SV), we identified putative binding sites for transcription factors associated with the regulation of detoxification genes. An inverse correlation was observed between the 4.3kb SV and malaria parasite infection, indicating that mosquitoes lacking the 4.3kb SV were more frequently infected compared to those possessing it. Our findings highlight the underexplored role and rapid spread of SVs in the evolution of insecticide resistance and provide additional tools for molecular surveillance of insecticide resistance.

Fine-scale adaptive divergence and population genetic structure of Aedes aegypti in Metropolitan Manila, Philippines

BACKGROUND

The adaptive divergence of Aedes aegypti populations to heterogeneous environments can be a driving force behind the recent expansion of their habitat distribution and outbreaks of dengue disease in urbanized areas. In this study, we investigated the population genomics of Ae. aegypti at a regional scale in Metropolitan Manila, Philippines.

METHODS

We used the Pool-Seq double digestion restriction-site association DNA sequencing (ddRAD-Seq) approach to generate a high number of single nucleotide polymorphisms (SNPs), with the aim to determine local adaptation and compare the population structure with 11 microsatellite markers. A total of 217 Ae. aegypti individuals from seven female and seven male populations collected from Metropolitan Manila were used in the assays.

RESULTS

We detected 65,473 SNPs across the populations, of which 76 were non-neutral SNPs. Of these non-neutral SNPs, the multivariate regression test associated 50 with eight landscape variables (e.g. open space, forest, etc.) and 29 with five climate variables (e.g. air temperature, humidity, etc.) (P-value range 0.005-0.045) in female and male populations separately. Male and female populations exhibited contrasting spatial divergence, with males exhibiting greater divergence than females, most likely reflecting the different dispersal abilities of male and female mosquitoes. In the comparative analysis of the same Ae. aegypti individuals, the pairwise FST values of 11 microsatellite markers were lower than those of the neutral SNPs, indicating that the neutral SNPs generated via pool ddRAD-Seq were more sensitive in terms of detecting genetic differences between populations at fine-spatial scales.

CONCLUSIONS

Overall, our study demonstrates the utility of pool ddRAD-Seq for examining genetic differences in Ae. aegypti populations in areas at fine-spatial scales that could inform vector control programs such as Wolbachia-infected mosquito mass-release programs. This in turn would provide information on mosquito population dispersal patterns and the potential barriers to mosquito movement within and around the release area. In addition, the potential of environmental adaptability observed in Ae. aegypti could help population control efforts.

A population genomics analysis of the Aotearoa New Zealand endemic rewarewa tree (Knightia excelsa)

Rewarewa (Knightia excelsa, Proteaceae) is a tree species endemic to Aotearoa New Zealand, with a natural distribution spanning Te Ika-a-Māui (North Island) and the top of Te Waipounamu (South Island). We used the pseudo-chromosome genome assembly of rewarewa as a reference and whole genome pooled sequencing from 35 populations sampled across Aotearoa New Zealand, including trees growing on Māori-owned land, to identify 1,443,255 single nucleotide polymorphisms (SNPs). Four genetic clusters located in the northern North Island (NNI), eastern North Island (NIE), western and southern North Island (NIWS), and the South Island (SI) were identified. Gene flow was revealed between the SI and NIE genetic clusters, plus bottleneck and contraction events within the genetic clusters since the mid-late Pleistocene, with divergence between North and South Island clusters estimated to have occurred ~115,000-230,000 years ago. Genotype environment analysis (GEA) was used to identify loci and genes linked with altitude, soil pH, soil carbon, slope, soil size, annual mean temperature, mean diurnal range, isothermality, annual precipitation, and precipitation seasonality. The location of the SNPs associated with these environmental variables was compared with the position of 52,192 gene-coding sequences that were predicted in the rewarewa genome using RNA sequencing. This new understanding of the genetic variation present in rewarewa and insights into the genetic control of adaptive traits will inform efforts to incorporate the species in restoration plantings and for marketing rewarewa honey based on provenance.

ViTAL: Vision TrAnsformer based Low coverage SARS-CoV-2 lineage assignment

MOTIVATION

Rapid spread of viral diseases such as Coronavirus disease 2019 (COVID-19) highlights an urgent need for efficient surveillance of virus mutation and transmission dynamics, which requires fast, inexpensive and accurate viral lineage assignment. The first two goals might be achieved through low-coverage whole-genome sequencing (LC-WGS) which enables rapid genome sequencing at scale and at reduced costs. Unfortunately, LC-WGS significantly diminishes the genomic details, rendering accurate lineage assignment very challenging.

RESULTS

We present ViTAL, a novel deep learning algorithm specifically designed to perform lineage assignment of low coverage-sequenced genomes. ViTAL utilizes a combination of MinHash for genomic feature extraction and Vision Transformer for fine-grain genome classification and lineage assignment. We show that ViTAL outperforms state-of-the-art tools across diverse coverage levels, reaching up to 87.7% lineage assignment accuracy at 1× coverage where state-of-the-art tools such as UShER and Kraken2 achieve the accuracy of 5.4% and 27.4% respectively. ViTAL achieves comparable accuracy results with up to 8× lower coverage than state-of-the-art tools. We explore ViTAL's ability to identify the lineages of novel genomes, i.e. genomes the Vision Transformer was not trained on. We show how ViTAL can be applied to preliminary phylogenetic placement of novel variants.

AVAILABILITY AND IMPLEMENTATION

The data underlying this article are available in https://github.com/zuherJahshan/vital and can be accessed with 10.5281/zenodo.10688110.

Experimental evolution of hybrid populations to identify Dobzhansky-Muller incompatibility loci

Epistatic interactions between loci that reduce fitness in interspecies hybrids, Dobzhansky-Muller incompatibilities (DMIs), contribute genetically to the inviability and infertility within hybrid populations. It remains a challenge, however, to identify the loci that contribute to DMIs as causes of reproductive isolation between species. Here, we assess through forward simulation the power of evolve-and-resequence (E&R) experimental evolution of hybrid populations to map DMI loci. We document conditions under which such a mapping strategy may be most feasible and demonstrate how mapping power is sensitive to biologically relevant parameters such as one-way versus two-way incompatibility type, selection strength, recombination rate, and dominance interactions. We also assess the influence of parameters under direct control of an experimenter, including duration of experimental evolution and number of replicate populations. We conclude that an E&R strategy for mapping DMI loci, and other cases of epistasis, can be a viable option under some circumstances for study systems with short generation times like Caenorhabditis nematodes.

Genomic Diversity Illuminates the Environmental Adaptation of Drosophila suzukii

Biological invasions carry substantial practical and scientific importance, and represent natural evolutionary experiments on contemporary timescales. Here, we investigated genomic diversity and environmental adaptation of the crop pest Drosophila suzukii using whole-genome sequencing data and environmental metadata for 29 population samples from its native and invasive range. Through a multifaceted analysis of this population genomic data, we increase our understanding of the D. suzukii genome, its diversity and its evolution, and we identify an appropriate genotype-environment association pipeline for our data set. Using this approach, we detect genetic signals of local adaptation associated with nine distinct environmental factors related to altitude, wind speed, precipitation, temperature, and human land use. We uncover unique functional signatures for each environmental variable, such as a prevalence of cuticular genes associated with annual precipitation. We also infer biological commonalities in the adaptation to diverse selective pressures, particularly in terms of the apparent contribution of nervous system evolution to enriched processes (ranging from neuron development to circadian behavior) and to top genes associated with all nine environmental variables. Our findings therefore depict a finer-scale adaptive landscape underlying the rapid invasion success of this agronomically important species.

Integrating Pool-seq uncertainties into demographic inference

Next-generation sequencing of pooled samples (Pool-seq) is a popular method to assess genome-wide diversity patterns in natural and experimental populations. However, Pool-seq is associated with specific sources of noise, such as unequal individual contributions. Consequently, using Pool-seq for the reconstruction of evolutionary history has remained underexplored. Here we describe a novel Approximate Bayesian Computation (ABC) method to infer demographic history, explicitly modelling Pool-seq sources of error. By jointly modelling Pool-seq data, demographic history and the effects of selection due to barrier loci, we obtain estimates of demographic history parameters accounting for technical errors associated with Pool-seq. Our ABC approach is computationally efficient as it relies on simulating subsets of loci (rather than the whole-genome) and on using relative summary statistics and relative model parameters. Our simulation study results indicate Pool-seq data allows distinction between general scenarios of ecotype formation (single versus parallel origin) and to infer relevant demographic parameters (e.g. effective sizes and split times). We exemplify the application of our method to Pool-seq data from the rocky-shore gastropod Littorina saxatilis, sampled on a narrow geographical scale at two Swedish locations where two ecotypes (Wave and Crab) are found. Our model choice and parameter estimates show that ecotypes formed before colonization of the two locations (i.e. single origin) and are maintained despite gene flow. These results indicate that demographic modelling and inference can be successful based on pool-sequencing using ABC, contributing to the development of suitable null models that allow for a better understanding of the genetic basis of divergent adaptation.

Genomes from historical Drosophila melanogaster specimens illuminate adaptive and demographic changes across more than 200 years of evolution

The ability to perform genomic sequencing on long-dead organisms is opening new frontiers in evolutionary research. These opportunities are especially notable in the case of museum collections, from which countless documented specimens may now be suitable for genomic analysis-if data of sufficient quality can be obtained. Here, we report 25 newly sequenced genomes from museum specimens of the model organism Drosophila melanogaster, including the oldest extant specimens of this species. By comparing historical samples ranging from the early 1800s to 1933 against modern-day genomes, we document evolution across thousands of generations, including time periods that encompass the species' initial occupation of northern Europe and an era of rapidly increasing human activity. We also find that the Lund, Sweden population underwent local genetic differentiation during the early 1800s to 1933 interval (potentially due to drift in a small population) but then became more similar to other European populations thereafter (potentially due to increased migration). Within each century-scale time period, our temporal sampling allows us to document compelling candidates for recent natural selection. In some cases, we gain insights regarding previously implicated selection candidates, such as ChKov1, for which our inferred timing of selection favors the hypothesis of antiviral resistance over insecticide resistance. Other candidates are novel, such as the circadian-related gene Ahcy, which yields a selection signal that rivals that of the DDT resistance gene Cyp6g1. These insights deepen our understanding of recent evolution in a model system, and highlight the potential of future museomic studies.

Novel insights into symbiont population structure: Globe-trotting avian feather mites contradict the specialist-generalist variation hypothesis

Researchers often examine symbiont host specificity as a species-level pattern, but it can also be key to understanding processes occurring at the population level, which are not as well understood. The specialist-generalist variation hypothesis (SGVH) attempts to explain how host specificity influences population-level processes, stating that single-host symbionts (specialists) exhibit stronger population genetic structure than multi-host symbionts (generalists) because of fewer opportunities for dispersal and more restricted gene flow between populations. However, this hypothesis has not been tested in systems with highly mobile hosts, in which population connectivity may vary temporally and spatially. To address this gap, we tested the SGVH on proctophyllodid feather mites found on migratory warblers (family Parulidae) with contrasting host specificities, Amerodectes protonotaria (a host specialist of Protonotaria citrea) and A. ischyros (a host generalist of 17 parulid species). We used a pooled-sequencing approach and a novel workflow to analyse genetic variants obtained from whole genome data. Both mite species exhibited fairly weak population structure overall, and contrary to predictions of the SGVH, the generalist was more strongly structured than the specialist. These results may suggest that specialists disperse more freely among conspecifics, whereas generalists sort according to geography. Furthermore, our results may reflect an unexpected period for mite transmission - during the nonbreeding season of migratory hosts - as mite population structure more closely reflects the distributions of hosts during the nonbreeding season. Our findings alter our current understanding of feather mite biology and highlight the potential for studies to explore factors driving symbiont diversification at multiple evolutionary scales.

Contrasting Evolutionary Patterns Between Sexual and Asexual Lineages in a Genomic Region Linked to Reproductive Mode Variation in the pea aphid

Although asexual lineages evolved from sexual lineages in many different taxa, the genetics of sex loss remains poorly understood. We addressed this issue in the pea aphid Acyrthosiphon pisum, whose natural populations encompass lineages performing cyclical parthenogenesis (CP) and producing one sexual generation per year, as well as obligate parthenogenetic (OP) lineages that can no longer produce sexual females but can still produce males. An SNP-based, whole-genome scan of CP and OP populations sequenced in pools (103 individuals from 6 populations) revealed that an X-linked region is associated with the variation in reproductive mode. This 840-kb region is highly divergent between CP and OP populations (FST = 34.9%), with >2,000 SNPs or short Indels showing a high degree of association with the phenotypic trait. In OP populations specifically, this region also shows reduced diversity and Tajima's D, consistent with the OP phenotype being a derived trait in aphids. Interestingly, the low genetic differentiation between CP and OP populations at the rest of the genome (FST = 2.5%) suggests gene flow between them. Males from OP lineages thus likely transmit their op allele to new genomic backgrounds. These genetic exchanges, combined with the selection of the OP and CP reproductive modes under different climates, probably contribute to the long-term persistence of the cp and op alleles.

Polygenic adaptation to overnutrition reveals a role for cholinergic signaling in longevity

Overnutrition by high-sugar (HS) feeding reduces both the lifespan and healthspan across taxa. Pressuring organisms to adapt to overnutrition can highlight genes and pathways important for the healthspan in stressful environments. We used an experimental evolution approach to adapt four replicate, outbred population pairs of Drosophila melanogaster to a HS or control diet. Sexes were separated and aged on either diet until mid-life, then mated to produce the next generation, allowing enrichment for protective alleles over time. All HS-selected populations increased their lifespan and were therefore used as a platform to compare allele frequencies and gene expression. Pathways functioning in the nervous system were overrepresented in the genomic data and showed evidence for parallel evolution, although very few genes were the same across replicates. Acetylcholine-related genes, including the muscarinic receptor mAChR-A, showed significant changes in allele frequency in multiple selected populations and differential expression on a HS diet. Using genetic and pharmacological approaches, we show that cholinergic signaling affects Drosophila feeding in a sugar-specific fashion. Together, these results suggest that adaptation produces changes in allele frequencies that benefit animals under conditions of overnutrition and that it is repeatable at the pathway level.

Deep haplotype analyses of target-site resistance locus ACCase in blackgrass enabled by pool-based amplicon sequencing

Rapid adaptation of weeds to herbicide applications in agriculture through resistance development is a widespread phenomenon. In particular, the grass Alopecurus myosuroides is an extremely problematic weed in cereal crops with the potential to manifest resistance in only a few generations. Target-site resistances (TSRs), with their strong phenotypic response, play an important role in this rapid adaptive response. Recently, using PacBio's long-read amplicon sequencing technology in hundreds of individuals, we were able to decipher the genomic context in which TSR mutations occur. However, sequencing individual amplicons are costly and time-consuming, thus impractical to implement for other resistance loci or applications. Alternatively, pool-based approaches overcome these limitations and provide reliable allele frequencies, although at the expense of not preserving haplotype information. In this proof-of-concept study, we sequenced with PacBio High Fidelity (HiFi) reads long-range amplicons (13.2 kb), encompassing the entire ACCase gene in pools of over 100 individuals, and resolved them into haplotypes using the clustering algorithm PacBio amplicon analysis (pbaa), a new application for pools in plants and other organisms. From these amplicon pools, we were able to recover most haplotypes from previously sequenced individuals of the same population. In addition, we analysed new pools from a Germany-wide collection of A. myosuroides populations and found that TSR mutations originating from soft sweeps of independent origin were common. Forward-in-time simulations indicate that TSR haplotypes will persist for decades even at relatively low frequencies and without selection, highlighting the importance of accurate measurement of TSR haplotype prevalence for weed management.

The redlegged earth mite draft genome provides new insights into pesticide resistance evolution and demography in its invasive Australian range

Genomic data provide valuable insights into pest management issues such as resistance evolution, historical patterns of pest invasions and ongoing population dynamics. We assembled the first reference genome for the redlegged earth mite, Halotydeus destructor (Tucker, 1925), to investigate adaptation to pesticide pressures and demography in its invasive Australian range using whole-genome pool-seq data from regionally distributed populations. Our reference genome comprises 132 autosomal contigs, with a total length of 48.90 Mb. We observed a large complex of ace genes, which has presumably evolved from a long history of organophosphate selection in H. destructor and may contribute towards organophosphate resistance through copy number variation, target-site mutations and structural variants. In the putative ancestral H. destructor ace gene, we identified three target-site mutations (G119S, A201S and F331Y) segregating in organophosphate-resistant populations. Additionally, we identified two new para sodium channel gene mutations (L925I and F1020Y) that may contribute to pyrethroid resistance. Regional structuring observed in population genomic analyses indicates that gene flow in H. destructor does not homogenize populations across large geographic distances. However, our demographic analyses were equivocal on the magnitude of gene flow; the short invasion history of H. destructor makes it difficult to distinguish scenarios of complete isolation vs. ongoing migration. Nonetheless, we identified clear signatures of reduced genetic diversity and smaller inferred effective population sizes in eastern vs. western populations, which is consistent with the stepping-stone invasion pathway of this pest in Australia. These new insights will inform development of diagnostic genetic markers of resistance, further investigation into the multifaceted organophosphate resistance mechanism and predictive modelling of resistance evolution and spread.

Genomic differentiation in Pacific cod using Pool-Seq

Patterns of genetic differentiation across the genome can provide insight into selective forces driving adaptation. We used pooled whole genome sequencing, gene annotation, and environmental covariates to evaluate patterns of genomic differentiation and to investigate mechanisms responsible for divergence among proximate Pacific cod (Gadus macrocephalus) populations from the Bering Sea and Aleutian Islands and more distant Washington Coast cod. Samples were taken from eight spawning locations, three of which were replicated to estimate consistency in allele frequency estimation. A kernel smoothing moving weighted average of relative divergence (F ST) identified 11 genomic islands of differentiation between the Aleutian Islands and Bering Sea samples. In some islands of differentiation, there was also elevated absolute divergence (d XY) and evidence for selection, despite proximity and potential for gene flow. Similar levels of absolute divergence (d XY) but roughly double the relative divergence (F ST) were observed between the distant Bering Sea and Washington Coast samples. Islands of differentiation were much smaller than the four large inversions among Atlantic cod ecotypes. Islands of differentiation between the Bering Sea and Aleutian Island were associated with SNPs from five vision system genes, which can be associated with feeding, predator avoidance, orientation, and socialization. We hypothesize that islands of differentiation between Pacific cod from the Bering Sea and Aleutian Islands provide evidence for adaptive differentiation despite gene flow in this commercially important marine species.

Genomic landscape of drug response reveals mediators of anthelmintic resistance

Like other pathogens, parasitic helminths can rapidly evolve resistance to drug treatment. Understanding the genetic basis of anthelmintic drug resistance in parasitic nematodes is key to tracking its spread and improving the efficacy and sustainability of parasite control. Here, we use an in vivo genetic cross between drug-susceptible and multi-drug-resistant strains of Haemonchus contortus in a natural host-parasite system to simultaneously map resistance loci for the three major classes of anthelmintics. This approach identifies new alleles for resistance to benzimidazoles and levamisole and implicates the transcription factor cky-1 in ivermectin resistance. This gene is within a locus under selection in ivermectin-resistant populations worldwide; expression analyses and functional validation using knockdown experiments support that cky-1 is associated with ivermectin survival. Our work demonstrates the feasibility of high-resolution forward genetics in a parasitic nematode and identifies variants for the development of molecular diagnostics to combat drug resistance in the field.

Genome resequencing clarifies phylogeny and reveals patterns of selection in the toxicogenomics model Pimephales promelas

Background

The fathead minnow (Pimephales promelas) is a model species for toxicological research. A high-quality genome reference sequence is available, and genomic methods are increasingly used in toxicological studies of the species. However, phylogenetic relationships within the genus remain incompletely known and little population-genomic data are available for fathead minnow despite the potential effects of genetic background on toxicological responses. On the other hand, a wealth of extant samples is stored in museum collections that in principle allow fine-scale analysis of contemporary and historical genetic variation.

Methods

Here we use short-read shotgun resequencing to investigate sequence variation among and within Pimephales species. At the genus level, our objectives were to resolve phylogenetic relationships and identify genes with signatures of positive diversifying selection. At the species level, our objective was to evaluate the utility of archived-sample resequencing for detecting selective sweeps within fathead minnow, applied to a population introduced to the San Juan River of the southwestern United States sometime prior to 1950.

Results

We recovered well-supported but discordant phylogenetic topologies for nuclear and mitochondrial sequences that we hypothesize arose from mitochondrial transfer among species. The nuclear tree supported bluntnose minnow (P. notatus) as sister to fathead minnow, with the slim minnow (P. tenellus) and bullhead minnow (P. vigilax) more closely related to each other. Using multiple methods, we identified 11 genes that have diversified under positive selection within the genus. Within the San Juan River population, we identified selective-sweep regions overlapping several sets of related genes, including both genes that encode the giant sarcomere protein titin and the two genes encoding the MTORC1 complex, a key metabolic regulator. We also observed elevated polymorphism and reduced differentation among populations (F_ST) in genomic regions containing certain immune-gene clusters, similar to what has been reported in other taxa. Collectively, our data clarify evolutionary relationships and selective pressures within the genus and establish museum archives as a fruitful resource for characterizing genomic variation. We anticipate that large-scale resequencing will enable the detection of genetic variants associated with environmental toxicants such as heavy metals, high salinity, estrogens, and agrichemicals, which could be exploited as efficient biomarkers of exposure in natural populations.

Identification of loci controlling timing of stem elongation in red clover using genotyping by sequencing of pooled phenotypic extremes

MAIN CONCLUSION

Through selective genotyping of pooled phenotypic extremes, we identified a number of loci and candidate genes putatively controlling timing of stem elongation in red clover. We have identified candidate genes controlling the timing of stem elongation prior to flowering in red clover (Trifolium pratense L.). This trait is of ecological and agronomic significance, as it affects fitness, competitivity, climate adaptation, forage and seed yield, and forage quality. We genotyped replicate pools of phenotypically extreme individuals (early and late-elongating) within cultivar Lea using genotyping-by-sequencing in pools (pool-GBS). After calling and filtering SNPs and GBS locus haplotype polymorphisms, we estimated allele frequencies and searched for markers with significantly different allele frequencies in the two phenotypic groups using BayeScan, an F_ST-based test utilizing replicate pools, and a test based on error variance of replicate pools. Of the three methods, BayeScan was the least stringent, and the error variance-based test the most stringent. Fifteen significant markers were identified in common by all three tests. The candidate genes flanking the markers include genes with potential roles in the vernalization, autonomous, and photoperiod regulation of floral transition, hormonal regulation of stem elongation, and cell growth. These results provide a first insight into the potential genes and mechanisms controlling transition to stem elongation in a perennial legume, which lays a foundation for further functional studies of the genetic determinants regulating this important trait.

Whole genome sequencing reveals footprints of adaptive genetic variation in populations of Eruca sativa

Populations of Eruca sativa (Brassicaceae) derived from arid and Mediterranean habitats exhibit ecotypic differentiation. Here, pooled DNA sequencing was used to assess adaptive genome differentiation in the two ecotypes. Differentiated SNP loci were scanned with the empirical FST outlier method and by correlating allele frequencies with environmental parameters. Genetic diversity values were relatively higher in the pooled arid genome, whereas the pooled Mediterranean genome exhibited stronger directional selection, indicating the impact of climatic conditions on genetic diversity. GO enrichment analysis categorized the annotated differentiated loci according to biological processes, revealing a large set of candidate genes related to abiotic and biotic stress responses. Allelic variation was detected in regulatory elements and coding regions (synonymous and non-synonymous mutations) of genes belonging to different transcription factors and phytohormone signaling, suggesting adaptation to both abiotic and biotic conditions. Furthermore, SNP mutations were also found in genic regions belonging to the synthesis of secondary metabolites, including aliphatic glucosinolates and their hydrolyzed bioactive compounds, among others. The results of this eco-genomic study demonstrate the role of divergent abiotic and biotic selection factors in evolutionary processes leading to adaptive ecotypic differentiation.

Population genomics of ancient and modern Trichuris trichiura

The neglected tropical disease trichuriasis is caused by the whipworm Trichuris trichiura, a soil-transmitted helminth that has infected humans for millennia. Today, T. trichiura infects as many as 500 million people, predominantly in communities with poor sanitary infrastructure enabling sustained faecal-oral transmission. Using whole-genome sequencing of geographically distributed worms collected from human and other primate hosts, together with ancient samples preserved in archaeologically-defined latrines and deposits dated up to one thousand years old, we present the first population genomics study of T. trichiura. We describe the continent-scale genetic structure between whipworms infecting humans and baboons relative to those infecting other primates. Admixture and population demographic analyses support a stepwise distribution of genetic variation that is highest in Uganda, consistent with an African origin and subsequent translocation with human migration. Finally, genome-wide analyses between human samples and between human and non-human primate samples reveal local regions of genetic differentiation between geographically distinct populations. These data provide insight into zoonotic reservoirs of human-infective T. trichiura and will support future efforts toward the implementation of genomic epidemiology of this globally important helminth.

Detecting Genetic Variation of Colonizing Streptococcus agalactiae Genomes in Humans: A Precision Protocol

Deciphering the genotypic diversity of within-individual pathogens and verifying the evolutionary model can help elucidate resistant genotypes, virulent subpopulations, and the mechanism of opportunistic pathogenicity. However, observed polymorphic mutations (PMs) are rare and difficult to be detected in the “dominant-lineage” model of bacterial infection due to the low frequency. The four pooled group B Streptococcus (GBS) samples were collected from the genital tracts of healthy pregnant women, and the pooled samples and the isogenic controls were genomically sequenced. Using the PMcalling program, we detected the PMs in samples and compared the results between two technical duplicates, GBS-M001T and GBS-M001C. Tested with simulated datasets, the PMcalling program showed high sensitivity especially in low-frequency PMs and reasonable specificity. The genomic sequence data from pooled samples of GBS colonizing carrier pregnant women were analyzed, and few high-frequency PMs and some low-frequency PMs were discovered, indicating a dominant-lineage evolution model. The PMs mainly were nonsynonymous and enriched in quorum sensing, glycolysis/gluconeogenesis, ATP-binding cassette (ABC) transporters, etc., suggesting antimicrobial or environmental selective pressure. The re-analysis of the published Burkholderia dolosa data showed a diverse-community model, and only a few low-frequency PMs were shared between different individuals. Genes of general control non-repressible 5-related N-acetyltransferases family, major facilitator superfamily (MFS) transporter, and ABC transporter were positive selection candidates. Our findings indicate an unreported nature of the dominant-lineage model of GBS colonization in healthy women, and a formerly not observed mutation pool in a colonized microbial community, possibly maintained by selection pressure.

The ectodysplasin-A receptor is a candidate gene for lateral plate number variation in stickleback fish

Variation in lateral plating in stickleback fish represents a classical example of rapid and parallel adaptation in morphology. The underlying genetic architecture involves polymorphism at the ectodysplasin-A gene (EDA). However, lateral plate number is influenced by additional loci that remain poorly characterized. Here, we search for such loci by performing genome-wide differentiation mapping based on pooled whole-genome sequence data from a European stickleback population variable in the extent of lateral plating, while tightly controlling for the phenotypic effect of EDA. This suggests a new candidate locus, the EDA receptor gene (EDAR), for which additional support is obtained by individual-level targeted Sanger sequencing and by comparing allele frequencies among natural populations. Overall, our study illustrates the power of pooled whole-genome sequencing for searching phenotypically relevant loci and opens opportunities for exploring the population genetics and ecological significance of a new candidate locus for stickleback armor evolution.

Genetic composition and diversity of Arabica coffee in the crop's center of origin and its impact on four major fungal diseases

Conventional wisdom states that genetic variation reduces disease levels in plant populations. Nevertheless, crop species have been subject to a gradual loss of genetic variation through selection for specific traits during breeding, thereby increasing their vulnerability to biotic stresses such as pathogens. We explored how genetic variation in Arabica coffee sites in southwestern Ethiopia was related to the incidence of four major fungal diseases. Sixty sites were selected along a gradient of management intensity, ranging from nearly wild to intensively managed coffee stands. We used genotyping-by-sequencing of pooled leaf samples (pool-GBS) derived from 16 individual coffee shrubs in each of the sixty sites to assess the variation in genetic composition (multivariate: reference allele frequency) and genetic diversity (univariate: mean expected heterozygosity) between sites. We found that genetic composition had a clear spatial pattern and that genetic diversity was higher in less managed sites. The incidence of the four fungal diseases was related to the genetic composition of the coffee stands, but in a specific way for each disease. In contrast, genetic diversity was only related to the within-site variation of coffee berry disease, but not to the mean incidence of any of the four diseases across sites. Given that fungal diseases are major challenges of Arabica coffee in its native range, our findings that genetic composition of coffee sites impacted the major fungal diseases may serve as baseline information to study the molecular basis of disease resistance in coffee. Overall, our study illustrates the need to consider both host genetic composition and genetic diversity when investigating the genetic basis for variation in disease levels.

Genome-wide analysis of the response to ivermectin treatment by a Swedish field population of Haemonchus contortus

Haemonchus contortus is a pathogenic gastrointestinal nematode of small ruminants and, in part due to its capacity to develop resistance to drugs, contributes to significant losses in the animal production sector worldwide. Despite decades of research, comparatively little is known about the specific mechanism(s) driving resistance to drugs such as ivermectin in this species. Here we describe a genome-wide approach to detect evidence of selection by ivermectin treatment in a field population of H. contortus from Sweden, using parasites sampled from the same animals before and seven days after ivermectin exposure followed by whole-genome sequencing. Despite an 89% reduction in parasites recovered after treatment measured by the fecal egg count reduction test, the surviving population was highly genetically similar to the population before treatment, suggesting that resistance has likely evolved over time and that resistance alleles are present on diverse haplotypes. Pairwise gene and SNP frequency comparisons indicated the highest degree of differentiation was found at the terminal end of chromosome 4, whereas the most striking difference in nucleotide diversity was observed in a region on chromosome 5 previously reported to harbor a major quantitative trait locus involved in ivermectin resistance. These data provide novel insight into the genome-wide effect of ivermectin selection in a field population as well as confirm the importance of the previously established quantitative trait locus in the development of resistance to ivermectin.

Genome-wide patterns of genetic diversity, population structure and demographic history in mānuka (Leptospermum scoparium) growing on indigenous Māori land

Leptospermum scoparium J. R. Forst et G. Forst, known as mānuka by Māori, the indigenous people of Aotearoa (New Zealand), is a culturally and economically significant shrub species, native to New Zealand and Australia. Chemical, morphological and phylogenetic studies have indicated geographical variation of mānuka across its range in New Zealand, and genetic differentiation between New Zealand and Australia. We used pooled whole genome re-sequencing of 76 L. scoparium and outgroup populations from New Zealand and Australia to compile a dataset totalling ~2.5 million SNPs. We explored the genetic structure and relatedness of L. scoparium across New Zealand, and between populations in New Zealand and Australia, as well as the complex demographic history of this species. Our population genomic investigation suggests there are five geographically distinct mānuka gene pools within New Zealand, with evidence of gene flow occurring between these pools. Demographic modelling suggests three of these gene pools have undergone expansion events, whilst the evolutionary histories of the remaining two have been subjected to contractions. Furthermore, mānuka populations in New Zealand are genetically distinct from populations in Australia, with coalescent modelling suggesting these two clades diverged ~9-12 million years ago. We discuss the evolutionary history of this species and the benefits of using pool-seq for such studies. Our research will support the management and conservation of mānuka by landowners, particularly Māori, and the development of a provenance story for the branding of mānuka based products.

Select and Sequence of a Segregating Sugar Beet Population Provides Genomic Perspective of Host Resistance to Seedling Rhizoctonia solani Infection

Understanding the genetic basis of polygenic traits is a major challenge in agricultural species, especially in non-model systems. Select and sequence (SnS) experiments carried out within existing breeding programs provide a means to simultaneously identify the genomic background of a trait while improving the mean phenotype for a population. Using pooled whole genome sequencing (WGS) of selected and unselected bulks derived from a synthetic outcrossing sugar beet population EL57 (PI 663212), which segregates for seedling rhizoctonia resistance, we identified a putative genomic background involved in conditioning a resistance phenotype. Population genomic parameters were estimated to measure fixation (He), genome divergence (F ST ), and allele frequency changes between bulks (DeltaAF). We report on the genome wide patterns of variation resulting from selection and highlight specific genomic features associated with resistance. Expected heterozygosity (He) showed an increased level of fixation in the resistant bulk, indicating a greater selection pressure was applied. In total, 1,311 biallelic loci were detected as significant FST outliers (p < 0.01) in comparisons between the resistant and susceptible bulks. These loci were detected in 206 regions along the chromosomes and contained 275 genes. We estimated changes in allele frequency between bulks resulting from selection for resistance by leveraging the allele frequencies of an unselected bulk. DeltaAF was a more stringent test of selection and recovered 186 significant loci, representing 32 genes, all of which were also detected using FST. Estimates of population genetic parameters and statistical significance were visualized with respect to the EL10.2 physical map and produced a candidate gene list that was enriched for function in cell wall metabolism and plant disease resistance, including pathogen perception, signal transduction, and pathogen response. Specific variation associated with these genes was also reported and represents genetic markers for validation and prediction of resistance to Rhizoctonia. Select and sequence experiments offer a means to characterize the genetic base of sugar beet, inform selection within breeding programs, and prioritize candidate variation for functional studies.

A Population Genomic Assessment of Three Decades of Evolution in a Natural Drosophila Population

Population genetics seeks to illuminate the forces shaping genetic variation, often based on a single snapshot of genomic variation. However, utilizing multiple sampling times to study changes in allele frequencies can help clarify the relative roles of neutral and non-neutral forces on short time scales. This study compares whole-genome sequence variation of recently collected natural population samples of Drosophila melanogaster against a collection made approximately 35 years prior from the same locality—encompassing roughly 500 generations of evolution. The allele frequency changes between these time points would suggest a relatively small local effective population size on the order of 10,000, significantly smaller than the global effective population size of the species. Some loci display stronger allele frequency changes than would be expected anywhere in the genome under neutrality—most notably the tandem paralogs Cyp6a17 and Cyp6a23, which are impacted by structural variation associated with resistance to pyrethroid insecticides. We find a genome-wide excess of outliers for high genetic differentiation between old and new samples, but a larger number of adaptation targets may have affected SNP-level differentiation versus window differentiation. We also find evidence for strengthening latitudinal allele frequency clines: northern-associated alleles have increased in frequency by an average of nearly 2.5% at SNPs previously identified as clinal outliers, but no such pattern is observed at random SNPs. This project underscores the scientific potential of using multiple sampling time points to investigate how evolution operates in natural populations, by quantifying how genetic variation has changed over ecologically relevant timescales.

Drosophila Evolution over Space and Time (DEST): A New Population Genomics Resource

Drosophila melanogaster is a leading model in population genetics and genomics, and a growing number of whole-genome data sets from natural populations of this species have been published over the last years. A major challenge is the integration of disparate data sets, often generated using different sequencing technologies and bioinformatic pipelines, which hampers our ability to address questions about the evolution of this species. Here we address these issues by developing a bioinformatics pipeline that maps pooled sequencing (Pool-Seq) reads from D. melanogaster to a hologenome consisting of fly and symbiont genomes and estimates allele frequencies using either a heuristic (PoolSNP) or a probabilistic variant caller (SNAPE-pooled). We use this pipeline to generate the largest data repository of genomic data available for D. melanogaster to date, encompassing 271 previously published and unpublished population samples from over 100 locations in >20 countries on four continents. Several of these locations have been sampled at different seasons across multiple years. This data set, which we call Drosophila Evolution over Space and Time (DEST), is coupled with sampling and environmental metadata. A web-based genome browser and web portal provide easy access to the SNP data set. We further provide guidelines on how to use Pool-Seq data for model-based demographic inference. Our aim is to provide this scalable platform as a community resource which can be easily extended via future efforts for an even more extensive cosmopolitan data set. Our resource will enable population geneticists to analyze spatiotemporal genetic patterns and evolutionary dynamics of D. melanogaster populations in unprecedented detail.

ƒ-statistics estimation and admixture graph construction with Pool-Seq or allele count data using the R package poolfstat

By capturing various patterns of the structuring of genetic variation across populations, f -statistics have proved highly effective for the inference of demographic history. Such statistics are defined as covariance of SNP allele frequency differences among sets of populations without requiring haplotype information and are hence particularly relevant for the analysis of pooled sequencing (Pool-Seq) data. We here propose a reinterpretation of the F (and D) parameters in terms of probability of gene identity and derive from this unified definition unbiased estimators for both Pool-Seq data and standard allele count data obtained from individual genotypes. We implemented these estimators in a new version of the R package poolfstat, which now includes a wide range of inference methods: (i) three-population test of admixture; (ii) four-population test of treeness; (iii) F4-ratio estimation of admixture rates; and (iv) fitting, visualization and (semi-automatic) construction of admixture graphs. A comprehensive evaluation of the methods implemented in poolfstat on both simulated Pool-Seq (with various sequencing coverages and error rates) and allele count data confirmed the accuracy of these approaches, even for the most cost-effective Pool-Seq design involving relatively low sequencing coverages. We further analyzed a real Pool-Seq data made of 14 populations of the invasive species Drosophila suzukii which allowed refining both the demographic history of native populations and the invasion routes followed by this emblematic pest. Our new package poolfstat provides the community with a user-friendly and efficient all-in-one tool to unravel complex population genetic histories from large-size Pool-Seq or allele count SNP data.

Chromosome-level genome assembly reveals genomic architecture of northern range expansion in the mountain pine beetle, Dendroctonus ponderosae Hopkins (Coleoptera: Curculionidae)

Genome sequencing methods and assembly tools have improved dramatically since the 2013 publication of draft genome assemblies for the mountain pine beetle, Dendroctonus ponderosae Hopkins (Coleoptera: Curculionidae). We conducted proximity ligation library sequencing and scaffolding to improve contiguity, and then used linkage mapping and recent bioinformatic tools for correction and further improvement. The new assemblies have dramatically improved contiguity and gaps compared to the originals: N50 values increased 26- to 36-fold, and the number of gaps were reduced by half. Ninety per cent of the content of the assemblies is now contained in 12 and 11 scaffolds for the female and male assemblies, respectively. Based on linkage mapping information, the 12 largest scaffolds in both assemblies represent all 11 autosomal chromosomes and the neo-X chromosome. These assemblies now have nearly chromosome-sized scaffolds and will be instrumental for studying genomic architecture, chromosome evolution, population genomics, functional genomics, and adaptation in this and other pest insects. We also identified regions in two chromosomes, including the ancestral-X portion of the neo-X chromosome, with elevated differentiation between northern and southern Canadian populations.

Unique Genomic and Phenotypic Responses to Extreme and Variable pH Conditions in Purple Urchin Larvae

Environmental variation experienced by a species across space and time can promote the maintenance of genetic diversity that may be adaptive in future global change conditions. Selection experiments have shown that purple sea urchin, Strongylocentrotus purpuratus, populations have adaptive genetic variation for surviving pH conditions at the "edge" (pH 7.5) of conditions experienced in nature. However, little is known about whether populations have genetic variation for surviving low-pH events beyond those currently experienced in nature or how variation in pH conditions affects organismal and genetic responses. Here, we quantified survival, growth, and allele frequency shifts in experimentally selected developing purple sea urchin larvae in static and variable conditions at three pH levels: pH 8.1 (control), pH 7.5 (edge-of-range), and pH 7.0 (extreme). Variable treatments recovered body size relative to static treatments, but resulted in higher mortality, suggesting a potential tradeoff between survival and growth under pH stress. However, within each pH level, allele frequency changes were overlapping between static and variable conditions, suggesting a shared genetic basis underlying survival to mean pH regardless of variability. In contrast, genetic responses to pH 7.5 (edge) versus pH 7.0 (extreme) conditions were distinct, indicating a unique genetic basis of survival. In addition, loci under selection were more likely to be in exonic regions than regulatory, indicating that selection targeted protein-coding variation. Loci under selection in variable pH 7.5 conditions, more similar to conditions periodically experienced in nature, performed functions related to lipid biosynthesis and metabolism, while loci under selection in static pH 7.0 conditions performed functions related to transmembrane and mitochondrial processes. While these results are promising in that purple sea urchin populations possess genetic variation for surviving extreme pH conditions not currently experienced in nature, they caution that increased acidification does not result in a linear response but elicits unique physiological stresses and survival mechanisms.

An overview of current population genomics methods for the analysis of whole-genome resequencing data in eukaryotes

Characterizing the population history of a species and identifying loci underlying local adaptation is crucial in functional ecology, evolutionary biology, conservation and agronomy. The constant improvement of high-throughput sequencing techniques has facilitated the production of whole genome data in a wide range of species. Population genomics now provides tools to better integrate selection into a historical framework, and take into account selection when reconstructing demographic history. However, this improvement has come with a profusion of analytical tools that can confuse and discourage users. Such confusion limits the amount of information effectively retrieved from complex genomic data sets, and impairs the diffusion of the most recent analytical tools into fields such as conservation biology. It may also lead to redundancy among methods. To address these isssues, we propose an overview of more than 100 state-of-the-art methods that can deal with whole genome data. We summarize the strategies they use to infer demographic history and selection, and discuss some of their limitations. A website listing these methods is available at www.methodspopgen.com.

Benchmarking the performance of Pool-seq SNP callers using simulated and real sequencing data

Population genomics is a fast-developing discipline with promising applications in a growing number of life sciences fields. Advances in sequencing technologies and bioinformatics tools allow population genomics to exploit genome-wide information to identify the molecular variants underlying traits of interest and the evolutionary forces that modulate these variants through space and time. However, the cost of genomic analyses of multiple populations is still too high to address them through individual genome sequencing. Pooling individuals for sequencing can be a more effective strategy in Single Nucleotide Polymorphism (SNP) detection and allele frequency estimation because of a higher total coverage. However, compared to individual sequencing, SNP calling from pools has the additional difficulty of distinguishing rare variants from sequencing errors, which is often avoided by establishing a minimum threshold allele frequency for the analysis. Finding an optimal balance between minimizing information loss and reducing sequencing costs is essential to ensure the success of population genomics studies. Here, we have benchmarked the performance of SNP callers for Pool-seq data, based on different approaches, under different conditions, and using computer simulations and real data. We found that SNP callers performance varied for allele frequencies up to 0.35. We also found that SNP callers based on Bayesian (SNAPE-pooled) or maximum likelihood (MAPGD) approaches outperform the two heuristic callers tested (VarScan and PoolSNP), in terms of the balance between sensitivity and FDR both in simulated and sequencing data. Our results will help inform the selection of the most appropriate SNP caller not only for large-scale population studies but also in cases where the Pool-seq strategy is the only option, such as in metagenomic or polyploid studies.

Growth genes are implicated in the evolutionary divergence of sympatric piscivorous and insectivorous rainbow trout (Oncorhynchus mykiss)

Background

Identifying ecologically significant phenotypic traits and the genomic mechanisms that underly them are crucial steps in understanding traits associated with population divergence. We used genome-wide data to identify genomic regions associated with key traits that distinguish two ecomorphs of rainbow trout (Oncorhynchus mykiss)—insectivores and piscivores—that coexist for the non-breeding portion of the year in Kootenay Lake, southeastern British Columbia. “Gerrards” are large-bodied, rapidly growing piscivores with high metabolic rates that spawn north of Kootenay Lake in the Lardeau River, in contrast to the insectivorous populations that are on average smaller in body size, with lower growth and metabolic rates, mainly forage on aquatic insects, and spawn in tributaries immediately surrounding Kootenay Lake. We used pool-seq data representing ~ 60% of the genome and 80 fish per population to assess the level of genomic divergence between ecomorphs and to identify and interrogate loci that may play functional or selective roles in their divergence.

Results

Genomic divergence was high between sympatric insectivores and piscivores (\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$F_{\text{ST}}$$\end{document}FST = 0.188), and in fact higher than between insectivorous populations from Kootenay Lake and the Blackwater River (\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$F_{\text{ST}}$$\end{document}FST = 0.159) that are > 500 km apart. A window-based \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$F_{\text{ST}}$$\end{document}FST analysis did not reveal “islands” of genomic differentiation; however, the window with highest \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$F_{\text{ST}}$$\end{document}FST estimate did include a gene associated with insulin secretion. Although we explored the use of the “Local score” approach to identify genomic outlier regions, this method was ultimately not used because simulations revealed a high false discovery rate (~ 20%). Gene ontology (GO) analysis identified several growth processes as enriched in genes occurring in the ~ 200 most divergent genomic windows, indicating many loci of small effect involved in growth and growth-related metabolic processes are associated with the divergence of these ecomorphs.

Conclusion

Our results reveal a high degree of genomic differentiation between piscivorous and insectivorous populations and indicate that the large body piscivorous phenotype is likely not due to one or a few loci of large effect. Rather, the piscivore phenotype may be controlled by several loci of small effect, thus highlighting the power of whole-genome resequencing in identifying genomic regions underlying population-level phenotypic divergences.

Supplementary Information

The online version contains supplementary material available at 10.1186/s12862-021-01795-9.

Drivers of linkage disequilibrium across a species' geographic range

While linkage disequilibrium (LD) is an important parameter in genetics and evolutionary biology, the drivers of LD remain elusive. Using whole-genome sequences from across a species’ range, we assessed the impact of demographic history and mating system on LD. Both range expansion and a shift from outcrossing to selfing in North American Arabidopsis lyrata were associated with increased average genome-wide LD. Our results indicate that range expansion increases short-distance LD at the farthest range edges by about the same amount as a shift to selfing. However, the extent over which LD in genic regions unfolds was shorter for range expansion compared to selfing. Linkage among putatively neutral variants and between neutral and deleterious variants increased to a similar degree with range expansion, providing support that genome-wide LD was positively associated with mutational load. As a consequence, LD combined with mutational load may decelerate range expansions and set range limits. Finally, a small number of genes were identified as LD outliers, suggesting that they experience selection by either of the two demographic processes. These included genes involved in flowering and photoperiod for range expansion, and the self-incompatibility locus for mating system.

Nearby genomic variants are often co-inherited because of limited recombination. The extent of non-random association of alleles at different loci is called linkage disequilibrium (LD) and is commonly used in genomic analyses, for example to detect regions under selection or to determine effective population size. Here we reversed testing and addressed how demographic history may affect LD within a species. Using genomic data from more than a thousand individuals of North American Arabidopsis lyrata from across the entire species’ range, we quantified the effect of postglacial range expansion and a shift in mating system from outcrossing to selfing on LD. We show that both factors lead to increased LD, and that the maximal effect of range expansion is comparable with a shift in mating system to selfing. Heightened LD involves deleterious mutations, and therefore, LD can also serve as an indicator of mutation accumulation. Furthermore, we provide evidence that some genes experienced stronger increases in LD possibly due to selection associated with the two demographic changes. Our results provide a novel and broad view on the evolutionary factors shaping LD that may also apply to the very many species that underwent postglacial range expansion.

Revealing environmentally driven population dynamics of an Arctic diatom using a novel microsatellite PoolSeq barcoding approach

Ecological stability under environmental change is determined by both interspecific and intraspecific processes. Particularly for planktonic microorganisms, it is challenging to follow intraspecific dynamics over space and time. We propose a new method, microsatellite PoolSeq barcoding (MPB), for tracing allele frequency changes in protist populations. We successfully applied this method to experimental community incubations and field samples of the diatom Thalassiosira hyalina from the Arctic, a rapidly changing ecosystem. Validation of the method found compelling accuracy in comparison with established genotyping approaches within different diversity contexts. In experimental and environmental samples, we show that MPB can detect meaningful patterns of population dynamics, resolving allelic stability and shifts within a key diatom species in response to experimental treatments as well as different bloom phases and years. Through our novel MPB approach, we produced a large dataset of populations at different time-points and locations with comparably little effort. Results like this can add insights into the roles of selection and plasticity in natural protist populations under stable experimental but also variable field conditions. Especially for organisms where genotype sampling remains challenging, MPB holds great potential to efficiently resolve eco-evolutionary dynamics and to assess the mechanisms and limits of resilience to environmental stressors.

Genomic Diversity and Evolution of Quasispecies in Newcastle Disease Virus Infections

Newcastle disease virus (NDV) infections are well known to harbour quasispecies, due to the error-prone nature of the RNA polymerase. Quasispecies variants in the fusion cleavage site of the virus are known to significantly change its virulence. However, little is known about the genomic patterns of diversity and selection in NDV viral swarms. We analyse deep sequencing data from in vitro and in vivo NDV infections to uncover the genomic patterns of diversity and the signatures of selection within NDV swarms. Variants in viruses from in vitro samples are mostly localised in non-coding regions and 3′ and 5′ untranslated regions (3′UTRs or 5′UTRs), while in vivo samples contain an order of magnitude more variants. We find different patterns of genomic divergence and diversity among NDV genotypes, as well as differences in the genomic distribution of intra-host variants among in vitro and in vivo infections of the same strain. The frequency spectrum shows clear signatures of intra-host purifying selection in vivo on the matrix protein (M) coding gene and positive or diversifying selection on nucleocapsid (NP) and haemagglutinin-neuraminidase (HN). The comparison between within-host polymorphisms and phylogenetic divergence reveals complex patterns of selective pressure on the NDV genome at between- and within-host level. The M sequence is strongly constrained both between and within hosts, fusion protein (F) coding gene is under intra-host positive selection, and NP and HN show contrasting patterns: HN RNA sequence is positively selected between hosts while its protein sequence is positively selected within hosts, and NP is under intra-host positive selection at the RNA level and negative selection at the protein level.

Demographic Processes Linked to Genetic Diversity and Positive Selection across a Species' Range

Demography determines the strength of genetic drift, which generally reduces genetic variation and the efficacy of selection. Here, we disentangled the importance of demographic processes at a local scale (census size and mating system) and at a species-range scale (old split between population clusters, recolonization after the last glaciation cycle, and admixture) in determining within-population genomic diversity and genomic signatures of positive selection. Analyses were based on re-sequence data from 52 populations of North American Arabidopsis lyrata collected across its entire distribution. The mating system and range dynamics since the last glaciation cycle explained around 60% of the variation in genomic diversity among populations and 52% of the variation in the signature of positive selection. Diversity was lowest in selfing compared with outcrossing populations and in areas further away from glacial refugia. In parallel, reduced positive selection was found in selfing populations and in populations with a longer route of postglacial range expansion. The signature of positive selection was also reduced in populations without admixture. We conclude that recent range expansion can have a profound influence on diversity in coding and non-coding DNA, similar in magnitude to the shift toward selfing. Distribution limits may in fact be caused by reduced effective population size and compromised positive selection in recently colonized parts of the range.

Genomics versus mtDNA for resolving stock structure in the silky shark (Carcharhinus falciformis)

Conservation genetic approaches for elasmobranchs have focused on regions of the mitochondrial genome or a handful of nuclear microsatellites. High-throughput sequencing offers a powerful alternative for examining population structure using many loci distributed across the nuclear and mitochondrial genomes. These single nucleotide polymorphisms are expected to provide finer scale and more accurate population level data; however, there have been few genomic studies applied to elasmobranch species. The desire to apply next-generation sequencing approaches is often tempered by the costs, which can be offset by pooling specimens prior to sequencing (pool-seq). In this study, we assess the utility of pool-seq by applying this method to the same individual silky sharks, Carcharhinus falciformis, previously surveyed with the mtDNA control region in the Atlantic and Indian Oceans. Pool-seq methods were able to recover the entire mitochondrial genome as well as thousands of nuclear markers. This volume of sequence data enabled the detection of population structure between regions of the Atlantic Ocean populations, undetected in the previous study (inter-Atlantic mitochondrial SNPs FST values comparison ranging from 0.029 to 0.135 and nuclear SNPs from 0.015 to 0.025). Our results reinforce the conclusion that sampling the mitochondrial control region alone may fail to detect fine-scale population structure, and additional sampling across the genome may increase resolution for some species. Additionally, this study shows that the costs of analyzing 4,988 loci using pool-seq methods are equivalent to the standard Sanger-sequenced markers and become less expensive when large numbers of individuals (>300) are analyzed.

Within-population sperm competition intensity does not predict asymmetry in conpopulation sperm precedence

Postcopulatory sexual selection can generate evolutionary arms races between the sexes resulting in the rapid coevolution of reproductive phenotypes. As traits affecting fertilization success diverge between populations, postmating prezygotic (PMPZ) barriers to gene flow may evolve. Conspecific sperm precedence is a form of PMPZ isolation thought to evolve early during speciation yet has mostly been studied between species. Here, we show conpopulation sperm precedence (CpSP) between Drosophila montana populations. Using Pool-seq genomic data we estimate divergence times and ask whether PMPZ isolation evolved in the face of gene flow. We find models incorporating gene flow fit the data best indicating populations experienced considerable gene flow during divergence. We find CpSP is asymmetric and mirrors asymmetry in non-competitive PMPZ isolation, suggesting these phenomena have a shared mechanism. However, we show asymmetry is unrelated to the strength of postcopulatory sexual selection acting within populations. We tested whether overlapping foreign and coevolved ejaculates within the female reproductive tract altered fertilization success but found no effect. Our results show that neither time since divergence nor sperm competitiveness predicts the strength of PMPZ isolation. We suggest that instead cryptic female choice or mutation-order divergence may drive divergence of postcopulatory phenotypes resulting in PMPZ isolation.

This article is part of the theme issue ‘Fifty years of sperm competition’.

An experimental test of temperature-dependent selection on mitochondrial haplotypes in Callosobruchus maculatus seed beetles

Mitochondrial DNA (mtDNA) consists of few but vital maternally inherited genes that interact closely with nuclear genes to produce cellular energy. How important mtDNA polymorphism is for adaptation is still unclear. The assumption in population genetic studies is often that segregating mtDNA variation is selectively neutral. This contrasts with empirical observations of mtDNA haplotypes affecting fitness-related traits and thermal sensitivity, and latitudinal clines in mtDNA haplotype frequencies. Here, we experimentally test whether ambient temperature affects selection on mtDNA variation, and whether such thermal effects are influenced by intergenomic epistasis due to interactions between mitochondrial and nuclear genes, using replicated experimental evolution in Callosobruchus maculatus seed beetle populations seeded with a mixture of different mtDNA haplotypes. We also test for sex-specific consequences of mtDNA evolution on reproductive success, given that mtDNA mutations can have sexually antagonistic fitness effects. Our results demonstrate natural selection on mtDNA haplotypes, with some support for thermal environment influencing mtDNA evolution through mitonuclear epistasis. The changes in male and female reproductive fitness were both aligned with changes in mtDNA haplotype frequencies, suggesting that natural selection on mtDNA is sexually concordant in stressful thermal environments. We discuss the implications of our findings for the evolution of mtDNA.

Genomic Signatures of Local Adaptation in Clam Shrimp (Eulimnadia texana) from Natural Vernal Pools

Vernal pools are unique in their isolation and the strong selection acting on their resident species. Vernal pool clam shrimp (Eulimnadia texana) are a promising model due to ease of culturing, short generation time, small genomes, and obligate desiccated diapaused eggs. Clam shrimp are also androdioecious (sexes include males and hermaphrodites), and here we use population-scaled recombination rates to support the hypothesis that the heterogametic sex is recombination free in these shrimp. We collected short-read sequence data from pooled samples from different vernal pools to gain insights into local adaptation. We identify genomic regions in which some populations have allele frequencies that differ significantly from the metapopulation. BayPass (Gautier M. 2015. Genome-wide scan for adaptive divergence and association with population-specific covariates. Genetics 201(4):1555-1579.) detected 19 such genomic regions showing an excess of population subdivision. These regions on average are 550 bp in size and had 2.5 genes within 5 kb of them. Genes located near these regions are involved in Malpighian tubule function and osmoregulation, an essential function in vernal pools. It is likely that salinity profiles vary between pools and over time, and variants at these genes are adapted to local salinity conditions.

Genome-Wide Changes in Genetic Diversity in a Population of Myotis lucifugus Affected by White-Nose Syndrome

Novel pathogens can cause massive declines in populations, and even extirpation of hosts. But disease can also act as a selective pressure on survivors, driving the evolution of resistance or tolerance. Bat white-nose syndrome (WNS) is a rapidly spreading wildlife disease in North America. The fungus causing the disease invades skin tissues of hibernating bats, resulting in disruption of hibernation behavior, premature energy depletion, and subsequent death. We used whole-genome sequencing to investigate changes in allele frequencies within a population of Myotis lucifugus in eastern North America to search for genetic resistance to WNS. Our results show low F_ST values within the population across time, i.e., prior to WNS (Pre-WNS) compared to the population that has survived WNS (Post-WNS). However, when dividing the population with a geographical cut-off between the states of Pennsylvania and New York, a sharp increase in values on scaffold GL429776 is evident in the Post-WNS samples. Genes present in the diverged area are associated with thermoregulation and promotion of brown fat production. Thus, although WNS may not have subjected the entire M. lucifugus population to selective pressure, it may have selected for specific alleles in Pennsylvania through decreased gene flow within the population. However, the persistence of remnant sub-populations in the aftermath of WNS is likely due to multiple factors in bat life history.

Sex determination in the GIFT strain of tilapia is controlled by a locus in linkage group 23

Background

Tilapias (Family Cichlidae) are the second most important group of aquaculture species in the world. They have been the subject of much research on sex determination due to problems caused by early maturation in culture and their complex sex-determining systems. Different sex-determining loci (linkage group 1, 20 and 23) have been detected in various tilapia stocks. The ‘genetically improved farmed tilapia’ (GIFT) stock, founded from multiple Nile tilapia (Oreochromis niloticus) populations, with some likely to have been introgressed with O. mossambicus, is a key resource for tilapia aquaculture. The sex-determining mechanism in the GIFT stock was unknown, but potentially complicated due to its multiple origins.

Results

A bulk segregant analysis (BSA) version of double-digest restriction-site associated DNA sequencing (BSA-ddRADseq) was developed and used to detect and position sex-linked single nucleotide polymorphism (SNP) markers in 19 families from the GIFT strain breeding nucleus and two Stirling families as controls (a single XY locus had been previously mapped to LG1 in the latter). About 1500 SNPs per family were detected across the genome. Phenotypic sex in Stirling families showed strong association with LG1, whereas only SNPs located in LG23 showed clear association with sex in the majority of the GIFT families. No other genomic regions linked to sex determination were apparent. This region was validated using a series of LG23-specific DNA markers (five SNPs with highest association to sex from this study, the LG23 sex-associated microsatellite UNH898 and ARO172, and the recently isolated amhy marker for individual fish (n = 284).

Conclusions

Perhaps surprisingly given its multiple origins, sex determination in the GIFT strain breeding nucleus was associated only with a locus in LG23. BSA-ddRADseq allowed cost-effective analysis of multiple families, strengthening this conclusion. This technique has potential to be applied to other complex traits. The sex-linked SNP markers identified will be useful for potential marker-assisted selection (MAS) to control sex-ratio in GIFT tilapia to suppress unwanted reproduction during growout.

Genomic release-recapture experiment in the wild reveals within-generation polygenic selection in stickleback fish

How rapidly natural selection sorts genome-wide standing genetic variation during adaptation remains largely unstudied experimentally. Here, we present a genomic release-recapture experiment using paired threespine stickleback fish populations adapted to selectively different lake and stream habitats. First, we use pooled whole-genome sequence data from the original populations to identify hundreds of candidate genome regions likely under divergent selection between these habitats. Next, we generate F2 hybrids from the same lake-stream population pair in the laboratory and release thousands of juveniles into a natural stream habitat. Comparing the individuals surviving one year of stream selection to a reference sample of F2 hybrids allows us to detect frequency shifts across the candidate regions toward the genetic variants typical of the stream population-an experimental outcome consistent with polygenic directional selection. Our study reveals that adaptation in nature can be detected as a genome-wide signal over just a single generation.