Comparative sex chromosome genomics in snakes: differentiation, evolutionary strata, and lack of global dosage compensation

Abstract

Snakes exhibit genetic sex determination, with female heterogametic sex chromosomes (ZZ males, ZW females). Extensive cytogenetic work has suggested that the level of sex chromosome heteromorphism varies among species, with Boidae having entirely homomorphic sex chromosomes, Viperidae having completely heteromorphic sex chromosomes, and Colubridae showing partial differentiation. Here, we take a genomic approach to compare sex chromosome differentiation in these three snake families. We identify homomorphic sex chromosomes in boas (Boidae), but completely heteromorphic sex chromosomes in both garter snakes (Colubridae) and pygmy rattlesnake (Viperidae). Detection of W-linked gametologs enables us to establish the presence of evolutionary strata on garter and pygmy rattlesnake sex chromosomes where recombination was abolished at different time points. Sequence analysis shows that all strata are shared between pygmy rattlesnake and garter snake, i.e., recombination was abolished between the sex chromosomes before the two lineages diverged. The sex-biased transmission of the Z and its hemizygosity in females can impact patterns of molecular evolution, and we show that rates of evolution for Z-linked genes are increased relative to their pseudoautosomal homologs, both at synonymous and amino acid sites (even after controlling for mutational biases). This demonstrates that mutation rates are male-biased in snakes (male-driven evolution), but also supports faster-Z evolution due to differential selective effects on the Z. Finally, we perform a transcriptome analysis in boa and pygmy rattlesnake to establish baseline levels of sex-biased expression in homomorphic sex chromosomes, and show that heteromorphic ZW chromosomes in rattlesnakes lack chromosome-wide dosage compensation. Our study provides the first full scale overview of the evolution of snake sex chromosomes at the genomic level, thus greatly expanding our knowledge of reptilian and vertebrate sex chromosomes evolution.

Figure 1. Normalized read coverage depth for female (red) and male (blue) scaffolds ordered along the Anolis genome for (A) boa, (B) garter snake, and (C) pygmy rattlesnake.

The points show the normalized log₂ coverage for each scaffold and the lines represent a smoothing spline drawn along the chromosome. Coverage was normalized by dividing the coverage for each scaffold-by-sex combination by the median coverage of all scaffolds in chromosomes 1–5 in that sex, resulting in a median log₂ coverage score for autosomes of 0. Under this normalization, hemizygous sequences are expected to have a median log₂ coverage of −1. Scaffolds were mapped to the Anolis macrochromosomes on the basis of the location of their gene content. The phylogenetic relationship between the species investigated is shown. Boids split from the other two groups about 100 MY ago, while colubrids and viperids diverged about 50 MY ago . The snake photographs are used under a Creative Commons Attribution 2.5 Generic license (CC BY 2.5). Credit for the photographs are as follows: Nick Turland for the western pygmy rattlesnake (http://www.flickr.com/photos/nturland/1436776818/); Guilherme Jófili for the Boa constrictor (http://www.flickr.com/photos/gjofili/5005623645/); Steve Jurvetson for the coastal garter snake (http://www.flickr.com/photos/jurvetson/825514494/). Additional permission to publish the western pygmy rattlesnake image was granted by Nick Turland.

Figure 2. Mapping of the best candidate female-specific W-candidate scaffolds of (A) boa, (B) garter snake, and (C) pygmy rattlesnake to the Anolis macrochromosomes.

The histograms (left) show the number of candidate W scaffolds mapped to the six major chromosomes of Anolis, with the green bar highlighting the Anolis homolog to the snake Z chromosome. The W candidates homologous to Z-linked scaffolds (which accounts for 7.6% of the genome) make up 49% and 35% of all female-biased scaffolds in pygmy rattlesnake and garter snake, respectively. This is a 6.4-fold excess in pygmy rattlesnake, and 4.5-fold excess in garter snake over random mapping based on chromosome size. In contrast, in Boa, only 6.5% of scaffolds map to chromosome 6, which does not differ from random mapping on the basis of chromosome size. The right panel shows color-coded mapping density of W-candidates along the Anolis macrochromosomes. The density of W-candidates is not uniform across the Z chromosome in both pygmy rattlesnake and garter snake (p<0.0001 for mean nearest neighbor distances). The data in this figure are from all female biased candidate W-linked scaffolds and will thus contain both non-coding scaffolds as well as scaffolds containing protein coding genes.

Figure 3. Evolutionary strata and sequence conservation between the pygmy rattlesnake and garter snake W-candidate scaffolds mapped along the Z chromosome.

The middle three tracks show the position of candidate W sequences along the Z chromosome in garter snake (top) and pygmy rattlesnake (bottom), and their overlap (center). The top and bottom plots show nucleotide identity between Z-W gametologs (grey dots), and the median (red line) inferred for each of the putative strata for garter snake and pygmy rattlesnake (blue boxes; the y range of the boxes represents the interquantile range of the identity values in each region). The gray shaded region represents identity below 30%, indicating low quality mapping. As in Figure 2, the data in this figure contain both non-coding and coding scaffolds.

Figure 4. Molecular evolution of snake Z chromosomes and autosomes at synonymous sites (K_s) and non-synonymous sites (K_a), and their ratio (K_a/K_s).

For each gene, we calculated the Anolis-boa, Anolis-pygmy rattlesnake, and Anolis-garter snake rates of synonymous and non-synonymous evolution. To detect branch-specific differences, we obtained for each gene the ratios of these evolutionary rates between the different snake species pairs (pygmy rattlesnake/boa, garter snake/boa, and garter snake/pygmy rattlesnake), and plotted them for each macrochromosome. Please note that in the figure, “Pygmy Rattlesnake/Boa” refers to the ratio: pygmy-Anolis divergence/boa-Anolis divergence, and so on. Significant differences between the Z-chromosome and the autosomes are marked with asterisks (***, p<0.001; NS, non-significant).

Figure 5. Log₂ of expression in female, log₂ of expression in male, and log₂ of female over male expression, for the different macrochromosomes of boa (A) and pygmy rattlesnake (B).

FPKM values were obtained for each gene using Cufflinks. Genes were assigned to different chromosomes according to their location in the Anolis genome. (C) log₂ of female over male expression along chromosome 6 of Anolis (Z of snakes), using a sliding window size of 30 genes.