Independently evolved pollution resistance in four killifish populations is largely explained by few variants of large effect

Abstract

The genetic architecture of phenotypic traits can affect the mode and tempo of trait evolution. Human-altered environments can impose strong natural selection, where successful evolutionary adaptation requires swift and large phenotypic shifts. In these scenarios, theory predicts that adaptation is due to a few adaptive variants of large effect, but empirical studies that have revealed the genetic architecture of rapidly evolved phenotypes are rare, especially for populations inhabiting polluted environments. Fundulus killifish have repeatedly evolved adaptive resistance to extreme pollution in urban estuaries. Prior studies, including genome scans for signatures of natural selection, have revealed some of the genes and pathways important for evolved pollution resistance, and provide context for the genotype-phenotype association studies reported here. We created multiple quantitative trait locus (QTL) mapping families using progenitors from four different resistant populations, and using RAD-seq genetically mapped variation in sensitivity (developmental perturbations) following embryonic exposure to a model toxicant PCB-126. We found that one to two large-effect QTL loci accounted for resistance to PCB-mediated developmental toxicity. QTLs harbored candidate genes that govern the regulation of aryl hydrocarbon receptor (AHR) signaling. One QTL locus was shared across all populations and another was shared across three populations. One QTL locus showed strong signatures of recent natural selection in the corresponding wild population but another QTL locus did not. Some candidate genes for PCB resistance inferred from genome scans in wild populations were identified as QTL, but some key candidate genes were not. We conclude that rapidly evolved resistance to the developmental defects normally caused by PCB-126 is governed by few genes of large effect. However, other aspects of resistance beyond developmental phenotypes may be governed by additional loci, such that comprehensive resistance to PCB-126, and to the mixtures of chemicals that distinguish urban estuaries more broadly, may be more genetically complex.

Keywords: adaptation; contemporary evolution; ecological genomics; genomics/proteomics; molecular evolution; quantitative genetics.

FIGURE 1

(a) Locations of source populations for founders of mapping families. Sites are characterized by complex mixtures of pollutants but dominated by either polychlorinated biphenyls (PCBs) or dioxins in the three Northern sites (New Bedford Harbor, MA [NBH], Bridgeport Harbor, CT [BRP], and Newark, NJ [NEW]), or polycyclic aromatic hydrocarbons (PAHs) in the Southern site (Elizabeth River, VA [ELR]). All polluted site populations (filled circles) have been previously shown to be resistant to PCB‐126. Open blue circle is the clean reference site (Block Island, RI [BLI]) and source for one of the progenitors for each of the mapping families. (b) Bar graphs show the phenotypic distribution (developmental deformity index, x‐axis) of embryos from each of the F2 intercross mapping families at 10 days after an exposure to a discriminating dose, (200 ng PCB‐126/L in sea water). Full bars show the total count of embryos collected and assigned a deformity score, which ranged from 0 (no deformities) to 5 (most severe deformities) (the BRP families were scored on a scale where a score of 4 was maximum). Colored portion of each bar indicates the subset of embryos that were selected for genotyping.

FIGURE 2

QTL and their contribution to variation in sensitivity to PCB‐induced toxicity. (a) LOD score results along each chromosome (x‐axis numbers) for single QTL scans of RAD‐Tag markers in each mapping family. Blue, green, orange, and red series correspond to mapping families from NBH, BRP, NEW, and ELR, respectively. Tall colored boxes (blue and gold) indicate the QTL that fall on the same chromosome/position in multiple families. Grey boxes indicate minor effect QTL that explain a relatively minor proportion of the resistance phenotype and are not shared between populations. Note that the vertical alignment of homologous chromosomes between families is not perfect, because estimated mapping distances for the same chromosome can vary between mapping families because of variation in recombination rates between families and genotyping error. (b) Summary of the contribution (drop one QTL analysis) of large‐effect QTL (blue and gold boxes) to the full QTL model (LOD and percent variance explained) in each mapping family. (c) Genotype‐by‐phenotype plot for the NBH family at the AIP candidate locus (QTL on chromosome 2) and AHR1b/2b candidate locus (QTL on chromosome 18). The Y‐axis groups individuals (points) as either resistant (R: phenotype malformation scores 0–1 following PCB‐126 exposure) or sensitive (S: phenotype malformation scores 4–5 following PCB‐126 exposure). The X‐axis distinguishes individuals based on their genotype at the AIP locus (QTL on chromosome 2), where individuals homozygous for the resistant allele, heterozygous, or homozygous for the sensitive allele, are represented by RR, RS, or SS, respectively. The color of points distinguishes individuals based on their genotype at the AHR1b/2b locus (QTL on chromosome 18), where individuals homozygous for the resistant allele, heterozygous, or homozygous for the sensitive allele, are represented by red, green, or blue, respectively. All individuals that carry the RR genotype at the AIP locus are phenotypically resistant (top left group). Those that carry one or no copies of the resistant allele at the AIP locus tend to be phenotypically sensitive (bottom middle and right groups), unless they are also homozygous for the resistant allele at the AHR1b/2b locus (red dots) in which case they tend to be phenotypically resistant (top middle and right groups).

FIGURE 3

Large‐effect QTL (first row) on CHR2 (left column) and CHR18 (right column) aligned with signatures of selection from genome scans for each resistant population (remaining rows). X‐axis tick marks represent Mb intervals along the entire length of the two chromosomes. QTL LOD scores (top row) are the results from a single‐QTL scan (marker regression) in each of the QTL mapping families (the blue series represents LOD scores from NBH mapping families, green from BRP, orange from NEW, and red from ELR). Panels in the remaining rows include F _ST and delta p _i between each of the resistant populations and their nearby reference populations (original genome scan data were collected by Reid et al., 2016). Elevated F _ST and reduced delta p _i are signatures of recent strong natural selection. Tall boxes highlight the position and rank of the top‐ranked signatures of selection that co‐localize to a large‐effect QTL. The large effect QTL on CHR2 (which contains AIP) coincides with a highly ranked selection signature region; although this region shows a strong signature of selection in all four populations, it includes a QTL for only three of the populations (NBH, BRP, and NEW). The large effect QTL on chromosome 18 (which contains AHR1b and AHR2b) is found in all four populations, but coincides with highly ranked selection signature regions in only two of those populations (ELR and NEW).