The genomic basis of cichlid fish adaptation within the deepwater "twilight zone" of Lake Malawi

Abstract

Deepwater environments are characterized by low levels of available light at narrow spectra, great hydrostatic pressure, and low levels of dissolved oxygen-conditions predicted to exert highly specific selection pressures. In Lake Malawi over 800 cichlid species have evolved, and this adaptive radiation extends into the "twilight zone" below 50 m. We use population-level RAD-seq data to investigate whether four endemic deepwater species (Diplotaxodon spp.) have experienced divergent selection within this environment. We identify candidate genes including regulators of photoreceptor function, photopigments, lens morphology, and haemoglobin, many not previously implicated in cichlid adaptive radiations. Colocalization of functionally linked genes suggests coadapted "supergene" complexes. Comparisons of Diplotaxodon to the broader Lake Malawi radiation using genome resequencing data revealed functional substitutions and signatures of positive selection in candidate genes. Our data provide unique insights into genomic adaptation within deepwater habitats, and suggest genome-level specialization for life at depth as an important process in cichlid radiation.

Keywords: Root effect; cichlid; hemoglobin; sensory drive; supergene.

Figure 1

Maps of (A) Africa (topographic) and (B) Lake Malawi (bathymetric) indicating the sampling location Nkhata Bay. (C) D. “limnothrissa black pelvic”; (D) D. “macrops offshore”; (E) D. “macrops black dorsal”; (F) D. “macrops ngulube.”

Figure 2

(A) Maximum‐likelihood tree (unrooted) of phylogenomic structure among Diplotaxodon species. Yellow—D. “macrops black dorsal”; red—D. “limnothrissa black pelvic”; black—D. “macrops offshore”; blue—D. “macrops ngulube.” Asterisks indicate >95% bootstrap branch support. Scale bar indicates genetic divergence (nucleotide divergence per site). (B) Genome‐wide pattern of pairwise divergence between populations. Dots represent pairwise F_ST values for SNPs as calculated by populations (Catchen et al. 2013). Only scaffolds with minimum length of 100 kb containing a minimum of 10 SNPs are displayed. Individual scaffold boundaries are indicated by alternate white/grey background. Dashed lines indicate the global pairwise F_ST average. Highlighted groups in the phylogenetic trees on the right‐hand side indicate the population pairs. Green arrow heads on top of the figure indicate locations of candidate regions supported by all three candidate outlier approaches (see Fig. S3).

Figure 3

(A) Diplotaxodon eye diameter (normalized by total length of the fish, TL) across the four species. Yellow—D. “macrops black dorsal”; red—D. “limnothrissa black pelvic”; black—D. “macrops offshore”; blue—D. “macrops ngulube.” (B) Pope plot, summarizing correlations of population allele frequency with vertical eye diameter, inferred from 20 independent Bayenv runs. Dots illustrate per locus average relative rank (ARR) versus relative rank standard deviation. The green vertical line delimits the 95th ARR percentile. The 42 loci most consistently correlated with eye size (ARR ≥ 0.95) are shown in green. (C) Pairwise F_ST divergence (top six panels) and global allele frequency correlations with interspecific vertical eye diameter differences (bottom panel). Displayed are the six scaffolds containing the most significantly correlated loci. Corresponding regions are highlighted in shades of green (light green—ARR ≥ 0.95; dark green—ARR ≥ 0.95, and smoothed ARR P < 0.001). Dots represent the kernel smoothed averages across 50 kb windows. Dashed lines indicate the genome wide average of F_ST/ARR. Population pairs are indicated by the highlighted populations in the phylogenetic trees on the right‐hand side. Black arrowheads on top indicate regions that were also supported as candidate outlier loci by at least two of three independent outlier detection approaches. Table S5 summarizes the gene complements in highlighted regions.