Next Article in Journal
Biomarker of Aflatoxin Ingestion: 1H NMR-Based Plasma Metabolomics of Dairy Cows Fed Aflatoxin B1 with or without Sequestering Agents
Previous Article in Journal
Occurrence and Levels of Aflatoxins in Fish Feeds and Their Potential Effects on Fish in Nyeri, Kenya
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Toxins
Volume 10
Issue 12
10.3390/toxins10120544
toxins-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Transcriptomic Characterization of the South American Freshwater Stingray Potamotrygon motoro Venom Apparatus

by
Filipe Silva
1,2,†,
Yu Huang
3,4,†,
Vítor Yang
1,2,
Xidong Mu
5,
Qiong Shi
3,4 and
Agostinho Antunes
1,2,*
1
CIIMAR/CIMAR, Interdisciplinary Centre of Marine and Environmental Research, University of Porto, Av. General Norton de Matos, s/n, 4450-208 Porto, Portugal
2
Department of Biology, Faculty of Sciences, University of Porto, Rua do Campo Alegre, s/n, 4169-007 Porto, Portugal
3
Shenzhen Key Lab of Marine Genomics, Guangdong Provincial Key Lab of Molecular Breeding in Marine Economic Animals, BGI Academy of Marine Sciences, Shenzhen 518083, China
4
BGI Education Center, University of Chinese Academy of Sciences, Shenzhen 518083, China
5
Pearl River Fisheries Research Institute, Chinese Academy of Fishery Sciences, Key Laboratory of Recreational Fisheries, Ministry of Agriculture, Guangdong Engineering Technology Research Center for Advanced Recreational Fisheries, Guangzhou 510380, China
*
Author to whom correspondence should be addressed.
†
These authors contributed equally to this work.
Toxins 2018, 10(12), 544; https://doi.org/10.3390/toxins10120544
Submission received: 15 November 2018 / Revised: 14 December 2018 / Accepted: 14 December 2018 / Published: 18 December 2018
(This article belongs to the Section Animal Venoms)
Browse Figures
Versions Notes

Abstract

:
Venomous animals are found through a wide taxonomic range including cartilaginous fish such as the freshwater stingray Potamotrygon motoro occurring in South America, which can injure people and cause venom-related symptoms. Ensuring the efficacy of drug development to treat stingray injuries can be assisted by the knowledge of the venom composition. Here we performed a detailed transcriptomic characterization of the venom gland of the South American freshwater stingray Potamotrygon motoro. The transcripts retrieved showed 418 hits to venom components (comparably to 426 and 396 hits in other two Potamotrygon species), with high expression levels of hyaluronidase, cystatin and calglandulin along with hits uniquely found in P. motoro such as DELTA-alicitoxin-Pse1b, Augerpeptide hhe53 and PI-actitoxin-Aeq3a. We also identified undescribed molecules with extremely high expression values with sequence similarity to the SE-cephalotoxin and Rapunzel genes. Comparative analyses showed that despite being closely related, there may be significant variation among the venoms of freshwater stingrays, highlighting the importance of considering elicit care in handling different envenomation cases. Since hyaluronidase represents a major component of fish venom, we have performed phylogenetic and selective pressure analyses of this gene/protein across all fish with the available information. Results indicated an independent recruitment of the hyaluronidase into the stingray venom relative to that of venomous bony fish. The hyaluronidase residues were found to be mostly under negative selection, but 18 sites showed evidence of diversifying positive selection (P < 0.05). Our data provides new insight into stingray venom variation, composition, and selective pressure in hyaluronidase.
Keywords:
venom; transcriptomics; stingray; hyaluronidase; next-generation sequencing; selective pressure
Key Contribution: The transcriptomic characterization of the venom gland of Potamotrygon motoro and its comparison to other Potamotrygon species’ transcriptomes, such as P. falkneri and P. amandae, allowed determining consistent components involved in the envenomation scenarios. Results suggest the recruitment of a different hyaluronidase (HYAL) gene into the venom cocktail of stingrays when compared to venomous bony fish (HYAL6 in stingrays and HYAL5 in Scorpaeniformes).

1. Introduction

Potamotrygon motoro, commonly known as the ocellate river stingray, is a freshwater stingray with a circular disk, eyes raised from the dorsal surface, typically beige or brown, with dark-ringed orange spots [1]. It belongs to the Potamotrygonidae family that represents the only group of cartilaginous fish living mostly in freshwater. Freshwater stingrays occur in South America, within some river systems, each with its own endemic stingrays. Recent classifications of this family included five genera, Heliotrygon, Paratrygon, Plesiotrygon, Potamotrygon and Styracura [2,3]. One of the first described freshwater stingray was P. motoro, which is the most widespread species of the whole family across South America, found in freshwater rivers in Uruguay, Paraná-Paraguay, Orinoco, and Amazon River basins [1,4,5].
Similar to other species of the genus, the P. motoro possesses a tail with a sting comprehending the venom apparatus consisting of a rigid structure of dentin, serrated barbs, and an enveloping layer of cells with venom-producing capabilities [6,7]. Stingray injuries usually have two components, an immediate physical trauma from the penetration of the spine and an envenomation in the created wound when the spine tears through the integumentary tissue.
The stings will be used by the animal through tail-flicking towards the desired target if frightened, such as when there are capture attempts or close proximity of bathers. Unsuspecting individuals stepping on hidden spines in the sand is another possibility. The result is tissue breakdown and necrosis, followed by heavy pain and potentially further infection at the wound site [7,8,9]. The physical trauma is described as incredibly painful and can provoke severe damage due to the serrations of the stingers, which may strike vessels or nerves [7,10]. Furthermore, mucus surrounding the stingers and bacteria present in the epithelium can drastically worsen the condition of the wound. Bacterial infections are reported to commonly be gram-negative species with wide antibiotic resistance [11,12]. In areas without appropriate treatment or facilities available, or when no treatment is procured, the injuries’ severity increases, and death may occur. Edema, erythema, inflammation, vomiting, and headaches are other commonly described symptoms. Tetanus is also a reported risk [7,13,14,15,16].
Wounds related to freshwater stingrays are common in some South American regions in which these species are abundant, and the locals openly interact with them, such as in bathing areas or within fishing communities. Mishandling captive stingrays can also result in envenomation [8,9]. Actual confirmation of the species involved is difficult, but it is likely that many recorded occurrences are due to the highly widespread Potamotrygon motoro [15].
Despite the recent advances in venom studies, most species that are suspected of being venomous or have shown such potential remain largely unexplored. In particular, when it comes to aquatic species, from invertebrates such as coleoids to vertebrates such as scorpionfish, the number of analyzed species compared to the prospective amount is considerably reduced [17,18,19,20,21,22]. Moreover, in snake’s species venom composition could vary depending on where the sting occurred (different environments) and whether it was captive or wild specimens being handled, which could lead to a bigger array of compounds and symptoms [23,24]. Within stingrays, there have been studies isolating some of the molecules involved in their venoms and assessing their potential impact [25,26,27]. More recently, three transcriptomes were generated from the freshwater stingrays P. amandae, P. falkneri [28] and the marine blue-spotted stingray Neotrygon kuhlii [29], illuminating the molecular diversity of venoms in stingrays relatively to other venomous species.
Here, we present the venom gland transcriptome characterization of the South American freshwater stingray P. motoro revealing its venom components and performing comparative analyses to identify the main differences across members of the same genus (i.e., Potamotrygon). We assessed the protein coding transcript variability in the venom of the closely related Potamotrygon species to elucidate the venom composition of the P. motoro stingray that can easily injure humans, providing information required for the development of safety and treatment measures. Findings show slight variations that should be taken into consideration when approaching freshwater stingray envenomation cases and components that require further proteomic clarification.
We also analyzed the hyaluronidase (HYAL) gene family among fish due to its presence in the P. motoro venom gland and its relevance in the venom of other animals. The hyaluronidase enzyme is capable of hydrolyzing hyaluronan, an extracellular matrix component, lowering its viscosity [30]. Our results support an independent venom recruitment of the hyaluronidase in cartilaginous freshwater stingrays relatively to other fish. We identified the HYAL1, HYAL2 and HYAL6 paralogs, the last being a candidate for the gene recruitment into the venom molecular network. Furthermore, hyaluronidase key sites were detected as being under diversifying selective pressure. We hypothesized that such residues likely do not impact the hyaluronidase enzymatic function but could instead affect the protein stability, affinity, and protein-protein interaction.

2. Results and Discussion

2.1. Denovo Assembly

The Illumina HiSeq sequencing resulted in 38,674,474 paired-end reads of 90 bp each. After filtering out the adapters and low-quality sequences, a total of 35,977,224 clean reads were obtained for further processing. De novo assembly using Trinity produced 140,078 contigs with a length average of ~498 bp, minimum of 174 bp and maximum of 15,040 bp. There were 14543 contigs longer than 1000 bp. Of the contigs with expression higher than one fragment per kilobase of transcript per million mapped reads (FPKM), 27,032 potentially coding transcripts with an average length of ~787 bp were identified by the Trinotate annotation pipeline (Table 1).
Our sequencing and assembly data provided a comparable number of contigs to previous studies [28], with both significant expression and predicted coding regions. However, the average length of contigs and length of a large amount of our dataset was considerably lower, which could have resulted from the sampling and sequencing methods procedures considered.

2.2. Transcriptome Functional and Pathway Annotation

For the full annotation, the filtered contigs were subjected to BLAST [31] runs against the SwissProt [32], ToxProt [33] and the NCBI non-redundant protein databases [34]. Overall, 84.77% of the contigs were annotated against NCBI, 80.26% to SwissProt, 1.69% to ToxProt and 13.40% were left unannotated (Table 2). The top five species hits in the NCBI database were Callorhinchus milii, Latimeria chalumnae, Lepisosteus oculatus, Xenopus tropicalis and Chrysemys picta bellii, matching previous data with hits to fish (even though very few fish species were represented in the databases).
Similar to previous studies, little sequence diversity and closeness was found in the BLAST hits to the venom gland contigs (Figure S1) [28]. This is due to the lack of fish sequences in the SwissProt and the NCBI databases, which increases the complexity in the proper identification and characterization of venom components in species not closely related. Our data will assist identification in future studies but more effort in a higher diversity of species would be welcome.

2.3. Gene Ontology and Metabolic Pathways

Gene Ontology (GO) [35] analysis was performed to understand the functional properties of the contigs. Of the 27,032 predicted coding transcripts, 19,194 (77.60%) were successfully mapped into the three major functional groups: cellular component, molecular function, and biological process (Figure 1A).
GO terms annotation in the cellular components presented peaks for cell and membrane, with organelles closely following, a slight difference in comparison to P. amandae and P. falkneri (Figure 2A). As with the other members of the genus, cellular, metabolic, and biological regulation processes were the highest matching categories in biological processes. Regarding molecular functions, binding and catalytic activity were also the most mapped (Figure 1A).
The consistency found in the results of the transcriptomic characterization of three different stingray species of the Potamotrygon genus suggest that despite the inter-specific variation, key processes and components appear overabundant in the venom apparatus.
The Kyoto Encyclopedia of Genes and Genomes (KEGG) [36] annotation can assist in the understanding of high-level functions and utilities of biological systems. Using the previous annotation retrieved from the SwissProt database it was possible to identify the KEGG pathways in which the transcripts were involved. There were 17,572 (71.04%) transcripts associated with KEGG pathways. Compared to the results obtained in [26], Cancers, Signal transduction, Immune system and Infectious diseases were still the most represented groups. However, the Endocrine system also reached values comparable to that of the Immune system (Figure 1B). It should also be noted that only around 6% of contigs were assigned to this pathway, compared to the 8% found in P. amandae and P. falkneri transcripts (Figure 2B). The distribution of the KEGG pathway annotation through the six major categories—Metabolism, Genetic Information Processing, Environmental Information Processing, Cellular Processes, Organismal Systems, and Human Diseases—differs from the other members of the genus, overall. Pathways involved in cancer were also the pathways with the highest matches, followed by the PI3K-Akt signaling pathway, HTLV-I infection, and Endocytosis. This suggests that while being closely related, the Potamotrygon genus possess significant diversity within its species and their venoms.
In respect to its metabolism, an overview of the KEGG pathway distribution showed 46 hits to the biosynthesis of amino acids and 138 to the biosynthesis of antibiotics. The amino acid metabolism is essential for the production of the venomous cocktail, while the antibiotic pathways is related with the production of bioactive components. There were also pathways identified in the carbon (78) and fatty acid (30) metabolisms.
The functional annotation of our dataset yielded a greater number of annotated contigs in KEGG pathways (17,572–71.0%) and in gene ontology (19,194–77.6%) compared to former works (KEGG: 14,131–56.3% and 13,147–59.5%; GO: 16,921–64.9% and 15,394–69.7%; P. amandae and P. falkneri, respectively) in the other two species of stingrays. This may be the underlying cause for the KEGG pathway and gene ontology to differ in relative values.

2.4. Venom Expression

The contigs of the P. motoro transcriptome were also searched for similarity to the UniProt ToxProt database, a curated list of toxins and venoms. High similarity hits included zinc metalloproteinases, hyaluronidases, venom prothrombin activator and alpha-latrocrustotoxin-Lt1a. These molecules play key roles in venom functionality. Several of the identified venoms have been previously described in stingrays. These include proteins involved in ion channels disruption and inhibition of muscle contraction, neuromuscular transmission and nervous system interference and affecting the circulatory system and its related processes [26]. The most abundant transcripts with hits to the database presented some common molecules with other species of the genus, such as cysteine-rich venom protein, cystatin-2, calglandulins, metalloproteinases, peroxiredoxin-4 and translationally controlled tumor protein homolog (Table 3). Interestingly, hyaluronidase is the second most abundant, but not the first as in the other two stingray species. Two particular transcripts that matched SE-cephalotoxin are highly expressed (Table 3). The fact that this molecule did not appear in the venom gland transcriptomes of both P. amandae and P. falkneri and possessed such high values led to a more in-depth identification of the coding region. BLAST searches with the sequence matched mostly hypothetical undescribed proteins, SE-cephalotoxin-like molecules and proteins encoded by the Rapunzel gene.
The reason the Rapunzel gene would be overrepresented in a venom gland transcriptome is not fully known. It is difficult to understand the underlying implications of this molecule primarily since it is not clearly identified through current databases. Assuming it is related to SE-cephalotoxin, then this result shows a very significant difference in venom composition within the Potamotrygon genus. It also represents a potential key component of the P. motoro venom. If similar to Rapunzel, a gene involved in bone structure formation, it will be necessary to determine first the function of the gene in stingrays and its implications as part of a venom cocktail. In the possibility of being another molecule, further studies would be required to understand its function and the possible evolutionary relationship with SE-cephalotoxin and Rapunzel. It should be noted that the hits to the NCBI database regarding SE-cephalotoxin-like genes were to non-venomous species and both these hits and those to Rapunzel, despite the low e-values (~1 × 10−70), showed only moderate sequence similarity (~40%).
Notably, we have not detected the presence of contigs in P. motoro with high similarities to phospholipase A2, previously described in other members of the genus, as well as in the common stingray Dasyatis pastinaca, and in envenomation scenarios causing symptoms such as hemorrhage and necrosis [37].

Unique Toxins Identified in P. motoro and Fish Venom Comparison

In addition to the SE-cephalotoxin/ Rapunzel transcripts, unexpected finds at high levels of expression included the DELTA-alicitoxin-Pse1b, Augerpeptide hhe53 and PI-actitoxin-Aeq3a. The transcripts for these molecules or similar have yet to be described in stingray venom. DELTA-alicitoxin-Pse1b is a pore-forming protein that can provoke hemolysis, making it a toxin focused on damaging the circulatory system [38]. Augerpeptide hhe53 is a protein expressed in venom ducts of the sea snail Hastula hectica and little is known about its function [39]. PI-actitoxin-Aeq3a is both a serine protease trypsin and potassium channel inhibitor and was originally described as a Kunitz-type serine protease inhibitor in the sea anemone Actinia equine [40], with consequences similar to those of the serine protease inhibitors, capable of affecting coagulation and the kinin system processes. Recently, P. amandae has been described as a distinct species from P. motoro that is believed to be a species complex [1]. Our results would sustain that hypothesis but could also be indicative of intraspecific variation.
It is worth mentioning that transcripts with similarity to venom phosphodiesterase were found, although at much lower expression values (<3.6 FPKM). Phosphodiesterase activity was previously found in the stingray Urolophus halleri [14] and in the bony fish Scatophagus argus, Gymnapistes marmoratus and Synanceia horrida [41]. However, there was no mention of these molecules in the most recent transcriptomic studies in stingrays, including those of the Potamotrygon genus. This protein is an enzyme with nuclease, pyrophosphatase and phosphatase activity and will hydrolyze nucleoside 5′-triphosphates and 5′-diphosphates, though not 5′-monophosphates. Together with 5′-nucleotidase, it has been reported as responsible for potential tissue necrosis and capable of inhibiting platelet aggregation induced by ADP [14,42,43,44].
Overall, within the high expression values, not including hyaluronidase, only phospholipase A2 of the key venom components found in bony fish (Table 4) was detected in the Potamotrygon species, although entirely absent from the Potamotrygon motoro transcriptome. Compared to the proteins of the marine stingray Neotrygon kuhlii extracted from the barb venom gland, only the cystatin and peroxiredoxin were found in common (Table 4).
The hyaluronidase transcripts showed elevated expression levels and this protein is believed to be a key component in the Potamotrygon genus venom. Thus, to better understand the hyaluronidase gene family evolution we performed phylogenetic analyses across all fish for which the hyaluronidase gene information is publicly available.

2.5. Hyaluronidase Phylogeny and Selective Pressure

We estimated a phylogenetic reconstruction of a total of 80 hyaluronidases sequences found across fish and mammal representatives of all HYAL functional groups using maximum likelihood methods. The produced trees strengthened the proximity previously inferred among hyaluronidase present in the venomous glands of Scorpaeniformes and the PH-20 isotype (HYAL5) [46], located on the sperm surface acting as a receptor for the zona pellucida surrounding the oocyte and responsible for the degradation of extracellular matrix (ECM) during oocyte fertilization [47]. It also reinforced the idea that the hyaluronidases found in the venoms of stingrays are not immediately related to the hyaluronidases of venomous bony fish, in particular those found in Scorpaeniformes. The ocellate river stingray P. motoro possessed three well-defined sequences of hyaluronidase. Interestingly, the phylogenetic analyses placed the most expressed transcript together with a zebrafish HYAL6 sequence in a clade, with no proximity to any potentially venomous fish (Figure 3). Another transcript was grouped in the clade containing the HYAL2 sequences and together with its low expression level indicates a likelihood of this particular protein being involved in house-keeping activities. Similarly, the other transcript was likely also part of house-keeping processes. We were unable to properly resolve the hyaluronidase group for this sequence. This may be due to the distance between cartilaginous and bony fish and the similarity of the HYAL1 and HYAL3 groups. BLAST identification indicated the transcript to be similar to the HYAL1 group sequences. These results do not create a monophyletic group with other venomous animals, supporting previous transcriptomic studies in stingrays [28,29] and the more precise phylogenetic study of venomous fish that showed up to 18 independent origins (convergent evolution) of venom in fish [48].
Interestingly, the fact that the most expressed type of hyaluronidase (HYAL6) found in the venom gland transcriptome of all three freshwater stingrays were all grouped together suggests a common origin of the venom specialization. The amino acid sequence similarity and the phylogenetic proximity to the zebrafish HYAL6 sequence compared to the other types of HYAL suggests this gene as a candidate for the hyaluronidase recruitment into the venom of the Potamotrygon lineage. However, data in fish regarding hyaluronidase is limited, especially in HYAL6. The sequences besides the zebrafish found in the group were mostly predicted hyaluronidase-like and with no group indicated. Given the support for both Danio rerio H6 and Mus musculus H6 phylogenetic position and the resolution of the hyaluronidase family members, it is likely that this group represents HYAL6. The Rhincodon typus and the Callorhinchus milii sequences were likely incorrectly annotated, potentially due to the absence of HYAL6 in the annotation database at the time, as it is a pseudogene in Homo sapiens (i.e., no protein sequence available). Performing a BLAST with Mus musculus sequence against the NCBI database returns only hits to itself and a HYAL6 sequence from Rattus norvegicus, while the remaining hits represent HYAL4, suggesting high sequence similarity and the likely reason for the annotation found in the cartilaginous fish. Regardless, this is the first identification of a potentially recruited gene that would also match with an independent evolution of venom in fish. The hyaluronidases of venomous bony fish were found to be closely related to the HYAL5 as opposed to HYAL6 (Figure 3), indicating that different groups in the hyaluronidase family were recruited into venoms. Currently, little is known about HYAL6 other than it is highly expressed in murine testicular tissue together with HYAL5. HYAL6 was determined to not have hyaluronidase activity when pH values were neutral [49].
Using the HyPhy software package [50], several selection analyses algorithms were used. As expected with coding sequences, most sites were found to be under negative selection ensuring that a functional molecule is encoded. Results from Single Likelihood Ancestor Counting (SLAC), Random Effects Likelihood (REL), Fixed Effects Likelihood (FEL) [51] and Fast Unconstrained Bayesian AppRoximation (FUBAR) [52] analyses did not reveal any sites with evidence of positive selection (317, 529, 374 and 415 negatively selected sites, respectively). The prominent negative selection profile may have masked possible sites under positive selection and thus the Mixed Effects Model of Evolution (MEME) [53] model was used to identify sites under episodic diversifying selection (codons: 20, 53, 58, 79, 103, 113, 114, 120, 159, 179, 190, 234, 245, 488, 501, 505, 523 and 526).
The three-dimensional protein structural modelling of hyaluronidase allowed the visualization of the potential location of the specific residues. To further improve the accuracy of the residue selective pressure, only sites still significant with P < 0.01 were marked (113, 114, 120 and 505; 523 could not be modelled). The residues under diversifying selection were found mostly in the protein surface. Interestingly, some of the residues under positive selection were found to be close to each other in the modelled three-dimensional structure (Figure 4). This pattern corroborates a relation between positive selection and the protein secondary structure elements studied in Drosophila. It was found that amino acids in disordered regions had a higher chance to be under positive selection in relation to their proportion in the proteins. In comparison, residues in helices and beta-structures had less sites under positive selection than expected. It is noteworthy to mention that sites under diversifying pressure occurred close to each other in the protein sequence more often than anticipated [54].
Inspecting the protein model revealed that the positively selected sites were in the proximity of the catalytic region. Although the proximity to the catalytic region could mean these sites affect the protein activity, it does not directly translate to a guaranteed effect. However, changes to these residues may result in slight structural shifts, altering interactions with other molecules, the molecule stability or affecting the orientation of the catalytic residues. Our results suggest that the enzymatic reaction is conserved, and the modifications required for the protein venom role would be optimizations on the surface recognition and related protein-protein interactions.
The fact that these changes are observed across an analysis of multiple family members of hyaluronidase involved in different processes, but all catalyzing the degradation of hyaluronic acid, shows that altering the tissue permeability is a key factor in a wide variety of processes, ranging from usual metabolic pathways to venomous activity.
This demonstrates that while hyaluronidase recruitment into venoms in fish was likely independent from other venomous animals, not all genes part of the venom arsenal had to undergo rapid mutation to retain relevance. For example, the king cobra genome showed that the hyaluronidase was not found under diversifying selection, suggesting that auxiliary genes not causing resistance in the targets are under lower selective pressure [57]. This would also be insightful to understand the multiple cases of convergent evolution seen in venoms across a wide taxonomic array.
Hyaluronidase has been recruited into venoms along a wide range of taxonomic groups and it has been suggested to use the protein as a therapeutic target [58]. This convergent evolution suggests that the molecule is well suited for an important role in venoms, whether as a facilitating agent or something else. The high expression values found in the venom glands transcriptomes of the freshwater stingrays and the small number of positively selected residues, along with their location in the protein structure, strongly suggests that hyaluronidase plays a more important role than previously recognized, which should be further assessed in future proteomic studies.

3. Conclusions

The transcriptomic characterization of the P. motoro venom gland revealed that despite the species being closely related to other freshwater stingrays (i.e., same genus) there might be significant variation in the venom composition among species. Notably, we found contigs with hits to sequences of DELTA-alicitoxin-Pse1b, Augerpeptide hhe53, PI-actitoxin-Aeq3a and SE-cephalotoxin with significant expression levels, when compared to the most expressed molecules of P. amandae and P. falkneri. The most notable absence in the venom of P. motoro was the phospholipase A2. The venom gland of P. motoro possessed multiple transcripts that were mapped to the pathways of antibiotic synthesis, such as terpenoid and polyketide metabolism, and biosynthesis of secondary metabolites, in particular the pathways of penicillin, cephalosporin, streptomycin, neomycin, kanamycin, and gentamicin. Also identified were pathways related to diseases, these being cancers, immune diseases, substance dependence, cardiovascular diseases, endocrine and metabolic diseases, neurodegenerative and infectious diseases.
Detailed analyses of the hyaluronidase showed that its recruitment into the venom of freshwater stingrays was an independent event relatively to other venomous bony fish species. The hyaluronidase genes have been the target of diversifying selection, indication that across duplications and gene mutations, these molecules were recruited and adjusted into several different processes. Interestingly, there is no evidence suggesting that the protein optimization required changes to the original enzymatic reaction. Rather, the molecule structural stability, substrate affinity and protein-protein interaction would be the elements most likely to have changed.
Our work expands the existing knowledge in cartilaginous fish venom composition. It also reinforced the importance of hyaluronidase in venoms, evidence of its independent recruitment and indication that not all venom components are subject to high diversifying pressures. The knowledge of venom composition of fish is still very narrow given the huge diversity of venomous species. In particular, to fully understand all the independent specializations occurring in this group more information on the involved molecules and their sequences is required. The more data becomes available the clearer the role and potential of the bioactive components used by venomous animals, such as freshwater stingrays.

4. Materials and Methods

4.1. Fish Collection and Sample Processing

Here we studied the venomous freshwater stingray P. motoro (PM, Taxonomy ID: 86373). The fish was bred by the Pearl River Basin Sub-center, National Sharing Service Infrastructure of Fishery Germplasm Resources, Guangzhou, China.
The species identification was supported by both morphology and DNA barcoding experiments (with identity > 99%). Overall, 10 spine samples, including the venom glands (one from each fish), were dissected, immediately snap-frozen and deposited in liquid nitrogen tanks for future processing. All the collection and processing procedures have been reviewed and approved by the Institutional Review Board on Bioethics and Biosafety of BGI (No. BGI-IRB 15139; 27 November 2015).

4.2. RNA Extraction and Sequencing

Total RNA was isolated from the 10 pooled venom glands using TRIzol® LS Reagent (Invitrogen, Carlsbad, CA, USA) and the quality of each RNA sample was assessed by Agilent 2100 Bioanalyzer (Agilent Technologies, Palo Alto, CA, USA). Afterwards, purification and isolation of mRNA with poly(A) tails was performed using Oligo-(dT)-attached magnetic beads, and the obtained mRNA was proceeded to Illumina cDNA library construction and sequencing through Illumina HiSeq2000 platform (Illumina, San Diego, CA, USA) at BGI-Tech (BGI, Shenzhen, Guangdong, China). The lack of a well-defined gland structure makes tissue extraction challenging and thus not all the molecules found may be representative of the stingray venom [59,60].

4.3. Data Filtering and De novo Assembly

Paired-end raw reads (China National Genebank Database, project accession CNP0000235) generated from the sequencing platform were filtered with SOAPnuke [61] to remove the junk reads with adaptors, more than 10% of N bases and more than 50% of low-quality bases (base quality score ≤ 10). The remaining clean data were then de novo assembled into contigs using the Trinity software (v2.1.1) with the specific parameter set to “--min_contig_length 150” and others set to default [35]. TGICL software (v2.0) was used to eliminate the redundancies in the assembly with given parameter “-v 25 -O ‘-repeat_stringency 0.95’” [62].

4.4. Transcript Expression Calculation

Expression values were calculated by mapping the raw reads against the assembled contigs using CLC’s Genomics Workbench (v9.5.2) RNA-Seq analysis option. Multiple statistics were obtained, including fragments per kilobase of transcript per million mapped reads.

4.5. Transcriptome, Functional and Pathway Annotations

The transcripts were subject to a filter dependent on their expression, in which only contigs with values higher than 1 FPKM were considered for further analyses. The remaining contigs were then processed through Trinity’s Trinotate software [63] to fully annotate and to determine the function and pathways of the molecules encoded by the transcripts retrieved from the venom gland. The annotation part of this process included the necessary SwissProt [32] and Pfam [64] databases, as well as NCBI’s non-redundant protein database [34] (restricted to vertebrates) and Uniprot’s ToxProt [33]. All of these were target databases for the contigs’ BLAST analysis [31]. Through this process, the TransDecoder [65] software determined coding regions in the contigs, retaining only opening reading frames (ORFs) longer than 100 amino acids. GO [35] and pathway data (KEGG) [36] was extracted from hits to the SwissProt and Pfam databases. GO term annotation results were exported using WEGO webserver [66].

4.6. Sequence Alignment and Phylogenetic Analyses

Nucleotide sequences corresponding to hyaluronidases were retrieved from the GenBank [34] and RefSeq [67]. The search was performed using BLAST. The search queries were the coding sequences (CDS) of hyaluronidases from Synanceia horrida already available in GenBank and the hyaluronidase transcripts from P. motoro. Representative sequences from different hyaluronidase functional groups were manually searched for and retrieved. The Homo sapiens sequences were used to more accurately identify functional HYAL groups. The Mus musculus HYAL6 was used instead of the human HYAL6 that is a pseudogene. All sequences used, and their accession or reference, are found in Table S1.
Nucleotide sequences were translated into amino acid sequences and aligned using the CLUSTAL [68] algorithm within SEAVIEW [69] (v4.5 4), which was also used for visualization and manual editing of the sequences [70]. Regions in which homology could not be guaranteed were removed from further analyses.
The most appropriate model of protein evolution for Maximum Likelihood phylogenetic trees was chosen using ProtTest [71,72], which returned WAG + I + G. Maximum Likelihood phylogenetic trees were reconstructed using IQ-TREE [73]. Maximum likelihood tree branches were supported by 1000 bootstraps.

4.7. Selective Pressure Analyses

In the event of recombination, single tree topologies cannot explain the evolutionary path of the recombined sequence. In this scenario, the likelihood-ratio test (LRT) viability may be doubtful. LRT was shown to be robust at low levels of recombination (less than three recombination events with 10 sequences) [74]. To account for possible interference, a genetic algorithm for recombination detection (GARD) and Single Breakpoints (SBP) algorithms [75,76,77,78,79] of the HyPhy [50] multiplatform were used. When relevant levels of recombination were found, sequences were partitioned prior to the selection analyses.
The HyPhy multiplatform was used to assess selective pressure, running the SLAC algorithm, a derivative of the Suzuki—Gojobori counting approach; the FEL algorithm, estimating non-synonymous and synonymous substitution rates at each site; the REL algorithm [51,78,79], which categorizes non-synonymous and synonymous rates variation across all sites and infers selective pressure in sites using an empirical Bayes approach; and the FUBAR algorithm [52,76,77,79], which uses a Markov chain Monte Carlo method and allows visualization of Bayesian inference for each site. To detect episodic diversifying selection, masked by heavy purifying pressure, the MEME algorithm [53,76,77,78,79] was also used.

4.8. Structure Modelling

Computational structural methods have been successfully used to identify relevant molecular interactions of protein involved in complex pathways such as toxin target and protein inhibition mechanisms [77,78,79,80,81]. The three-dimensional (3D) molecular structure of hyaluronidase was predicted using the Phyre2 webserver (normal modelling mode) [55,77,78,79]. The obtained models were viewed using VMD—Visual Molecular Dynamics [56]. Amino acids under positive selection were marked in the structure for posterior analysis. Hyaluronidase catalytic residues were annotated accordingly to The Catalytic Site Atlas 2.0 [82].

Supplementary Materials

The following data are available online at https://www.mdpi.com/2072-6651/10/12/544/s1, Figure S1: Hits distribution representation per species following BLAST searches of the Potamotrygon motoro transcriptome against the SwissProt database; Table S1: Accession numbers for the sequences used in phylogenetic and selective pressure analyses.

Author Contributions

Conceptualization, F.S., Y.H. and A.A.; methodology and software, F.S. and Y.H.; data curation, F.S. and Y.H.; validation and formal analysis, F.S., Y.H. and V.Y.; supervision, A.A.; project administration, A.A.; funding acquisition, A.A.; resources, X.M., Q.S., and A.A.; writing—original draft preparation, F.S. and A.A.; writing—review and editing, Y.H., Q.S., and A.A.

Funding

Filipe Silva was supported by a PhD grant (SFRH/BD/126560/2016) from the Portuguese Foundation for Science and Technology (FCT—Fundação para a Ciência e a Tecnologia, Portugal). Xidong Mu was funded by the National Sharing Service Infrastructure of Fishery Germplasm Resources (2018DKA30470). Agostinho Antunes was funded in part by the Strategic Funding UID/Multi/04423/2013 through national funds provided by FCT and the European Regional Development Fund (ERDF) in the framework of the program PT2020, by the European Structural and Investment Funds (ESIF) through the Competitiveness and Internationalization Operational Program—COMPETE 2020 and by National Funds through the FCT under the project PTDC/AAG-GLO/6887/2014 (POCI-01-0124-FEDER-016845).

Acknowledgments

We thank Dany Domínguez Pérez, Emanuel Maldonado, Daniela Almeida and Anoop Alex (CIIMAR, University of Porto) for their advices and comments. Special thanks go to all members of the Transcriptomes of 1000 Fish (Fish-T1K) project.

Conflicts of Interest

The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

References

  1. Loboda, T.S.; Carvalho, M.R. Systematic revision of the Potamotrygon motoro (Müller & Henle, 1841) species complex in the Paraná-Paraguay basin, with description of two new ocellated species (Chondrichthyes: Myliobatiformes; Potamotrygonidae). Neotrop. Ichthyol. 2013, 11, 693–737. [Google Scholar]
  2. Carvalho, M.; Lovejoy, N. Morphology and phylogenetic relationships of a remarkable new genus and two new species of Neotropical freshwater stingrays from the Amazon basin (Chondrichthyes: Potamotrygonidae). Zootaxa 2011, 2776, 13–48. [Google Scholar]
  3. De Carvalho, M.R.; Loboda, T.S.; da Silva, J.P.C.B. A new subfamily, Styracurinae, and new genus, Styracura, for Himantura schmardae (Werner, 1904) and Himantura pacifica (Beebe & Tee-Van, 1941) (Chondrichthyes: Myliobatiformes). Zootaxa 2016, 4175, 21. [Google Scholar] [CrossRef]
  4. Júlio Júnior, H.F.; Tós, C.D.; Agostinho, Â.A.; Pavanelli, C.S. A massive invasion of fish species after eliminating a natural barrier in the upper rio Paraná basin. Neotrop. Ichthyol. 2009, 7, 709–718. [Google Scholar] [CrossRef] [Green Version]
  5. Oddone, M.C.; Velasco, G.; Charvet, P. Record of the freshwater stingrays Potamotrygon brachyura and P. motoro (Chondrichthyes, Potamotrygonidae) in the lower Uruguay river, South America. Acta Amaz. 2012, 42, 299–304. [Google Scholar] [CrossRef]
  6. Barbaro, K.C.; Lira, M.S.; Malta, M.B.; Soares, S.L.; Garrone Neto, D.; Cardoso, J.L.; Santoro, M.L.; Haddad Junior, V. Comparative study on extracts from the tissue covering the stingers of freshwater (Potamotrygon falkneri) and marine (Dasyatis guttata) stingrays. J. Int. Soc. Toxinol. 2007, 50, 676–687. [Google Scholar] [CrossRef] [PubMed]
  7. Haddad, V., Jr.; Neto, D.G.; de Paula Neto, J.B.; de Luna Marques, F.P.; Barbaro, K.C. Freshwater stingrays: Study of epidemiologic, clinic and therapeutic aspects based on 84 envenomings in humans and some enzymatic activities of the venom. J. Int. Soc. Toxinol. 2004, 43, 287–294. [Google Scholar] [CrossRef] [PubMed]
  8. Pedroso, C.M.; Jared, C.; Charvet-Almeida, P.; Almeida, M.P.; Garrone Neto, D.; Lira, M.S.; Haddad, V., Jr.; Barbaro, K.C.; Antoniazzi, M.M. Morphological characterization of the venom secretory epidermal cells in the stinger of marine and freshwater stingrays. J. Int. Soc. Toxinol. 2007, 50, 688–697. [Google Scholar] [CrossRef] [PubMed]
  9. Williamson, J.A.; Burnett, J.W.; Fenner, P.J.; Rifkin, J.F. Venomous and Poisonous Marine Animals: A Medical and Biological Handbook; University of New South Wales Press: Kensington, Australia, 1996. [Google Scholar]
  10. Junghanss, T.; Bodio, M. Medically important venomous animals: Biology, prevention, first aid, and clinical management. Clin. Infect. Dis. 2006, 43, 1309–1317. [Google Scholar] [CrossRef] [PubMed]
  11. Barber, G.R.; Swygert, J.S. Necrotizing fasciitis due to Photobacterium damsela in a man lashed by a stingray. N. Engl. J. Med. 2000, 342, 824. [Google Scholar] [CrossRef] [PubMed]
  12. Domingos, M.O.; Franzolin, M.R.; dos Anjos, M.T.; Franzolin, T.M.; Barbosa Albes, R.C.; de Andrade, G.R.; Lopes, R.J.; Barbaro, K.C. The influence of environmental bacteria in freshwater stingray wound-healing. J. Int. Soc. Toxinol. 2011, 58, 147–153. [Google Scholar] [CrossRef] [PubMed]
  13. Clark, R.F.; Girard, R.H.; Rao, D.; Ly, B.T.; Davis, D.P. Stingray envenomation: A retrospective review of clinical presentation and treatment in 119 cases. J. Emerg. Med. 2007, 33, 33–37. [Google Scholar] [CrossRef] [PubMed]
  14. Diaz, J.H. The evaluation, management, and prevention of stingray injuries in travelers. J. Travel Med. 2008, 15, 102–109. [Google Scholar] [CrossRef] [PubMed]
  15. Junior, V.H.; Cardoso, J.L.; Neto, D.G. Injuries by marine and freshwater stingrays: History, clinical aspects of the envenomations and current status of a neglected problem in Brazil. J. Venom. Anim. Toxins Incl. Trop. Dis. 2013, 19, 16. [Google Scholar] [CrossRef] [PubMed]
  16. Torrez, P.P.; Quiroga, M.M.; Said, R.; Abati, P.A.; Franca, F.O. Tetanus after envenomations caused by freshwater stingrays. J. Int. Soc. Toxinol. 2015, 97, 32–35. [Google Scholar] [CrossRef] [PubMed]
  17. Fry, B.G.; Roelants, K.; Norman, J.A. Tentacles of venom: Toxic protein convergence in the Kingdom Animalia. J. Mol. Evol. 2009, 68, 311–321. [Google Scholar] [CrossRef] [PubMed]
  18. Fry, B.G.; Wuster, W.; Kini, R.M.; Brusic, V.; Khan, A.; Venkataraman, D.; Rooney, A.P. Molecular evolution and phylogeny of elapid snake venom three-finger toxins. J. Mol. Evol. 2003, 57, 110–129. [Google Scholar] [CrossRef] [PubMed]
  19. Ruder, T.; Sunagar, K.; Undheim, E.A.; Ali, S.A.; Wai, T.C.; Low, D.H.; Jackson, T.N.; King, G.F.; Antunes, A.; Fry, B.G. Molecular phylogeny and evolution of the proteins encoded by coleoid (cuttlefish, octopus, and squid) posterior venom glands. J. Mol. Evol. 2013, 76, 192–204. [Google Scholar] [CrossRef] [PubMed]
  20. Smith, W.L.; Wheeler, W.C. Venom Evolution Widespread in Fishes: A Phylogenetic Road Map for the Bioprospecting of Piscine Venoms. J. Hered. 2006, 97, 206–217. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  21. Undheim, E.A.; Sunagar, K.; Herzig, V.; Kely, L.; Low, D.H.; Jackson, T.N.; Jones, A.; Kurniawan, N.; King, G.F.; Ali, S.A.; et al. A proteomics and transcriptomics investigation of the venom from the barychelid spider Trittame loki (brush-foot trapdoor). Toxins 2013, 5, 2488–2503. [Google Scholar] [CrossRef] [PubMed]
  22. Sunagar, K.; Undheim, E.A.; Chan, A.H.; Koludarov, I.; Munoz-Gomez, S.A.; Antunes, A.; Fry, B.G. Evolution stings: The origin and diversification of scorpion toxin peptide scaffolds. Toxins 2013, 5, 2456–2487. [Google Scholar] [CrossRef] [PubMed]
  23. Freitas-de-Sousa, L.A.; Amazonas, D.R.; Sousa, L.F.; Sant’Anna, S.S.; Nishiyama, M.Y.; Serrano, S.M.T.; Junqueira-de-Azevedo, I.L.M.; Chalkidis, H.M.; Moura-da-Silva, A.M.; Mourão, R.H.V. Comparison of venoms from wild and long-term captive Bothrops atrox snakes and characterization of Batroxrhagin, the predominant class PIII metalloproteinase from the venom of this species. Biochimie 2015, 118, 60–70. [Google Scholar] [CrossRef] [PubMed]
  24. Sunagar, K.; Undheim, E.A.B.; Scheib, H.; Gren, E.C.K.; Cochran, C.; Person, C.E.; Koludarov, I.; Kelln, W.; Hayes, W.K.; King, G.F.; et al. Intraspecific venom variation in the medically significant Southern Pacific Rattlesnake (Crotalus oreganus helleri): Biodiscovery, clinical and evolutionary implications. J. Proteom. 2014, 99, 68–83. [Google Scholar] [CrossRef] [PubMed]
  25. Conceição, K.; Konno, K.; Melo, R.L.; Marques, E.E.; Hiruma-Lima, C.A.; Lima, C.; Richardson, M.; Pimenta, D.C.; Lopes-Ferreira, M. Orpotrin: A novel vasoconstrictor peptide from the venom of the Brazilian Stingray Potamotrygon gr. orbignyi. Peptides 2006, 27, 3039–3046. [Google Scholar] [CrossRef] [PubMed]
  26. Conceição, K.; Santos, J.M.; Bruni, F.M.; Klitzke, C.F.; Marques, E.E.; Borges, M.H.; Melo, R.L.; Fernandez, J.H.; Lopes-Ferreira, M. Characterization of a new bioactive peptide from Potamotrygon gr. orbignyi freshwater stingray venom. Peptides 2009, 30, 2191–2199. [Google Scholar] [CrossRef] [PubMed]
  27. Magalhães, M.R.; da Silva, N.J.; Ulhoa, C.J. A hyaluronidase from Potamotrygon motoro (freshwater stingrays) venom: Isolation and characterization. J. Int. Soc. Toxinol. 2008, 51, 1060–1067. [Google Scholar] [CrossRef] [PubMed]
  28. Júnior, N.G.d.O.; Fernandes, G.d.R.; Cardoso, M.H.; Costa, F.F.; Cândido, E.d.S.; Neto, D.G.; Mortari, M.R.; Schwartz, E.F.; Franco, O.L.; de Alencar, S.A. Venom gland transcriptome analyses of two freshwater stingrays (Myliobatiformes: Potamotrygonidae) from Brazil. Sci. Rep. 2016, 6, 21935. [Google Scholar] [CrossRef] [PubMed]
  29. Baumann, K.; Casewell, N.R.; Ali, S.A.; Jackson, T.N.; Vetter, I.; Dobson, J.S.; Cutmore, S.C.; Nouwens, A.; Lavergne, V.; Fry, B.G. A ray of venom: Combined proteomic and transcriptomic investigation of fish venom composition using barb tissue from the blue-spotted stingray (Crotalus adamanteus). J. Proteom. 2014, 109, 188–198. [Google Scholar] [CrossRef] [PubMed]
  30. Stern, R.; Kogan, G.; Jedrzejas, M.J.; Šoltés, L. The many ways to cleave hyaluronan. Biotechnol. Adv. 2007, 25, 537–557. [Google Scholar] [CrossRef] [PubMed]
  31. Camacho, C.; Coulouris, G.; Avagyan, V.; Ma, N.; Papadopoulos, J.; Bealer, K.; Madden, T.L. BLAST+: Architecture and applications. BMC Bioinform. 2009, 10, 421. [Google Scholar] [CrossRef] [PubMed]
  32. UniProt: The universal protein knowledgebase. Nucleic Acids Res. 2017, 45, D158–D169. [CrossRef] [PubMed]
  33. Jungo, F.; Bougueleret, L.; Xenarios, I.; Poux, S. The UniProtKB/Swiss-Prot Tox-Prot program: A central hub of integrated venom protein data. J. Int. Soc. Toxinol. 2012, 60, 551–557. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  34. Clark, K.; Karsch-Mizrachi, I.; Lipman, D.J.; Ostell, J.; Sayers, E.W. GenBank. Nucleic Acids Res. 2016, 44, D67–D72. [Google Scholar] [CrossRef] [PubMed]
  35. Ashburner, M.; Ball, C.A.; Blake, J.A.; Botstein, D.; Butler, H.; Cherry, J.M.; Davis, A.P.; Dolinski, K.; Dwight, S.S.; Eppig, J.T.; et al. Gene ontology: Tool for the unification of biology. The Gene Ontology Consortium. Nat. Genet. 2000, 25, 25–29. [Google Scholar] [CrossRef] [PubMed]
  36. Kanehisa, M.; Goto, S.; Sato, Y.; Furumichi, M.; Tanabe, M. KEGG for integration and interpretation of large-scale molecular data sets. Nucleic Acids Res. 2012, 40, D109–114. [Google Scholar] [CrossRef] [PubMed]
  37. Ben Bacha, A.; Abid, I.; Horchani, H.; Mejdoub, H. Enzymatic properties of stingray Dasyatis pastinaca group V, IIA and IB phospholipases A(2): A comparative study. Int. J. Biol. Macromol. 2013, 62, 537–542. [Google Scholar] [CrossRef] [PubMed]
  38. Uechi, G.; Toma, H.; Arakawa, T.; Sato, Y. Molecular characterization on the genome structure of hemolysin toxin isoforms isolated from sea anemone Actineria villosa and Phyllodiscus semoni. J. Int. Soc. Toxinol. 2010, 56, 1470–1476. [Google Scholar] [CrossRef] [PubMed]
  39. Imperial, J.S.; Kantor, Y.; Watkins, M.; Heralde, F.M., 3rd; Stevenson, B.; Chen, P.; Hansson, K.; Stenflo, J.; Ownby, J.P.; Bouchet, P.; et al. Venomous auger snail Hastula (Impages) hectica (Linnaeus, 1758): Molecular phylogeny, foregut anatomy and comparative toxinology. J. Exp. Zool. Part B Mol. Dev. Evol. 2007, 308, 744–756. [Google Scholar] [CrossRef] [PubMed]
  40. Ishida, M.; Minagawa, S.; Miyauchi, K.; Shimakura, K.; Nagashima, Y.; Shiomi, K. Amino Acid Sequences of Kunitz-type Protease Inhibitors from the Sea Anemone Actinia Equina. Fish. Sci. 1997, 63, 794–798. [Google Scholar] [CrossRef]
  41. Ziegman, R.; Alewood, P. Bioactive Components in Fish Venoms. Toxins 2015, 7, 1497. [Google Scholar] [CrossRef] [PubMed]
  42. Mitra, J.; Bhattacharyya, D. Phosphodiesterase from Daboia russelli russelli venom: Purification, partial characterization and inhibition of platelet aggregation. J. Int. Soc. Toxinol. 2014, 88, 1–10. [Google Scholar] [CrossRef] [PubMed]
  43. Belli, S.I.; Goding, J.W. Biochemical characterization of human PC-1, an enzyme possessing alkaline phosphodiesterase I and nucleotide pyrophosphatase activities. Eur. J. Biochem. 1994, 226, 433–443. [Google Scholar] [CrossRef] [PubMed]
  44. Margres, M.J.; McGivern, J.J.; Wray, K.P.; Seavy, M.; Calvin, K.; Rokyta, D.R. Linking the transcriptome and proteome to characterize the venom of the eastern diamondback rattlesnake (Crotalus adamanteus). J. Proteom. 2014, 96, 145–158. [Google Scholar] [CrossRef] [PubMed]
  45. Casewell, N.R.; Visser, J.C.; Baumann, K.; Dobson, J.; Han, H.; Kuruppu, S.; Morgan, M.; Romilio, A.; Weisbecker, V.; Mardon, K.; et al. The Evolution of Fangs, Venom, and Mimicry Systems in Blenny Fishes. Curr. Biol. 2017, 27, 1184–1191. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  46. Ng, H.C.; Ranganathan, S.; Chua, K.L.; Khoo, H.E. Cloning and molecular characterization of the first aquatic hyaluronidase, SFHYA1, from the venom of stonefish (Synanceja horrida). Gene 2005, 346, 71–81. [Google Scholar] [CrossRef] [PubMed]
  47. Cherr, G.N.; Yudin, A.I.; Overstreet, J.W. The dual functions of GPI-anchored PH-20: Hyaluronidase and intracellular signaling. J. Int. Soc. Matrix Biol. 2001, 20, 515–525. [Google Scholar] [CrossRef]
  48. Smith, W.L.; Stern, J.H.; Girard, M.G.; Davis, M.P. Evolution of Venomous Cartilaginous and Ray-Finned Fishes. Integr. Comp. Biol. 2016, 56, 950–961. [Google Scholar] [CrossRef] [PubMed]
  49. Reitinger, S.; Laschober, G.T.; Fehrer, C.; Greiderer, B.; Lepperdinger, G. Mouse testicular hyaluronidase-like proteins SPAM1 and HYAL5 but not HYALP1 degrade hyaluronan. Biochem. J. 2007, 401, 79–85. [Google Scholar] [CrossRef] [PubMed]
  50. Pond, S.L.K.; Frost, S.D.W.; Muse, S.V. HyPhy: Hypothesis testing using phylogenies. Bioinformatics 2005, 21, 676–679. [Google Scholar] [CrossRef] [PubMed]
  51. Kosakovsky Pond, S.L.; Frost, S.D.W. Not So Different After All: A Comparison of Methods for Detecting Amino Acid Sites Under Selection. Mol. Biol. Evol. 2005, 22, 1208–1222. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  52. Murrell, B.; Moola, S.; Mabona, A.; Weighill, T.; Sheward, D.; Kosakovsky Pond, S.L.; Scheffler, K. FUBAR: A Fast, Unconstrained Bayesian AppRoximation for Inferring Selection. Mol. Biol. Evol. 2013, 30, 1196–1205. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  53. Murrell, B.; Wertheim, J.O.; Moola, S.; Weighill, T.; Scheffler, K.; Kosakovsky Pond, S.L. Detecting Individual Sites Subject to Episodic Diversifying Selection. PLoS Genet. 2012, 8, e1002764. [Google Scholar] [CrossRef] [PubMed]
  54. Ridout, K.E.; Dixon, C.J.; Filatov, D.A. Positive Selection Differs between Protein Secondary Structure Elements in Drosophila. Genome Biol. Evol. 2010, 2, 166–179. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  55. Kelley, L.A.; Mezulis, S.; Yates, C.M.; Wass, M.N.; Sternberg, M.J.E. The Phyre2 web portal for protein modeling, prediction and analysis. Nat. Protoc. 2015, 10, 845–858. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  56. Humphrey, W.; Dalke, A.; Schulten, K. VMD—Visual Molecular Dynamics. J. Mol. Gr. 1996, 14, 33–38. [Google Scholar] [CrossRef]
  57. Vonk, F.J.; Casewell, N.R.; Henkel, C.V.; Heimberg, A.M.; Jansen, H.J.; McCleary, R.J.R.; Kerkkamp, H.M.E.; Vos, R.A.; Guerreiro, I.; Calvete, J.J.; et al. The king cobra genome reveals dynamic gene evolution and adaptation in the snake venom system. Proc. Natl. Acad. Sci. USA 2013, 110, 20651–20656. [Google Scholar] [CrossRef] [PubMed]
  58. Kemparaju, K.; Girish, K.S. Snake venom hyaluronidase: A therapeutic target. Cell Biochem. Funct. 2006, 24, 7–12. [Google Scholar] [CrossRef] [PubMed]
  59. Da Silva, N.J.; Clementino Ferreira, K.R.; Leite Pinto, R.N.; Aird, S.D. A Severe Accident Caused by an Ocellate River Stingray (Potamotrygon motoro) in Central Brazil: How Well Do We Really Understand Stingray Venom Chemistry, Envenomation, and Therapeutics? Toxins 2015, 7, 2272–2288. [Google Scholar] [CrossRef] [PubMed]
  60. Monteiro-dos-Santos, J.; Conceicao, K.; Seibert, C.S.; Marques, E.E.; Silva, P.I., Jr.; Soares, A.B.; Lima, C.; Lopes-Ferreira, M. Studies on pharmacological properties of mucus and sting venom of Potamotrygon cf. henlei. Int. Immunopharmacol. 2011, 11, 1368–1377. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  61. Li, R.; Li, Y.; Kristiansen, K.; Wang, J. SOAP: Short oligonucleotide alignment program. Bioinformatics 2008, 24, 713–714. [Google Scholar] [CrossRef] [PubMed]
  62. Pertea, G.; Huang, X.; Liang, F.; Antonescu, V.; Sultana, R.; Karamycheva, S.; Lee, Y.; White, J.; Cheung, F.; Parvizi, B.; et al. TIGR Gene Indices clustering tools (TGICL): A software system for fast clustering of large EST datasets. Bioinformatics 2003, 19, 651–652. [Google Scholar] [CrossRef] [PubMed]
  63. Grabherr, M.G.; Haas, B.J.; Yassour, M.; Levin, J.Z.; Thompson, D.A.; Amit, I.; Adiconis, X.; Fan, L.; Raychowdhury, R.; Zeng, Q.; et al. Trinity: Reconstructing a full-length transcriptome without a genome from RNA-Seq data. Nat. Biotechnol. 2011, 29, 644–652. [Google Scholar] [CrossRef] [PubMed]
  64. Punta, M.; Coggill, P.C.; Eberhardt, R.Y.; Mistry, J.; Tate, J.; Boursnell, C.; Pang, N.; Forslund, K.; Ceric, G.; Clements, J.; et al. The Pfam protein families database. Nucleic Acids Res. 2012, 40, D290–D301. [Google Scholar] [CrossRef] [PubMed]
  65. Haas, B.J.; Papanicolaou, A.; Yassour, M.; Grabherr, M.; Blood, P.D.; Bowden, J.; Couger, M.B.; Eccles, D.; Li, B.; Lieber, M.; et al. De novo transcript sequence reconstruction from RNA-seq using the Trinity platform for reference generation and analysis. Nat. Protoc. 2013, 8, 1494–1512. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  66. Ye, J.; Fang, L.; Zheng, H.; Zhang, Y.; Chen, J.; Zhang, Z.; Wang, J.; Li, S.; Li, R.; Bolund, L.; et al. WEGO: A web tool for plotting GO annotations. Nucleic Acids Res. 2006, 34, W293–W297. [Google Scholar] [CrossRef] [PubMed]
  67. O’Leary, N.A.; Wright, M.W.; Brister, J.R.; Ciufo, S.; Haddad, D.; McVeigh, R.; Rajput, B.; Robbertse, B.; Smith-White, B.; Ako-Adjei, D.; et al. Reference sequence (RefSeq) database at NCBI: Current status, taxonomic expansion, and functional annotation. Nucleic Acids Res. 2016, 44, D733–D745. [Google Scholar] [CrossRef] [PubMed]
  68. Sievers, F.; Wilm, A.; Dineen, D.; Gibson, T.J.; Karplus, K.; Li, W.; Lopez, R.; McWilliam, H.; Remmert, M.; Söding, J.; et al. Fast, scalable generation of high-quality protein multiple sequence alignments using Clustal Omega. Mol. Syst. Biol. 2011, 7. [Google Scholar] [CrossRef] [PubMed]
  69. Gouy, M.; Guindon, S.; Gascuel, O. SeaView Version 4: A Multiplatform Graphical User Interface for Sequence Alignment and Phylogenetic Tree Building. Mol. Biol. Evol. 2010, 27, 221–224. [Google Scholar] [CrossRef] [PubMed]
  70. Da Fonseca, R.R.; Antunes, A.; Melo, A.; Ramos, M.J. Structural divergence and adaptive evolution in mammalian cytochromes P450 2C. Gene 2007, 387, 58–66. [Google Scholar] [CrossRef] [PubMed]
  71. Abascal, F.; Zardoya, R.; Posada, D. ProtTest: Selection of best-fit models of protein evolution. Bioinformatics 2005, 21, 2104–2105. [Google Scholar] [CrossRef] [PubMed]
  72. Da Fonseca, R.R.; Johnson, W.E.; O’Brien, S.J.; Vasconcelos, V.; Antunes, A. Molecular evolution and the role of oxidative stress in the expansion and functional diversification of cytosolic glutathione transferases. BMC Evol. Biol. 2010, 10, 281. [Google Scholar] [CrossRef] [PubMed]
  73. Nguyen, L.-T.; Schmidt, H.A.; von Haeseler, A.; Minh, B.Q. IQ-TREE: A fast and effective stochastic algorithm for estimating maximum likelihood phylogenies. Mol. Biol. Evol. 2014. [Google Scholar] [CrossRef] [PubMed]
  74. Anisimova, M.; Nielsen, R.; Yang, Z. Effect of recombination on the accuracy of the likelihood method for detecting positive selection at amino acid sites. Genetics 2003, 164, 1229–1236. [Google Scholar] [PubMed]
  75. Kosakovsky Pond, S.L.; Posada, D.; Gravenor, M.B.; Woelk, C.H.; Frost, S.D.W. Automated Phylogenetic Detection of Recombination Using a Genetic Algorithm. Mol. Biol. Evol. 2006, 23, 1891–1901. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  76. Pineda, S.S.; Sollod, B.L.; Wilson, D.; Darling, A.; Sunagar, K.; Undheim, E.A.B.; Kely, L.; Antunes, A.; Fry, B.G.; King, G.F. Diversification of a single ancestral gene into a successful toxin superfamily in highly venomous Australian funnel-web spiders. BMC Genom. 2014, 15, 177. [Google Scholar] [CrossRef] [PubMed]
  77. Low, D.H.; Sunagar, K.; Undheim, E.A.; Ali, S.A.; Alagon, A.C.; Ruder, T.; Jackson, T.N.; Pineda Gonzalez, S.; King, G.F.; Jones, A.; et al. Dracula’s children: Molecular evolution of vampire bat venom. J. Proteom. 2013, 89, 95–111. [Google Scholar] [CrossRef] [PubMed]
  78. Fry, B.G.; Undheim, E.A.B.; Ali, S.A.; Jackson, T.N.W.; Debono, J.; Scheib, H.; Ruder, T.; Morgenstern, D.; Cadwallader, L.; Whitehead, D.; et al. Squeezers and leaf-cutters: Differential diversification and degeneration of the venom system in toxicoferan reptiles. Mol. Cell. Proteom. 2013, 12, 1881–1899. [Google Scholar] [CrossRef] [PubMed]
  79. Sunagar, K.; Fry, B.G.; Jackson, T.N.W.; Casewell, N.R.; Undheim, E.A.B.; Vidal, N.; Ali, S.A.; King, G.F.; Vasudevan, K.; Vasconcelos, V.; et al. Molecular Evolution of Vertebrate Neurotrophins: Co-Option of the Highly Conserved Nerve Growth Factor Gene into the Advanced Snake Venom Arsenalf. PLoS ONE 2013, 8, e81827. [Google Scholar] [CrossRef]
  80. Pereira, S.R.; Vasconcelos, V.M.; Antunes, A. Computational study of the covalent bonding of microcystins to cysteine residues—A reaction involved in the inhibition of the PPP family of protein phosphatases. FEBS J. 2013, 280, 674–680. [Google Scholar] [CrossRef] [PubMed]
  81. Pereira, S.R.; Vasconcelos, V.M.; Antunes, A. The phosphoprotein phosphatase family of Ser/Thr phosphatases as principal targets of naturally occurring toxins. Crit. Rev. Toxicol. 2011, 41, 83–110. [Google Scholar] [CrossRef] [PubMed]
  82. Furnham, N.; Holliday, G.L.; de Beer, T.A.; Jacobsen, J.O.; Pearson, W.R.; Thornton, J.M. The Catalytic Site Atlas 2.0: Cataloging catalytic sites and residues identified in enzymes. Nucleic Acids Res. 2014, 42, D485–D489. [Google Scholar] [CrossRef] [PubMed]
Toxins 10 00544 g001 550
Figure 1. (A) Gene Ontology distribution of the coding region containing contigs retrieved from the freshwater stingray Potamotrygon motoro venom gland transcriptome. (B) Kyoto Encyclopedia of Genes and Genomes (KEGG) Classification of the coding region containing contigs of the freshwater stingray Potamotrygon motoro venom gland transcriptome.
Figure 1. (A) Gene Ontology distribution of the coding region containing contigs retrieved from the freshwater stingray Potamotrygon motoro venom gland transcriptome. (B) Kyoto Encyclopedia of Genes and Genomes (KEGG) Classification of the coding region containing contigs of the freshwater stingray Potamotrygon motoro venom gland transcriptome.
Toxins 10 00544 g001
Toxins 10 00544 g002 550
Figure 2. Graphical representation of the comparison between the most represented (A) gene ontology categories and (B) KEGG pathways found in Potamotrygon motoro, Potamotrygon amandae and Potamotrygon falkneri (data from the last two species reproduced from [26] 2009, Peptides). The percentage values correspond to the annotation values in each work for GO term and KEGG pathway annotation. The images of P. motoro (by Steven G. Johnson) and P. falkneri (by Andrew Kuchling) are used under the Creative Commons Attribution-Share Alike 3.0 Unported license. The image of P. amandae was reproduced from [1] 2013, Neotropical Ichtyhology.
Figure 2. Graphical representation of the comparison between the most represented (A) gene ontology categories and (B) KEGG pathways found in Potamotrygon motoro, Potamotrygon amandae and Potamotrygon falkneri (data from the last two species reproduced from [26] 2009, Peptides). The percentage values correspond to the annotation values in each work for GO term and KEGG pathway annotation. The images of P. motoro (by Steven G. Johnson) and P. falkneri (by Andrew Kuchling) are used under the Creative Commons Attribution-Share Alike 3.0 Unported license. The image of P. amandae was reproduced from [1] 2013, Neotropical Ichtyhology.
Toxins 10 00544 g002
Toxins 10 00544 g003 550
Figure 3. Maximum likelihood tree reconstruction of the hyaluronidase protein family (80 amino acid sequences of 529 residues) with 1000 bootstrap replications. Only branch supports higher than 80 are shown. Identified enzyme functional groups are labelled after each species. Highly expressed hyaluronidase sequences are indicated by green colored names. Known venomous fish hyaluronidase sequences are indicated by red colored names. P. motoro sequence with similarity to HYAL2 is highlighted with a blue colored name. HYAL1 sequence from P. motoro is indicated in bold. In the clade containing HYAL6, given the clade resolution and the phylogenetic proximity of Danio rerio H6 and Mus musculus H6, the HYAL4 designation for the Rhincodon typus and the Callorhinchus milii are considered mislabels.
Figure 3. Maximum likelihood tree reconstruction of the hyaluronidase protein family (80 amino acid sequences of 529 residues) with 1000 bootstrap replications. Only branch supports higher than 80 are shown. Identified enzyme functional groups are labelled after each species. Highly expressed hyaluronidase sequences are indicated by green colored names. Known venomous fish hyaluronidase sequences are indicated by red colored names. P. motoro sequence with similarity to HYAL2 is highlighted with a blue colored name. HYAL1 sequence from P. motoro is indicated in bold. In the clade containing HYAL6, given the clade resolution and the phylogenetic proximity of Danio rerio H6 and Mus musculus H6, the HYAL4 designation for the Rhincodon typus and the Callorhinchus milii are considered mislabels.
Toxins 10 00544 g003
Toxins 10 00544 g004 550
Figure 4. The three-dimensional protein model of the Potamotrygon motoro hyaluronidase inferred by Phyre2 [55]. Positively selected amino acids (P < 0.01) are colored in red and the catalytic residues are represented in white (green label). Images were generated using the Visual Molecular Dynamics software [56]. (A) Side-view of the catalytic cleft and protein overview. (B) Catalytic cleft close look.
Figure 4. The three-dimensional protein model of the Potamotrygon motoro hyaluronidase inferred by Phyre2 [55]. Positively selected amino acids (P < 0.01) are colored in red and the catalytic residues are represented in white (green label). Images were generated using the Visual Molecular Dynamics software [56]. (A) Side-view of the catalytic cleft and protein overview. (B) Catalytic cleft close look.
Toxins 10 00544 g004
Table
Table 1. Statistics of the freshwater stingray Potamotrygon motoro venom gland transcriptome sequencing and assembly (retrieved from the Trinotate annotation pipeline).
Table 1. Statistics of the freshwater stingray Potamotrygon motoro venom gland transcriptome sequencing and assembly (retrieved from the Trinotate annotation pipeline).
ProcessValue
Number of raw reads38,674,474
Raw data (bp)3,480,702,660
Read length (bp)90
Number of high-quality reads35,977,224
Average high-quality read length (bp)80
Number of contigs140,078
Number of contigs  ≥  1 FPKM107,129
Number of contigs  ≥  1 FPKM and containing coding sequences27,032
Contigs (bp)69,861,618
N50690
Average contig length (bp)498
Min. contig length (bp)174
Max. contig length (bp)15,040
Table
Table 2. Statistics of the coding region containing contigs annotated in the venom gland transcriptome of the freshwater stingray Potamotrygon motoro.
Table 2. Statistics of the coding region containing contigs annotated in the venom gland transcriptome of the freshwater stingray Potamotrygon motoro.
DatabaseNumber of Annotated Contigs
NCBI20,967 (84.77%)
SwissProt19,851 (80.26%)
NCBI and SwissProt19,398 (78.43%)
NCBI or SwissProt21,420 (86.60%)
ToxProt418 (1.69%)
PFAM16,823 (68.02%)
GO19,194 (77.60%)
KEGG17,572 (71.04%)
Table
Table 3. Top 25 most expressed transcripts from the Potamotrygon motoro venom gland transcriptome with hits to ToxProt. The corresponding Top 25 transcripts from P. amandae and P. falkneri are also shown (retrieved from [26]). Rows marked in red indicate differences in the most expressed molecules for each species.
Table 3. Top 25 most expressed transcripts from the Potamotrygon motoro venom gland transcriptome with hits to ToxProt. The corresponding Top 25 transcripts from P. amandae and P. falkneri are also shown (retrieved from [26]). Rows marked in red indicate differences in the most expressed molecules for each species.
Potamotrygon motoroPotamotrygon amandaePotamotrygon falkneri
TranscriptProteinUniprot AccessionFPKMProteinUniprot AccessionFPKMProteinUniprot AccessionFPKM
TR1637SE-cephalotoxinB2DCR842139.08HyaluronidaseJ3S82022201.33HyaluronidaseJ3S8209488.18
TR45956HyaluronidaseJ3S8207271.53Translationally controlled tumor protein homologU3EQ601225.86Hemolytic toxin Avt-1Q5R231757.19
TR7202SE-cephalotoxinB2DCR8489.97Cysteine-rich venom protein latiseminQ8JI38891.77Translationally controlled tumor protein homologU3EQ60734.38
TR53378DELTA-alicitoxin-Pse1bP0DL56399.42Venom allergen 5P81656854.22CalglandulinQ3SB11644.06
TR15069Cystatin-2J3SE80373.92CalglandulinQ8AY75835.98Putative Kunitz-type serine protease inhibitorB2BS84450.57
TR10238Augerpeptide hhe53P0CI21324.45Cystatin-2J3SE80629.15Peroxiredoxin-4P0CV91397.49
TR3474Translationally controlled tumor protein homologU3EQ60315.67Hemolytic toxin Avt-1Q5R231401.11Cysteine-rich venom protein 1Q8T0W5386.31
TR112682CalglandulinQ3SB11202.31Putative Kunitz-type serine protease inhibitorB2BS84354.6Venom allergen 5 P81656357.2
TR24275Putative Kunitz-type serine protease inhibitorB2BS84198.6Peroxiredoxin-4P0CV91304.79CalglandulinQ3SB11318.87
TR20805Venom allergen 5P81656174.48Cysteine-rich venom protein 1Q8T0W5211.07Cysteine-rich venom protein latiseminQ8JI38299.5
TR20804Venom allergen 5.02P35782165.85Alpha-latroinsectotoxin-Lt1aQ02989206.6VesprynQ2XXL4172.92
TR14068PI-actitoxin-Aeq3aP0DMW6158.3Kunitz-type serine protease inhibitor bitisilin-3Q6T269141.02Venom serine protease 34 Q8MQS8156.83
TR15070Cystatin-2J3SE80146.47Analgesic polypeptide HC3C0HJF3121.14CalglandulinQ8AY75150.86
TR86455Insulin-like growth factor-binding protein-related protein 1G4V4G1142.04Kunitz-type serine protease inhibitor kappa-theraphotoxin-Hh1aP68425117.78Alpha-latrocrustotoxin-Lt1aQ9XZC0150.43
TR110388CalglandulinQ3SB11128.78Venom prothrombin activator porpharin-DQ58L93104.65Alpha-latrotoxin-Lh1aG0LXV8150.2
TR112679CalglandulinQ3SB11126.9Zinc metalloproteinase-disintegrin-like BmMPA8QL49101.89Zinc metalloproteinase-disintegrin-like BmMPA8QL49143.04
TR113284Cysteine-rich venom protein 1Q8T0W582.72Acidic phospholipase A2Q9DF56100.95Cystatin-2J3SE80140.7
TR67254Putative Kunitz-type serine protease inhibitorB2BS8473.96VesprynQ2XXL495.12Venom protease Q7M4I3118.09
TR9752Zinc metalloproteinase-disintegrin-like BmMPA8QL4969.87Insulin-like growth factor-binding protein-related protein 1 G4V4G193.55Alpha-latroinsectotoxin-Lt1aQ02989114.49
TR1967Alpha-latrocrustotoxin-Lt1aQ9XZC068.33Delta-latroinsectotoxin-Lt1aQ2533887.69Kunitz-type serine protease inhibitor HNTX-852P0DJ69108.68
TR119403Alpha-latrocrustotoxin-Lt1aQ9XZC068.3OhaninP8323486.84Venom prothrombin activator porpharin-DQ58L93107.52
TR110386CalglandulinQ3SB1167.47Venom prothrombin activator vestarin-D2A6MFK871.33Analgesic polypeptide HC3C0HJF3106.37
TR7292DELTA-thalatoxin-Avl1aQ5R23166.92Blarina toxinQ76B4570.84Snake venom metalloprotease inhibitorA8YPR986.64
TR53095Peroxiredoxin-4P0CV9164.91Venom proteaseQ7M4I359.88Kunitz-type serine protease inhibitor bitisilin-3Q6T26981.19
Table
Table 4. Relevant protein venom components found in bony fish (22 species) and the marine stingray Neotrygon kuhlii. Molecules found in the freshwater stingrays (this study and [26]) are indicated in red.
Table 4. Relevant protein venom components found in bony fish (22 species) and the marine stingray Neotrygon kuhlii. Molecules found in the freshwater stingrays (this study and [26]) are indicated in red.
Venom Components
Bony Fish [41]Marine Stingray [29]
Trachynilysin (TLY)60S acidic ribosomal protein
Stonustoxin (SNTX)ATP synthase
Verrucotoxin (VTX)Coronin
Neoverrucotoxin (neoVTX)Cystatin
CardioleputinCytochrome C
NocitoxinFerritin
KaratoxinGalectin
Sp-CTxGanglioside GM2 activator
PlumieribetinGlutathione S-transferase mu
SP-CL 1-5Hemoglobin subunit alpha
DracotoxinLeukocyte elastase inhibitor
TrachinineNucleoside diphosphate kinase
SA-HTPeroxiredoxin 6
TmC4-47.2Transaldolase
NattectinType III intermediate filament
Toxin-PCVoltage-dependent anion channel
Wap65-
Natterin-
Hyaluronidase-
Phospholipase A2 [45]-
Proenkephalin [45]-
Neuropeptide Y [45]-

Share and Cite

MDPI and ACS Style

Silva, F.; Huang, Y.; Yang, V.; Mu, X.; Shi, Q.; Antunes, A. Transcriptomic Characterization of the South American Freshwater Stingray Potamotrygon motoro Venom Apparatus. Toxins 2018, 10, 544. https://doi.org/10.3390/toxins10120544

AMA Style

Silva F, Huang Y, Yang V, Mu X, Shi Q, Antunes A. Transcriptomic Characterization of the South American Freshwater Stingray Potamotrygon motoro Venom Apparatus. Toxins. 2018; 10(12):544. https://doi.org/10.3390/toxins10120544

Chicago/Turabian Style

Silva, Filipe, Yu Huang, Vítor Yang, Xidong Mu, Qiong Shi, and Agostinho Antunes. 2018. "Transcriptomic Characterization of the South American Freshwater Stingray Potamotrygon motoro Venom Apparatus" Toxins 10, no. 12: 544. https://doi.org/10.3390/toxins10120544

APA Style

Silva, F., Huang, Y., Yang, V., Mu, X., Shi, Q., & Antunes, A. (2018). Transcriptomic Characterization of the South American Freshwater Stingray Potamotrygon motoro Venom Apparatus. Toxins, 10(12), 544. https://doi.org/10.3390/toxins10120544

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop