CROSSTALK OF CONSERVATION OF SEQUENCES OF miRNAs AND ENZYMATIC MACHINERY OF miRNAs PRODUCTION IN AMPULLARIIDAE

REGISTRO DOI: 10.69849/revistaft/ar10202408192041


Raphael Rodrigues Porto1*, Flávio César Thiago1*, Carlos Bruno de Araújo1, Tamires Caixeta Alves1, Laurence Rodrigues do Amaral1, Matheus de Souza Gomes1**.


ABSTRACT:  

BACKGROUND: The Ampullariidae family of molluscs is an emerging model for evolutionary studies due to its high diversity, ancient history and wide geographic distribution. miRNAs are essential for the development of the organism, as they help and act in the control of gene expression. The analysis of miRNAs in molluscs, especially in the family Ampullariidae, is especially important due to the low representativeness of analyzed species and the possibility of analyzing the conservation of miRNAs among ampullariids species. OBJECTIVE: Identification and characterization of miRNAs (precursors and matures) and their processing pathway genes. METHODS: Computational prediction and characterization of miRNAs and genes involved in the miRNA pathway were performed using public database of Lanistes nyassanus, Marisa cornuarietis, Pomacea canaliculata and Pomacea maculata. The in silico analysis was performed using a robust algorithm to identify and characterize miRNAs and their precursors in genome of ampullariids species. To search for the putative proteins involved in the miRNA biogenesis the putative proteomes from 4 Ampullariids species and Blastp tool were used. Characterization of conserved protein domains was performed using the PFAM and CDD. Phylogenetic analyzes were performed for the ampullariids miRNA precursors and their orthologs and also for the putative ampullariids proteins involved in the miRNA pathway and their orthologs using MEGA program. FINDINGS: 141 pre-miRNAs and 162 mature miRNAs were identified in the genome of L. nyassanus, 279 pre-miRNAs and 297 mature miRNAs in the genome of M. cornuarietis, 269 pre-miRNAs and 296 mature miRNAs in the genome of P. canaliculata and 299 pre-miRNAs and 316 mature miRNAs in the genome of P. maculata. We identified and characterized 24 putative key proteins involved in the miRNA pathway including Argonaute, DICER, DROSHA and EXPORTIN protein families in the predicted proteome of the 4 ampullariids. The data obtained in this work will support studies of phylogeny, population divergence, speciation and patterns of diversity in the Ampullariidae Family. MAIN CONCLUSIONS: The searching for novel miRNAs and their processing pathway genes in 4 species of ampullariids was able to predict new structures expanding the study of miRNAs in molluscs and in Ampullariidae family, as well as open an avenue to study the roles of miRNAs in the organisms. 

Keyword: miRNAs. Computational Analysis. Snails. Ampullariidae. Genome 

INTRODUCTION 

Molluscs of the Ampullariidae family are freshwater snails found in tropical and subtropical regions of Africa, Asia and the Americas. The animals belong to the subclass Caenogastropoda, which comprises about 50,000 species and approximately 60% of Gastropod species (Hayes et al., 2009; Strong et al., 2008). As they are a family that derived from the subclass Caenogastropoda in early evolution (COWIE, 2015) ampullariids are an excellent model to study speciation and adaptation (Hayes et al., 2009).  

The scarcity of sequenced molluscs genomes was identified by the Global Invertebrate Genomics Alliance as a problem for species characterization (Voolstra et al., 2017). There are four ampullariids with sequenced genomes including the Old World (Lanistes nyassanus) and New World (Pomacea canaliculata, Pomacea maculata and Marisa cornuarietis) clades that can be the object of comparative genomic studies. The division between New and Old World lineages is consistent with the time of separation from Gondwana, around 120 million years ago (Jokat et al., 2003; Sun et al., 2019). The availability of the genome and transcriptome of the cited species makes it promising to study gene regulation, intrinsically the posttranscriptional regulation mediated by miRNAs. Almost all biological processes, including cell differentiation, tissue development and cell metabolism can be modulated by miRNAs. This regulation of gene expression is carried out through a fine control, where miRNAs and their silencing pathways carry out a specific regulation of gene expression and maintenance of genome integrity (Bartel, 2009). This may help to better understand the biology of the ampullariids studied and their adaptive properties. miRNAs have a low homoplasy, that is, their presence in two different species is more likely to be linked to a common ancestry, thus correlating the distribution of miRNAs in evolutionary relationships of the species (Kenny et al., 2015). This will allow a better understanding of the biology of the Ampullariidae family in molluscs, particularly its medical and economic importance for the humanity of ampullariids species. 

P. canaliculata is involved in the transmission of the helminth Angiostrongylus cantonensis that causes eosinophilic meningitis in humans. Molluscs are the intermediate hosts when they are contaminated by ingestion of Larva L1. The infecting larva lodges in the Central Nervous System in humans, causing eosinophilic meningitis or ocular angiostrongyliasis in the eye ball (Luessi et al., 2009; Wang et al., 2012, 2008). Despite eosinophilic meningitis usually course in humans with low severity there may be some severe cases of the disease with severe neurological sequelae, coma and death (especially in children) (Chen et al., 2014; Lindo et al., 2002; Tsai et al., 2001). miRNAs are able to modulate the immune system of molluscs in response to stress caused by infectious processes (Picone et al., 2017).  

P. canaliculata is used as a source of dietary protein in Southeast Asia. Ingestion of undercooked animals can transmit A. cantonensis causing eosinophilic meningitis in humans. In the Nueva Ecija region of the Philippines, the prevalence of P. canaliculata contaminated with nematodes was 17.75%  (Cawas et al., 2020). The lack of awareness of the population is notorious, being considered a neglected disease and an important public health problem in the region (Dalton et al., 2017). 

In mainland China, eosinophilic meningitis caused by A. cantonensis is considered an emerging infectious disease (Lv et al., 2008). Epidemiological evidence indicates that P. canaliculata is becoming the most prevalent intermediate host of A. cantonensis in this country, due to its high susceptibility to the parasite and its wide environmental tolerance (Yang et al., 2013). A major outbreak occurred in Beijing in 2006 with 160 cases of human angiostrongyliasis causing eosinophilic meningitis. The fact occurred due to the consumption of contaminated P. canaliculata imported from regions in southern China (Wang et al., 2008). Thus, we note the extreme importance of more epidemiological studies in the country regarding this serious public health problem. The study of the biology and adaptability of P. canaliculata has the same degree of importance for better management and control of eosinophilic meningitis. 

P. canaliculata is described as an invasive species and was ranked among the top 100 pests worldwide. The species is resistant to elevated temperatures and drought (Matsukura et al., 2009; Yusa et al., 2006) . P. maculata has physiological and adaptive characteristics similar to P. canaliculata, both species being recognized for their rapid propagation and population increase, constituting a serious threat to the balance of an ecosystem. The economic costs of damage, control and repair can be immense (Pimentel et al., 2005). 

P. canaliculata is considered an agricultural pest in many parts of Asia (Wang et al., 2020). In Indonesia, when introduced into rice plantations, it causes a 15% drop in production (Novarino, 2011). The dispersion of P. canaliculata in crops planted in humid regions causes serious damage in Southeast Asia  (Carlsson et al., 2016) in Japan (Yusa, 1999) and also in China (Yang et al., 2013). 

In 2003, China’s State Department of the Environment ranked P. canaliculata in a list of 16 invasive pests. Since the mollusc was introduced to China around 1980, there is no accurate data on the spread of P. canaliculata in agriculture. Over the next twenty years, the molluscs spread widely in the country and as the snails moved into agriculture, control measures against P. canaliculata were implemented (Yang et al., 2013). More epidemiological studies are needed regarding the impact of P. canaliculata on agriculture in China. Work on the biology of the snail such as adaptability and environmental plasticity is also needed. The study of miRNAs in the species is an important step towards understanding part of the physiology of P. canaliculata, for possible mechanisms of population control of the species in agriculture. 

The inherent skills of P. canaliculata and P. maculata referring to the invasive aspect could be confirmed by mapping characteristics by gene locus. The ability of P. canaliculata and P. maculata species to colonize and occupy new territories is attributed in part to genes that encode the G protein family (GRL101) (Sun et al., 2019), highly expressed in tentacles and cephalic palps (chemoreception in aquatic snails) (Schultz and Adema, 2017) . Still referring to invasiveness, studies have suggested that the voracious appetite of P. canaliculata and P. maculata species are attributed to the expansion of gene families related to cellulose digestion (Sun et al., 2019). 

L. nyassanus and M. cornuarietis are snails related to freshwater (Strong et al., 2008) and have divergences regarding their adaptability. While M. cornuarietis is better able to resist factors such as low temperature and hypoxia (Matsukura et al., 2016; Mu et al., 2018) L. nyassanus has a low adaptability and is endemic to Lake Nyasa, or Malawi, geographically located between Malawi, Tanzania and Mozambique (COWIE, 2015; Dohrn, 1865; Van Bocxlaer, 2017). 

Since miRNA-mediated post-transcriptional silencing performs several regulatory aspects in animals, it is necessary to understand the processing mechanisms of the miRNA pathway and their respective genes. The identification of major miRNA biogenesis factors including, AGO2, DGCR8, DICER, DROSHA, RAN and XPO5 was elucidated based on homologous sequence search in several molluscs (Huang et al., 2021) and in Biomphalaria Glabrata by researchers from our group (Queiroz et al., 2017). 

Sun et al., published in 2019 the sequencing of the total genome of four species of Ampullariidae, L. nyassanus, P. canaliculata, P. maculata and M. cornuarietis. From then on, it was putative to carry out the mapping of characteristics by gene locus. Given the relevant regulatory aspect of non-protein coding RNAs in gene expression, specifically miRNAs, it is necessary to study mature and precursor miRNAs, as well as the genes involved in the miRNA processing pathway. In this work, mature and precursor miRNAs and the genes involved in miRNA pathway were identified and characterized using the genome of L. nyassanus, P. canaliculata, P. maculata and M. cornuarietis

MATERIALS AND METHODS 

PREDICTION AND CHARACTERIZATION OF MATURE AND PRECURSOR MIRNAS 

The genomes of L. nyassanus, M. cornuarietis, P. canaliculata and P. maculata were retrieved from the NCBI database (https://www.ncbi.nlm.nih.gov/genome). A robust algorithm was used to predict miRNAs, precursor and mature forms. The sequences were processed for identification and characterization of predicted mature miRNAs and their precursors, following the methodology established by Gomes et al. (2011). 

Genome sequences with potential tendency to formation of hairpin structures or similarity to precursor miRNAs were obtained using the Einverted program (EMBOSS tool) and the BLASTn program (NCBI), for identification of homologous miRNAs from miRBase (http: //www.mirbase.org/). Then, the sequences were submitted through a series of filters, selecting, according to the desired parameters, those that corresponded to putative-putative precursors of miRNAs. These filters were based on conserved characteristics of precursor miRNAs and unwanted sequences were discarded. The filters used were the Minimum Free Energy (MFE) (Gruber et al., 2008; Hofacker, 2009), GC content (Zhou et al., 2009), homology with mature miRNAs already identified. Finally, sequences with high similarity to genes encoding proteins, repetitive elements and non-coding RNAs, except miRNAs, were discarded. 

The putative precursors of miRNAs were compared to their orthologs at various levels. The sequences of pre-miRNAs and their respective miRNAs were submitted to multiple alignment, using the tools ClustalX 2.1(Larkin et al., 2007) and RNAalifold (Bernhart et al., 2008). Alignment patterns for ClustalX 2.1 were based on adjusted parameters (gap opening: 22.50; gap extension: 0.83). Secondary structures were obtained using the RNAfold platform (Gruber et al., 2008; Hofacker, 2009). The phylogenetic analysis of the putative pre-miRNAs of L. nyassanus, M. cornuarietis, P. canaliculata and P. maculata was carried out using MEGA program version X (Tamura et al., 2007) and the Neighbor-joining method applying the Kimura two-parameters model (“Kimura-twoparameters”), to stimulate the divergence between the sequences (N Saitou and Nei, 1987). Consensus trees were obtained using a bootstrap of 5000 replicates. 

PREDICTION AND CHARACTERIZATION OF PROTEINS INVOLVED IN MIRNAS BIOGENESIS 

In silico identification and characterization of putative miRNA pathway proteins were performed using the amino acid sequences of the proteins from the model organisms Biomphalaria glabrata, Drosophila melanogaster and C. elegans and the Blastp tool on the ampubase website (https://www.comp.hkbu.edu.hk/~db/AmpuBase/index.php#&panel1-2). Using PFAM (31.0) (Mistry et al., 2021) and conserved Domains CDD (Marchler-Bauer et al., 2017, 2015, 2011; Marchler-Bauer and Bryant, 2004; Shennan Lu et al., 2020) the putative domain architectures of the ampullariids proteins were established. Analysis of the active sites were used for the proteins Argonaut, Dicer and Drosha using Clustal X 2.1 and Weblogo. 

Phylogenetic analysis was performed using the MEGA version X program (Kumar et al., 2018), after alignment in default parameters, using the Neighbor-Joining method (Naruya Saitou and Nei, 1987) and the JTT model. The statistical reliability of each branch in the generated tree was evaluated using a 2000-replica bootstrap. 

RESULTS

IDENTIFICATION AND CHARACTERIZATION OF PUTATIVE PROTEINS INVOLVED IN MIRNA PATHWAY IN AMPULLARIIDS SPECIES

The predicted proteome of the species L. nyassanus. M. cornuarietis, P. canaliculata and P. maculata were retrieved from the NCBI database (National Center for Biotechnology Information – http://.ncbi.nlm.nih.gov/) to identify and characterize the putative proteins responsible for the biogenesis of the miRNAs. The identification was based on conservation, using as parameters protein length, conserved domains, active sites conservation and phylogeny analysis. Lny23421, Lny25187,Lny28936, Lny28171, Lny18581, Lny14710, Lny12572, Lny27692, Lny22624, Lny12572, Lny29643, Lny27624, Lny9931, Mco26370, Mco27042, Mco25093, Mco22837, Mco24030, Mco81614, Mco51286, Mco13761, Mco21592, Mco12258, Pca61914, Pca59652, Pca69814, Pca67531, Pca67887,Pca62316, Pca65193, Pca60969, Pca65193, Pca65324, Pca59979, Pca69254, Pma35144, Pma60501, Pma63715, Pma65850, Pma51393, Pma59467, Pma63331, Pma32975, Pma62419, Pma51302 and Pma61707 putative proteins involved in miRNA pathway were identified in predicted proteome of L. nyassanus, M. cornuarietis, P. canaliculate and P. maculate, respectively [Supplementary data I – (Supplementary Tables 1,2,3 and 4)]. The key proteins of the miRNA pathway were selected, Argonaute, DICER, DROSHA and Exportin, and thus used for further analysis.

Ten putative Argonaute proteins were identified in the predicted proteome of the ampullariids species, being 3 sequences for L. nyassanus (Lny23421_c0_g1, Lny25187_c0_g1 and Lny28936_c7_g7), 3 sequences for M. cornuarietis (Mco25093_c0_g2, Mco26370_c0_g1, Mco27042_c1_g1), 2 sequences for P. canaliculata (Pca69814_c1_g2, Pca61914_c0_g1) and 2 sequences for P. maculata (Pma63715_c0_g1_, Pma60501_c2_g16). The sequences Lny23421_c0_g1, Mco26370_c0_g, Pca69814_c1_g2 and Pma63715_c0_g1 exhibited six conserved domains, while the protein sequences Lny25187_c0_g1, Lny28936_c7_g7, Mco25093_c0_g2, Mco27042_c1_g1, Pca61914_c0_g1 e Pma60501_c2_g16 displayed only the PAZ and PIWI domains (Figure 1).

Figure 1: Argonaute’s Conserved Domains, L. nyassanus. M. cornuarietis, P. canaliculata, P. maculata and orthologs.  

The PIWI domain of these proteins presented a catalytic triad formed by DDH (aspartic acid, aspartic acid and histidine). The global alignment of the PIWI domains of the ampullariids Argonaute proteins and their orthologs showed a conservation of amino acid sequences, being accentuated in the surroundings of the catalytic sites (Figure 2).

Figure 2: Global alignment analysis of the PIWI domain of Argonaute proteins de L. nyassanus, M. cornuarietis, P. canaliculata and P. maculata.

The two subfamilies, AGO and PIWI compose the AGO protein family. The predicted sequences of L. nyassanus, M. cornuarietis, P. canaliculata and P. maculata were used in the phylogenetic tree and their orthologs, in order to emphasize the conservation of these putative proteins identified in the Ampullariids species and wide distribution among deuterostomes and protostomes species. The distribution of Argonaute and PIWI proteins corroborated the distribution of species in the tree of life. As expected, the ampullariids species grouped together with the phylum Molluscs. The evolutionary proximity of species ampullariids, L. nyassanus, M. cornuarietis, P. canaliculata and P. maculata is represented in the tree (Figure 3).

Figure 3: Phylogenetic analysis of AGO proteins, L. nyassanus, M. cornuarietis, P. canaliculata e P. maculata and their orthologs.

Five putative Drosha and two putative Dicer proteins were identified and characterized in the ampullariids species. For Drosha, 1 sequence for L. nyassanus (Lny18581_c0_g2), 1 sequence for M. cornuarietis (Mco22837_c2_g1), 1 sequence for P. canaliculata (Pca67887_c3_g1) and 2 sequences for P. maculata (Pma51393_c0-g1 and Pma65850_c2g4) were identified. The amino acid sequences of Pca67887_c3_g1 and Mco22837_c2_g1 had three conserved domains, Lny18581_c0_g2 and Pma65850_c2_g4 had two conserved domains, while the protein sequence Pma51393_c0_g1 had only one conserved domain. For Dicer 1 amino acid sequence were found for P. canaliculata (Pca67531_c1_g2), showing 4 conserved domains and for L. nyassanus (Lny28171_c0_g1) showing only 3 conserved domains (Figure 4 and Figure 5).

Figure 4: Analysis of conserved domains Dicer for L. nyassanus, M. cornuarietis, P. canaliculata, P. maculata and their orthologs.

Figure 5: Analysis of conserved domains Drosha for L. nyassanus, M. cornuarietis, P. canaliculata, P. maculata and their orthologs.

In the analysis of the Drosha/Dicer proteins, only one conserved domain of Ribonuclease III was identified in the ampullariids studied. The Riboc I domains presented the catalytic amino acids in the conserved site regions (EDDE): glutamic acid (E), aspartate (D), aspartate (D) and glutamic acid (E). This was verified in the sequences Drosha, Lny18581_c0_g2, Mco22837_c2_g1, Pca67887_c3_g1 and Dicer, Lny28171_c0_g1. Note the conservation of amino acids in the active site itself and in the surroundings of the region (Figure 6).

Figure 6: Analysis of the conserved domains of Ribonuclease III (Riboc) of the Drosha and Dicer proteins from ampullariids species.

Phylogenetic analysis of the predicted proteins DROSHA and DICER showed the evolutionary relationship with their orthologs. Proteins from L. nyassanus, M. cornuarietis, P. canaliculata, P. maculata grouped in the Molluscs clade, revealing close similarity among the ampullariids studied. The disposition of the orthologs of the putative DICER and DROSHA in the clades in the tree corroborated with the tree of life, with a visible distinction between Deuterostomes and Protostomes for the two classes of RNase (Figure 7).

Figure 7: Phylogenetic analysis of DICER and DROSHA proteins from ampullariids and their orthologs. 

Similarity of the predicted protein Drosha from L. nyassanus, M. cornuarietis and Pcanaliculata with model organism D. melanogaster and evolutionarily close organisms. The percent identity, query coverage, length of the predicted protein Drosha from P. canaliculata are similar among closely related organisms (Table 1).  

Code Query coverage E Value Identity % Length of string 
NP_477436.1 Drosha [D. melanogaster] 100% 0.0 100,00 1327 
XP_013069067.1 (B. glabrata73% 0.0 54,91 1469 
XP_005107264.1 (Aplysia californica72% 0.0 54,81 1528 
XP_041368618.1 (G. aegis)  72% 0,0 59,09 1412 
Pca67887_c3_g1 (P.canaliculata) 66% 0.0 57,75 1525 
Mco22837_c2_g1 (M. cornuarietis56% 0.0 60,96 760 
Lny18581_c0_g2 (L. yassanus)  27% 7e-150 69,57 317 
Table 1: Putative DROSHA of L. nyassanus, M. cornuarietis and P. canaliculata cover against NP_477436.1 DROSHA [D. melanogaster], E-Value, % of String Identity and Length. 

Comparison of the Dicer protein from the model organism D. melanogaster with the putative Dicer protein in L. nyassanus and P. canaliculata and evolutionarily close organisms. Once again, the predicted protein Dicer of P. canaliculata showed greater similarity with its orthologs when compared to the other Ampullariidae (Table 2). 

Código Cobertura da query E ValueIdentidade %Comprimento da sequência
NP_524453.1 Dicer-1 [D. melanogaster] 100% 0.0 100% 2249 
XP_012942349.1 (A. californica)  77% 2e-14150% 2503 
XP_013067888.1 (B. glabrata)  74% 2e-148 51,90% 2332 
XP_021347347.1 (M. yessoensis72% 7e-152 52,87% 2284 
XP_041353083.1 (G. aegis70% 2e-150 50,26% 2063 
Pca67531_c1_g2 (P. canaliculata71% 7e-140 45.10% 1985 
Lny28171_c0_g1 (L. nyassanus38% 2e-162 53,32% 1434 
Table 2: Putative Dicer L. nyassanus e P. canaliculata, coverage against NP_524453.1 Dicer1 [D. melanogaster], E-Value, % identity and string length.

Nine putative protein sequences belonging to the Exportin protein family (XPO-5, XPO1, XPO-T) were identified and characterized in the predicted proteome of ampullariid species; XPO-5: L. nyassanus (Lny14170_c0_g1), P. canaliculata (Pca62316_c0_g2); XPO-1: L. nyassanus (Lny12572_c0_g1), M. cornuarietis (Mco24030_c0_g2), P. canaliculata (Pca65193_c2_g1) and P. maculata (Pma59467_c0_g1), XPO-T: L. nyassanus (Lny27692_c0_g2), M. cornuarietis (Mco81614_c0_g2) P. maculata (Pma63331_c0_g1). The sequences of XPO-5, Lny14170_c0_g1 and Pca62316_c0_g2, exhibited 2 domains. For the XPO-1 protein sequences, the putative proteins Lny12572_c0_g1, Mco24030_c0_g2 and Pma59467_c0_g1 displayed six conserved domains, except for the Pca65193_c2_g1, which exhibited only five conserved domains. Likewise, XPO-T putative proteins Lny27692_c0_g2, Mco81614_c0_g2 and Pma63331_c0_g1 showed one conserved domain (Figure 8).

Figure 8: Distribution of conserved domains in putative proteins of the Exportin family in ampullariids species and orthologs.

In the phylogenetic analysis of putative proteins of the Exportin family (XPO-5, XPO-1 and XPO-T) the evolutionary relationship with their orthologs was verified. Protein sequences from L. nyassanus. M. cornuarietis, P. canaliculata and P. maculata were grouped in the phylum Molluscs, confirming once again a close proximity between the ampullariids species, as expected by the degree of relatedness of the species. The phylogenetic tree presented three distinct clades (XPO-5, XPO-1 and XPO-T clades), as well as species belonging to the Protostome and Deuterostome groups, corroborating with the literature (Figure 9).

Figure 9: Phylogenetic distribution of putative proteins, XPO-5 XPO-1 and XPO-T identified in L. nyassanus. M. cornuarietis, P. canaliculata e P. maculata and their orthologs. 

PREDICTION AND CHARACTERIZATION OF PERCURSOR AND MATURE MIRNAS.  

A robust algorithm was applied to the genome sequences of the four ampullariids species to identify mature miRNAs and their precursors. 141 pre-miRNAs and 162 miRNAs of L. nyassanus, 279 pre-miRNAs and 297 miRNAs of M. cornuarietis, 269 pre-miRNAs and 296 miRNAs of P. canaliculata; 299 pre-miRNAs and 316 miRNAs of P. maculata were predicted [Supplementary data II – (Supplementary Tables 5,6,7,8 and 9)]. 

An accurate and stringent structural and thermodynamic analysis were used to identify and characterize conserved miRNAs(Y. Q. Zhang et al., 2009). All pre-miRNAs identified were analyzed for these particular characteristics [Supplementary data II – (Supplementary Tables 5,6,7 and 8)].The pre-miRNAs for L. nyassanus exhibited MFE (minimum free energy) with a mean of -30.30 kcal/mol, with values between -18.50 and -50.90 kcal/mol; for M. cornuarietis the mean MFE of -28.05 kcal/mol, with values between -18.50 and -67.10 kcal/mol; for P. canaliculata the mean of MFE -28.70 kcal/mol, with values between -18.70 and -48.4 kcal/mol; for P. maculata the mean of -27.55 kcal/mol, with values between -18.50 and -52.7 kcal/mol (Table 3). 

Species Average size of precursors (nucleotides)  Precursor size variation (nucleotides) Mean value MFE – Kcal/mol Variation MFE – Kcal/mol 
L. nyassanus  87.04  66 to 100 -30,3 (-50.90 to -18.50) 
M. cornuarietis  88.16  65 to 100 -28,05 (-67.10 to -18.50) 
P. canaliculata 89,55 63 to 144 -28,7 (-48.40 to -18.70) 
P. maculata  89.27 65 to 100 -27,55 (-52.70 to 18.50) 
Table 3: Average size and size variation of precursors, average minimum free energy (MFE) and MFE variation of precursors in 4 species Ampullariidae 

The precursors showed mean values for the guanine-cytosine (GC) content for L. nyassanus, M. cornuarietis, P. canaliculata and P. maculata; 41.26%, 41.64%, 41.70%, 40.74%, respectively. The maximum and minimum values were respectively 69.05 and 20.25% for L. nyassanus, 67.03% and 20.20% for M. cornuarietis, 68.29% and 20.83% for P. canaliculata, 68, 75% and 20.40% for P. maculata. In this study, the average size of miRNA precursors for L. nyassanus, M. cornuarietis, P. canaliculata, P. maculata were respectively; 87.04 nucleotides, 88.16 nucleotides, 89.55 nucleotides and 89.27 nucleotides. Size range of pre-miRNAs for L. nyassanus (66 to 100 nucleotides), M. cornuarietis (65 to 100 nucleotides), P. canaliculata (63 to 144 nucleotides), P. maculata (65 to 100 nucleotides).   

The sequences of miRNAs have been shown to be evolutionarily conserved (De Wit et al., 2009). This conservation is accentuated in the seed region of the mature miRNAs, being considered important for recognition of their mRNA targets (Lewis et al., 2005). All mature miRNAs identified in this work showed 100% identity in the seed regions (position 2 to 8 nucleotides) with the structure of their respective ortholog deposited in miRBase, on which the homology comparison was based. 

The mature miRNAs identified in this work had an average of 22.22 nucleotides in L. nyassanus, 21.92 nucleotides in M. cornuarietis, 21.86 nucleotides in P. canaliculata and 21.77 nucleotides for P. maculata, ranging in size from 17 to 25 nucleotides for all species analyzed. 

The Argonaute protein has a preference for the uracil nucleotide in the first position of the miRNA sequence at the 5′ end to be loaded, to be inserted into the RISC complex (Seitz et al., 2011). In this work mature miRNAs from L. nyassanus, M. cornuarietis. P. canaliculata and P. maculata respectively presented uracil as the main nucleotides in the first position; 41.97%, 44.44%, 41.2% and 42.40% (uracil), 17.28%, 17.84%, 20.9% and 17.8% (cytosine), 26.54%, 22.56%, 20.4% and 24.5% (adenine) and 14.19%, 15.15%, 17.5% and 15.3% (guanine) [Supplementary data II – (Supplementary table 9)].  

The presence of miRNAs in two different species is more likely to be linked to a common ancestry, thus correlating the distribution of miRNAs in evolutionary relationships of species (Kenny et al., 2015). In order to emphasize the conservation of these putative miRNAs identified in the ampullariid species, phylogenetic trees were constructed. The global alignment and secondary structure characterization of the precursors of miRNAs was performed, to show the similarity of the sequences with their orthologs already identified. The taxa that were used as a parameter for the selection of specific are Bilateria, Protostomes, Lophotrocozoa and Molluscs. The selection of pre-miRNAs based on these taxons was performed having as parameter the miRNAs prevalent simultaneously in the 4 ampullariids species (Table 4). The characterization of miRNAs by taxons was shown in figures 1 to 11 and Supplementary data III – supplementary figures 1 to 34 

TaxSpeciesPre-miRNAs
BilateriaL. nyassanus lny-mir-33, lny-mir-190, lny-mir-252a, lny-mir-281, lny-mir-981  
M. cornuarietis mcr-mir-33, mcr-mir-190, mcr-mir-252a, mcr-mir-281, mcrmir-981  
P. canaliculata pcn-mir-33, pcn-mir-190, pcn -mir-252a, pcn-mir-281, pcn-mir981 
P. maculata pmc -mir-33, pmc-mir-190, pmc-mir-252a, pmc-mir-281, pmcmir-981 
ProtostomiesL. nyassanus lny-bantam, lny-mir-2a, lny-mir-67, lny-mir-750 
M. cornuarietis mcr-bantam, mcr-mir-2a, mcr-mir-67, mcr-mir-750 
P. canaliculata pcn-bantam,  pcn-mir-750 
P. maculata pcr-bantam, pcr-mir-2a, pcr-mir-67, pcr-mir-750 
LophotrocozoaL. nyassanus lny-mir-1990, lny-mir-745a, lny-mir-96b, lny-mir-1994 
M. cornuarietis mcr-mir-1990, mcr-mir-96b, mcr-mir-1994 
P. canaliculatapcn-mir-1990, pcn-mir-745a, pcn-mir-1992, pcn-mir-96b, pcnmir-1994 
P. maculata pmc-mir-1990, pmc-mir-745a, pmc-mir-1992, pmc-mir-96b, pmc-mir-1994 
MolluscsL. nyassanus lny-mir-1985, lny-mir-12096b 
M. cornuarietis mcr-mir-1985, mcr-mir-12096b 
P. canaliculata pcn-mir-1985, pcn-mir-12096b 
P. maculata Pmc-mir-1985, Pmc-mir-12096b 
Table 4: miRNAs selected by taxon 

The precursors of miR-33, Bilateria-specific-miRNA, were found in this study:  lnymiR-33, mcr-miR-33, pcn-miR-33, pmc-miR-33. All precursors found the 4 species showed two mature miRNAs, 3p and 5p. The sequences of the precursors showed high conservation in the primary and secondary structures. The similarity of the ampullariids miRNA structures were found to be highly conserved each other (Figure 10 and 11).

Figure 10: Alignment of mir-33 in L. nyassanus, M. cornuarietis, P. canaliculata and P. maculata and their orthologs; pcn – P. canaliculata, pmc – P. maculata, mcr – M. cornuarietis, lny- L. nyassanus, lgi – L. gigantea, tcf – Triops cancriformis, cte – Capitella teleta, sko – Saccoglossus kowalevskii, pny – Pundamilia nyererei.

Figure 11: Secondary structure of mir-33 in L. Nyassanus, M. cornuarietis, P. canaliculata, P. maculata and their orthologs; pcn – P. canaliculata, pmc – P. maculata, mcr – M. cornuarietis, lny – L. nyassanus, lgi – L. gigantea, tcf – Triops cancriformis, cte – C. teleta, sko – Saccoglossus kowalevskii, pny – Pundamilia nyererei. 

The Protostome and Deuterostome clades were observed in the phylogenetic. The first clade was subdivided into Lophotrochozoa and Ecdysozoa, being represented respectively by molluscs and arthropods. In the phylum Molluscs, the presence of gastropods was observed, where the ampullariids species were grouped. Deuterostome organisms were represented by clades; mammals, birds and amphibians, verifying wide distribution in phylogenetic clades. This distribution corroborated with the tree of life (Figure 12).

Figure 12: mir-33 phylogenetic distribution for L. Nyassanus, M. cornuarietis P. canaliculata, P. maculata  and their orthologs.

The precursors of miRNAs (protostomes specific) lny-miR-67, mcr-miR-67, pcn-miR67, pmc-miR-67 and a mature miRNA were identified for each Ampullariidae species analyzed. The primary structures of the mi-67 family in L. nyassanus, M. cornuarietis. P. canaliculata and P. maculata were highly conserved between themselves and their orthologs, especially in the region of mature miRNAs. Note the extreme similarity in the secondary structures of the precursors in M. cornuarietis and P. maculata (Figure 13 andFigure 14).

Figure 13: mir-67 alignment to L. nyassanus, M. cornuarietis. P. canaliculata e P. maculata and their orthologs; asu – Ascaris suum, cbn – Caenorhabditis brenneri, cbr – Caenorhabditis briggsae, cel – C. elegans, lny – L. nyassanus, lgi – L. gigantea, mcr – M. cornuarietis, pcn – P. canaliculata, pmc – P. maculata.

Figure 14: Secondary structure of mir-67 in L. nyassanus, M. cornuarietis, P. maculata and their orthologs; asu – Ascaris suum, cbn – Caenorhabditis brenneri, lny – L. nyassanus, lgi – L. gigantea, mcr-M. cornuarietis, pmc – P. maculata.

The Lophotrochozoa clade was subdivided into molluscs, annelids and flatworms. In the phylum Molluscs the presence of bivalves and gastropods was observed. For the first class of molluscs it was represented by the species: Crassostrea virginica, Crassostrea hongkongensis, Pinctada martensii, Pinctada Fucata, Patinopecten yessoensis and Argopecten purpuratus. The gastropod animals were represented exclusively by the analyzed Ampullariidae species: L. nyassanus, M. cornuarietis. P. canaliculata and P. maculata. This distribution is consistent with the tree of life, where a wide distribution was found in the phylogenetic clades (Figure 15).

Figure 15: Phylogenetic distribution of mir-67 in L. nyassanus, M. cornuarietis. P. canaliculata e P. maculata and their orthologs.

Alignment of mir-96b precursors from Ampullariidae species and their orthologs demonstrated the presence of two mature miRNAs (5p and 3p strands) in L. nyassanus, M. cornuarietis. P. canaliculata and P. maculata. The 3p strand of mature miRNAs from the analyzed Ampullariidae species showed 100% identity to each other and high similarity with orthologs. The sequence of the precursors of miR-96b, specific Lophochotrozoa, demonstrated high similarity in secondary structure in the gastropod organisms L.nyassanus, M. cornuarietis and B.glabrata (Figure 16 and Figure 17).

Figure 16: Alignment of miRNA-96b precursors of L. nyassanus. M. cornuarietis, P. canaliculata e P. maculata and their orthologs. bgl: B. glabrata, csi: Cyclina sinensis, lny: Lanistes nyassanus, mcr: Marisa cornuarietis, pfu: Pinctada fucata

Figure 17: Secondary structure of ampullariids miRNA-96b precursors studied with their orthologs; bgl: B. glabrata, csi: Cyclina sinensis, lny: Lanistes nyassanus, mcr: Marisa cornuarietis, pfu: Pinctada fucata

The phylogenetic analysis of the precursors of mir-96b (Lophotrochozoa- specific) was verified the distribution among the species in three phyla (Molluscs, Anelideos and Platelmiltos). A wide distribution was observed in the phylum Molluscs, where the presence of cephalopods, bivalves and gastropods was verified. The first class of molluscs was represented by Architeuthis dux. The bivalve class was represented by Crassostrea gigas, C. virginica, C. hongkongensis, P. Fucata and P. martensii. The presence of snails of the Ampullariidae family analyzed once again appeared grouped together, as expected. In addition, ampullariids appear segregated in the gastropod class with molluscs of the genus Haliotis and Radix auricularia species. A wide distribution was observed in the Lophotrochozoa clade, being constructed a phylogenetic tree consistent with the tree of life. (Figure 18).

Figure 18: Phylogenetic analysis of miR-96b (Lophotrochozoa-specific) for L. nyassanus. M. cornuarietis, P. canaliculata e P. maculata and their orthologs. 

The mir-1990 family predicted for species L. nyassanus, M. cornuarietis, P. canaliculata and P. maculata is composed of the precursors lny-mir-1990, mcr-mir-1990, pcnmir-1990 pcm-mir-1990 and the mature miRNAs, 3p and 5p, in each species ampullariids. For this family, the conservation and great similarity of the primary structure was observed, as well as secondary structures demonstrate conservation among orthologs (Figure 19 and Figure 20).

Figure 19: mir-1990 alignment for L. nyassanus, M. cornuarietis. P. canaliculata e P. maculata and their orthologs; pcn – P. canaliculata, lny – L. nyassanus, mcr – M. cornuarietis, pmc – P. maculata, lgi- L. gigantea, mle– M. leonina, cgi – C. gigas, bgl – B. glabrata, hru- H. rufescens

Figure 20: Secondary structure of the mir-1990 in L. nyassanus, M. cornuarietis. P. canaliculata e P. maculata and their orthologs 

DISCUSSION OF RESULTS 

In this study, three putative Argonaute proteins were found in L. nyassanus (Lny23421_c0_g1, Lny25187_c0_g1 and Lny28936_c7_g7). The same number encountered in M. cornuarietis (Mco25093_c0_g2, Mco26370_c0_g1, Mco27042_c1_g1). For each specie of the genus pomacea, two probable piwi protein were found, in orders P. canaliculata (Pca69814_c1_g2, Pca61914_c0_g1) and P. maculata (Pma63715_c0_g1_, Pma60501_c2_g16). The Argonaute protein is characterized by the presence of the PIWI and PAZ domains (PIWI-Argonaute-Zwille) (Song et al., 2004). The phylogenetic tree presented these two clades in a distinct way (figure 14). The argonaute proteins of eukaryotes can present in addition to these the following domains, N (N-terminal), and MID, together with two domain ligands L1 and L2 (Swarts et al., 2014). The PIWI domain in Argonaute contains three conserved catalytic residues composed of two aspartates and a histidine, called ‘DDH’ (Yang and Steitz, 1995). It was confirmed the presence of these conserved domains in Argonaute putative proteins identified. The analysis of the amino acid sequence of the PIWI domains confirmed the maintenance of the DDH amino acids (aspartic acid, aspartic acid and histidine) in the proper positions for the activity performed. The PIWI domains have great functional homology with RNase H, which is known for enabling the fragmentation of RNA sequences (Parker et al., 2004). The size of the predicted proteins was compatible with the average size of evolutionarily close organisms. 

In this study, Lny28171_c0_g1 and Pca67531_c1_g2 DICER proteins were found in L. nyassanus and P. canaliculata respectively. Lny18581_c0_g2, Mco22837_c2_g1 and Pca67887_c3_g1 DROSHA proteins were found in L. nyassanus, M. cornuarietis and P. canaliculata respectively. The DROSHA and DICER proteins belong to the family of endoribonucleases, ribonuclease III (RNase III) that have specificity for cleavage of doublestranded RNA (dsRNA). The DROSHA protein is a class II RNAase III containing two endonuclease domains and a dsRBD. Dicer is a class III RNAase III that has two endonuclease domains, a dsRBD, an N-terminal helicase domain and a PAZ domain. Functionally class III enzymes produce single-stranded RNA (ssRNA) products of approximately 22 nucleotides in length from long dsRNA substrates (Bernstein et al., 2001; Blaszczyk et al., 2001; Filippov et al., 2000). The canonical pathway of miRNA biogenesis is the dominant pathway in which miRNAs are processed by class II (Drosha) and class III (Dicer) RNAases consecutively. These RNases III contain two catalytic centers with amino acids important for function. The Riboc I domain is composed of the amino acids; E (Glutamic acid), E (Glutamic acid), D (Aspartic acid) D (Aspartic acid); and the Riboc II domain is formed by E (Glutamic acid), Q (Glutamine), E (Glutamic acid), D (Aspartic acid), D (Aspartic acid), E (Glutamic acid); responsible for the cleavage of diester bonds (Blaszczyk et al., 2001).  

The DROSHA protein predicted from P. canaliculata (Pca67887_c3_g1) showed compatible length and similarity with its orthologs (Table 1). The analysis of conserved domains showed the two endonuclease domains and one dsrm, characteristic of class II RNAases III. In the search for catalytic sites, he presented the first set composed of amino acids E-E-D-D in the Riboc I domain (figure 16). The 1525 amino acid length was compatible with the 1523 amino acid P. canaliculata protein (ID: XP_025111754.1) registered at the NCBI (National Center for Biotechnology Information ― http://.ncbi.nlm.nih.gov/). These data corroborated the identification of the putative Drosha protein in P. canaliculata. On the other hand, in the DROSHA protein of M. cornuarietis (Mco22837 c2 g1), despite having a relatively smaller size than its orthologs (figure 15), the lack of coverage in the prediction of the protein did not affect the prediction of conserved domains, since it occurred, probably in the N-terminal portion of the protein, which, for DROSHA, functional conserved domains have not yet been identified. This predicted Drosha protein had two RIBOc domains and a double-stranded RNA binding domain and the presence of complete catalytic amino acid sequence in the Riboc I domain. These data support the identification for the putative Drosha protein in M. cornuarietis (Mco22837 c2 g1). Regarding another DROSHA candidate, L. nyassanus (Lny18581_c0_g2) presented a Riboc domain, a fragment of the second Riboc domain and the presence of a complete catalytic amino acid sequence in the Riboc I domain. The regions in the genome predicted for this gene are probably in error of annotation leading to the truncation of the protein sequences, consequently the size of the amino acid sequences for these genes were underestimated in size and structure. Further versions of the genome may solve this problem. 

The putative DICER protein of P. canaliculata (Pca67531_c1_g2) was shorter in length (1984 aa) than the evolutionarily close organisms A. californica (2503 aa) and B. glabrata (2332 aa) (Table 2). Only conserved domains were identified; Ribonucleas_3_3, PAZ, Dicer dimer and Helicase_C; not showing the second Riboc domain. When verifying the starting position of the second Riboc domain of closely related organisms, A. californica (2273 aa) and B. glabrata (2100 aa), corroborated the fact that the second ribonuclease was not identified in the probable putative sequence. In Dicer, probably due to the lack of coverage in the sequencing of this region of the genome of P. maculata, it was not putative to evidence the presence of the second catalytic domain Riboc. The P. canaliculata Dicer protein (XP_025104500) identified at the NCBI with a length 2249 amino acids was orthologous to Pca67531_c1_g2, with a query coverage of 99% and a percentage of identification of 98.83%. The protein sequence Pca67531_c1_g2 using the BLASTp tool against the genome of P. canaliculata in the AmbuBase database (Department of Computer Science HKBU – https://www.comp.hkbu.edu.hk/~db/AmpuBase/#&panel1- 4) presented a percentage of identification and coverage of 100% (1984/1984 amino acids). Part of the first catalytic site of amino acids (ED) was observed in Riboc in the candidate protein cited (Pca67531_c1_g2) When performing the search for catalytic residues having as query NP_524453.1 Dicer-1 D. melanogaster the orthologous molluscs presented only the amino acids in the first domain Riboc catalysts (ED), with the exception of Pectem maximus (EDDE) NCBI (National Center for Biotechnology Information ― http://.ncbi.nlm.nih.gov/). The data presented corroborated the identification of the putative protein. On the other hand, the putative DICER protein from L. nyassanus (Lny28171_c0_g1), despite having a relatively smaller size than its orthologs (figure 15), the lack of coverage in the protein prediction did not affect the prediction of the two conserved Riboc domains and a PAZ domain, since it was probably in the N-terminal portion of the protein, which for DICER the catalytic domains lie at the C-terminal end. This predicted DICER protein showed the presence of the complete catalytic amino acid sequence in the Riboc I domain. The data support the identification of the putative DICER protein in L. nyassanus. The regions in the genome predicted for this gene are likely to have annotation error leading to the truncation of the protein sequences. Consequently, the size of the amino acid sequences for these genes are underestimated in size and structure (Queiroz et al., 2017).  

In this study, Lny14170_c0_g1, Lny12572_c0_g1, Lny27692_c0_g2, Mco24030_c0_g2, Mco81614_c0_g2, Pca62316_c0_g2, Pca65193_c2_g1, Pma59467_c0_g1, and Pma63331_c0_g1 Exportin protein family (XPO-5, XPO1 and XPO-T) were found L. nyassanus, M. cornuarietis, P. canaliculata and P. maculata respectively. Intrinsically, the transport of miRNA precursors is performed by the XPO-5 protein. Occasionally the transport of pre-miRNAs can be accomplished by XPO-1 entering the non-canonical miRNA processing pathway without Drosha cleavage. Displacement of pre-miRNAs is also performed by exportin T (XPOT) in a Ran-GTP-dependent manner (Ruby et al., 2007). The search for Exportins proteins (XPO-5, XPO-1 and XPO-T) using the methodology of this work, resulted in the putative proteins found. 

In order to emphasize the conservation of these putative proteins identified in the ampullariids species, a phylogenetic tree and studies of conserved domains were constructed. The phylogenetic analysis revealed a distribution of the studied amino acid sequences, in the genes of the Exportin protein family (XPO-5, XPO1 and XPO-T), grouped close to the Deuterostomes clade. Putative proteins from L. nyassanus. M. cornuarietis, P. canaliculata and P. maculata aggregated with their respective orthologs close to the phylum Molluscs and the clade Deuterostomes. The similarity of the identified Exportins of the Ampullariidae species was verified together with their orthologs, observing the level of preservation of the conserved domains and protein length. For Exportin 5 the putative proteins (Lny14170_c0_g1, Pca62316_c0_g2) presented the two conserved domains (Xpo1 and IBN_N) characteristic of the XPO-5 protein and a size compatible with the average dimension of its orthologs, confirming the identification of the putative Exportin 5 proteins. Exportin 1 putative proteins (Lny12572_c0_g1, Mco24030_c0_g2, Pca65193_c2_g1, Pma59467_c0_g1) showed six conserved domains characteristic of XPO-1 (IBN_N, Xpo1, CRM1_repeat, CRM1_repeat_2, CRM1_repeat_3 and only five domains except for the CRM1_cc1. The Pca65193_c2 did not show the conserved IBN_N domain. The lack of protein sequence prediction in the N-terminal region was already reported in this work for P. canaliculata, resulting in a relatively smaller size when compared to the analyzed species. Even the length of the protein sequences analyzed for XPO-1 (Lny12572_c0_g1, Mco24030_c0_g2, Pca65193_c2_g1, Pma59467_c0_g1) were strictly the same (1070 aa). This dimension also matched to the average size of the Exportin 1 protein for orthologous organisms. All these factors corroborated the identification of putative XPO-1 proteins in the analyzed species. Likewise, Exportin T in genes (Lny27692_c0_g2, Mco81614_c0_g2 and Pma63331_c0_g1) showed a characteristic conserved domain for all protein sequences (Xpo1). This time, lack of protein sequence prediction was verified in the Nterminal region for M. cornuarietis (Mco81614_c0_g2), presenting a smaller size (331 aa) than the analyzed species. Interestingly, the protein sequences (Lny27692_c0_g2 and Pma63331_c0_g1) also had the same size (966 aa) and matched the average size of XOP-T proteins from orthologous organisms. They were convincing evidence for the identification of putative Exportin T in the analyzed species. 

The genomes of four ampullariids (Old World – Lanistes nyassanus and New World Pomacea canaliculata, P. maculata, and Marisa cornuarietis) were published by Sun et al., 2019, in order to understand the genomic basis underlying the diversity of the Ampullariidae species and also their behavioral, morphological, and physiological adaptations. In our study we were able to apply a robust analysis to identify and characterize miRNAs and the putative proteins involved in miRNA pathway using the genome and predicted proteome from four Ampullariidae species. Either miRNAs or proteins of their biosynthetic pathway were described in several organisms, including several kingdoms such as plants and animals. In model species such as D. melanogaster and C. elegans mature and precursor miRNAs and also their pathway proteins have been extensively studied. Even though Mollusca is a phylum containing many species, little has been studied about miRNAs and their processing pathway in this phylum. 

Sun et al., 2019, sequenced and analyzed the four genomes of ampullariids species, P. canaliculata, P. maculata, and L. nyassanus. The genomes were assembled in 535.5 Mb million base pairs (Mb) in total length or smaller and produced thousands of scaffolds anchored or not to chromosomes. Automatic annotation analyses of the four genomes displayed around 20 thousand gene models (Sun et al., 2019). So far, no study on miRNAs and their processing pathway in ampullariids species has been reported in the literature. In order to improve and increase the annotation of non-protein coding genes, such as miRNAs, and also of proteincoding genes of small RNA processing pathways, this work proposed an in-depth study of mature and precursory miRNAs and their pathway processing in the genome of the four species of ampullariids. 

In our study, we were able to identify conserved miRNAs in four species, being 162 mature miRNAs and 141 miRNA precursors in L. nyassanus, 297 mature miRNAs and 279 miRNA precursors in M. cornuarietis, 296 mature miRNAs and 269 miRNA precursors in P. canaliculata and 316 mature miRNAs and 299 miRNA precursors in P. maculata. Two miRNAs found in all species was Mollusca-specific. 

The minimum free energy value, -18 kcal/mol, in general, is the value necessary for a miRNA precursor to be stable and generates miRNAs (de Souza Gomes et al., 2011; Hofacker, 2009). The sequences of pre-miRNAs presented in this work MFE mean of (minimum free energy) and variation of values (Table 3). The mean value of the studied precursors was similar to the mean value of pre-miRNAs of lophotrochozoan species (-31.27 kcal/mol) (de Souza Gomes et al., 2011). Of the pre-miRNAs deposited in miRBase  approximately 75% of those had minimum free energy values below -18.50 kcal/mol (Leclercq et al., 2013). 

The size of the precursors of the lophotrochozoan species is between 57 and 153 nucleotides, with an average of 90 nucleotides per sequence(de Souza Gomes et al., 2011). In this study, the average size of the precursors and the length variation were showed in Table 3. The size of pre-miRNAs in animal species ranges from 47 to 177 nucleotides with a mean of 87 ± 29 nucleotides (Zhou et al., 2009). The data from this study corroborated with the average size of pre-miRNAs and the reference range for the size of lophotrochozoan species and metazoans in general. 

GC content is one of the main parameters for the identification of pre-miRNAs and is very important for the stability of the secondary structure of precursors. The percentage of premiRNAs present in miRBase that had GC content above 20% and below 65% was 80% (Leclercq et al., 2013). In this work, the variation of GC content and the percentage of miRNAs in the range between 20% and 65% were respectively; L. nyassanus (20.25% to 69.05%) and 98.5%, M. cornuarietis (20.20% to 67.03%) and 98%, P. canaliculata (20.83% and 68.29 %) and 99.25%, P. maculata (20.40% and 68.75%) and 98.66%. The precursors presented mean values for the guanine-cytosine (GC) content for L. nyassanus, M. cornuarietis, P. canaliculata and P. maculata; 41.26%, 41.64%, 41.70%, 40.74% respectively. The mean value of the GC content of the studied ampullariids precursors was similar to the mean of the lophotrochozoan species (40.49%) (de Souza Gomes et al., 2011). Comparatively the GC content of the miRNA precursors were not observed statistical differences (p > 0.05) between the Ampullariidae analyzed and lophotrochozoan species. 

Mature miRNAs were identified in the genome of 4 ampullariids species: 162 miRNAs for L. nyassanus, 297 miRNAs for M. cornuarietis, 296 miRNAs for P. canaliculata and 316 miRNAs for P. maculata. The number of miRNAs identified in L. nyassanus genome was similar to the gastropod species Conus tribblei (159 miRNAs) (Huang et al., 2021). The number of miRNAs found in M. cornuarietis (297 miRNAs), P. canaliculata (296 miRNAs), P. maculata (316 miRNAs) corroborated with the number of miRNAs identified for Lymnaea stagnalis (freshwater gastropod molluscs), 264 conserved mature miRNAs (Walker et al., 2018). There are plausible justifications for the smaller quantitative number of miRNAs for L. nyassanus compared to other ampullariids under study. L. nyassanus has fewer adaptive and functional properties when compared to other ampullariids under study (Sun et al., 2019). L. nyassanus is endemic to Lake Nyasa in Africa (COWIE, 2015; Dohrn, 1865; Van Bocxlaer, 2017), while P. canaliculata and P. maculata are geographically distributed on almost every continent, fast growing, high reproduction rate, tolerance and adaptation to environmental stress (Liu et al., 2018; Pimentel et al., 2005). Species of the genus Pomacea have the adaptive property of terrestrial egg laying, while L. nyassanus does egg laying in the aquatic environment (Sun et al., 2019). In addition, M. cornuarietis is better able to resist factors such as low temperature and hypoxia (Matsukura et al., 2016; Mu et al., 2018). As reported in previous studies, there is a significant increase in the number of miRNAs as the morphological complexity of the species increases (Heimberg et al., 2008), thus, justifying a lower number of miRNAs in L. nyassanus, when compared to other ampullariids in study. 

The data presented corroborated the identification of miRNAs with the literature (Ha and Kim, 2014; Rachagani et al., 2010). The number of nucleotides in the mature miRNAs from ampullariids species had an arithmetic mean of 21.94 nucleotides/miRNA, with a size range from 17 to 25 nucleotides. 

The Argonaute protein prefers the uracil nucleotide in the first position of the miRNA sequence, to be inserted in the RISC complex (Seitz et al., 2011). In this work mature miRNAs from L. nyassanus, M. cornuarietis. P. canaliculata and P. maculata presented, respectively, quantitative uracil nucleotides in the first position; 41.97%, 44.44%, 41.2% and 42.40%. Mostly presented this nucleotide in the first position in the miRNA sequence.  

The structural and thermodynamic characteristics of pre-miRNAs were more similar among M. cornuarietis, P. canaliculata and P. maculata species. L. nyassanus had more dissimilar values, however it is the organism that most diverges with reference to the three species mentioned. The justifications for this probably in phylogeny. The division between the New World (M. cornuarietis, P. canaliculata and P. maculata) and Old World (L. nyassanus) lineages may have occurred at the time of the separation from Gondwana, around 120 million years ago (Jokat et al., 2003; Sun et al., 2019). 

miR-33, Bilateria-specific, was found in the genomes of L. nyassanus, M. cornuarietis, P. canaliculata and P. maculata and had two mature miRNAs (3p and 5p chains), corroborating studies observed in bivalve molluscs and D. melanogaster (Bao et al., 2014; Clerbaux et al., 2021). Mir-33-5p supposedly regulates genes involved in the heavy metal-induced stress response in the mollusc Tegillarca granosa (Bao et al., 2014). Likewise, P. canaliculata is able to survive in environments with heavy metals, where there is a high concentration of these metals in the gills (Kruatrachue et al., 2011), probably regulated by this micro-RNA. In this study, miR-33-5p showed 100% similarity between the ampullariid and orthologous species (figure 1). Furthermore, miR -33 is a pleiotropic regulator of metabolic and developmental processes in D. melanogaster. It is suggested that mir-33 is a conserved regulator of lipid homeostasis (Clerbaux et al., 2021). 

The family of precursors of mir-67, protostomes specific, was found in the genome of the 4 species Ampullariidae with the presence of mature miR-67-3p for each species, corroborating the wide distribution of this miRNA in molluscs and C. elegans (Huang et al., 2021; Ma et al., 2017). Like miRNA-33, miRNA-67 is potentially associated with gene regulation when exposed to toxic levels of cadmium (Cd) in the bivalve T. granosa (Bao et al., 2014). Regulation of the immune system is also performed by mi-67 in C. elegans when exposed to pathogenic bacteria Pseudomonas aeruginosa. miRNA-67 inhibits the expression of the SAX-7 protein in the worm, that induces pathogen avoidance behavior (Ma et al., 2017). Furthermore, mir-67 was identified in 34 species of molluscs (Huang et al., 2021). 

The sequences of the precursors of miR-96b presented two mature miRNAs, 3p and 5p, in the genomes of the species: L. nyassanus, M. cornuarietis. P. canaliculata and P. maculata. The mir-96b, Lophotrochozoa-specific, regulates genes post-transcriptionally in flatworms and molluscs. In the freshwater worm: Schmidtea mediterranea may have functions in wound healing, neoblast proliferation (mainly responsible for regeneration in planarians) and blastema differentiation (Sasidharan et al., 2013). In molluscs, miR-96b is involved together with other miRNAs in the pigmentation of the shells of C. gigas, inducing the synthesis of melanin, carotenoid or tetrapyrrole (Feng et al., 2020). The miR-96b was also found in the genome of gastropods evolutionarily close to the 4 species Ampullariidae, B. glabrata (Queiroz et al., 2020) and A. californica (Huang et al., 2021). Findings in the scientific literature corroborate the data of this work. 

The family of mir-1990, Lophotrochozoa-specific miRNA, presented two mature miRNAs, 3p and 5p. A wide distribution of mir-1990 in molluscs was verified, mainly in 15 bivalve animals, 10 gastropods and 1 cephalopod (Huang et al., 2021). Few studies describe the role of miRNAs in molluscs development. However, mir-1990-3p may play a vital role in the shell biomineralization process in P. fucata, as it is highly expressed in mantle tissues (Huang et al., 2021). mir-1990 was found only in molluscs, however the mir-1990 family is Lophotrochozoa-specific as cte-miR-1990c-3p, cte-miR-1990a and cte-miR-1990b have been identified in Capitella teleta (http: //www.mirbase.org/). The wide distribution of mir-1990 in molluscs and the likely mechanism of gene control in the species corroborated the data in this study. 

The functional form of a miRNA is generally 22 nucleotides in length, with a size range from 17 to 25 nucleotides (Ha and Kim, 2014; Rachagani et al., 2010). These small RNAs are generated by two cleavage reactions in which miRNAs are processed by RNAases III Drosha and Dicer consecutively (Lee et al., 2002; Tanzer and Stadler, 2004; Y. F. Zhang et al., 2009), in canonical miRNA processing pathway. The microprocessor complex, composed of a protein (DGCR8) and a ribonuclease III enzyme (Drosha), performs the cleavage of pri-miRNA in premiRNA in the cell nucleus. DGCR8 binds to a structure in pri-miRNA while Drosha performs the duplex scission. After the pre-miRNAs are generated, they are transported to the cytoplasm through the exportin 5 (XPO5) / RanGTP complex (Alarcón et al., 2015; Denli et al., 2004; Han et al., 2004; Okada et al., 2009). In the cytoplasm, pre-miRNAs are cleaved by Dicer into mature. The miRNAs can effect their action through the miRISC complex, in which a protein from the AGO family is loaded with a miRNA guide strand that directs the complex to base pairing, which by different mechanisms induces the silencing of gene expression (Ipsaro and Joshua-Tor, 2015; Jo et al., 2015). The AGO-2 protein has endonucleolytic action and the level of complementarity formed between the miRISC complex and the target mRNA is indicative of the complex’s action (Jo et al., 2015). All these factors in our analysis were of great relevance.  

CONCLUSION 

The search for miRNAs in L. nyassanus, M. cornuarietis, P. canaliculata and P. maculata was able to predict 162 mature forms of miRNAs for L. nyassanus, 297 miRNAs for M. cornuarietis, 296 miRNAs for P. canaliculata and 316 miRNAs for P. maculata

It was possible to identify the biogenesis pathway of miRNAs in species of the Ampullariidae family, highlighting several putative key pathway proteins; Argonaute, DROSHA, DICER, XPOs, among others. The prediction of some proteins of the miRNA pathway was not verified in some species, probably due to annotation error in the part of the sequenced gene specific to each organism. 

Our results expand the study of miRNAs in molluscs, bringing new challenges to the understanding of the essential processes related to the invasiveness of P. canaliculata and P. maculata species and supposed advances in the control of these molluscs considered agricultural pests for humanity. 

Thus, the search for miRNAs and their processing pathways in the species L. nyassanus, M. cornuarietis, P. canaliculata and P. maculata were able to predict structures that found the established criteria, also complementing previous results from our study group as well, as creating the need to establish putative new target genes for the discovered miRNAs. 

REFERENCES 

Alarcón, C.R., Lee, H., Goodarzi, H., Halberg, N., Tavazoie, S.F., 2015. N6-methyladenosine marks primary microRNAs for processing. Nature 519, 482–485. https://doi.org/10.1038/nature14281 

Bao, Y., Zhang, L., Dong, Y., Lin, Z., 2014. Identification and comparative analysis of the tegillarca granosa haemocytes MicroRNA transcriptome in response to Cd using a deep sequencing approach. PLoS One 9. https://doi.org/10.1371/journal.pone.0093619 

Bartel, D.P., 2009. MicroRNAs: Target Recognition and Regulatory Functions. Cell 136, 215–233. https://doi.org/10.1016/j.cell.2009.01.002 

Bernhart, S.H., Hofacker, I.L., Will, S., Gruber, A.R., Stadler, P.F., 2008. RNAalifold: Improved consensus structure prediction for RNA alignments. BMC Bioinformatics 9. https://doi.org/10.1186/1471-2105-9-474 

Bernstein, E., Caudy, A.A., Hammond, S.M., Hannon, G.J., 2001. Role for a bidentate ribonuclease in the initiation step of RNA interference. Nature 409, 363–366. https://doi.org/10.1038/35053110 

Blaszczyk, J., Tropea, J.E., Bubunenko, M., Routzahn, K.M., Waugh, D.S., Court, D.L., Ji, X., 2001. Crystallographic and modeling studies of RNase III suggest a mechanism for double-stranded RNA cleavage. Structure 9, 1225–1236. https://doi.org/10.1016/S09692126(01)00685-2 

Carlsson, N.O.L., Brönmark, C., Hansson, L., Carlsson, N.O.L., Bronmark, C., Hansson, L., 2016. Invading Herbivory : The Golden Apple Snail Alters Ecosystem Functioning in Asian Wetlands Stable URL : http://www.jstor.org/stable/3450580 REFERENCES Linked references are available on JSTOR for this article : You may need to log in to JSTOR to access th 85, 1575–1580. 

Cawas, J.R., Quisao, C.J.T., Castillo, D.S.C., Pornobi, K.O., 2020. Prevalence of Angiostrongylus cantonensis among different species of snails in the village of Bagong Sikat Muñoz, Nueva Ecija, Philippines and its associated risk factors for zoonotic transmission. J. Parasit. Dis. 44, 388–394. https://doi.org/10.1007/s12639-020-01200-0 

Chen, C.Y., Kuo, C.C., Lo, C.P., Huang, M.Y., Wang, Y.M., Wang, W.Y., 2014. Eosinophilic meningoencephalitis caused by Angiostrongylus cantonesis. Qjm 107, 573–575. https://doi.org/10.1093/qjmed/hcr261 

Clerbaux, L.-A., Schultz, H., Roman-Holba, S., Ruan, D., Yu, R., Lamb, A., Bommer, G., Kennell, J., 2021. The microRNA miR-33 is a pleiotropic regulator of metabolic and developmental processes in Drosophila melanogaster. Dev. Dyn. https://doi.org/10.1002/dvdy.344

COWIE, R.H., 2015. The recent apple snails of Africa and Asia (Mollusca: Gastropoda: Ampullariidae: Afropomus, Forbesopomus, Lanistes, Pila, Saulea): a nomenclatural and type catalogue. The apple snails of the Americas: addenda and corr. Zootaxa 3940, 1. https://doi.org/10.11646/zootaxa.3940.1.1 

Dalton, M.F., Fenton, H., Cleveland, C.A., Elsmo, E.J., Yabsley, M.J., 2017. Eosinophilic meningoencephalitis associated with rat lungworm (Angiostrongylus cantonensis) migration in two nine-banded armadillos (Dasypus novemcinctus) and an opossum (Didelphis virginiana) in the southeastern United States. Int. J. Parasitol. Parasites Wildl. 6, 131–134. https://doi.org/10.1016/j.ijppaw.2017.05.004 

de Souza Gomes, M., Muniyappa, M.K., Carvalho, S.G., Guerra-Sá, R., Spillane, C., 2011. Genome-wide identification of novel microRNAs and their target genes in the human parasite Schistosoma mansoni. Genomics 98, 96–111. https://doi.org/10.1016/j.ygeno.2011.05.007 

De Wit, E., Linsen, S.E.V., Cuppen, E., Berezikov, E., 2009. Repertoire and evolution of miRNA genes in four divergent nematode species. Genome Res. 19, 2064–2074. https://doi.org/10.1101/gr.093781.109 

Denli, A.M., Tops, B.B.J., Plasterk, R.H.A., Ketting, R.F., Hannon, G.J., 2004. Processing of primary microRNAs by the Microprocessor complex. Nature 432, 231–235. https://doi.org/10.1038/nature03049 

Dohrn, H., 1865. List of the land and freshwater shells of the Zambesi and lake Nyasa, eastern tropical Africa, collected by John Kirk. Proc. Zool. Soc. London 1865, 231–234. 

Feng, D., Li, Q., Yu, H., Liu, S., Kong, L., Du, S., 2020. Integrated analysis of microRNA and mRNA expression profiles in Crassostrea gigas to reveal functional miRNA and miRNA-targets regulating shell pigmentation. https://doi.org/10.1038/s41598-02077181-0 

Filippov, V., Solovyev, V., Filippova, M., Gill, S.S., 2000. A novel type of RNase III family proteins in eukaryotes. Gene 245, 213–221. https://doi.org/10.1016/S03781119(99)00571-5 

Gruber, A.R., Lorenz, R., Bernhart, S.H., Neuböck, R., Hofacker, I.L., 2008. The Vienna RNA websuite. Nucleic Acids Res. 36. https://doi.org/10.1093/nar/gkn188 

Ha, M., Kim, V.N., 2014. Regulation of microRNA biogenesis. Nat. Rev. Mol. Cell Biol. https://doi.org/10.1038/nrm3838 

Han, J., Lee, Y., Yeom, K.H., Kim, Y.K., Jin, H., Kim, V.N., 2004. The Drosha-DGCR8 complex in primary microRNA processing. Genes Dev. 18, 3016–3027. https://doi.org/10.1101/gad.1262504 

Hayes, K.A., Cowie, R.H., Jørgensen, A., Schultheiß, R., Albrecht, C., Thiengo, S.C., 2009. Molluscan Models in Evolutionary Biology: Apple Snails (Gastropoda: Ampullariidae) as a System for Addressing Fundamental Questions*. Am. Malacol. Bull. 27, 47–58. https://doi.org/10.4003/006.027.0204

Heimberg, A.M., Sempere, L.F., Moy, V.N., Donoghue, P.C.J., Peterson, K.J., 2008. MicroRNAs and the advent of vertebrate morphological complexity. Proc. Natl. Acad. Sci. U. S. A. 105, 2946–2950. https://doi.org/10.1073/pnas.0712259105 

Hofacker, I.L., 2009. Secondary Structure Analysis Using the Vienna Package. Curr. Protoc. Bioinforma. 26, 12.2.1-12.2.16. https://doi.org/10.1002/0471250953.bi1202s26 

Huang, S., Yoshitake, K., Asaduzzaman, M., Kinoshita, S., Watabe, S., Asakawa, S., 2021. Discovery and functional understanding of MiRNAs in molluscs: a genome-wide profiling approach. RNA Biol. https://doi.org/10.1080/15476286.2020.1867798 

Ipsaro, J.J., Joshua-Tor, L., 2015. From guide to target: Molecular insights into eukaryotic RNA-interference machinery. Nat. Struct. Mol. Biol. https://doi.org/10.1038/nsmb.2931 

Jo, M.H., Shin, S., Jung, S.R., Kim, E., Song, J.J., Hohng, S., 2015. Human Argonaute 2 Has Diverse Reaction Pathways on Target RNAs. Mol. Cell 59, 117–124. https://doi.org/10.1016/j.molcel.2015.04.027 

Jokat, W., Boebel, T., König, M., Meyer, U., 2003. Timing and geometry of early Gondwana breakup. J. Geophys. Res. Solid Earth 108. https://doi.org/10.1029/2002jb001802 

Kenny, N.J., Namigai, E.K.O., Marlétaz, F., Hui, J.H.L., Shimeld, S.M., 2015. Draft genome assemblies and predicted microRNA complements of the intertidal lophotrochozoans Patella vulgata (Mollusca, Patellogastropoda) and Spirobranchus (Pomatoceros) lamarcki (Annelida, Serpulida). Mar. Genomics 24, 139–146. https://doi.org/10.1016/j.margen.2015.07.004 

Kruatrachue, M., Sumritdee, C., Pokethitiyook, P., Singhakaew, S., 2011. Histopathological effects of contaminated sediments on golden apple snail (Pomacea canaliculata, Lamarck 1822). Bull. Environ. Contam. Toxicol. 86, 610–614. https://doi.org/10.1007/s00128011-0265-4 

Kumar, S., Stecher, G., Li, M., Knyaz, C., Tamura, K., 2018. MEGA X: Molecular evolutionary genetics analysis across computing platforms. Mol. Biol. Evol. 35, 1547– 1549. https://doi.org/10.1093/molbev/msy096 

Larkin, M.A., Blackshields, G., Brown, N.P., Chenna, R., Mcgettigan, P.A., Mcwilliam, H., Valentin, F., Wallace, I.M., Wilm, A., Lopez, R., Thompson, J.D., Gibson, T.J., Higgins, D.G., Bateman, A., 2007. Clustal W and Clustal X version 2.0 23, 2947–2948. https://doi.org/10.1093/bioinformatics/btm404 

Leclercq, M., Diallo, A.B., Blanchette, M., 2013. Computational prediction of the localization of microRNAs within their pre-miRNA. Nucleic Acids Res. 41, 7200–7211. https://doi.org/10.1093/nar/gkt466 

Lee, Y., Jeon, K., Lee, J.T., Kim, S., Kim, V.N., 2002. MicroRNA maturation: Stepwise processing and subcellular localization. EMBO J. 21, 4663–4670. https://doi.org/10.1093/emboj/cdf476 

Lewis, B.P., Burge, C.B., Bartel, D.P., 2005. Conserved seed pairing, often flanked by adenosines, indicates that thousands of human genes are microRNA targets. Cell 120, 15–20. https://doi.org/10.1016/j.cell.2004.12.035 

Lindo, J.F., Waugh, C., Hall, J., Cunningham-Myrie, C., Ashley, D., Eberhard James J Sullivan, M.L., Bishop, H.S., Robinson, D.G., Holtz, T., Robinson, R.D., 2002. Enzootic Angiostrongylus cantonensis in rats and snails after an outbreak of human eosinophilic meningitis, Jamaica. Emerg. Infect. Dis. 8, 324–326. https://doi.org/10.3201/eid0803.010316 

Liu, N., Zhang, Y., Ren, Y., Wang, H., Li, S., Jiang, F., Yin, L., Qiao, X., Zhang, G., Qian, W., Liu, B., Fan, W., 2018. The genome of the golden apple snail Pomacea canaliculata provides insight into stress tolerance and invasive adaptation. Gigascience 7, 1–13. https://doi.org/10.1093/gigascience/giy101 

Luessi, F., Sollors, J., Torzewski, M., Müller, H.D., Siegel, E., Blum, J., Sommer, C., Vogt, T., Thömke, F., 2009. Eosinophilic meningitis due to angiostrongylus cantonensis in Germany. J. Travel Med. 16, 292–294. https://doi.org/10.1111/j.17088305.2009.00337.x 

Lv, S., Zhang, Y., Steinmann, P., Zhou, X.N., 2008. Emerging angiostrongyliasis in mainland China. Emerg. Infect. Dis. 14, 161–164. https://doi.org/10.3201/eid1401.061529 

Ma, Y.C., Zhang, L., Dai, L.L., Khan, R.U., Zou, C.G., 2017. mir-67 regulates P. aeruginosa avoidance behavior in C. elegans. Biochem. Biophys. Res. Commun. 494, 120–125. https://doi.org/10.1016/j.bbrc.2017.10.069 

Marchler-Bauer, A., Bo, Y., Han, L., He, J., Lanczycki, C.J., Lu, S., Chitsaz, F., Derbyshire, M.K., Geer, R.C., Gonzales, N.R., Gwadz, M., Hurwitz, D.I., Lu, F., Marchler, G.H., Song, J.S., Thanki, N., Wang, Z., Yamashita, R.A., Zhang, D., Zheng, C., Geer, L.Y., Bryant, S.H., 2017. CDD/SPARCLE: Functional classification of proteins via subfamily domain architectures. Nucleic Acids Res. 45, D200–D203. https://doi.org/10.1093/nar/gkw1129 

Marchler-Bauer, A., Bryant, S.H., 2004. CD-Search: Protein domain annotations on the fly. Nucleic Acids Res. 32. https://doi.org/10.1093/nar/gkh454 

Marchler-Bauer, A., Derbyshire, M.K., Gonzales, N.R., Lu, S., Chitsaz, F., Geer, L.Y., Geer, R.C., He, J., Gwadz, M., Hurwitz, D.I., Lanczycki, C.J., Lu, F., Marchler, G.H., Song, J.S., Thanki, N., Wang, Z., Yamashita, R.A., Zhang, D., Zheng, C., Bryant, S.H., 2015. CDD: NCBI’s conserved domain database. Nucleic Acids Res. 43, D222–D226. https://doi.org/10.1093/nar/gku1221 

Marchler-Bauer, A., Lu, S., Anderson, J.B., Chitsaz, F., Derbyshire, M.K., DeWeese-Scott, C., Fong, J.H., Geer, L.Y., Geer, R.C., Gonzales, N.R., Gwadz, M., Hurwitz, D.I., Jackson, J.D., Ke, Z., Lanczycki, C.J., Lu, F., Marchler, G.H., Mullokandov, M., Omelchenko, M. V., Robertson, C.L., Song, J.S., Thanki, N., Yamashita, R.A., Zhang, D., Zhang, N., Zheng, C., Bryant, S.H., 2011. CDD: A Conserved Domain Database for the functional annotation of proteins. Nucleic Acids Res. 39. https://doi.org/10.1093/nar/gkq1189

Matsukura, K., Izumi, Y., Yoshida, K., Wada, T., 2016. Cold tolerance of invasive freshwater snails, Pomacea canaliculata, P. maculata, and their hybrids helps explain their different distributions. Freshw. Biol. 61, 80–87. https://doi.org/10.1111/fwb.12681 

Matsukura, K., Tsumuki, H., Izumi, Y., Wada, T., 2009. Physiological response to low temperature in the freshwater apple snail, Pomacea canaliculata (Gastropoda: Ampullariidae). J. Exp. Biol. 212, 2558–2563. https://doi.org/10.1242/jeb.031500 

Mistry, J., Chuguransky, S., Williams, L., Qureshi, M., Salazar, G.A., Sonnhammer, E.L.L., Tosatto, S.C.E., Paladin, L., Raj, S., Richardson, L.J., Finn, R.D., Bateman, A., 2021. Pfam: The protein families database in 2021. Nucleic Acids Res. 49, D412–D419. https://doi.org/10.1093/nar/gkaa913 

Mu, H., Sun, J., Cheung, S.G., Fang, L., Zhou, H., Luan, T., Zhang, H., Wong, C.K.C., Qiu, J.W., 2018. Comparative proteomics and codon substitution analysis reveal mechanisms of differential resistance to hypoxia in congeneric snails. J. Proteomics 172, 36–48. https://doi.org/10.1016/j.jprot.2017.11.002 

Novarino, A.A., 2011. Record of Small-Clawed Otters (Aonyx cinereus) Foraging on an Invasive Pest Species, Golden Apple Snails (Pomacea canaliculata) in a West Sumatra Rice Field. IUCN/SCC Otter Spec. Gr. Bull. 28, 34–38. 

Okada, C., Yamashita, E., Lee, S.J., Shibata, S., Katahira, J., Nakagawa, A., Yoneda, Y., Tsukihara, T., 2009. A high-Resolution structure of the pre-microrna nuclear export machinery. Science (80-. ). 326, 1275–1279. https://doi.org/10.1126/science.1178705 

Parker, J.S., Roe, S.M., Barford, D., 2004. Crystal structure of a PIWI protein suggests mechanisms for siRNA recognition and slicer activity. EMBO J. 23, 4727–4737. https://doi.org/10.1038/sj.emboj.7600488 

Picone, B., Rhode, C., Roodt-Wilding, R., 2017. Identification and characterization of miRNAs transcriptome in the South African abalone, Haliotis midae. Mar. Genomics 31, 9–12. https://doi.org/10.1016/j.margen.2016.10.005 

Pimentel, D., Zuniga, R., Morrison, D., 2005. Update on the environmental and economic costs associated with alien-invasive species in the United States. Ecol. Econ. 52, 273–288. https://doi.org/10.1016/j.ecolecon.2004.10.002 

Queiroz, F.R., Portilho, L.G., Jeremias, W. de J., Babá, É.H., do Amaral, L.R., Silva, L.M., Coelho, P.M.Z., Caldeira, R.L., Gomes, M. de S., 2020. Deep sequencing of small RNAs reveals the repertoire of miRNAs and piRNAs in Biomphalaria glabrata. Mem. Inst. Oswaldo Cruz 115, 1–13. https://doi.org/10.1590/0074-02760190498 

Queiroz, F.R., Silva, L.M., De Jesus Jeremias, W., Hideo Babá, É., Caldeira, R.L., Coelho, P.M.Z., De Souza Gomes, M., 2017. Differential expression of small RNA pathway genes associated with the Biomphalaria glabrata/Schistosoma mansoni interaction. PLoS One 12, 1–20. https://doi.org/10.1371/journal.pone.0181483 

Rachagani, S., Kumar, S., Batra, S.K., 2010. MicroRNA in pancreatic cancer: Pathological, diagnostic and therapeutic implications. Cancer Lett. 292, 8–16. https://doi.org/10.1016/j.canlet.2009.11.010 

Ruby, J.G., Jan, C.H., Bartel, D.P., 2007. Intronic microRNA precursors that bypass Drosha processing. Nature 448, 83–86. https://doi.org/10.1038/nature05983 

Saitou, N, Nei, M., 1987. The neighbor-joining method: a new method for reconstructing phylogenetic trees. Mol. Biol. Evol. 4, 406–25. https://doi.org/10.1093/oxfordjournals.molbev.a040454 

Saitou, Naruya, Nei, M., 1987. ESCALA CIWA-AR Escala CIWA-Ar(Clinical Institute Withdrawal Assesment for Alcohol) Evaluación del Síndrome de Abstinencia Alcohólica. Mol. Biol. Evol 4, 406–425. https://doi.org/citeulike-article-id:93683 

Sasidharan, V., Lu, Y.C., Bansal, D., Dasari, P., Poduval, D., Seshasayee, A., Resch, A.M., Graveley, B.R., Palakodeti, D., 2013. Identification of neoblast- And regenerationspecific miRNAs in the planarian Schmidtea mediterranea. Rna 19, 1394–1404. https://doi.org/10.1261/rna.038653.113 

Schultz, J.H., Adema, C.M., 2017. Comparative immunogenomics of molluscs. Dev. Comp. Immunol. 75, 3–15. https://doi.org/10.1016/j.dci.2017.03.013 

Seitz, H., Tushir, J.S., Zamore, P.D., 2011. A 5′-uridine amplifies miRNA/miRNA* asymmetry in Drosophila by promoting RNA-induced silencing complex formation. Silence 2, 1–10. https://doi.org/10.1186/1758-907X-2-4 

Shennan Lu et al., 2020. CDD: NCBI’s conserved domain database [WWW Document]. Nucleic Acids Res. URL https://www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi (accessed 6.7.20). 

Song, J.J., Smith, S.K., Hannon, G.J., Joshua-Tor, L., 2004. Crystal structure of argonaute and its implications for RISC slicer activity. Science (80-. ). 305, 1434–1437. https://doi.org/10.1126/science.1102514 

Strong, E.E., Gargominy, O., Ponder, W.F., Bouchet, P., 2008. Global diversity of gastropods (Gastropoda; Mollusca) in freshwater, in: Freshwater Animal Diversity Assessment. Springer Netherlands, Dordrecht, pp. 149–166. https://doi.org/10.1007/978-1-40208259-7_17 

Sun, J., Mu, H., Ip, J.C.H.H., Li, R., Xu, T., Accorsi, A., Sánchez Alvarado, A., Ross, E., Lan, Y., Sun, Y., Castro-Vazquez, A., Vega, I.A., Heras, H., Ituarte, S., Van Bocxlaer, B., Hayes, K.A., Cowie, R.H., Zhao, Z., Zhang, Y., Qian, P.-Y.Y., Qiu, J.-W.W., Alvarado, A.S., Ross, E., Lan, Y., Sun, Y., Castro-Vazquez, A., Vega, I.A., Heras, H., Ituarte, S., Van Bocxlaer, B., Hayes, K.A., Cowie, R.H., Zhao, Z., Zhang, Y., Qian, P.-Y.Y., Qiu, J.-W.W., 2019. Signatures of Divergence, Invasiveness, and Terrestrialization Revealed by Four Apple Snail Genomes. Mol. Biol. Evol. 36, 1507–1520. https://doi.org/10.1093/molbev/msz084 

Swarts, D.C., Makarova, K., Wang, Y., Nakanishi, K., Ketting, R.F., Koonin, E. V., Patel, D.J., Van Der Oost, J., 2014. The evolutionary journey of Argonaute proteins. Nat. Struct. Mol. Biol. 21, 743–753. https://doi.org/10.1038/nsmb.2879 

Tamura, K., Dudley, J., Nei, M., Kumar, S., 2007. MEGA4: Molecular Evolutionary Genetics Analysis (MEGA) software version 4.0. Mol. Biol. Evol. 24, 1596–1599. https://doi.org/10.1093/molbev/msm092 

Tanzer, A., Stadler, P.F., 2004. Molecular evolution of a microRNA cluster. J. Mol. Biol. 339, 327–335. https://doi.org/10.1016/j.jmb.2004.03.065 

Tsai, H.C., Liu, Y.C., Kunin, C.M., Lee, S.S.J., Chen, Y.S., Lin, H.H., Tsai, T.H., Lin, W.R., Huang, C.K., Yen, M.Y., Yen, C.M., 2001. Eosinophilic meningitis caused by Angiostrongylus cantonensis: Report of 17 cases. Am. J. Med. 111, 109–114. https://doi.org/10.1016/S0002-9343(01)00766-5 

Van Bocxlaer, B., 2017. Hierarchical structure of ecological and non-ecological processes of differentiation shaped ongoing gastropod radiation in the malawi basin. Proc. R. Soc. B Biol. Sci. 284. https://doi.org/10.1098/rspb.2017.1494 

Voolstra, C.R., GIGA Community of Scientists (COS), Wörheide, G., Lopez, J. V., 2017. Corrigendum to: Advancing genomics through the Global Invertebrate Genomics Alliance (GIGA). Invertebr. Syst. 31, 231. https://doi.org/10.1071/is16059_co 

Walker, S.E., Spencer, G.E., Necakov, A., Carlone, R.L., 2018. Identification and characterization of microRNAs during retinoic acid-induced regeneration of a molluscan central nervous system. Int. J. Mol. Sci. 19. https://doi.org/10.3390/ijms19092741 

Wang, J., Lu, X., Zhang, J., Ouyang, Y., Qin, Z., Zhao, B., 2020. Using golden apple snail to mitigate its invasion and improve soil quality: a biocontrol approach. Environ. Sci. Pollut. Res. 27, 14903–14914. https://doi.org/10.1007/s11356-020-07998-9 

Wang, Q.P., Lai, D.H., Zhu, X.Q., Chen, X.G., Lun, Z.R., 2008. Human angiostrongyliasis. Lancet Infect. Dis. 8, 621–630. https://doi.org/10.1016/S1473-3099(08)70229-9 

Wang, Q.P., Wu, Z.D., Wei, J., Owen, R.L., Lun, Z.R., 2012. Human angiostrongylus cantonensis: An update. Eur. J. Clin. Microbiol. Infect. Dis. 31, 389–395. https://doi.org/10.1007/s10096-011-1328-5 

Yang, T.B., Wu, Z.D., Lun, Z.R., 2013. The apple snail Pomacea canaliculata, a novel vector of the rat lungworm, Angiostrongylus cantonensis: its introduction, spread, and control in China. Hawaii. J. Med. Public Health 72, 23–25. 

Yang, W., Steitz, T.A., 1995. Recombining the structures of HIV integrase, RuvC and RNase H. Structure 3, 131–134. https://doi.org/10.1016/S0969-2126(01)00142-3 

Yusa, Y. and T.W., 1999. Impact of the introduction of Apple Snails and their control in Japan. Naga, ICLARM Q. 22, 9–13. 

Yusa, Y., Wada, T., Takahashi, S., 2006. Effects of dormant duration, body size, self-burial and water condition on the long-term survival of the apple snail, Pomacea canaliculata (Gastropoda: Ampullariidae). Appl. Entomol. Zool. 41, 627–632. https://doi.org/10.1303/aez.2006.627

Zhang, Y.F., Zhang, R., Su, B., 2009. Diversity and evolution of MicroRNA gene clusters. Sci. China, Ser. C Life Sci. 52, 261–266. https://doi.org/10.1007/s11427-009-0032-5 

Zhang, Y.Q., Chen, D.L., Tian, H.F., Zhang, B.H., Wen, J.F., 2009. Genome-wide computational identification of microRNAs and their targets in the deep-branching eukaryote Giardia lamblia. Comput. Biol. Chem. 33, 391–396. https://doi.org/10.1016/j.compbiolchem.2009.07.013 

Zhou, M., Wang, Q., Sun, J., Li, X., Xu, L., Yang, H., Shi, H., Ning, S., Chen, L., Li, Y., He, T., Zheng, Y., 2009. In silico detection and characteristics of novel microRNA genes in the Equus caballus genome using an integrated ab initio and comparative genomic approach. Genomics 94, 125–131. https://doi.org/10.1016/j.ygeno.2009.04.006


1Federal University of Uberlândia, Laboratory of Bioinformatics and Molecular Analysis, Patos de Minas, MG, Brazil 
*These authors contributed equally to this work.  
**Corresponding Author: Matheus de Souza Gomes, Federal University of Uberlandia, Institute of Genetics and Biochemistry, Laboratory of Bioinformatics and Molecular Analysis, GetúlioVargas Avenue, 230, Patos de Minas, MG, Brazil.  Phone: +55 34 38233714   matheusgomes@ufu.br.