Independent Recruitment of Different Types of Phospholipases A2 to the Venoms of Caenophidian Snakes: The Rise of PLA2-IIE within Pseudoboini (Dipsadidae).

Snake venoms harbor a wide and diverse array of enzymatic and nonenzymatic toxic components, allowing them to exert myriad effects on their prey. However, they appear to trend toward a few optimal compositional scaffolds, dominated by four major toxin classes: SVMPs, SVSPs, 3FTxs, and PLA2s. Nevertheless, the latter appears to be restricted to vipers and elapids, as it has never been reported as a major venom component in rear-fanged species. Here, by investigating the original transcriptomes from 19 species distributed in eight genera from the Pseudoboini tribe (Dipsadidae: Xenodontinae) and screening among seven additional tribes of Dipsadidae and three additional families of advanced snakes, we discovered that a novel type of venom PLA2, resembling a PLA2-IIE, has been recruited to the venom of some species of the Pseudoboini tribe, where it is a major component. Proteomic and functional analyses of these venoms further indicate that these PLA2s play a relevant role in the venoms from this tribe. Moreover, we reconstructed the phylogeny of PLA2s across different snake groups and show that different types of these toxins have been recruited in at least five independent events in caenophidian snakes. Additionally, we present the first compositional profiling of Pseudoboini venoms. Our results demonstrate how relevant phenotypic traits are convergently recruited by different means and from homologous and nonhomologous genes in phylogenetically and ecologically divergent snake groups, possibly optimizing venom composition to overcome diverse adaptative landscapes.

Whole snake genomes from eighteen families of snakes (Serpentes: Caenophidia) and their applications to systematics.

We present genome assemblies for 18 snake species representing 18 families (Serpentes: Caenophidia): Acrochordus granulatus, Aparallactus werneri, Boaedon fuliginosus, Calamaria suluensis, Cerberus rynchops, Grayia smithii, Imantodes cenchoa, Mimophis mahfalensis, Oxyrhabdium leporinum, Pareas carinatus, Psammodynastes pulverulentus, Pseudoxenodon macrops, Pseudoxyrhopus heterurus, Sibynophis collaris, Stegonotus admiraltiensis, Toxicocalamus goodenoughensis, Trimeresurus albolabris, and Tropidonophis doriae. From these new genome assemblies, we extracted thousands of loci commonly used in systematic and phylogenomic studies on snakes, including target-capture datasets composed of ultraconserved elements (UCEs) and anchored hybrid enriched loci (AHEs), as well as traditional Sanger loci. Phylogenies inferred from the two target-capture loci datasets were identical with each other and strongly congruent with previously published snake phylogenies. To show the additional utility of these non-model genomes for investigative evolutionary research, we mined the genome assemblies of two New Guinea island endemics in our dataset (S. admiraltiensis and T. doriae) for the ATP1a3 gene, a thoroughly researched indicator of resistance to toad toxin ingestion by squamates. We find that both these snakes possess the genotype for toad toxin resistance despite their endemism to New Guinea, a region absent of any toads until the human-mediated introduction of Cane Toads in the 1930s. These species possess identical substitutions that suggest the same bufotoxin resistance as their Australian congenerics (Stegonotus australis and Tropidonophis mairii) which forage on invasive Cane Toads. Herein, we show the utility of short-read high-coverage genomes, as well as improving the deficit of available squamate genomes with associated voucher specimens.

Proliferative strongyloidiasis in a colony of colubrid snakes.

Strongyloides are small rhabditid nematodes primarily associated with enteric disease in a variety of animal species, including reptiles. Strongyloides spp life stages were associated with a disease outbreak in a large breeding colony of snakes. Multiple Pantherophis and Lampropeltis colubrids exhibited respiratory distress, anorexia, stomatitis, facial deformation, and waning body condition that resulted in death or necessitated euthanasia. Postmortem examinations of 13 snakes revealed epithelial hyperplasia and inflammation of the alimentary and respiratory tracts associated with varying numbers of adult and larval nematodes and embryonated or larvated ova. In a subset of snakes, aberrant nematode migration was also observed in the eye, genitourinary system, coelom, and vasculature. Histomorphology and gross examination of parasitic adult female nematodes from host tissues were consistent with a Strongyloides spp. Sedimented fecal material from 101/160 (63%) snakes housed in the affected facility was positive for nematodes and/or larvated ova. Polymerase chain reaction amplification and sequencing of portions of the 18S and 28S ribosomal ribonucleic acid (RNA) genes and the internal transcribed spacer region of adult female parasites and positive fecal samples supported the diagnosis of strongyloidiasis. Strongyloides spp possess a unique life cycle capable of alternating between parasitic (homogonic) and free-living (heterogonic) stages, resulting in the production of directly infective larvae. Commonly utilized husbandry practices in reptile collections can amplify the numbers of infective larvae generated in the captive environment, increasing the risk for rhabditid hyperinfections. This report documents morbidity, mortality, and non-enteric disease manifestations due to Strongyloides hyperinfections in a captive colubrid snake colony.

Conspicuous coloration of toxin-resistant predators implicates additional trophic interactions in a predator-prey arms race.

Antagonistic coevolution between natural enemies can produce highly exaggerated traits, such as prey toxins and predator resistance. This reciprocal process of adaptation and counter-adaptation may also open doors to other evolutionary novelties not directly involved in the phenotypic interface of coevolution. We tested the hypothesis that predator-prey coevolution coincided with the evolution of conspicuous coloration on resistant predators that retain prey toxins. In western North America, common garter snakes (Thamnophis sirtalis) have evolved extreme resistance to tetrodotoxin (TTX) in the coevolutionary arms race with their deadly prey, Pacific newts (Taricha spp.). TTX-resistant snakes can retain large amounts of ingested TTX, which could serve as a deterrent against the snakes' own predators if TTX toxicity and resistance are coupled with a conspicuous warning signal. We evaluated whether arms race escalation covaries with bright red coloration in snake populations across the geographic mosaic of coevolution. Snake colour variation departs from the neutral expectations of population genetic structure and covaries with escalating clines of newt TTX and snake resistance at two coevolutionary hotspots. In the Pacific Northwest, bright red coloration fits an expected pattern of an aposematic warning to avian predators: TTX-resistant snakes that consume highly toxic newts also have relatively large, reddish-orange dorsal blotches. Snake coloration also seems to have evolved with the arms race in California, but overall patterns are less intuitively consistent with aposematism. These results suggest that interactions with additional trophic levels can generate novel traits as a cascading consequence of arms race coevolution across the geographic mosaic.

A reference genome assembly for the continentally distributed ring-necked snake, Diadophis punctatus.

Snakes in the family Colubridae include more than 2,000 currently recognized species, and comprise roughly 75% of the global snake species diversity on Earth. For such a spectacular radiation, colubrid snakes remain poorly understood ecologically and genetically. Two subfamilies, Colubrinae (788 species) and Dipsadinae (833 species), comprise the bulk of colubrid species richness. Dipsadines are a speciose and diverse group of snakes that largely inhabit Central and South America, with a handful of small-body-size genera that have invaded North America. Among them, the ring-necked snake, Diadophis punctatus, has an incredibly broad distribution with 14 subspecies. Given its continental distribution and high degree of variation in coloration, diet, feeding ecology, and behavior, the ring-necked snake is an excellent species for the study of genetic diversity and trait evolution. Within California, six subspecies form a continuously distributed "ring species" around the Central Valley, while a seventh, the regal ring-necked snake, Diadophis punctatus regalis is a disjunct outlier and Species of Special Concern in the state. Here, we report a new reference genome assembly for the San Diego ring-necked snake, D. p. similis, as part of the California Conservation Genomics Project. This assembly comprises a total of 444 scaffolds spanning 1,783 Mb and has a contig N50 of 8.0 Mb, scaffold N50 of 83 Mb, and BUSCO completeness score of 94.5%. This reference genome will be a valuable resource for studies of the taxonomy, conservation, and evolution of the ring-necked snake across its broad, continental distribution.

Genome mapping of a LYST mutation in corn snakes indicates that vertebrate chromatophore vesicles are lysosome-related organelles.

Reptiles exhibit a spectacular diversity of skin colors and patterns brought about by the interactions among three chromatophore types: black melanophores with melanin-packed melanosomes, red and yellow xanthophores with pteridine- and/or carotenoid-containing vesicles, and iridophores filled with light-reflecting platelets generating structural colors. Whereas the melanosome, the only color-producing endosome in mammals and birds, has been documented as a lysosome-related organelle, the maturation paths of xanthosomes and iridosomes are unknown. Here, we first use 10x Genomics linked-reads and optical mapping to assemble and annotate a nearly chromosome-quality genome of the corn snake Pantherophis guttatus The assembly is 1.71 Gb long, with an N50 of 16.8 Mb and L50 of 24. Second, we perform mapping-by-sequencing analyses and identify a 3.9-Mb genomic interval where the lavender variant resides. The lavender color morph in corn snakes is characterized by gray, rather than red, blotches on a pink, instead of orange, background. Third, our sequencing analyses reveal a single nucleotide polymorphism introducing a premature stop codon in the lysosomal trafficking regulator gene (LYST) that shortens the corresponding protein by 603 amino acids and removes evolutionary-conserved domains. Fourth, we use light and transmission electron microscopy comparative analyses of wild type versus lavender corn snakes and show that the color-producing endosomes of all chromatophores are substantially affected in the LYST mutant. Our work provides evidence characterizing xanthosomes in xanthophores and iridosomes in iridophores as lysosome-related organelles.

Hepatozoon Miller, 1908 parasites in the Colubridae snakes Clelia clelia (Daudin, 1803) and Drymarchon corais (Boie, 1827) from the Eastern Amazonia.

Based on the genetic, morphological, and morphometric data of blood gamonts, we identified Hepatozoon parasites in colubrid snakes sampled in the Eastern Amazon region. Hepatozoon trigeminum was detected in the mussurana snake Clelia clelia and exhibited wide and elongated gamonts (mean dimensions: 14.25±0.65 × 4.31±0.43 µm) with an evident parasitophorous vacuole. Hepatozoon odwyerae sp. nov. was described in the indigo snake Drymarchon corais, whose gamonts have elongated and thin bodies (mean dimensions: 13.41±0.79 × 3.72±0.35 µm) with one end more tapered than the other. Phylogenetic analyses, based on the amplification of a 441 bp fragment of the 18S rRNA gene, revealed that the novel sequences of Hepatozoon spp. from our study were closely related to hemogregarine lineages found in lizards and snakes from Brazil, forming a well-supported monophyletic clade with them. The present study provides the first species description of Hepatoozon in D. corais and a new record of a host species for C. clelia using the integrated taxonomic data. We also highlight the importance of further investigations into the diversity of Hepatozoon in snakes, a rich but underestimated group of parasites, especially in the Amazonian biome.

Uncovering Cryptic Coevolution.

AbstractStudies of coevolution in the wild have largely focused on reciprocally specialized species pairs with striking and exaggerated phenotypes. Textbook examples include interactions between toxic newts and their garter snake predators, long-tongued flies and the flowers they pollinate, and weevils with elongated rostra used to bore through the defensive pericarp of their host plants. Although these studies have laid a foundation for understanding coevolution in the wild, they have also contributed to the widespread impression that coevolution is a rare and quirky sideshow to the day-to-day grind of ecology and evolution. In this perspective, we argue that the focus of coevolution has been biased toward the obvious and ignored the cryptic. We have focused on the obvious-studies of reciprocally specialized species pairs with exaggerated phenotypes-mainly because we have lacked the statistical tools required to study coevolution in more generalized and phenotypically mundane systems. Building from well-established coevolutionary theory, we illustrate how model-based approaches can be used to remove this barrier and begin estimating the strength of coevolutionary selection indirectly using routinely collected data, thus uncovering cryptic coevolution in more typical communities. By allowing the distribution of coevolutionary selection to be estimated across genomes, phylogenies, and communities and over deep timescales, these novel approaches have the potential to revolutionize the way we study coevolution. As we develop a road map to these next-generation approaches, we highlight recent studies making notable progress in this direction.

The road not taken: Evolution of tetrodotoxin resistance in the Sierra garter snake (Thamnophis couchii) by a path less travelled.

The repeated evolution of tetrodotoxin (TTX) resistance provides a model for testing hypotheses about the mechanisms of convergent evolution. This poison is broadly employed as a potent antipredator defence, blocking voltage-gated sodium channels (Nav ) in muscles and nerves, paralysing and sometimes killing predators. Resistance in taxa bearing this neurotoxin and a few predators appears to come from convergent replacements in specific Nav residues that interact with TTX. This stereotyped genetic response suggests molecular and phenotypic evolution may be constrained and predictable. Here, we investigate the extent of mechanistic convergence in garter snakes (Thamnophis) that prey on TTX-bearing newts (Taricha) by examining the physiological and genetic basis of TTX resistance in the Sierra garter snake (Th. couchii). We characterize variation in this predatory adaptation across populations at several biological scales: whole-animal TTX resistance; skeletal muscle resistance; functional genetic variation in three Nav encoding loci; and levels of gene expression for one of these loci. We found Th. couchii possess extensive geographical variation in resistance at the whole-animal and skeletal muscle levels. As in other Thamnophis, resistance at both levels is highly correlated, suggesting convergence across the biological levels linking organism to organ. However, Th. couchii shows no functional variation in Nav loci among populations or difference in candidate gene expression. Local variation in TTX resistance in Th. couchii cannot be explained by the same relationship between genotype and phenotype seen in other taxa. Thus, historical contingencies may lead different species of Thamnophis down alternative routes to local adaptation.

Phylogenomic analyses resolve relationships among garter snakes (Thamnophis: Natricinae: Colubridae) and elucidate biogeographic history and morphological evolution.

Garter snakes (Thamnophis) are a successful group of natricines endemic to North America. They have become important natural models for ecological and evolutionary research, yet prior efforts to resolve phylogenetic relationships have resulted in conflicting topologies and weak support for certain relationships. Here, we use genomic data generated with a reduced representation double-digest RADseq approach to reassess evolutionary relationships across Thamnophis. We then use the resulting phylogeny to better understand how biogeography and feeding ecology have influenced lineage diversification and morphological evolution. We recovered highly congruent and strongly supported topologies from maximum likelihood and Bayesian analyses, but some discordance with a multispecies coalescent approach. All phylogenomic estimates split Thamnophis into two clades largely defined by northern and southern North American species. Divergence time estimates and biogeographic analyses indicate a mid-Miocene origin of Thamnophis in Mexico. In addition, historic vicariant events thought to explain biogeographic patterns in other lineages (e.g., Isthmus of Tehuantepec, Rocky Mountain Range, and Trans-Mexican Volcanic Belt) appear to have influenced patterns of diversification in Thamnophis as well. Analyses of morphological traits associated with feeding ecology showed moderate to strong phylogenetic signal. Nevertheless, phylogenetic ANOVA suggested significant differences in certain cranial morphologies between aquatic specialists and garter snakes that are terrestrial-aquatic generalists, independent of evolutionary history. Our new estimate of Thamnophis phylogeny yields an improved understanding of the biogeographic history and morphological evolution of garter snakes, and provides a robust framework for future research on these snakes.

Cryptic diversity and phylogeography of the Rhabdophis nuchalis group (Squamata: Colubridae).

Previous studies, have found that the rapid uplift of the Tibetan plateau accelerated the diversification of species. However, there are few relevant biogeographic data for the Colubridae in this region. We conducted a comprehensive study of the Rhabdophis nuchalis Group, which presently contains four nominal species, R. nuchalis, R. pentasupralabialis, R. leonardi, and R. chiwen. Building upon previous studies with specimens we have recently examined, greater interspecific and intraspecific diversity has been revealed. Here we address three questions: (1) Do the intraspecific differences represent only geographic variation within lineages, or are there cryptic species? (2) What are the interspecific relationships among members of the R. nuchalis Group? (3) What has been the biogeographic history of this species group? To resolve these questions we used four mitochondrial gene sequences and one nuclear sequence to investigate the molecular phylogenetic and geographic relationships among populations. Our molecular analysis reveals cryptic species diversity within the R. nuchalis Group, and seven clades were identified in the analysis. Ancestral area estimation suggests that the R. nuchalis Group originated in the Hengduan Mountains approximately 6.24 Mya and expanded its range northward to the Qinling-Daba Mountains. The Sichuan Basin appears to have been a barrier to migration. Species divergence seems to have been related to the rapid uplift of the Qinghai-Tibet Plateau.

Reference genome of the California glossy snake, Arizona elegans occidentalis: A declining California Species of Special Concern.

The glossy snake (Arizona elegans) is a polytypic species broadly distributed across southwestern North America. The species occupies habitats ranging from California's coastal chaparral to the shortgrass prairies of Texas and southeastern Nebraska, to the extensive arid scrublands of central México. Three subspecies are currently recognized in California, one of which is afforded state-level protection based on the extensive loss and modification of its preferred alluvial coastal scrub and inland desert habitat. We report the first genome assembly of A. elegans occidentalis as part of the California Conservation Genomics Project (CCGP). Consistent with the reference genome strategy of the CCGP, we used Pacific Biosciences HiFi long reads and Hi-C chromatin-proximity sequencing technologies to produce a de novo assembled genome. The assembly comprises a total of 140 scaffolds spanning 1,842,602,218 base pairs, has a contig NG50 of 61 Mb, a scaffold NG50 of 136 Mb, and a BUSCO complete score of 95.9%, and is one of the most complete snake genome assemblies. The A. e. occidentalis genome will be a key tool for understanding the genomic diversity and the basis of adaptations within this species and close relatives within the hyperdiverse snake family Colubridae.

Integration of transcriptomic and proteomic approaches for snake venom profiling.

INTRODUCTION: Snake venoms contain many protein and peptide isoforms with high levels of sequence variation, even within a single species. AREAS COVERED: In this review, we highlight several examples, from both published and unpublished work in our lab, demonstrating how a combined venom gland transcriptome and proteome methodology allows for comprehensive characterization of venoms, including those from understudied rear-fanged snake species, and we provide recommendations for using these approaches. EXPERT OPINION: When characterizing venoms, peptide mass fingerprinting using databases built predominately from protein sequences originating from model organisms can be disadvantageous, especially when the intention is to document protein diversity. Therefore, the use of species-specific venom gland transcriptomes corrects for the absence of these unique peptide sequences in databases. The integration of transcriptomics and proteomics improves the accuracy of either approach alone for venom profiling.

Molecular phylogenetics of sub-Saharan African natricine snakes, and the biogeographic origins of the Seychelles endemic Lycognathophis seychellensis.

Phylogenetic relationships of sub-Saharan African natricine snakes are understudied and poorly understood, which in turn has precluded analyses of the historical biogeography of the Seychelles endemic Lycognathophis seychellensis. We inferred the phylogenetic relationships of Seychelles and mainland sub-Saharan natricines by analysing a multilocus DNA sequence dataset for three mitochondrial (mt) and four nuclear (nu) genes. The mainland sub-Saharan natricines and L. seychellensis comprise a well-supported clade. Two maximally supported sets of relationships within this clade are (Limnophis,Natriciteres) and (Afronatrix,(Hydraethiops,Helophis)). The relationships of L. seychellensis with respect to these two lineages are not clearly resolved by analysing concatenated mt and nu data. Analysed separately, nu data best support a sister relationship of L. seychellensis with (Afronatrix,(Hydraethiops,Helophis)) and mt data best support a sister relationship with all mainland sub-Saharan natricines. Methods designed to cope with incomplete lineage sorting strongly favour the former hypothesis. Genetic variation among up to 33 L. seychellensis from five Seychelles islands is low. Fossil calibrated divergence time estimates support an overseas dispersal of the L. seychellensis lineage to the Seychelles from mainland Africa ca. 43-25 million years before present (Ma), rather than this taxon being a Gondwanan relic.

Molecular phylogeny of Neotropical Parrot Snakes (Serpentes: Colubrinae: Leptophis) supports underestimated species richness.

Tetrapod taxa with broad geographic distributions across the Neotropics are often composed of multiple evolutionary lineages. In this paper, we present the most complete phylogeny of Leptophis to date and assess morphology-based species limits within the broadly distributed green parrot snake Leptophis ahaetulla sensu lato, which occurs from Mexico to Argentina. Although L. ahaetulla sensu stricto, L. nigromarginatus and L. occidentalis were recovered as paraphyletic, tree topology tests failed to reject their monophyly. Monophyly of L. bocourti, L. coeruleodorsus, L. cupreus, L. depressirostris, L. marginatus, L. riveti and L. sp. nov. was strongly supported. Our phylogenetic trees support recognition of multiple species within Leptophis ahaetulla sensu lato and suggest that color evolution and the uplift of the Andes played an important role in the diversification of parrot snakes.

Biogeographic barriers, Pleistocene refugia, and climatic gradients in the southeastern Nearctic drive diversification in cornsnakes (Pantherophis guttatus complex).

The southeastern Nearctic is a biodiversity hotspot that is also rich in cryptic species. Numerous hypotheses (e.g., vicariance, local adaptation, and Pleistocene speciation in glacial refugia) have been tested in an attempt to explain diversification and the observed pattern of extant biodiversity. However, previous phylogeographic studies have both supported and refuted these hypotheses. Therefore, while data support one or more of these diversification hypotheses, it is likely that taxa are forming within this region in species-specific ways. Here, we generate a genomic data set for the cornsnakes (Pantherophis guttatus complex), which are widespread across this region, spanning both biogeographic barriers and climatic gradients. We use phylogeographic model selection combined with hindcast ecological niche models to determine regions of habitat stability through time. This combined approach suggests that numerous drivers of population differentiation explain the current diversity of this group of snakes. The Mississippi River caused initial speciation in this species complex, with more recent divergence events linked to adaptations to ecological heterogeneity and allopatric Pleistocene refugia. Lastly, we discuss the taxonomy of this group and suggest there may be additional cryptic species in need of formal recognition.

Out of the Hengduan Mountains: Molecular phylogeny and historical biogeography of the Asian water snake genus Trimerodytes (Squamata: Colubridae).

The Asian water snake genus Trimerodytes (formerly Sinonatrix) is endemic to East and Southeast Asia. Although several species have been included in various phylogenetic studies previously, the evolution and relationships among members of this genus as a whole remain unexplored. In this study, we report the sequencing two protein-coding mitochondrial gene fragments (MTCYB and ND2) and three nuclear genes (c-mos, NT3, and Rag1), reconstruct interspecific phylogeny, and explore biogeography for the genus Trimerodytes. Both Bayesian inference and maximum likelihood analyses consistently recover the monophyly of Trimerodytes with strong support, with T. yapingi the sister-group to the remaining species. The divergence date and ancestral area estimation suggest that Trimerodytes likely originated in Hengduan Mountains (eastern Tibetan Plateau) in western China at 23.93 Ma (95% HPD: 17.09-31.30), and intraspecific divergence began at about 4.23 Ma (95% HPD: 2.74-6.10). Analyses support the validity of T. yunnanensis.

Sex linkage of the skeletal muscle sodium channel gene (SCN4A) explains apparent deviations from Hardy-Weinberg equilibrium of tetrodotoxin-resistance alleles in garter snakes (Thamnophis sirtalis).

The arms race between tetrodotoxin-bearing Pacific newts (Taricha) and their garter snake predators (Thamnophis) in western North America has become a classic example of coevolution, shedding light on predator-prey dynamics, the molecular basis of adaptation, and patterns of convergent evolution. Newts are defended by tetrodotoxin (TTX), a neurotoxin that binds to voltage-gated sodium channels (Nav proteins), arresting electrical activity in nerves and muscles and paralyzing would-be predators. However, populations of the common garter snake (T. sirtalis) have overcome this defense, largely through polymorphism at the locus SCN4A, which renders the encoded protein (Nav1.4) less vulnerable to TTX. Previous work suggests that SCN4A commonly shows extreme deviations from Hardy-Weinberg equilibrium (HWE) in these populations, which has been interpreted as the result of intense selection imposed by newts. Here we show that much of this apparent deviation can be attributed to sex linkage of SCN4A. Using genomic data and quantitative PCR, we show that SCN4A is on the Z chromosome in Thamnophis and other advanced snakes. Taking Z-linkage into account, we find that most apparent deviations from HWE can be explained by female hemizygosity rather than low heterozygosity. Sex linkage can affect mutation rates, selection, and drift, and our results suggest that Z-linkage of SCN4A may make significant contributions to the overall dynamics of the coevolutionary arms race between newts and snakes.

The Harderian gland transcriptomes of Caraiba andreae, Cubophis cantherigerus and Tretanorhinus variabilis, three colubroid snakes from Cuba.

The Harderian gland is a cephalic structure, widely distributed among vertebrates. In snakes, the Harderian gland is anatomically connected to the vomeronasal organ via the nasolacrimal duct, and in some species can be larger than the eyes. The function of the Harderian gland remains elusive, but it has been proposed to play a role in the production of saliva, pheromones, thermoregulatory lipids and growth factors, among others. Here, we have profiled the transcriptomes of the Harderian glands of three non-front-fanged colubroid snakes from Cuba: Caraiba andreae (Cuban Lesser Racer); Cubophis cantherigerus (Cuban Racer); and Tretanorhinus variabilis (Caribbean Water Snake), using Illumina HiSeq2000 100 bp paired-end. In addition to ribosomal and non-characterized proteins, the most abundant transcripts encode putative transport/binding, lipocalin/lipocalin-like, and bactericidal/permeability-increasing-like proteins. Transcripts coding for putative canonical toxins described in venomous snakes were also identified. This transcriptional profile suggests a more complex function than previously recognized for this enigmatic organ.

Shifts in Selective Pressures on Snake Phototransduction Genes Associated with Photoreceptor Transmutation and Dim-Light Ancestry.

The visual systems of snakes are heavily modified relative to other squamates, a condition often thought to reflect their fossorial origins. Further modifications are seen in caenophidian snakes, where evolutionary transitions between rod and cone photoreceptors, termed photoreceptor transmutations, have occurred in many lineages. Little previous work, however, has focused on the molecular evolutionary underpinnings of these morphological changes. To address this, we sequenced seven snake eye transcriptomes and utilized new whole-genome and targeted capture sequencing data. We used these data to analyze gene loss and shifts in selection pressures in phototransduction genes that may be associated with snake evolutionary origins and photoreceptor transmutation. We identified the surprising loss of rhodopsin kinase (GRK1), despite a low degree of gene loss overall and a lack of relaxed selection early during snake evolution. These results provide some of the first evolutionary genomic corroboration for a dim-light ancestor that lacks strong fossorial adaptations. Our results also indicate that snakes with photoreceptor transmutation experienced significantly different selection pressures from other reptiles. Significant positive selection was found primarily in cone-specific genes, but not rod-specific genes, contrary to our expectations. These results reveal potential molecular adaptations associated with photoreceptor transmutation and also highlight unappreciated functional differences between rod- and cone-specific phototransduction proteins. This intriguing example of snake visual system evolution illustrates how the underlying molecular components of a complex system can be reshaped in response to changing selection pressures.