Our laboratories have carried out an extensive analysis of the transcriptomes of venom glands of cone snails (family Conidae) (Safavi-Hemami et al., 2016; Li et al., 2017; Robinson et al., 2017; Li et al., 2018; Lu et al., 2020). This information has greatly facilitated the characterization of bioactive venom components. We have also analyzed a number of taxa in the superfamily Conoidea which are not cone snails; many of these belong to the family Turridae (Figures 1–3, Supplementary Tables S1, S2). However, species in some other conoidean families, including Terebridae and Drilliidae, have also been analyzed (see Figures 2, 3).

Using these datasets, we have developed a method for phylogenetic analysis based on common transcripts recovered from transcriptomic data; a detailed description of the procedures and protocols is provided in the Methods section. This approach should be applicable to any molluscan group. The results of this analysis, shown in Figure 3, include nine different species conventionally assigned to the genus Turris. The results clearly demonstrate that Turris, as conventionally defined, is polyphyletic, and that a subset of species is only distantly related to the major group. We describe a new genus, Purpuraturris gen. nov., to include this subset of outlier species (see Appendix 1). The phylogenetic tree in Figure 3 demonstrates that the family Turridae has two major divisions. The four species (cryptorraphe, cristata, undosa, nadaensis) formerly assigned to Turris are not in the same branch as the other five species in Turris (guidopoppei, hidalgoi, spectabilis, dollyae, normandavidsoni).

The proposed new genus includes these four species conventionally assigned to Turris, but for which the genetic data are incompatible with that assignment. The molecular results suggest that these species are unrelated to their former congeners in the genus Turris. Additionally, as a group these species are morphologically distinct and share a number of characters that justify creation of a new genus. Because they have some morphological characters similar to Turris (i.e., the position of the slit, the thin linear cleft in the margin of the aperture), but also have differentiating characters, notably a distinctive purplish color, particularly within the aperture, we have chosen the name Purpuraturris gen. nov. The type species of the new genus is Purpuraturris cryptorrhaphe (formerly Turris cryptorrhaphe). The transcriptomic data indicate that Purpuraturris nadaensis, Purpuraturris undosa and Purpuraturris cristata are congeners. This taxonomic placement is supported by the phylogenetic tree constructed using standard genetic markers (Figure 4). A conventional taxonomic description of the new genus is provided in Appendix 1.

The specific family of venom peptides that is the focus of this study are venom peptides from species formerly assigned to the genus Turris; these belong to a group known as the P-like turripeptides (Watkins et al., 2006). This nomenclature is derived from the number and arrangement of cysteine residues in these peptides, which are shared with the P-superfamily of conotoxins (Robinson and Norton, 2014). A description of the nomenclature adopted for turrid venom peptides is provided in Appendix 2. Based on their signal sequences, propeptides, and cysteine loop sizes, these can either be regarded as a single superfamily of venom peptides, or as two closely-related superfamilies (see PII and PIII superfamilies, Supplementary Table S1). This superfamily was chosen due to the presence of multiple venom peptides from both Turris and Purpuraturris. The predicted mature peptide sequences are shown in Table 1, while a tree that shows how the sequences are related, along with the alignment of the precursor sequences used to generate the tree are presented in Figure 5. The precursor sequences are listed in Supplementary Table S1, divided into P-like turripeptide superfamily PII, consisting of Clades I and II in Figure 5, and P-like turripeptide superfamily PIII, consisting of Clades III and IV in Figure 5.

It is clear from this tree that the peptide sequences fall into discrete branches, which we have labeled Clades I–IV. It is notable that each branch of the tree either has species from Turris s.s. (as redefined above), or only species that have been removed from the genus Turris and reassigned to the proposed new genus, Purpuraturris gen. nov. None of the branches have a mixture of species from both genera. Thus, the grouping of several species in a given branch is consistent with and supportive of the phylogenetic trees shown in Figures 3, 4. The peptide tree reflects the phylogenetic relationships between the various taxa analyzed: the species in Turris are not closely-related to species in Purpuraturris, leading to the pattern of divergence in their venom peptides observed in Figure 5. The tree of peptide sequences has two major branches; in one major branch, there are two clades, with each subgroup containing sequences only from Turris s.s., or Purpuraturris species; in contrast, the second major branch has two clades, both having only species assigned to Turris whereas peptides from Purpuraturris are absent. Thus, a juxtaposition of the trees based on evolutionary relationships between species in Figures 3, 4, and on peptide sequence similarities in Figure 5 suggests that the divergence between the two major divisions of the P-like turripeptides shown in Figure 5 may have predated the divergence of Purpuraturris from Turris. Alternatively, it is possible that the divergence of Clades III and IV may have predated the divergence of Clade II (Purpuraturris) from Clade I (Turris) with a subsequent loss of these two lineages from Purpuraturris.

We note that the sequences of the mature peptides in the two major branches have strikingly divergent features. The two clades in the major branch found only in Turris species (Clades III and IV) have P-like turripeptides each with a single amino acid between the first two cysteine residues; in Clades I and II there are six amino acids between the first two cysteine residues, with the consensus sequence of: (E/P) EN (E/L) (A/E) X.

The venom peptide family that is analyzed in Figure 5 and Table 1 has over 30 related sequences. In Figure 5, we have subdivided the major branches into four defined clades; the phylogenetic framework indicates that the peptides in Clade II should be considered separately from the peptides in Clade I, III and IV, since the divergence established by the phylogenetic trees in Figures 3, 4 indicate that these peptides may have evolved as part of an entirely different biological strategy, from two lineages of venomous animals that are long diverged. Nevertheless, Clades I, III and IV, still raise the issue of how to efficiently analyze this library of diverse peptides. Using phylogenetic methodology to compare sequences resulted in the discovery that peptides from Turris venoms that belong to the P-like turripeptide group fall into three clades; if sufficient resources were available, at least one peptide from each clade should be characterized, since it is likely that these may act through divergent mechanisms. Even within each clade of peptides, there are too many to comprehensively characterize. In the most complex group, Clade III, the mature peptide sequences can be clustered into two groups; the peptide sequences from Turris hidalgoi have representatives in each of these Clade III subgroups, i.e., Thd9.2 defining one group, and Thdg9.1 defining the second.

The use of a transcriptomic approach to building the phylogenetic tree automatically reveals relative expression levels (Figure 5C). If the expression levels of the two T. hidalgoi peptides in Clade III are compared, by far the most highly-expressed is Thd9.1. The high-expression level of Thd9.1 would more easily facilitate characterization and thus favor prioritization of this peptide for further study.

In the superfamily Conoidea, the increasing number of comprehensive transcriptomic analyses of venom glands has made it possible to routinely generate native venom peptide libraries. This is particularly facile for conoidean venoms because of how venom peptide precursors are organized. The typical conoidean venom gland has a repertoire of approximately 200 different peptides encoded and expressed by epithelial cells lining the gland in each species. The accelerated evolution of these venom peptide genes results in the mature peptide regions being hypervariable, except for the cysteine residues which provide a conserved framework for genetically related peptides. Because the signal sequences are also very highly conserved, it is straightforward to identify related peptides in each conoidean lineage. In this study, we demonstrate how the transcriptomic analysis of one lineage can be used to identify a peptide family, and in addition, to reveal phylogenetic relationships. Until recently, this type of analysis was significantly more difficult, since obtaining a comprehensive transcriptome was much more cumbersome and expensive.

The peptide family we analyzed in this study was obtained from the venom glands of species conventionally assigned to a single genus, Turris. One previously characterized group of venom peptides in turrids are the P-like turripeptides (Heralde et al., 2008), which have six cysteine residues in a pattern characteristic of Conus P-superfamily peptides (Robinson and Norton, 2014) but which are otherwise unrelated. There is more than one gene superfamily that encodes peptides of this type in these species (Watkins et al., 2006) — in this study, we focused on a subset of P-like turripeptides with related signal sequences. The precursor sequences, which are the data used to predict sequences of the mature bioactive venom peptides, are shown in Figure 5 and Supplementary Table S1. In some cases, P-like turripeptide precursor sequences were obtained using PCR or sequencing of cDNA libraries rather than from transcriptomes, indicated in Supplementary Table S1 by asterisks. In a few cases, P-like turritoxin precursor sequences were obtained from a transcriptome of a specimen other than that used for the neutral marker analysis in Figure 3. Those sequences are included in Supplementary Table S1, indicated by double asterisks. Expression levels of all P-like turripeptides identified in transcriptomes relative to the mean expression of HKG in the same transcriptome are shown in Figure 5C and provide potential insight into the functional role of each peptide within a species. The species assigned to the genus Turris sensu lato for which a transcriptome was obtained and analyzed in this work are shown in Figure 1.

In addition to venom peptide genes, transcriptomic analysis also provides sequences for all of the other genes expressed in the venom glands, such as housekeeping genes and genes involved in the normal cellular processes required to produce the venom components (i.e., proteins involved in secretion and membrane trafficking). This sequence information was used to carry out a “neutral” phylogenetic analysis of the species that are simultaneously being analyzed for the identity of venom components. Thus, transcriptomic analysis leads to the easy identification of peptide families from a particular lineage of venomous animals, provides neutral estimates of the phylogenetic relationships between the species being analyzed, and establishes a robust basis for the relationship of that species group to other taxa for which similar transcriptomic data are available. The results of the phylogenetic analyses of the species shown in Figures 1, 2 are discussed in the following section; the relationship of the species in Figure 1 to other conoideans was assessed, using transcriptomic information from the genera illustrated in Figure 2.

The phylogenetic tree in Figure 3, which was constructed using the transcriptomic approach described under Methods, includes four different conoidean families; all family assignments are firmly supported. The family Terebridae includes two genera, Terebra, represented by Terebra subulata and Terebra guttata, and Hastula, including Hastula nitida and Hastula matheroniana, that define two subbranches on the tree. Similarly, the family Drilliidae is represented by Drillia and Clavus, and the three genera in Conidae; Californiconus (Conus californicus), Conasprella and Conus, are well-supported in this phylogenetic analysis. Thus, the approach used is consistent with present generic assignments in the families Clavidae, Terebridae and Conidae (Bouchet et al., 2011).

In the family Turridae, species conventionally assigned to the genus Turris (shown in Figure 1) are split between two major branches of the family, with most species analyzed remaining in Turris, but some conventionally assigned to Turris have proven to be only distantly related. We have proposed a new genus, Purpuraturris gen. nov. (see the Appendix 1 and Figures 3, 4). The problem of proper generic assignments does not apply only to Turris, but also to other major genera in Turridae, such as Gemmula which could be addressed using similar methodologies in the future.

Linking primary amino acid sequences of venom peptides to the phylogeny of the animals in which they evolved is a prerequisite for integrating the biochemical and molecular information obtained from the standard approach used in venom peptide research with the broader biology (and even ecology) of the venomous animals that evolved that library of peptides. The phylogenetic correlation provided for the P-like turripeptides in Figure 5 illustrates a strategy for the further definition and characterization of these peptides. We have shown that there are four major groups of peptides based on sequence similarities, one found in those species that we have moved to the new genus, Purpuraturris gen. nov., and three found in species that remain within the genus Turris. The presence of the venom peptides found in Purpuraturris species but not in any Turris may suggest that these play a role in aspects of Purpuraturris biology divergent from the biology of true Turris species. The three divergent branches of the tree found in species of Turris sensu stricto in turn, suggests that these may play a functional role exclusive to the biology of members of that genus.

While we cannot definitively conclude that clade II P-like turripeptides are without exception absent from Turris species, nor that clade III and IV P-like turripeptides are universally absent from Purpuraturris species, by differentiating the bioactivity of a P-like turripeptide representative from each branch of the phylogenetic tree in Figure 5, an understanding of the differences in underlying molecular mechanisms can be achieved more efficiently. A correlation of the bioactivity of the various peptides with the biology of the snails that evolved these venom peptides helps to build a foundation for understanding how a peptide in a particular venom contributes to the biology of that species. Furthermore, details of the effects of the peptide on a target animal (presumably prey, predator or competitor) can connect molecular mechanisms to the systems level and behavior. The distribution of homologous peptides across a particular lineage, in this case, the species in the genus Turris, will suggest whether that specific venom peptide is involved in a biological strategy shared widely across the genus Turris by many species, or whether a particular type of peptide represents an innovative facet of a subset of species in that lineage. Behavioral observations of the species in a particular lineage may then be correlated with particular venom peptides. For example, phylogenetic analysis of A-superfamily toxins from fish-hunting cone snails revealed that distinct clades of toxins correlate with distinct subtype selectivity profiles at the nicotinic acetylcholine receptor (nAChR) (Puillandre et al., 2010).

The assessment of which peptides in a large venom peptide library should be selected for further analysis can be viewed as a three-step process. The first step is to obtain the library of venom peptide sequences from a lineage of related venomous animals, in this study, the species conventionally assigned to a single genus, Turris. The phylogenetic analysis revealed that a group of four species were not congeneric, but in fact, were only distantly related to species correctly assigned to Turris. By carrying out the analysis of both groups of species and analyzing the peptide sequences obtained, peptides more closely related to each other could be grouped into discrete clusters. Each cluster likely represents a different underlying mechanism for bioactivity, and analyzing at least one representative of each cluster would be a priority goal.

Some clusters however, are represented by multiple peptides even within the venom of a single species; these may have related molecular targets. If the potential macromolecular targets are of high interest, it may be desirable to investigate the subgroups within each cluster of related peptides. Investigation of dissimilar peptides based on branch lengths in the venom peptide phylogenetic tree may also help prioritize the selection peptides likely to have differences in their bioactivities. This initial phase of analysis is based purely on the peptide sequences initially obtained. Additional input data for developing an experimental strategy are the level of expression of peptides that define a particular cluster and the ability to carry out a thorough proteomic analysis so that the native venom peptide can be compared as comprehensively as possible with a peptide that is either synthesized chemically or produced through recombinant expression. In the case of conoidean venom peptide libraries, this is particularly important because of the frequent presence of posttranslational modifications.

An additional consideration for which peptides to characterize is which venom can be obtained in greater quantity and without negative ecological impact; in addition to the phylogenetic data and the level of expression derived from the transcriptomic data, a proteomic analysis confers significant advantages, and this generally requires a larger venom sample than is required for the other assessments. Thus, in the Thd9.1 cluster, whether it should be the peptide from T. hidalgoi or from T. spectabilis that is investigated depends on which venom is more easily available. The more abundant venom can be fractionated using the standard HPLC methods to purify individual peptides, and a mass spectrometric analysis with appropriate proteomic bioinformatics should identify peptides in the crude venom related to the Thd9.1 cluster. This analysis is important because it would uncover any posttranslational modifications, which is essential information required before initiating a comprehensive biochemical characterization. Furthermore, if a sufficient amount of venom were available to purify the native venom peptide in large amounts, after the peptide has been chemically synthesized (or produced by recombinant methods), a thorough comparison between the synthetic and the native peptides can be carried out, both with respect to biochemical identity, as well as to demonstrate similar bioactivity. These experiments are greatly facilitated by having larger amounts of native material readily available.

Thus, even before any biochemical work is initiated, the venom peptide library can be better assessed if the following were available: 1. A phylogeny of species in the lineage from which the peptide library was obtained. 2. A tree clustering individual members of the library into groups in order to assess their relatedness, 3. The level of expression of each venom peptide in the library. 4. A proteomic analysis of the venom to identify post-translational modifications and proteolytic cleavage sites, 5. Purified native peptide from venom for subsequent structural and functional comparison to a synthetic analog. Together, these can be used to prioritize which specific peptides will be the leads for chemical synthesis and biological characterization. Assembling these data leads to an informed decision of which peptides in the venom peptide library should be investigated. Since the characterization of each individual peptide is a major commitment of time and resources, this is a key decision point.

Once a significant quantity of the authentic bioactive peptide becomes available, either through chemical synthesis or recombinant expression, assessing potential biomedical applications can then be initiated. Several major factors influence how this assessment should be carried out: which biomedical applications are of the highest priority to the researchers (e.g., cancer applications, pain, neurodegenerative disease), what range of relevant assays are accessible, and the quantity of peptide available for the initial characterization. Since these factors vary considerably, which peptides to focus on for biomedical applications such as drug development will clearly be dictated by the results of the initial characterization. Venom peptides that have a novel activity and that score on an assay that has relevance for specific biomedical applications could then be prioritized for further comprehensive characterization and assessment.

RNA isolation and cDNA synthesis were performed as previously described for Iotyrris olangoensis (Watkins et al., 2006). First strand cDNA from Turris normandavidsoni was used as template to amplify the genes encoding PII and PIII turripeptides. Oligonucleotide primers were designed from untranslated region sequences flanking the open reading frame of PII turripeptide-encoding clones identified from Iotyrris olangoensis (Watkins et al., 2006) cDNA to amplify PII turripeptide sequences from Turris normandavidsoni (Tnr9.2; Table 1, Supplementary Table S1): 5’ (GAA CCR GCC AGC RAG ATG GG) and 3’ (GCC ATC AGT TGA CAA CAC GTG). Primers to the signal sequence and 3′UTR regions of PIII turripeptide-encoding genes identified from Iotyrris olangoensis (Watkins et al., 2006) were used to amplify PIII turritoxin sequences from Turris normandavidsoni (Tn9.3 and Tnr9.4; Table 1, Supplementary Table S1): 5’ (GAA CCR GCC AGC RAG) and 3’ (GCC ATC AGT TGA CAA CAC GTG). Polymerase chain reactions were performed with Advantage 2 polymerase Mix as described for mitochondrial gene products. The PCR cycling profiles are as follows: Initial denaturation (95°C, 60s); followed by 40 cycles of denaturation (95°C, 20s); annealing (58°C, 20s) and extension (72°C, 30s).

The resulting PCR products were purified, annealed to pNEB206A vector and transformed into competent cells, as described for mitochondrial gene products.

The nucleic acid sequences of these PII and PIII turritoxin-encoding clones were determined by automated sequencing (Core Sequencing Facility, University of Utah, United States).

Venom gland transcriptomes were generated as previously described (Li et al., 2018). Briefly, adapter clipping and quality trimming of raw reads were performed using fqtrim software (Version 0.9.4, http://ccb.jhu.edu/software/fqtrim/) and PRINSEQ (Version 0.20.4) (Schmieder and Edwards, 2011). After processing, sequences shorter than 70 bps and those containing more than 5% ambiguous bases (Ns) were discarded. De novo transcriptome assembly was performed using Trinity Version 2.0.5 (Grabherr et al., 2011) with a kmer size for building De Bruijn Graphs of 31, a minimum kmer coverage of 10, and a minimum glue of 10. Assembled transcripts were annotated using Blastx [(NCBI-Blast-2.2.28+, (Altschul et al., 1990)] against conotoxin sequences extracted from the ConoServer (Kaas et al., 2008) and (UniProt Consortium, 2015). The species, number of transcripts, toxins identified, toxin transcript expression values (tpm) and the number of toxin superfamilies identified in each of the 35 transcriptomes used here are listed in Supplementary Table S3. These transcriptomes were assessed using BUSCO analysis (busco:v5.2.2_cv2) (Manni et al., 2021) with the database mollusca_odb10 (5295 entries) and found to be of comparable quality (Supplementary Table S3). Common genes that were shared between all samples were identified using the blast identities of assembled contigs (e < 10^–4). This set of 66 housekeeping genes (HKG), half of which are ribosomal genes, is listed in Supplementary Table S4, along with Uniprot accession numbers. For each HKG, we aligned the nucleotide sequences from all species using mafft v7.407 (Katoh et al., 2005) (with the “auto” flag). For each alignment file we used the RAxML v8.2.12 (Stamatakis, 2006) program with a GTRGAMMA model (unpartitioned) to determine the “best” tree. Finally, we used Astral (version 5.7.3) (Zhang et al., 2018) to generate a consensus tree from all the individual single-gene trees. The HKG dataset generated for this study is included in the DRYAD repository (DOI): doi:10.5061/dryad.ksn02v751) as fasta files, aligned nt mafft fasta files and newick formatted “best” trees.

The identities of the turrid specimens used for transcriptome-based phylogeny were confirmed using the COI barcode sequence, in the same manner as described in Methods for standard phylogeny. These GenBank accession numbers are listed in Supplementary Table S2 and are consistent with the COI barcode sequences generated for Figure 4 and listed in Supplementary Table S5.

Phylogenetic analysis was performed using concatenated 12S and 16S rRNA and the mitochondrial cytochrome oxidase I (COI) barcode sequences for each species as previously described (Nam et al., 2009, Puillandre et al., 2014). GenBank accession numbers for these standard phylogenetic markers, which were used to generate Figure 4, are listed in Supplementary Table S5.

Multiple sequence alignment was performed using MAFFT version 7 (RRID:SCR_011811) (Katoh et al., 2017), maximum likelihood tree reconstruction was performed using IQ tree (RRID:SCR_021163) (Nguyen et al., 2015), consisting of best-fit model TVM+F+I+G4 using ModelFinder (Kalyaanamoorthy et al., 2017) and ultrafast bootstrap approximation (Minh et al., 2013) and the tree was visualized using iTOL version 6.4 (RRID:SCR_018174) (Letunic and Bork, 2021). The alignment file has been included in the supplementary data as Supplementary Data Sheet S1.

The sequences of the P-type turripeptide precursors were aligned using the MAFFT multiple sequence alignment program (RRID:SCR_011811) (wasabi.org) (Veidenberg et al., 2015) and refined by eye. Phylogenetic analysis was performed using IQ tree (RRID:SCR_021163) (Nguyen et al., 2015), consisting of best-fit model JTTDCMut+G4 using ModelFinder (Kalyaanamoorthy et al., 2017) and ultrafast bootstrap approximation (Minh et al., 2013) and the tree was visualized using iTOL version 6.4 (RRID:SCR_018174) (Letunic and Bork, 2021). The phylogenetic tree, which is rooted at the midpoint, is shown in Figure 5A. The alignment has been included in Figure 5B. Complete peptide precursors are shown in Supplementary Table S1, predicted mature peptides are shown in Table 1, and relative expression levels, where available, are shown in Figure 5C. The sequences were deposited into the GenBank nucleotide database, accession numbers: OK247621-OK247652.