Integrative multi-omics analysis reveals the contribution of neoVTX genes to venom diversity of Synanceia verrucosa

doi:10.21203/rs.3.rs-5243476/v1

Download PDF

Research Article

Integrative multi-omics analysis reveals the contribution of neoVTX genes to venom diversity of Synanceia verrucosa

https://doi.org/10.21203/rs.3.rs-5243476/v1

This work is licensed under a CC BY 4.0 License

You are reading this latest preprint version

Background

Animal venom systems have been considered as valuable model for investigating molecular mechanisms underlying phenotypic evolution. The stonefish were reported as the most venomous and dangerous fish due to sever human envenomation and occasionally fatality, whereas the genomic background of their venom remained under investigated and poorly explored compared with other venomous animals.

Results

In this study, we followed modern venomic pipelines to decode the Synanceia verrucosa venom components. A catalog of 478 toxin genes were annotated based on our assembled chromosome-level genome. Integrative analysis of the high-quality genome, transcriptome of venom gland and proteome of crude venom revealed a mechanism underlined the venom complexity in S. verrucosa. Six tandem-duplicated neoVTX subunit genes were evidenced as the major source for the neoVTX protein production. Further isoform sequencing enabled us to uncover massive alternative splicing events with a total of 411 isoforms demonstrated by the six genes, further contributing to the venom diversity. We then characterized 12 dominantly expressed toxin genes in the venom gland, and 11 of them were evidenced to produce the venom protein components, with the neoVTX proteins as the most abundant for granted. Other major venom proteins included a presumed CRVP, Kuntiz-type serine protease inhibitor, calglandulin protein, and hyaluronidase. Besides, a few of highly abundant non-toxin proteins were also characterized and they were hypothesized to imply housekeeping or hemostasis maintaining roles in the venom gland. Notably, a gastrotropin like non-toxin proteins ranked as the second highest abundant proteins in the venom, which had never been reported in other venomous animals, contributing to the unique venom property of S. verrucosa.

Conclusions

The results decoded the major venom composition of S. verrucosa, and highlighted the contribution of neoVTX genes to venom composition diversity by demonstrating tandem-duplication and alternative splicing. The diverse neoVTX proteins in the venom as lethal particles are hypothesized to be pivotal to understand adaptive evolution of S. verrucosa. Further functional studies are encouraged to exploit venom components of S. verrucosa for pharmaceutical innovation.

Multi-omics

venom diversity

Synanceia verrucosa

tandem duplication

alternative splicing

Venom is one of the most successful convergent evolutionary novelties across all major animal lineages for defense, predation, and competition[1]. Animal venom systems have been considered as valuable model for investigating molecular mechanisms underlying phenotypic evolution, and widely used for exploring and engineering biologically active molecules for pharmaceutical innovation[2].

Modern venomics are integrative research involves genetics, evolution, ecology as well as translational application, facilitating our understanding and exploiting of the venomous animals including scorpions, snakes and spiders[3, 4]. Most of the animal venom are highly complexed mixture of molecules with modifications at the genomic, transcriptomic and protein level during evolution. Recently the emerging chromosome level genomes greatly motivated venom toxins screening and gene identification in snake[5–7], spider[8, 9], and other venomous animals[10–13]. Large-scale comparative genomic analysis provided evolutionary views of toxin genes involved co-option of pre-existing nontoxic genes, massive tandem duplication and domain loss[14–17]. Despite different venom composition and apparatus spanning venomous metazoans, a comparative analysis of 20 species’ venom gland transcriptomes from wasps, spiders to fish and snakes revealed that a conserved genetic toolkit was recruited to evolve the ability to secret toxins[18]. On the other hand, full-length transcriptomic analysis found that more than half of the venom genes generated multiple isoforms within the parasitoid wasp venom gland, suggesting that alternative splicing highly contributed to venom diversity and evolution[19]. In the future, the integrative genomic, transcriptomic and proteomic approaches will greatly motivate novel toxin discovery and venom evolution deciphering.

It has been reported that nearly 2500 fish species are venomous across the world, whereas most of their venom were left as under investigated and unexplored resource[20]. The stonefish (Scorpaeniformes, Synanceiidae, Synanceia) were recognized as the most venomous and dangerous fish due to sever human envenomation and occasionally fatality[21, 22]. There are five species of stonefish, referring to S.verrucosa, S.horrida, S.alula, S.nana, and S.platyrhyncha. Encompassing the development of venomic approaches, attempts have been denotated to decode the stonefish venom composition and functional activities[23] The inceptive research in stonefish venom was predominantly reported in S.horrida, with the first fish-derived protein toxin characterized as SNTX, inspiring the venom research in the teleost[24–26].

As the most widely distributed stonefish species across the world, S. verrucosa demonstrated great potential in understanding venom pathology and toxin gene evolution in fish[27, 28]. The lethal particles have been purified in the venom at early 1990s, including a tetramer protein VTX and a dimer neoVTX, which both demonstrated hemolytic, hypotensive and cardiotoxic activities[29–31]. Besides, two widely existed components in the animal venom, hyaluronidase[32] and lectin[33], have also been characterized in S. verrucosa venom. The above work were limited to single particle-level purification and molecular cloning, leaving most of the venom components remained uncharacterized, which hindered our understanding that how the venom took poisonous actions and how the different components added to the physiological effects of envenomation.

Two criteria were widely adopted to identify major functional venom components across animals: (1) The genes should exhibit highly reliable expression in the venom gland; (2) The highly expressed genes should be further characterized in the venom proteome. In this work, we assembled chromosomal level genome of S. verrucosa with a set of toxin genes preliminary annotated. The major venom protein components were then deciphered and characterized by comprehensive analysis with genomic, transcriptomic and proteomic data. Notably, we revealed a mechanism of tandem-duplicated neoVTX genes underwent massive alternative splicing, which underlined the venom diversity and evolution in S. verrucosa.

Animals

Adult S.verrucosa fish were purchased from local fishery market and kept in the laboratory aquaria system for at least 2 weeks before sampling. The fish were anesthetized using MS222 (General reagent, China) before dissection. Muscle, skin and venom gland were dissected and rinsed by autoclaved artificial sea water, followed by flash frozen in liquid nitrogen. The samples were then stored at -80 ℃ freezer. All the experiments were conducted in accordance with guidelines approved by the Animal Ethical and Welfare Committee of Southern Marine Science and Engineering Guangdong Laboratory (Guangzhou).

Genome extraction, sequencing

Genomic DNA was extracted from the muscle using the DNeasy Blood and Tissue Kit (Qiagen, USA). The integrity, purity and concentration of extracted DNA was analyzed by gel electrophoresis and NanoDrop™ (Thermo Fisher Scientific, USA). DNA concentration was further quantified using the Qubit® 4.0 Fluorometer (Invitrogen, USA). A total of 15 µg DNA were then fragmented for HiFi sequencing library construction and the library was then sequenced using the PacBio Sequel II series sequencer (Pacific Biosciences, CA, USA). At the same time, whole genome sequencing (WGS) was included to facilitate genome assembly. Briefly, a library with 300 bp inserted size was prepared using genomic DNA from muscle followed the protocol of the MGIEasy FS DNA Library Prep Kit (MGI-Tech, China) and sequenced on the DNBSEQ-T7 platform (MGI Tech, China).

For Hi-C sequencing, the fresh muscle tissue were fixed with 2% formaldehyde to enable cross-linking between genomic DNA and chromatin. The crosslinked DNA were then extracted, and a library with 150-bp inserted size was constructed. The sequencing were performed using the DNBSEQ T7 platform (MGI Tech, China).

RNA extraction and sequencing

Total RNA was isolated from muscle, skin and venom gland (n = 7) using TRIzol reagents (Thermo Fisher Scientific, USA). The RNA was analyzed by NanoDrop™ for concentration, and Bioanalyzer 2000 (Agilent, USA) for integrity. For next generation sequencing (NGS), the cDNA libraries were constructed following the manufacturer’s recommendations and sequenced using the MGI DNBSEQ-T7 sequencing platform. Further, a venom gland full-length transcriptome sequencing was performed according to the pipeline of Iso-Seq PacBio platform (Pacific Biosciences, CA, USA).

Genome assembly and annotation

To generate a chromosomal-level genome of S.verrucosa, a high-quality contig genome assembly was generated using Hifiasm v0.16.1[34], followed by redundant sequences removal using Purge_dups v1.2.5[35]. The assembled contigs were then aligned to the Hi-C sequencing data of S.verrucosa with default parameters using chromap[36], followed by genome assembly with default parameters using the yahs pipeline[37]. Finally, the completeness of the genome assembly was assessed by BUSCO v5.4.7 with genome mode referring to the database of actinopterygii_odb10[38].

The repetitive sequences in the genome were screened using RepeatModeler v2.0[39]. The S.verrucosa genome sequences were then mapped against de novo repeat libraries and Repbase TE library v16.02 to identify transposable elements(TEs) using RepeatMasker v4.1.1[40], and TE-related proteins were further explored using RepeatProteinMask v.4.0.6[41]. Besides, Tandem Repeat Finder (TRF) v.4.07 was employed to search for tandem repeat sequences[42].

Homology-based prediction, transcriptome-based prediction and de novo prediction were employed to annotate genes in the genome, with the combined use of EvidenceModeler (EVM) v1.1.1 and PASA v2.4.0[43]. For homology-based predictions, the protein sequences of Danio rerio, Bufo gargarizans, Cyclopterus lumpus, Pseudoliparis swirei, Sebastes umbrosus, Gasterosteus aculeatus, Scatophagus argus, Thalassophryne amazonica, Siniperca chuatsi, Oncorhynchus mykiss, Oryzias latipes, Salmo salar, Oreochromis niloticus, Electrophorus electricus, Ostracion cubicus, and Acanthochromis polyacanthus obtained from NCBI were aligned to the assembled S. verrucosa genome[44]. For transcriptome-based prediction, RNA-seq data were aligned to the genome with default parameters using Hisat v2.1.0[45]. Subsequently, transcripts were reconstructed with default parameters using StringTie v2.0[46]. For de novo gene predictions, Augustus v3.4.0, Genscan v3.1, and GlimmerHMM v3.0.1 were utilized to analyze the repeat-masked genome[47–49].

The functional annotation of genes were conducted using Blastp and Blastn (-e 1e-5) to blast protein or nucleotide sequences against various databases including the NCBI nonredundant protein (NR), the NCBI nonredundant nucleotide sequence (NT), as well as SwissProt databases. Gene Ontology (GO) terms were then retrieved and assigned to S. verrucosa query sequences, and enzyme codes (EC) corresponding to GO terms were mapped to KEGG pathway annotations for further functional insights.

Long-read RNA alignment and isoform identification

The raw reads obtained by the full-length transcriptome sequence were analyzed following the Iso-Seq3 pipeline (v3.1). Briefly, circular consensus sequences (CCSs) were created from raw subreads, and the cDNA primers were removed using lima. Full-length reads were then filtered, clustered, and polished (Iso-Seq3, v3.1). The polished sequences were aligned to the S.verrucosa genome using pbmm2 (v1.13.0). Isoforms identification was conducted using TAMA[50]. Redundant isoforms were removed by comparing with the raw data of NGS transcriptome data.

An in-depth characterization of isoforms and removal of artifacts was completed using SQANTI3[51]. Briefly, isoform characterization was performed using sqanti_qc.py, and filtering was completed using sqanti_filter.py. The filtered data set was rerun through sqanti_qc.py for the final confirmation. The SQANTI filtered isoforms were used for the remaining analyses.

Expression analysis of the genes in the venom gland

HISAT2 v2.1.0 was used to align the quality-controlled RNA-seq reads to the genome of S. verrucosa. Aligned reads were assembled into transcripts, and expression was quantified using featureCounts[52]. Differential expression analysis was conducted using DESeq2 v1.26.0 to identify genes that were differentially expressed in each tissue type[53]. Genes with an adjusted p-value < 0.05 and a log₂(fold change) > 1 were considered to be differentially expressed. The visualization of differentially expressed toxin genes across all samples was performed using the pheatmap package v1.0.12 in R. Transcript per million (TPM) values were used to assess the expression level among the venom glands, muscles, and skin of S. verrucosa.

Toxin genes annotation

Due to the propensity of toxin gene families being organized in extensive tandem arrays, the number of paralogs in particular gene families may be compromised, and thereby impeding the accuracy of gene annotation[54]. Consequently, additional annotation steps were included in this study. Combining the ToxProt based annotation with the isoforms generated from the third-generation sequencing (TGS) transcriptome of venom gland, we performed empirical annotation using FGENESH (http://www.softberry.com), as well as manual annotation using IGV-GSAman (https://gitee.com/CJchen/IGV-sRNA). Finally, each gene re-annotated was searched against ToxProt database back to confirm its identity and further confirmed its homology through phylogenetic analysis. We also downloaded particular HMM models from Pfam database v35.0[55], and then blasted the annotated toxin genes set using hmmsearch in HMMER v3.1b1 (http://hmmer.org/). Sequences with E-values below the default inclusion threshold (E-value < 1e-05) were retrieved and identified as toxin homologs.

Proteomic analysis of the Venom

100 µg of lyophilized venom was dissolved in a Triethylammonium bicarbonate buffer. 10 mM of tris (2-carboxyethyl) phosphine were then added and the venom solution were kept at 37°C for 60 minutes to reduce the disulfide bonds in the proteins. Subsequently, 40 mM of iodoacetamide was supplemented to block the sulfhydryl groups of cysteine residues, preventing the reformation of disulfide bonds. Followed by the addition of pancreatic trypsin at an enzyme-to-protein ratio of 1:50, the venom was digested overnight at 37°C. The obtained peptides were then desalted using HLB and quantified (Thermo Fisher Scientific). The samples were then dissolved in the sample loading buffer (2% acetonitrile, 0.1% formic acid) at a concentration of 0.25 µg/µL for mass spectrometry analysis. The samples separated by high-performance liquid chromatography were analyzed using the timsTOF Pro2 mass spectrometer (Bruker, Germany) in data-dependent acquisition (DDA) mode. The ionization mode is set to positive, with an ion source voltage of 1.5 kV, and both MS and MS/MS data were acquired using time-of-flight (TOF) analysis. The mass spectrometry scanning range was set to 100–1700 m/z. The data acquisition mode utilized parallel accumulation serial fragmentation (PASEF), where 10 rounds of PASEF were performed after each full MS scan. Second-level spectra with charge states ranging from 0 to 5 were acquired, and a dynamic exclusion time of 24 seconds was set for tandem mass spectrometry scans to avoid repeated ion selection. The Scorpaeniformes protein database was searched using MaxQuant v2.4.2.0, and the results were filtered based on the parameter Peptide FDR ≤ 0.01.

For native-PAGE electrophorese, the lyophilized venom was dissolved in a phosphate buffer solution (10 mg/mL) and centrifuged at 5,000g for 15 minutes at room temperature to remove insoluble material and primarily cellular debris. The equivalent of 50 µg of venom was then subjected to preparative native-PAGE (6%) gels after dilution with Tris-glycine buffer (0.05 M, pH 8.0). The gels were run at a constant voltage of 120 V at 4°C. Gel strips were stained with Coomassie brilliant blue R-250 for 30 minutes and subsequently destained using water. The protein bands with interested molecular weight were cut and sent for mass spectrometry analysis (Sangon Biotech, China).

Chromosome-scale de novo assembly and toxin genes annotation

High quality genome greatly facilitated the genetic decoding of animal venom components. In this work, a total of 19.80 Gb HiFi and 154.59 Gb Hi-C library sequencing data were generated respectively, resulting in an assembled genome size of 677.90 Mb for S.verricosa (Fig. 1A). The assembled genome was further evaluated, and the results showed that high quality of S.verricosa genome were obtained in this study (Fig. 1B and Supplementary Fig. 1).

The toxin gene identification was performed with a combined analysis with the assembled genome and venom gland transcriptomic data as described in the Methodology. A total of 478 potential toxin genes were preliminary identified at the whole-genome level. Most of them clustered to 56 gene families, with the peptidase S1 family (102), venom metalloproteinase family (36), the SNTX/VTX toxin family (32), the type-B carboxylesterase/lipase family (19), the ficolin lectin family (18), the true venom lectin family (18), peptidase A1 family (17), venom Kunitz-type family (17) as the most abundant toxin gene family (Fig. 1C, Supplementary Data 1). Further GO annotation revealed that two thirds (321 out of 476) of the genes showed explicit toxin activity (GO:0090729). The other 155 annotated genes were assigned to other functional categories, including calcium channel regulator activity (GO:0005246), killing of cells of another organism (GO:0031640), hyalurononglucosaminidase activity (GO:0004415), serine-type endopeptidase inhibitor activity (GO:0004867) (Supplementary Fig. 2).

Toxin genes usually demonstrated remarkable tandem duplication events during evolution. In the present study, a total of 13 toxin genes were characterized as tandem repeated genes in the genome of S. verricosa, including the SNTX/VTX genes, true venom lectin genes, peptidase S1 genes, glycosyl hydrolase 56 genes and etc. Especially, 8 annotated SNTX/VTX genes were clustered on chromosome 14 (Fig. 1D). The front 6 SNTX/VTX genes showed high sequence similarity (> 85%) with the reported neovtxa (GenBank: AB262392.1) and neovtxb (GenBank: AB262393.1), (Fig. 1D and Supplementary Fig. 3). The six genes were arranged in the following order and named as: neovtx-b1, neovtx-b2, neovtx-a1, neovtx-a2, neovtx-b3, neovtx-a3. Interestingly, the other two SNTX/VTX genes (neovtx-like1 and neovtx-like2) clustered in this region shared a much lower sequence identity (< 50%) and clustered into a separate branch in the phylogenetic tree, whereas InterPro scan showed that they all demonstrated similar conserved domains, including Stonustoxin helical domain (IPR048997), Thioredoxin_11 domian (IPR040581), PRY domain (IPR006574) and SPRY domain (IPR003877) (Fig. 1D).

Toxin genes expression in the venom gland

Next-generation sequencing (NGS) of transcriptome was performed using venom gland (n = 7) to further assist the toxin genes characterization. A total of 14, 380 genes were considered as actively expressed in the venom gland with mean transcripts per million (TPM) > 1. Further analysis revealed that 240 of them corresponded to the previous annotated toxin genes set (Supplementary Data 1).

We then characterized a number of dominantly expressed toxin genes (DETGs) with the TPM > 500 in the venom gland. A total of 12 toxin genes were designated as DETGs, including one cysteine rich venom protein (CRVP) gene (sv_crvp1), the six tandem-repeated neoVTX subunit genes on chromosome 14, one hyaluronidase gene (sv_hyal1), one translationally-controlled tumor protein gene (sv_tpt1), two calglandulin protein genes (sv_caglp1 and sv_caglp2), and one Kunitz-type serine protease inhibitor gene (sv_kspi1) (Fig. 2A and B). The 12 DETGs accounted for 93.70% of the totally expressed toxin genes. Meanwhile, we also looked into the highly expressed non-toxin genes in the venom gland, which were suggested to be pivotal to regulate the toxin gene expression and/or produce important non-toxic components in the venom. The highest 10 non-toxin genes expressed in the venom gland were shown in Fig. 2B, including sv_fabp6, sv_tmsb4, sv_eef1a, sv_gpx1, sv_rpl39, sv_clu, sv_rpl26, sv_cstb, sv_ybx1 and sv_krt14. Interestingly, 4 unannotated genes with TPM > 500 were identified in the venom gland, suggesting them as potential unique genes in S. verrucosa not previously found in other venomous animals (Fig. 2B).

Comparative analysis was further performed between the venom gland (VG) and its closely related anatomic tissues, including skin (SK) and muscle (MU). The results showed that the VG demonstrated distinct expression profiles, and a higher transcriptomic similarity with the SK than MU (Supplementary Fig. 4). The results indicated that skin may serve as the potential evolutionary origin of venom gland in S. verrucosa. Differentially expressed genes (DEGs) in the venom gland compared with other tissues usually evidenced their preferred contributing roles to the venom. A total of 2,700 genes showed significantly higher expression in the VG (log₂ Fold Change > 1 and padj < 0.05) compared to the SK. By blasting with the actively expressed 240 toxin genes, we obtained a total of 87 VG biasedly expressed toxin genes including 10 of the 12 previously characterized DETGs, with the 6 neoVTX subunit genes demonstrated the highest differential expression level (Fig. 2. C and Supplementary Data 2). The other 2 DTEGs, the sv_tpt1 and sv_caglp2 (Sver.16G.00747), exhibited high expression levels in VG, as well as in SK, indicating their potential non-toxic roles outside venom gland. We further investigated the DEGs of non-toxin genes set (with TPM > 500 in VG) between the VG and SK. A total of 84 genes were characterized and predominantly showed molecular function of enzyme inhibitor, antioxidant and peroxidase activities, and involved in biological process of detoxification, metabolism and biosynthesis (Fig. 2D and Supplementary Data 3). They were suggested to be partially involved in the resistance of S. verrucosa to its own venom.

Proteomic evidence of the toxin gene in venom

As the major components of the venom, the proteins produced by the previously identified genes should be further characterized in the venom proteome. A 4D-label free proteomic approach was employed to analyze the venom protein components of S. verrucosa. A total of 662 proteins were recognized, and 33 of them hit the previously annotated toxin genes (Supplementary Data 4). We characterized the major venom proteins as two categories, toxin proteins and non-toxin proteins. Quantitative analysis showed that the toxin proteins accounted for 32.24% of the total venom proteins (Fig. 3A). We then analyzed the proteomic evidence of the 12 DETGs. As shown in Fig. 3B, the most abundant toxin proteins were neoVTX proteins, accounting for about 18.69% of the total venom protein. Corresponding to the venom gland transcriptomic data, almost all the neoVTX proteins came from the 6 tandem-repeated neoVTX subunit genes. Among them, neoVTX a2 (47.15%) and neoVTX b3 (46.98%) provided the most of the neoVTX proteins, and only a small proportion came from the other 4 genes, indicating their dominant roles in protein production (Fig. 3C). Interestingly, the native-PAGE analysis showed that two major protein bands with molecular weights around 160 and 420 kDa were both identified as containing neoVTX proteins by mass spectrometry (Fig. 3D). Proteins from the other DETGs, including helothermine like protein (SvCRVP) from sv_crvp1 (9.19%), Kunitz type serine protease inhibitors (SvKSPI) from sv_kspi1 (1.43%), the calglandulin like protein (SvCAGL, 0.60%), and hyaluronidase (0.49%) were also evidenced in the venom proteome (Fig. 3B). The above 11 DTEGs accounted 93.75% of all the designated toxin proteins, served as the major toxin genes in S. verricosa (Fig. 3E). The only DETGs with no proteomic evidence in the venom was tpt1, indicating its potential major role in venom gland otherwise as the venom components.

Other notable toxin proteins in the venom included homologues of putative endothelial lipase (Uniprot ID: J3RZ81), lactotransferrin (K9IMD0), snaclec cosgulation factor IX/X binding protein subunit A (SLXA, Q9DG39), snaclec B1 (B4XT00), peroxiredoxin-4 (P0CV91), and cobra venom factor (Q91132), whose abundance were higher than the lowest expressed DETG (neoVTX-b2) in the VG. As for the non-toxin protein, gastrotropin protein (13.14%) ranked the highest abundance in the venom of S. verricosa. The other highly expressed (abundance > 1%) non-toxin proteins were EI24, GPx1, PTGS1, cystatin-B, TCAF1, BHMT1, Actin and CSRP1. (Fig. 3B).

Alternative splicing contributed to venom diversity

Alternative splicing (AS) of toxin genes was proposed as the important molecular mechanism underlined venom diversity and evolution. To further explore the complexity of the venom components at the transcriptional level, an integrative analysis pipeline of ISO-Seq and NGS and was performed with the same RNA samples from the venom gland (Fig. 4A).

A total of 379,860 ISO-Seq reads were obtained, with a N50 length of 2,408 bp. All the reads were then aligned to the above assembled genome and 374,821 of them hit corresponding target genome block. A total of 27,600 unique isoforms were retained, corresponding to 10,945 gene loci. These unique isoforms were further blasted with the 476 annotated toxin genes, and a total of 825 cis-splicing isoforms were mapped to 191 of them (Fig. 4B, Supplementary Data 5). Further analysis showed that the 12 DETGs accounted for 52.48% of the isoforms among all the toxin genes (Fig. 4C). The results indicated that AS predominantly occurred in the highly expressed toxin genes in the venom. Among all the toxins genes, the 6 neoVTX genes were found to encompass a total of 411 AS isoforms, accounting for about half of the total number of toxin gene isoforms (Fig. 4C).

Figure 4 Massive alternative isoforms from the neovtx genes in the venom gland a) The work flow for sequencing and analysis for full length transcriptome construction. b) Distribution of the number of alternative splicing isoforms per annotated toxin gene.

c) Percent of the AS isoforms from each DETGs in the total isoforms (856) mapped to the annotated toxin genes (476). d) The 6 tandem-repeated neovtx genes demonstrated the most isoforms. The isoforms were mapped to each gene with the exon usage shown. 5 trans-splicing isoforms were listed at the lower panel, showing a single isoform embraced exons from different genes. e) Expression level of 27 cis-splicing and 5 trans-splicing isoforms from the 6 neovtx genes, the value was illustrated as mean ± SD, n = 7.

Mapping the isoforms to the 6 neoVTX genes showed that neovtx-a2 (154) and neovtx-b3 (101) demonstrated the highest number of AS isoforms, followed by neovtx-a3 (71), neovtx-a1 (35), neovtx-b1 (17) and neovtx-b2 (4) (Fig. 4D). Interestingly, a total of 29 trans-splicing isoforms were found among the 6 neoVTX genes, which utilized exons across these genes (Fig. 4D). The expression of all these AS isoforms was further analyzed referring to NGS data for quantification analysis. A total of 27 isoforms accounted for about 90% expression of the 6 neoVTX genes predominantly coming from the neovtx-a2, neovtx-b3 and neovtx-a3 genes. Additionally, the heterogenous expression pattern of the AS isoforms among all the 7 samples revealed the diversity of venom genetic regulation at individual level (Supplementary Fig. 5). By contrast, all the trans-splicing isoforms showed much lower expression (Fig. 4E and Supplementary Data 5). These results indicated that the alternative splicing of neoVTX genes contributes to the venom diversity.

High-throughput sequencing and rapid development of bioinformatic analysis have greatly facilitated the modern venomic research. In this work, we investigated the venom components of S. verrucosa by integrative analysis of genomic, transcriptomic and proteomic data. The chromosome level genome provides the fundamental genetic resource for venom components characterization and in-depth biomedical study. We first assembled the high-quality S. verrucosa genome, which enabled us to preliminary annotate a total of 476 toxin genes. Further transcriptomic analysis of the venom gland characterized 12 DTEGs as the potential major expressed toxic genes in the venom gland and 11 of them were eventually evidenced by the venom proteome.

Unique venom composition of S. verrucosa underlined its evolution adaption

A total of 467 toxin genes were annotated in the S. verrucosa genome, with peptidase S1 family, venom metalloproteinase family, the SNTX/VTX toxin family listed as the most abundant (Fig. 1C). Peptidase S1 family genes encode serine proteases, which are responsible for coordinating various physiological functions, including digestion, immune response, blood coagulation and reproduction. Snake venom serine protease (SVSP) was characterized as one of the four dominant protein families, together with phospholipase A₂ (PLA₂), snake venom metalloprotease (SVMP), and three-finger toxins (3FXs)[56]. In S. verrucosa genome, peptidase S1 family genes (102) were listed as the most abundant toxin genes in S. verrucosa genome (Fig. 1C). Interesting only two of them showed expression level with TPM > 100, and peptidase S1 genes as a whole accounted for only a small proportion of the VG transcripts. This was similar with the cases in scorpion, spider and centipede[57–59]. Correspondingly, only around 0.13% of peptidase S1 proteins were found in the crude venom (Supplementary Data 3). A total of 36 venom metalloproteinase family genes were annotated in our S. verrucosa genome, all of them showed low level expression with TPM < 50, and no metalloproteinase protein identified in the crude venom. Similarly, 6 PLA₂ genes (TPM < 1) were barely detected in the VG transcriptome. Instead, neoVTX proteins from the 6 tandem-repeated genes were the most abundant toxin proteins in the venom (Fig. 3B). Though we characterized 12 DETGs in the VG transcripts, only 3 toxin proteins demonstrated abundance more than 1% in the crude venom, including a pair of neoVTX subunits, a CRVP like protein, and a Kuntiz-type serine protease inhibitor. The results indicated the venom is relatively less complex compared to other venomous animals, which was also evidenced in a similar study in S.horrida[25]. The last but not the least, we found 4 unannotated genes with high expression (TPM > 500) in the VG but no proteins in the venom (Fig. 2B), indicating their potential roles as unique genes involved in the VG function of S. verrucosa. All the above results suggested that the unique venom composition underlined the venom phenotype evolution of S. verrucosa.

Integrative analysis revealed contribution of neoVTX genes to venom diversity of S. verrucosa

The major lethal molecules of S. verrucosa venom has been purified and molecular characterized as early as 1990s. The VTX and neoVTX were purified by two different research group and characterized as different proteins and remained unclarified so far. The VTX was first characterized as a tetrameric protein composing of VTXa and VTXb subunits[30, 31], whereas the neoVTX was a dimeric protein composing of neoVTXa and neoVTXb subunits[60]. They were considered as different proteins mainly based on the difference of their deduced amino acid sequence, neoVTXb showed 90% identity compared with the reported VTXb, lower than that with the SNTXb from another species, S.horrida. Further analysis found that VTXb demonstrated an additional nucleotide at the 3’ end of the CDS, leading to a frame shit, which contributed to the major difference between them. It was suggested that the geographical variation of the fish in the two studies contributed to major loci difference between the two toxin genes[60]. In our assembled genome, the discriminative loci in the 3 tandem-repeated beta subunit genes were all consistent with the neoVTXb (GenBank: AB262393.1) (Supplementary Fig. 6). Interestingly, we did find two protein bands with different molecular weights close to the reported dimeric and tetrameric neoVTX proteins by native PAGE electrophoresis, which were both further confirmed to contain by mass spectrometry sequencing (Fig. 3D). The results indicated that neoVTX proteins were the major lethal toxins in the fish strain we used and they may present as either dimeric or tetrameric, even other unknown polymeric proteins.

Gene duplication have been proved to play pivotal roles in the evolution and diversification of toxin genes[61]. The major toxins in the snake venom, including metalloproteinases, phospholipases A₂, and three-finger toxins, all experienced gene duplication events during the evolution[15, 62]. A total of 32 SNTX/VTX family genes were annotated in our assembled S. verrucosa genome. Integrative analysis with the venom gland transcriptome and venom proteome enabled us to evidenced that most of neoVTX proteins in the venom came from the six tandem-repeated neoVTX subunits, predominantly a pair of a2 and b3 genes (Fig. 3C). Whereas another two SNTX-like genes were also located adjacent to this region. They all showed the 4 conserved domains but with much lower sequence identity and were clustered to a separate branch in the phylogenetic tree (Fig. 1D). Additionally, the transcriptomic data showed that these two SNTX-like genes barely expressed in the venom gland, which was further evidenced that no protein was detected in the venom proteome. The results indicated that they may share the same origin with the neoNTX genes but showed no functional roles in S. verrucosa.

Alternative splicing of venom genes was considered as a complementary mechanism for the generation of venom complexity across animal lineages[19, 63, 64]. We found a total of 856 isoforms expressed in the venom gland from 198 toxin genes, and 411 of them including 29 trans-splicing isoforms were mapped to the 6 tandem-repeated neoVTX subunit genes, greatly increased the potential to produce diverse neoVTX protein. Expression analysis of these AS isoforms revealed that a pair of neoVTX a2 and b3 genes transcribed the several highest expressed isoforms (Fig. 4E). Interestingly, though we identified a total of 29 trans-splicing isoform among the six tandem-duplicated neoVTX genes, all of them showed much lower expression level (Fig. 4D and E), leaving their functional roles and biological significance unknown. Within the crude venom, 94.13% of the SNTX/VTX proteins were mapped to a pair of a2 and b3 genes, whereas transcriptomic data showed that the expression of them only accounted for 77.50% of all the 6 tandem-repeated genes (Fig. 3C), indicating an unknown regulation mechanism to biasedly promote translation of neoVTX a2 and b3 existed in the venom gland. Unfortunately, our proteomic sequencing approach failed to distinguish these proteins with various sizes and charges. A 2D gel electrophoresis of the crude venom from the closely related species S. horrida revealed 15 SNTXa isomers and 26 SNTXb isomers with varied length of ORF, partially supported the alternative splicing of neoVTX genes in this study. The results indicated that the AS of the toxin genes were conserved mechanism in stonefish underlined venom diversity[25]. Further work is required to investigate the biological significance of these neoVTX isomers.

The other major toxin proteins characterized in the venom of S. verrucosa

The most abundant toxin protein except the neoVTX proteins was characterized as CRVP in this study, whose homologue in S.horrida was defined as Golgi-associated plant pathogenesis related protein (ShGAPR) instead. Though ShGAPR showed the higher sequence similarity with the human GAPR1, it was suggested to act as a venom toxin similar to CVRPs[25]. The transcriptomic analysis of S. verrucosa venom gland revealed that, sv_crvp1 was found to be the highest expressed toxin gene in the VG and categorized into DEG compared to SK (Fig. 2B and C), further supporting its potential role as toxin gene in S. verrucosa. The most homologous of SvCRVP annotated in the ToxProt was helothermine, which was first isolated from the venom of the Mexican beaded lizard (Heloderma horridum horridum)[65]. Lizard helothermine was a peptide toxin could block ryanodine receptors, which were responsible for Ca²⁺ release within skeletal, cardiac and neuronal cells[66, 67]. Further work is encouraged to illustrate the potential role of SvCRVP as the most abundant venom component except neoVTX proteins in S. verrucosa.

Protease inhibitors in the venom were widely utilized by snakes as weapons to disrupt the homeostasis of prey’s physical biochemical reactions, such as the blood coagulation and blood pressure regulation, which resulting immobilization or death of the prey[68]. Besides, protease inhibitors were hypothesized to play roles in protecting peptide/protein toxins in the venom from degradation by proteases from the prey or predators[69]. Kuntiz-type serine protease inhibitors, C-type lectins and Phospholipase A₂ were the most extensively distributed protease inhibitors in the snake venom[70]. We annotated 17 Kunitz-type family genes in S. verrucosa genome and only one of them (sv_kspi1/SvKSPI) were characterized as DETGs in the venom gland (Fig. 2B and 3B) with an abundance of 1.84% in total toxin proteins. Besides, homologues of Snaclec B1 and SLXA belonging to C-type lectins were also found in the crude venom (Fig. 3B), together with SvKSPI served as the major protease inhibitors in S. verrucosa.

Calglandulin, previously identified from the venom gland of Bothrops insularis, was characterized as a putative Ca²⁺ binding protein with four EF-hand motifs[71]. Using antiserum against the recombinant calglandulin, it was found that calglandulin was only detectable in the venom gland of Bothrops insularis, but not in the crude venom or other tissue, indicating that Calglandulin may be primarily involved in the cellular control mechanism of the secretion of toxins from the gland into the venom[71, 72]. Whereas Calglandulin was screened out as a major component in the crude venom of fire ant Solenopsis invicta[73] and stingray[74], indicating its unknown functions beyond the toxins secretion. In the S. verrucosa venom, a small proportion of Calglandulin proteins (0.60%) were also evidenced, indistinguishably corresponding to the two sv_caglp1 and sv_caglp2 due to their same deduced amino acid sequence (Supplementary Fig. 7). Notably, transcriptomic analysis revealed that only sv_caglp1 gene was characterized as DEGs between venom gland and the skin, indicating a transcriptional regulation mechanism underlined driving sv_caglp1 as the major one to exert the function in toxins secretion (Fig. 2B and C).

Hyaluronidase is one of the enzymes commonly seen in the venom of snake, scorpions, spiders and leeches[75]. Hyaluronidase activity is considered critical for the spreading of toxins and interrupting the integrity of target’s extracellular matrix through the degradation of hyaluronic acid, which have been considered as a therapeutic target against snake venom[76, 77]. The cDNA of S. verrucosa hyaluronidase has been previously cloned and characterized with a conserved active site composing of one catalytic residue and four substrate positioning residues[32]. In the present study, we characterized hyaluronidase gene as one of DETGs in the venom gland and further supported by the proteomic data (Fig. 2B and 3B).

Among all the 12 DETGs defined in the S. verrucosa venom gland, only the sv_tpt1 gene demonstrated no proteomic evidence in the crude venom (Fig. 3B). Sv_tpt1 gene encodes a Translationally Controlled Tumor Protein (TCTP), which was predominantly described as venom toxins from spiders, contributing the allergic and inflammatory activity of the venom[78, 79]. Few studies paid attention to the TCTP toxins in other venomous animals, leaving its biological function under investigated. The expression of TCTP gene was reported in the venom gland from three scorpionfishes (Scorpaenidae.spp), indicating its widely distribution in the venomous teleost[80]. In the S. verrucosa, sv_tpt1 gene showed comparable high expression level both in the VG and SK (Fig. 2C), together with the absence in the crude venom, indicating its potential role beyond as toxins in the S. verrucosa.

A group of toxin proteins beyond the 12 DETGs were also identified based on the abundance in the crude venom (Fig. 3B). A putative endothelial lipase first found in venom gland transcriptome of the eastern diamondback rattlesnake, were suggested to have phospholipase and triglyceride lipase activities[81]. Lactotransferrin protein was a major iron-binding protein usually found in exocrine fluids including mucosal secretions and mammal breast milk. Lactotransferrin was annotated as toxins as its presence in the proteome of venom accessory glands in the vampire bat (Desmodus rotundus) and acted as antimicrobials[82]. Similarly, peroxiredoxin-4[83], Snaclec B1[84], Cobra venom factor like protein[85] were also evidenced as venom proteins in other venomous animals. The presence of these proteins in the crude venom of S. verrucosa further supported their characterization as a toxin protein, whereas the molecular mechanism underlined need to be further investigated.

A notable non-toxin protein in the venom

A set of non-toxin proteins were shown in Fig. 3B. Though some of them were never reported to be present in the crude venom, they were hypothesized to imply housekeeping or hemostasis maintaining roles in the VG. Among them, a gastrotropin-like protein was found to be the second most enriched followed neoVTX proteins, encoded by the FABP6 gene, which also showed high level expression in VG. (Fig. 2B and 3B). Gastrotropin/FABP6 is a fatty acid binding protein and therapeutic target related to immune infiltration[86]. Gastrotropin was first characterized as venom protein in a S. horrida, but not any other venomous animals. It was suggested to contribute the unique activities of its venom, by maintaining venom hemostasis and protecting venom components from degradation [25]. We did find that the extracted crude venom of S. verrucosa remained high activity at room temperature after hours. Further work is encouraged to investigate the biological significance of gastrotropin-like protein in S. verrucosa venom.

Overall, the major venom components of S. verrucosa were thoroughly investigated by integrative multi-omics approaches in this work, for the first attempt in the venomous teleost. Six tandem-duplicated neoVTX subunit genes were explored in the genome and discovered to embrace massive alternative splicing isoforms and diverse neoVTX proteins consequently. The diverse neoVTX proteins in the venom as lethal particles are hypothesized to be pivotal to understand adaptive evolution of S. verrucosa, especially when facing different nature predators. Further functional studies are necessary to elucidate the mechanism of alternative isoforms recruitment in venom. Additionally, more biomedical studies are encouraged to exploit venom components of S. verrucosa for pharmaceutical innovation.

Ethics approval and consent to participate

All the experiments were conducted in accordance with guidelines approved by the Animal Ethical and Welfare Committee of Southern Marine Science and Engineering Guangdong Laboratory (Guangzhou).

Consent for publication

Not applicable

Competing interests

The authors declare no competing interests

Funding

This work was supported by National Natural Science Foundation of China (32222014), Ministry of Science and Technology of China (2023YFF1304900, 2021YFF0502804), Science and Technology Department of Guangdong Province (2021QN02H103, 2024A1515013196), Guangdong Forestry Administration (SLYJ2023B4004), and the PI Project of Southern Marine Science and Engineering Guangdong Laboratory (Guangzhou) (GML2020GD0804, GML2022GD0804).

Author Contribution

Conceptualization, Z.Z. and W.Z.; Methodology, Q.L., H.L., W.Y., Z.P., and Z.Z.; Investigation, Z.Z., W.Z., and F.W.; Writing-Original Draft, Q.L. and Z.Z.; Writing-Review & Editing, Z.Z., S.W., and W.Z.; Funding Acquisition, Z.Z. and W.Z.; Resources, W.Z. and F.W.; Supervision, W.Z. and F.W.

Acknowledgements

Not applicable

Data Availability

The datasets used and/or analysed during the current study are available from the corresponding author on reasonable request

Zancolli G, Casewell NR. Venom systems as models for studying the origin and regulation of evolutionary novelties. Mol Biol Evol. 2020;37(10):2777–90.
Bordon KCF, Cologna CT, Fornari-Baldo EC, Pinheiro-Júnior EL, Cerni FA, Amorim FG, Anjolette FAP, Cordeiro FA, Wiezel GA, Cardoso IA, et al. From animal poisons and venoms to medicines: achievements, challenges and perspectives in drug discovery. Front Pharmacol. 2020;11:1132.
von Reumont BM, Anderluh G, Antunes A, Ayvazyan N, Beis D, Caliskan F, Crnković A, Damm M, Dutertre S, Ellgaard L et al. Modern venomics-Current insights, novel methods, and future perspectives in biological and applied animal venom research. GigaScience 2022, 11.
Suranse V, Srikanthan A, Sunagar K. Animal Venoms: Origin, Diversity and Evolution. In: eLS. John Wiley & Sons, Ltd; 2018: 1–20.
Suryamohan K, Krishnankutty SP, Guillory J, Jevit M, Schröder MS, Wu M, Kuriakose B, Mathew OK, Perumal RC, Koludarov I, et al. The Indian cobra reference genome and transcriptome enables comprehensive identification of venom toxins. Nat Genet. 2020;52(1):106–17.
Ludington AJ, Hammond JM, Breen J, Deveson IW, Sanders KL. New chromosome-scale genomes provide insights into marine adaptations of sea snakes (Hydrophis: Elapidae). BMC Biol 2023, 21(1).
Margres MJ, Rautsaw RM, Strickland JL, Mason AJ, Schramer TD, Hofmann EP, Stiers E, Ellsworth SA, Nystrom GS, Hogan MP et al. The Tiger Rattlesnake genome reveals a complex genotype underlying a simple venom phenotype. Proc Natl Acad Sci U S A 2021, 118(4).
Zhu B, Jin P, Hou Z, Li J, Wei S, Li S. Chromosomal-level genome of a sheet-web spider provides insight into the composition and evolution of venom. Mol Ecol Resour. 2022;22(6):2333–48.
You Y, Tang Y, Yin W, Liu X, Gao P, Zhang C, Tembrock LR, Zhao Y, Yang Z. From genome to proteome: Comprehensive identification of venom toxins from the Chinese funnel-web spider (Macrothelidae: Macrothele yani). Int J Biol Macromol 2024, 268(Pt 2).
Herráez-Pérez A, Pardos-Blas JR, Afonso CML, Tenorio MJ, Zardoya R. Chromosome-level genome of the venomous snail Kalloconus canariensis: a valuable model for venomics and comparative genomics. Gigascience 2022, 12.
Peng C, Huang Y, Bian C, Li J, Liu J, Zhang K, You X, Lin Z, He Y, Chen J, et al. The first Conus genome assembly reveals a primary genetic central dogma of conopeptides in C. betulinus. Cell Discov. 2021;7(1):11.
Farhat S, Modica MV, Puillandre N. Whole genome duplication and gene evolution in the hyperdiverse venomous gastropods. Mol Biol Evol. 2023;40(8):msad171.
Wei S, Zhou W, Fan H, Zhang Z, Guo W, Peng Z, Wei F. Chromosome-level genome assembly of the yellow boxfish (Ostracion cubicus) provides insights into the evolution of bone plates and ostracitoxin secretion. Front Mar Sci. 2023;10:1170704.
Martinson EO, Mrinalini, Kelkar YD, Chang C-H, Werren JH. The evolution of venom by co-option of single-copy genes. Curr Biol. 2017;27(13):2007–13. e2008.
Zhang Z-Y, Lv Y, Wu W, Yan C, Tang C-Y, Peng C, Li J-T. The structural and functional divergence of a neglected three-finger toxin subfamily in lethal elapids. Cell Rep. 2022;40(2):111079.
Almeida DD, Viala VL, Nachtigall PG, Broe M, Gibbs HL, Serrano SMT, Moura-da-Silva AM, Ho PL, Nishiyama-Jr MY, Junqueira-de-Azevedo ILM. Tracking the recruitment and evolution of snake toxins using the evolutionary context provided by the Bothrops jararaca genome. Proc Natl Acad Sci U S A. 2021;118(20):e2015159118.
Koludarov I, Senoner T, Jackson TNW, Dashevsky D, Heinzinger M, Aird SD, Rost B. Domain loss enabled evolution of novel functions in the snake three-finger toxin gene superfamily. Nat Commun. 2023;14(1):4861.
Zancolli G, Reijnders M, Waterhouse RM, Robinson-Rechavi M. Convergent evolution of venom gland transcriptomes across Metazoa. Proc Natl Acad Sci U S A. 2022;119(1):e2111392119.
Ye X, He C, Yang Y, Sun YH, Xiong S, Chan KC, Si Y, Xiao S, Zhao X, Lin H, et al. Comprehensive isoform-level analysis reveals the contribution of alternative isoforms to venom evolution and repertoire diversity. Genome Res. 2023;33(9):1554–67.
Harris RJ, Jenner RA. Evolutionary ecology of fish venom: Adaptations and consequences of evolving a venom system. Toxins. 2019;11(2):60.
Maillaud C, Hoang-Oppermann T, Hoang-Oppermann V, Rigot H, Girardot S, Nour M. Is stonefish Synanceia verrucosa envenomation potentially lethal? Toxicon 2020, 184:78–82.
Lopez CAJ, Magee CAJ, Belyea CCM, Gumboc L. Finger flexor tenosynovitis from Stonefish envenomation injury. J Am Acad Orthop Surg Glob Res Rev. 2019;3(5):e024.
Khoo HE. Bioactive proteins from stonefish venom. Clin Exp Pharmacol Physiol. 2002;29(9):802–6.
Poh CH, Yuen R, Khoo HE, Chung M, Gwee M, Gopalakrishnakone P. Purification and partial characterization of stonustoxin (lethal factor) from Synanceja horrida venom. Comp Biochem Physiol B. 1991;99(4):793–8.
Ziegman R, Undheim EAB, Baillie G, Jones A, Alewood PF. Investigation of the estuarine stonefish (Synanceia horrida) venom composition. J Proteom. 2019;201:12–26.
Ellisdon AM, Reboul CF, Panjikar S, Huynh K, Oellig CA, Winter KL, Dunstone MA, Hodgson WC, Seymour J, Dearden PK, et al. Stonefish toxin defines an ancient branch of the perforin-like superfamily. Proc Natl Acad Sci U S A. 2015;112(50):15360–5.
Saggiomo SL, Firth C, Wilson DT, Seymour J, Miles JJ, Wong Y. The geographic distribution, venom components, pathology and treatments of stonefish (Synanceia spp.) venom. Mar Drugs. 2021;19(6):302.
Tang T, Huang Y, Peng C, Liao Y, Lv Y, Shi Q, Gao B. A chromosome-level genome assembly of the reef stonefish (Synanceia verrucosa) provides novel insights into Stonustoxin (sntx) Genes. Mol Biol Evol. 2023;40(10):msad215.
Garnier P, Sauviat M-p, Goudey-perriere F, Perriere C. Cardiotoxicity of verrucotoxin, a protein isolated from the venom of Synanceia verrucosa. Toxicon. 1997;35(1):47–55.
Garnier P, Goudey-Perrière F, Breton P, Dewulf C, Petek F, Perrière C. Enzymatic properties of the stonefish (Synanceia verrucosa Bloch and Schneider, 1801) venom and purification of a lethal, hypotensive and cytolytic factor. Toxicon. 1995;33(2):143–55.
Garnier P, Ducancel F, Ogawa T, Boulain J-C, Goudey-Perrière F, Perrière C, Ménez A. Complete amino-acid sequence of the β-subunit of VTX from venom of the stonefish (Synanceia verrucosa) as identified from cDNA cloning experiments. Biochim Biophys Acta. 1997;1337(1):1–5.
Madokoro M, Ueda A, Kiriake A, Shiomi K. Properties and cDNA cloning of a hyaluronidase from the stonefish Synanceia verrucosa venom. Toxicon. 2011;58(4):285–92.
Kato K, Nakagawa H, Shinohara M, Ohura K. Purification of a novel lectin from the dorsal spines of the stonefish, Synanceia verrucosa. J Osaka Dent Univ. 2016;50(2):55–61.
Cheng H, Concepcion GT, Feng X, Zhang H, Li H. Haplotype-resolved de novo assembly using phased assembly graphs with hifiasm. Nat Methods. 2021;18(2):170–5.
Guan D, McCarthy SA, Wood J, Howe K, Wang Y, Durbin R. Identifying and removing haplotypic duplication in primary genome assemblies. Bioinformatics. 2020;36(9):2896–8.
Zhang H, Song L, Wang X, Cheng H, Wang C, Meyer CA, Liu T, Tang M, Aluru S, Yue F, et al. Fast alignment and preprocessing of chromatin profiles with Chromap. Nat Commun. 2021;12(1):6566.
Zhou C, McCarthy SA, Durbin R. YaHS: yet another Hi-C scaffolding tool. Bioinformatics 2023, 39(1).
Manni M, Berkeley MR, Seppey M, Simão FA, Zdobnov EM. BUSCO update: novel and streamlined workflows along with broader and deeper phylogenetic coverage for scoring of eukaryotic, prokaryotic, and viral genomes. Mol Biol Evol. 2021;38(10):4647–54.
Flynn JM, Hubley R, Goubert C, Rosen J, Clark AG, Feschotte C, Smit AF. RepeatModeler2 for automated genomic discovery of transposable element families. Proc Natl Acad Sci U S A. 2020;117(17):9451–7.
Bao W, Kojima KK, Kohany O. Repbase Update, a database of repetitive elements in eukaryotic genomes. Mob DNA. 2015;6:11.
Tarailo-Graovac M, Chen N. Using RepeatMasker to identify repetitive elements in genomic sequences. Curr Protoc Bioinf 2009, Chap. 4:Unit 4.10.
Benson G. Tandem repeats finder: a program to analyze DNA sequences. Nucleic Acids Res. 1999;27(2):573–80.
Haas BJ, Salzberg SL, Zhu W, Pertea M, Allen JE, Orvis J, White O, Buell CR, Wortman JR. Automated eukaryotic gene structure annotation using EVidenceModeler and the program to assemble spliced alignments. Genome Biol. 2008;9(1):R7.
Chen L, Qiu Q, Jiang Y, Wang K, Lin Z, Li Z, Bibi F, Yang Y, Wang J, Nie W, et al. Large-scale ruminant genome sequencing provides insights into their evolution and distinct traits. Science. 2019;364(6446):eaav6202.
Kim D, Paggi JM, Park C, Bennett C, Salzberg SL. Graph-based genome alignment and genotyping with HISAT2 and HISAT-genotype. Nat Biotechnol. 2019;37(8):907–15.
Pertea M, Pertea GM, Antonescu CM, Chang TC, Mendell JT, Salzberg SL. StringTie enables improved reconstruction of a transcriptome from RNA-seq reads. Nat Biotechnol. 2015;33(3):290–5.
Stanke M, Diekhans M, Baertsch R, Haussler D. Using native and syntenically mapped cDNA alignments to improve de novo gene finding. Bioinformatics. 2008;24(5):637–44.
Burge CB, Karlin S. Finding the genes in genomic DNA. Curr Opin Struct Biol. 1998;8(3):346–54.
Majoros WH, Salzberg SL. An empirical analysis of training protocols for probabilistic gene finders. BMC Bioinformatics. 2004;5:206.
Kuo RI, Cheng Y, Zhang R, Brown JWS, Smith J, Archibald AL, Burt DW. Illuminating the dark side of the human transcriptome with long read transcript sequencing. BMC Genomics. 2020;21(1):751.
Tardaguila M, de la Fuente L, Marti C, Pereira C, Pardo-Palacios FJ, Del Risco H, Ferrell M, Mellado M, Macchietto M, Verheggen K, et al. SQANTI: extensive characterization of long-read transcript sequences for quality control in full-length transcriptome identification and quantification. Genome Res. 2018;28(3):396–411.
Liao Y, Smyth GK, Shi W. featureCounts: an efficient general purpose program for assigning sequence reads to genomic features. Bioinformatics. 2014;30(7):923–30.
Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15(12):550.
Margres MJ, Rautsaw RM, Strickland JL, Mason AJ, Schramer TD, Hofmann EP, Stiers E, Ellsworth SA, Nystrom GS, Hogan MP, et al. The Tiger Rattlesnake genome reveals a complex genotype underlying a simple venom phenotype. Proc Natl Acad Sci U S A. 2021;118(4):e2014634118.
Mistry J, Chuguransky S, Williams L, Qureshi M, Salazar Gustavo A, Sonnhammer ELL, Tosatto SCE, Paladin L, Raj S, Richardson LJ, et al. Pfam: The protein families database in 2021. Nucleic Acids Res. 2020;49(D1):D412–9.
Tasoulis T, Isbister GK. A current perspective on snake venom composition and constituent protein families. Arch Toxicol. 2023;97(1):133–53.
Morgenstern D, Rohde BH, King GF, Tal T, Sher D, Zlotkin E. The tale of a resting gland: transcriptome of a replete venom gland from the scorpion Hottentotta judaicus. Toxicon. 2011;57(5):695–703.
Almeida DD, Scortecci KC, Kobashi LS, Agnez-Lima LF, Medeiros SR, Silva-Junior AA, Junqueira-de-Azevedo Ide L. Fernandes-Pedrosa Mde F: Profiling the resting venom gland of the scorpion Tityus stigmurus through a transcriptomic survey. BMC Genomics. 2012;13:362.
Liu ZC, Zhang R, Zhao F, Chen ZM, Liu HW, Wang YJ, Jiang P, Zhang Y, Wu Y, Ding JP, et al. Venomic and transcriptomic analysis of centipede Scolopendra subspinipes dehaani. J Proteome Res. 2012;11(12):6197–212.
Ueda A, Suzuki M, Honma T, Nagai H, Nagashima Y, Shiomi K. Purification, properties and cDNA cloning of neoverrucotoxin (neoVTX), a hemolytic lethal factor from the stonefish Synanceia verrucosa venom. Biochim Biophys Acta. 2006;1760(11):1713–22.
Wong ESW, Belov K. Venom evolution through gene duplications. Gene. 2012;496(1):1–7.
Almeida DD, Viala VL, Nachtigall PG, Broe M, Gibbs HL, Serrano SMT, Moura-da-Silva AM, Ho PL, Nishiyama-Jr MY, Junqueira-de-Azevedo ILM. Tracking the recruitment and evolution of snake toxins using the evolutionary context provided by the Bothrops jararaca genome. Proc Natl Acad Sci U S A. 2021;118(20):e2015159118.
Haney RA, Matte T, Forsyth FS, Garb JE. Alternative transcription at venom genes and its role as a complementary mechanism for the generation of venom complexity in the common house spider. Front Ecol Evol. 2019;7:85.
Ogawa T, Oda-Ueda N, Hisata K, Nakamura H, Chijiwa T, Hattori S, Isomoto A, Yugeta H, Yamasaki S, Fukumaki Y, et al. Alternative mRNA splicing in three venom families underlying a possible production of divergent venom proteins of the habu snake, Protobothrops flavoviridis. Toxins (Basel). 2019;11(10):581.
Mochca-Morales J, Martin BM, Possani LD. Isolation and characterization of Helothermine, a novel toxin from Heloderma horridum horridum (Mexican beaded lizard) venom. Toxicon. 1990;28(3):299–309.
Van Petegem F. Ryanodine receptors: structure and function. J Biol Chem. 2012;287(38):31624–32.
Morrissette J, Krätzschmar J, Haendler B, el-Hayek R, Mochca-Morales J, Martin BM, Patel JR, Moss RL, Schleuning WD, Coronado R, et al. Primary structure and properties of helothermine, a peptide toxin that blocks ryanodine receptors. Biophys J. 1995;68(6):2280–8.
Inagaki H. Snake Venom Protease Inhibitors: Enhanced Identification, Expanding Biological Function, and Promising Future. In: Snake Venoms. Edited by Inagaki H, Vogel C-W, Mukherjee AK, Rahmy TR, Gopalakrishnakone P. Dordrecht: Springer Netherlands; 2017: 161–186.
Ma H, Xiao-Peng T, Yang S-L, Lu Q-M, Lai R. Protease inhibitor in scorpion (Mesobuthus eupeus) venom prolongs the biological activities of the crude venom. Chin J Nat Med. 2016;14(8):607–14.
Thakur R, Mukherjee AK. Pathophysiological significance and therapeutic applications of snake venom protease inhibitors. Toxicon. 2017;131:37–47.
Junqueira-de-Azevedo ILM, Pertinhez T, Spisni A, Carreño FR, Farah CS, Ho PL. Cloning and expression of calglandulin, a new EF-hand protein from the venom glands of Bothrops insularis snake in E. coli. Biochim Biophys Acta. 2003;1648(1):90–8.
Pierre LS, Woods R, Earl S, Masci PP, Lavin MF. Identification and analysis of venom gland-specific genes from the coastal taipan (Oxyuranus scutellatus) and related species. Cell Mol Life Sci. 2005;62(22):2679–93.
Cai L, Yang F, Wang Y, Yang J, Zhu Y, Ma X, Höfer J, Wang Y, Ma Y, Xiao L. A combined protein toxin screening based on the transcriptome and proteome of Solenopsis invicta. Proteome Sci. 2022;20(1):15.
Kirchhoff KN, Billion A, Voolstra CR, Kremb S, Wilke T, Vilcinskas A. Stingray venom proteins: Mechanisms of action revealed using a novel network pharmacology approach. Mar Drugs. 2021;20(1):27.
Bala E, Hazarika R, Singh P, Yasir M, Shrivastava R. A biological overview of Hyaluronidase: A venom enzyme and its inhibition with plants materials. Materials Today 2018, 5(2, Part 1):6406–6412.
Silva de França F, Tambourgi DV. Hyaluronan breakdown by snake venom hyaluronidases: From toxins delivery to immunopathology. Front Immunol. 2023;14:1125899.
Bittenbinder MA, van Thiel J, Cardoso FC, Casewell NR, Gutiérrez J-M, Kool J, Vonk FJ. Tissue damaging toxins in snake venoms: mechanisms of action, pathophysiology and treatment strategies. Commun Biol. 2024;7(1):358.
Boia-Ferreira M, Moreno KG, Basílio ABC, da Silva LP, Vuitika L, Soley B, Wille ACM, Donatti L, Barbaro KC, Chaim OM, et al. TCTP from Loxosceles Intermedia (Brown Spider) venom contributes to the allergic and inflammatory response of cutaneous loxoscelism. Cells. 2019;8(12):1489.
Senff-Ribeiro A. Translationally Controlled Tumor Protein (TCTP/HRF) in Animal Venoms. In: TCTP/tpt1 - Remodeling Signaling from Stem Cell to Disease. Edited by Telerman A, Amson R. Cham: Springer International Publishing; 2017: 193–200.
Xie B, Yu H, Kerkkamp H, Wang M, Richardson M, Shi Q. Comparative transcriptome analyses of venom glands from three scorpionfishes. Genomics. 2019;111(3):231–41.
Rokyta DR, Lemmon AR, Margres MJ, Aronow K. The venom-gland transcriptome of the eastern diamondback rattlesnake (Crotalus adamanteus). BMC Genomics. 2012;13(1):312.
Francischetti IMB, Assumpção TCF, Ma D, Li Y, Vicente EC, Uieda W, Ribeiro JMC. The Vampirome: Transcriptome and proteome analysis of the principal and accessory submaxillary glands of the vampire bat Desmodus rotundus, a vector of human rabies. J Proteom. 2013;82:288–319.
Calvete JJ, Fasoli E, Sanz L, Boschetti E, Righetti PG. Exploring the venom proteome of the western diamondback rattlesnake, Crotalus atrox, via snake venomics and combinatorial peptide ligand library approaches. J Proteome Res. 2009;8(6):3055–67.
Jebali J, Bazaa A, Sarray S, Benhaj K, Karboul A, El Ayeb M, Marrakchi N, Gargouri A. C-type lectin protein isoforms of Macrovipera lebetina: cDNA cloning and genetic diversity. Toxicon. 2009;53(2):228–37.
Fritzinger DC, Bredehorst R, Vogel CW. Molecular cloning and derived primary structure of cobra venom factor. Proc Natl Acad Sci U S A. 1994;91(26):12775–9.
McKillop IH, Girardi CA, Thompson KJ. Role of fatty acid binding proteins (FABPs) in cancer development and progression. Cell Signal. 2019;62:109336.

No competing interests reported.

Download PDF

Editorial decision: Revision requested
24 Oct, 2024
Editor assigned by journal
12 Oct, 2024
Submission checks completed at journal
12 Oct, 2024
First submitted to journal
11 Oct, 2024

You are reading this latest preprint version

Integrative multi-omics analysis reveals the contribution of neoVTX genes to venom diversity of Synanceia verrucosa

Status:

Version 1

Abstract

Background

Results

Conclusions

Figures

Background

Methods

Animals

Genome extraction, sequencing

RNA extraction and sequencing

Genome assembly and annotation

Long-read RNA alignment and isoform identification

Expression analysis of the genes in the venom gland

Toxin genes annotation

Proteomic analysis of the Venom

Results

Toxin genes expression in the venom gland

Proteomic evidence of the toxin gene in venom

Alternative splicing contributed to venom diversity

Discussion

A notable non-toxin protein in the venom

Conclusions

Declarations

Competing interests

Funding

Author Contribution

Acknowledgements

Data Availability

References

Additional Declarations

Supplementary Files

Status:

Version 1