Non-parasite genome encoded virus-like RNAs reprogram the pathogenicity of human blood flukes

Domain: biology
Published: 2025-12-26
Editorial outcome: minor_revision

Abstract
Non-parasite genome encoded RNAs (ngRNAs) are often overlooked in their roles, especially their functions in shaping the biology of organisms. Here we show four non-parasite genome encoded virus-like RNAs that present in Schistosoma japonicum, a human blood fluke responsible for significant morbidity and mortality. These ngRNAs are predominantly localized in the germ cells that play crucial roles in oogenesis and egg production, pivotal processes contributing to host pathology. The ngRNA encode RNA-dependent RNA polymerase (RdRp) that synthesize complementary RNA strands, governing embryo viability and egg production. Notably, similar virus-like ngRNAs are also found in the planarian Dugesia japonica, where RdRp sustains worm survival. In this work, we integrate metatranscriptomics, single-cell RNA sequencing, RNA interference and various molecular methods to comprehensively identify ngRNAs and then characterize their functions. We further demonstrate the mechanism by which ngRNAs shape critical biological traits in parasitic helminths, providing novel insights into potential therapeutic targets for schistosomiasis control and advancing our understanding of ngRNA functions in certain organisms.

Introduction
Schistosomiasis is a zoonotic disease caused by schistosome parasites, a serious zoonotic disease with global distribution1. Currently, the primary treatment for schistosomiasis relies heavily on a single drug, praziquantel. However, the long-term usage of a single-drug strategy raises concerns about the development of drug resistance, necessitating further research into alternative treatment approaches. Deeply investigating biological characteristics of schistosomes and identifying novel and effective drug targets is essential. Previous studies among unicellular protozoan parasites such as Trichomonas vaginali2,3, Giardia lamblia4, Leptomonas spp5,6, and Leishmania spp7–9. have demonstrated that these organisms are associated with certain viruses and/or virus-like RNA segments. More importantly, these exogenous RNAs play critical roles in the pathogenicity of pathogens and the immune responses of infected hosts10–13, suggesting that the targeting of functional virus-like RNAs could become an effective strategy against certain diseases.

In this work, we find four non-parasite genome encoded virus-like RNAs (ngRNAs) associated with Schistosoma japonicum. Through in-depth analyzes employing in situ hybridization and single-cell RNA sequencing (scRNA-seq), we demonstrate that these ngRNAs are predominantly localized in the germ cells. Functional analyses reveal that the ngRNAs play crucial roles in oogenesis and egg production, which are pivotal processes contributing to host pathology and disease dissemination. Additionally, we find that ngRNAs encoding RNA-dependent RNA polymerases (RdRp) can synthesize complementary ngRNAs, which are involved in egg production and embryo viability. Similar virus-like ngRNAs are also found in the planarian Dugesia japonica, where RdRp ngRNA contributes to worm survival. Our results indicate the key biological roles of these ngRNAs in certain organisms, which may represent unique targets for controlling certain diseases caused by some pathogens.

Results
Identification of ngRNAs in S. japonicum
We conducted high-throughput sequencing to analyze the libraries generated from total RNA isolated from S. japonicum (Fig. 1a). Briefly, the parasites were collected from the rabbits infected with S. japonicum cercariae at 16-, 22-, and 28-days post-infection (dpi), respectively. Then, the parasites at different timepoints were subject to independently lysis and treated with DNase and RNase. Subsequently, total RNAs were isolated, and libraries were prepared for sequencing. After removing reads mapped to S. japonicum genomes (NCBI accession: GCA_006368765.1 and GCA_021461655.1) as well as the rabbit host, assembly of the remaining reads yielded 8355, 7246, 7069 contigs, respectively (Fig. 1b). BLASTx analyses indicated that the majority of these assemble sequences are related to four contigs (Fig. 1c). Using the rapid amplification of cDNA end (RACE) technique, we successfully amplified full length sequences for each contig (Supplementary Fig. 1). Further analyses revealed that two long contigs exhibited high similarity to the RNA-dependent RNA polymerase protein (RdRp) of Bunyaviricetes, while two shorter contigs showed high homology to the nucleocapsid protein of Phenuiviridae in Bunyaviricetes. These findings led us to designate these four contigs as rp1/2 and np1/2.

To exclude the possibility these RNAs being detected under opportunistic conditions, we shifted our experimental model from rabbits to mice. Mice were infected with S. japonicum cercariae, and the resulting parasites were collected at 16, 22, and 28 dpi, respectively. These samples were pooled for metatranscriptomic analyses. The results confirmed the consistent presence of these RNAs in the S. japonicum-infected mouse samples (Fig. 1d). Additionally, direct metatranscriptomic analyses of the parasite lysates further corroborated their presence in schistosomes (Fig. 1e). Furthermore, RT-PCR analyses indicated that these RNAs were absent in the livers of uninfected mice and intermediate host Oncomelania hupensis lacking S. japonicum infection (Fig. 1f). To further validate our findings, we conducted a comprehensive search of Schistosoma RNA-seq data deposited in the NCBI database. Among 177 available datasets for various species and geographic regions of Schistosoma, the results indicated that these RNAs are likely specific for S. japonicum (Fig. 1g and Supplementary Fig. 2) (Supplementary data 1). Despite efforts to expand our understanding by exploring potential presence in other Schistosoma species like S. mansoni and S. haematobium, the rules governing biological safety and species availability leave uncertainty regarding the presence of similar RNA in these genera.

Through bioinformatical analyses of metatranscriptomic results (Fig. 1b, d), we removed reads mapped to both S. japonicum and host genomes. The remaining contigs were not theoretically encoded by the parasitic genome, as confirmed by BLAST analyses showing no coding regions for these RNAs in the S. japonicum genome. To further validate the results, we used genomic DNA or cDNA synthesized from reverse transcript of total RNA (both treated and untreated with DNase or RNase) as templates for PCR amplification and also followed sequencing confirmation. These results indicated that these RNAs were not encoded by the S. japonicum genome (Fig. 1h). Therefore, we named them as non-parasite genome encoded RNAs (ngRNAs). Phylogenetic analyses revealed that either RP1 and RP2 or NP1 and NP2 can cluster together, belonging to the Phenuiviridae family (Fig. 1i, j). Generally, they also showed to be closely related to Trichobiphenuilivirus (Fig. 1i, j). Protein domain and motif analyses indicated that RP1/2 contained an RNA-dependent RNA polymerase motif related to Bunyaviricetes and NP1/2 contained a Phenuiviridae nucleocapsid protein motif (Fig. 1k). Additionally, the NP1 shows 51.02% identity to NP2 while RP1 displayed 57.84% identity to RP2 (Supplementary Fig. 3).

ngRNAs are enriched in the sexual organs of S. japonicum
Schistosomes undergo multiple developmental stages characterized by distinct morphology. We further evaluated the abundance of these ngRNAs across several key developmental stages, including eggs, schistosomula (7d, 14d), and adult worms (21d, 28d). The results indicated that the ngRNAs were present in all investigated stages, with the highest abundance in eggs and adult worms (Fig. 2a). This finding aligns with analyses of previous RNA-seq datasets from different development stages14, which confirmed the consistently high abundance of ngRNAs in adult worms (Supplementary Fig. 4a). Given the crucial role of eggs produced by adult worms in disease transmission and host pathology, we conducted Whole mount in situ hybridization (WISH) to localize these ngRNAs within adult females (AF) and males (AM). Since schistosome is diecious, the distribution patterns were analyzed separately for males and females (Fig. 2b). Our findings revealed that the majority of detected ngRNAs were localized in the gonads, including ovaries and testes (Fig. 2b). RT-qPCR analyses of the isolated testes and ovaries from adult worms exhibited consistent results (Fig. 2c).

To further characterize ngRNA roles, we performed scRNA-seq on males and females (16 and 26 dpi). Since S. japonicum genome lacks the coding regions for these ngRNAs, we manually supplemented the genome with the corresponding ngRNA sequences prior to scRNA-seq analyses. After removing potential doublets and low-quality cells, we obtained 12,076 and 8964 cell populations for AM and AF, respectively. Unsupervised clustering of all integrated datasets resulted in a total of 47 molecularly distinct cell populations (Fig. 2d and Supplementary Fig. 4b, c). This included clusters enriched with specific developmental cell types, such as rp1/2 and np1/2, which showed significant abundance in mature oocytes (female late stage), GSC progeny, neoblast, S1, S1 progeny, and vitellocyte cells in females; and spermatocyte (male late stage), neoblast, GSC, GSC progeny, and tegument progenitors in males (Fig. 2e, f).

ngRNAs are essential to maintain ovarian development and egg production
Interestingly, these ngRNAs were predominantly localized within the gonads of schistosomes, particularly in AF that contribute to egg production. Consequently, we mainly focused on AF in the following study. Next, we employed RNA interference (RNAi) to investigate their functional roles. Optimal siRNA was carefully screened and optimized to achieve maximal knockdown efficiency (Fig. 3a and Supplementary Fig. 5). Following siRNA treatment, significant reductions of ngRNAs were observed (Fig. 3b) alongside severe morphological defects of oocytes and impaired oogenesis (Fig. 3c). Morphological defects were further corroborated by Z-stack scanning, which vividly demonstrated ovarian abnormalities in treated individuals (Supplementary Movie 1–5).

To assess the impact on egg production, we monitored egg production at various time points of post-siRNA treatments (24, 48, 72, and 96 h). The results indicated a significantly decreased number of egg production following ngRNA suppression, especially for rp1 ngRNA (Fig. 3d). Subsequent RT-qPCR analysis of the eggs also found significantly reduced ngRNA for rp1/2 and np1 (Fig. 3e). Notably, acridine orange staining revealed substantial reductions in egg viability following ngRNA suppression (Fig. 3f, g). The finding was further validated through cell viability assays using dissociated eggs (Fig. 3g). To further corroborate the results, we conducted in silico docking studies of rp1/2 encoding RdRp proteins from commercially available inhibitors. Among these, Favipiravir, an antiviral agent targeting influenza viruses, demonstrated potential inhibition of activities towards RP1/2 RdRp protein (Supplementary Fig. 6). In vivo validation showed that treatment with Favipiravir significantly reduced both egg production and the viability of eggs in the worms compared to untreated controls (Fig. 3h, i).

ngRNA encoding proteins interacting S. japonicum proteins regulate egg production
To investigate the regulatory roles of ngRNAs in egg production, we employed EdU staining to assess cell proliferation in ovaries and Fast Blue BB staining to evaluate vitellaria morphology following rp1 and np1 suppression. The results revealed that suppression of rp1 and np1 led to decreased cell proliferation in ovaries (Fig. 4a), accompanied by reduced Fast Blue BB staining intensity in the vitellarium, suggesting a potential role for these ngRNAs in governing cell proliferation in the gonad. RNA-seq analysis of rp1- and np1-inhibited parasites revealed 56 and 27 down-regulated genes, respectively, of which 5 genes were shown in both (Fig. 4b, c and Supplementary data 2–3), including tyrosinase (EWB00_006115), a known regulator of egg production15. These findings were further corroborated by scRNA-seq data, which demonstrated co-expression of these target genes with RP1-encoded RdRp in specific cell clusters such as S1/S1 progeny, vitellocytes, and female late-stage cells (Fig. 4c).

To elucidate the molecular basis underlying these regulatory effects, we utilized the Yeast two-hybrid method to identify the interaction partners of RP1 RdRp. The result indicated RP1 may interact with ribosomal protein LP1 (EWB00_007383), Lin-10 family protein (EWB00_001962), proteasome beta subunit, Peptidase C1A subfamily and Ribosomal protein S5a (Fig. 4d and Supplementary Fig. 7). Notably, LP1 and Lin-10 have been previously implicated in regulating development in schistosomes and Caenorhabditis elegans16,17. In addition, ribosomal protein RPLP1 is a known host factor in flavivirus infection18. Single cell RNA-seq results further supported the co-expressions of these interactors with RP1 in specific cell clusters, such as vitellocytes and female late-stage cells (Fig. 4e). To integrate these findings into a coherent molecular framework, we constructed an interaction map based on the STRING database and our Yeast two-hybrid results (Fig. 4f). This map revealed protein–protein interactions mediated by RP1, including tyrosinase, RpS5a, which are known to play roles in egg production within schistosomes15 and Drosophila melanogaster19. Next, we verified the functional roles of these identified interactors using RNAi. The results indicated that suppression of Lin-10 and Lp-1 (Fig. 4g) led to abnormal morphology of vitellocytes (Fig. 4h) and reduced number of egg production (Fig. 4i). Collectively, these results indicated the molecular mechanisms by which RP1 RdRp regulates egg production and ovarian development.

ngRNA encoding RdRp involves the synthesis of complementary ngRNAs
Since these ngRNAs are consistently present in multiple stages of schistosomes, we want to understand the mechanisms of ngRNA replication during cell division. RdRp, an RNA-dependent RNA polymerase found in viral replication machinery, is implied to play a role in the replication of ngRNAs in schistosomes20. To validate this hypothesis, we designed the experiments to determine whether RdRp can generate complementary (positive-sense) RNA and double-stranded RNAs (dsRNAs), which are essential for producing progeny RNA copies of these ngRNAs. We used a strand-specific RT-PCR approach to specifically detect the complementary RNAs (Supplementary Table 1). The PCR results confirmed the presence of complementary RNAs for these ngRNAs (Fig. 5a), suggesting RdRp is capable of synthesizing the complementary strand necessary for ngRNA replication.

To further assess the functional role of RP1 RdRp in ngRNA replication, we suppressed rp1 by siRNA and observed significantly reduced levels of other ngRNAs (Fig. 5b) (Supplementary data 2–3), indicating that the RdRp is critical for maintaining the abundance of ngRNAs. Next, we determined whether that RP1 RdRp is responsible of ngRNA replication and its activity can be inhibited by Favipiravir. Upon co-transfection of recombinant plasmids for expressing RP1 RdRp domain and reversed complementary Gaussia luciferase reporter gene embedding into 5′UTR and 3′UTR of rp1 or reversed complementary np1 into 293T cells, luciferase assay and RT-qPCR analyses indicated that RP1 RdRp can drive the expression of either Gaussia reporter gene or np1 gene (Fig. 5c, d). When applying Favipiravir in transfected cells, their expressions were significantly suppressed (Fig. 5c, d). Recognizing that dsRNA often serves as an intermediate product during RNA replication, we further performed J2 antibody-based Dot blot analyses to detect its presence. The results indicated that S. japonicum not only contains dsRNA but also shows a reduction in dsRNA levels upon suppression of Rp1 RdRp (Fig. 5e). Collectively, these findings suggest that RP1 RdRp has activity involving ngRNA replication. To ensure that our observations were biologically relevant, we performed Whole-mount immunofluorescence using J2 antibodies in the ovary of female worms. Control animals exhibited robust signals for dsRNA (Fig. 5f), while rp1 RdRp suppressed worms showed undetectable signal intensity, further supporting the important role of RP1 RdRp in maintaining ovarian development. Moreover, scRNA-seq data revealed that rp1/2 and np1/2 are primarily co-expressed in germline cells and vitellocyte (Fig. 5g). These results suggest that RdRp-mediated replication of these ngRNAs is highly active in cell clusters related to egg production, aligning with our experimental findings.

Physical characterization of ngRNA formation in S. japonicum
Since the ngRNAs can be successfully isolated even after DNase and RNase treatment, we hypothesized that these ngRNAs may be physically protected from ribozyme degradation in a micro-sized particle encapsulated form. To this end, we used sucrose density gradient centrifugation to separate the components of worm lysates. The most abundant fraction was collected for RNA isolation and quantification by RT-qPCR, indicating. Electron microscopy (EM) analysis of the fraction with the highest abundance of ngRNAs revealed significant enrichment of micro-sized particles (Fig. 6a), suggesting that these ngRNAs may be present in a particle-like form.

To corroborate this result, transmission electron microscopy (TEM) was employed to examine thin sections from the schistosome ovaries. Given the high enrichment of ngRNAs in this tissue region, TEM analysis revealed that the micro-sized particles were indeed abundant in the ovary-associated areas (Fig. 6b). To identify whether these particles are associated with ngRNAs, we cloned the full-length open reading frame (ORF) encoding np1 into the pET-28 vector and expressed the recombinant protein. The purified protein was then used to generate a polyclonal antiserum (Supplementary Fig. 8). When applying the antibodies in an immune TEM assay, some of the micro-sized particles were clearly recognized (Fig. 6c–e), further confirming that these particles are indeed associated with ngRNAs.

Targeting rp1 ngRNA in animals infected with S. japonicum reduces egg production and host pathology
Given the absence of ngRNA coding genes in mammalian genomes21 and certain mammals can serve as definitive hosts for disease dissemination, targeting these ngRNAs may not only represent a unique strategy for disease control with minimal side effects, but also mitigate host pathology. Therefore, we administered rp1 siRNA into the mice infected with S. japonicum to access the effect of ngRNA inhibition on egg production and the pathology of mice livers (Fig. 7a). Upon five injections, at 28 days of post-infection, worms were perfused from the mice. There was a statistically significant reduction in worm recovery [control (RNAi)=86% versus rp1 (RNAi)=50%, P<0.0001, Student’s t test] (Fig. 7b). Additionally, total RNA was extracted from the collected worms, and the transcript levels of rp1 were analyzed by RT-qPCR. The result indicated that the transcript level of rp1 was significantly decreased in remaining worms collected from the mice administered with rp1 siRNA (Fig. 7c). Moreover, mice receiving rp1 siRNA (RNAi) targeting parasites exhibited less eggs and non-significant granuloma in livers, in contrast to the abundance of egg-induced granulomata in livers of control mice (Fig. 7d, f). Subsequently, hatching experiments analyzing the eggs collected from mice livers revealed a significantly lower number of viable eggs [control (RNAi)=267.8 versus rp1 (RNAi)=37.17, P<0.0001, Student’s t test] (Fig. 7e). Overall, these results clearly indicated that targeting rp1 ngRNA in infected animals could reduce egg production and egg viability as well as host pathology, implying the alleviation of disease dissemination.

RdRp ngRNA is also present in planarian Dugesia japonica that contributes to worm survival
We investigated whether similar features of ngRNAs are also present in other flatworms, such as planarian D. japonica. Additionally, recent studies indicated that RNA viruses are present in Schmidtea mediterranea22,23. Through metatranscriptomic and bioinformatic analyses of RNA libraries from D. japonica, we identified that one ngRNA also encoding RdRp in D. japonica (Fig. 8a), named as rp1 RdRp. Similarly, rp1 ngRNA is not encoded by the D. japonica genome (Fig. 8b). Localization studies by WISH and FISH confirmed that the rp1 ngRNA is intracellularly localized within the worms (Fig. 8c, e). Furthermore, treatment of D. japonica with dsRNA targeting rp1 RdRp resulted in the reduced expression of target ngRNA (Fig. 8d) and also led to compromised epidermis integrity (Fig. 8f) and increased worm mortality (Fig. 8g), highlighting its critical role in maintaining worm survival.

Methods
Animal and parasites
BALB/c mice (6–8 weeks, male) and New Zealand Rabbits (male, three months) were purchased from the Shanghai Experimental Animal Centre of Chinese Academy of Sciences. Mice were housed under specific-pathogen-free, climate-controlled conditions on a 12-h light/dark cycle, with free access to standard chow and water. S. japonicum cercariae were provided by the National Institute of Parasitic Diseases of Chinese Center for Disease Control and Prevention, China. Each mouse and rabbit were percutaneously infected with approximately 200 and 2000 cercariae, respectively. All animal experiments were approved by the Institutional Animal Care and Use Committee (IACUC) of the Shanghai Veterinary Research Institute of Chinese Academy of Agricultural Sciences (Permit number: SV-20170929-01). Dugesia japonica (D. japonica) was obtained from Boshan fountain of China. Animals were reared in spring water at 20 ± 1 °C and fed raw chicken liver once or twice a week. Asexual planarians (a length between 2 and 4 mm) were starved for 1 week before any experiments.

Metatranscriptomics analyses
For metatranscriptomics analyses of parasites collected from rabbits, schistosomes [16, 22, 28 days post-infection (dpi), 1 g in weight for each stage] were homogenized with steel balls in 2 mL phosphate-buffered saline (PBS). The lysates were subjected to three freeze-thaw cycles on dry ice. After centrifugation (4°C, 10 min, 15,000 × g), the supernatants were collected carefully. To remove cell-sized particles, 300 mL supernatant was filtered through a 0.45 μm filter (Millipore), of which the filtrates were digested using a mixture of DNase (Turbo DNase from Thermo Fisher Scientific, Baseline-ZERO from Epicenter, and benzonase from Novagen) and RNase (Thermo Fisher Scientific, USA) to remove free nucleic acid at 37°C for 1.5 h33. Nucleic acids were isolated using the QIAamp RNA Mini Kit (Qiagen, Germany) according to the manufacturer’s protocols. Three libraries were generated using Nextera XT DNA Sample Preparation Kit (Illumina, USA), and then subjected to paired-end (PE) sequenced on the MiSeq Illumina platform (Majorbio, Shanghai)34. For metatranscriptomic analyses of parasites from mice, worms were collected from mice at 16 dpi, 22 dpi, and 26 dpi, respectively. The parasites were pooled and then subjected to metatranscriptomic analyses as described above.

For metatranscriptomic analyses of parasite lysates, adult schistosomes (28 dpi) were amputated and then centrifuged to remove the tissue and cell debris. The supernatants were subject to total RNA isolation for library preparation. The libraries were constructed using the Illumina TruSeq Stranded mRNA LT Sample Prep Kit (Illumina, USA). After screening by an Agilent Bioanalyzer 2100 (Agilent Technologies, USA), the selected libraries were quantified with a Promega QuantiFluor RNA System (Promega, USA), and the samples with a concentration >2 nM were further subjected to PE sequencing. Following serial dilution and mixing, the library was denatured into single strands for sequencing. The Illumina HiSeq X-Ten platform (Illumina, USA) was used for sequencing, along with the PE150 strategy.

Upon the identification of the species of planarian D. japonica by amplifying the Cox gene, worms were carefully washed using PBS (Thermo Fisher Scientific) at least three times, and then were subject to total RNA isolation using TRIzol (Thermo Fisher Scientific). The quality of isolated RNAs was evaluated by Agilent 2100. The libraries were prepared using the Illumina Stranded Total RNA Prep Kit (Illumina) and then sequenced as described above.

Bioinformatical analyses
Sequencing reads were first adapter- and quality-trimmed using the Trim-galore program (v0.6.10) (-phred33 -length 100 -stringency 3 -paired, https://zenodo.org/records/7598955). The remaining reads were mapped to the S. japonicum reference genomes (NCBI accession: GCA_006368765.1 and GCA_021461655.1), rabbit (NCBI accession: GCA_964237555.1), and mice (NCBI accession: GCF_000001635.27) genomes using the Bowtie2 program (v2.5.0)35. The remaining reads were assembled de novo using MEGAHIT (v1.2.9) with default parameters (-min-contig-len 200) to generate the longest contigs and singlets36. Next, the contigs and singlet reads were analyzed by DIAMOND (v2.1.12) with default parameters. Each sequence was aligned against the NCBI virus reference proteome (ftp://ftp.ncbi.nih.gov/refseq/release/viral/) in addition to viral protein sequences from the NCBI nr Fasta file (based on annotation taxonomy in Virus Kingdom)34,37, and an E-value threshold of <10⁻⁵ was set to screen potential viral sequences. To further eliminate potential false-positive viral sequences, the viral sequences obtained through the above steps were re-aligned against the NCBI Internal Non-Protein Non-Redundant (NVNR) Protein Database. As to planarian D. japonica, the remaining reads were mapped to the planarian reference genome (GCA_001938525.1), and the following analyses were performed for the assembly as described above.

As to metatranscriptomic analyses, the raw sequencing data were trimmed to remove adapters and filtered to remove low-quality reads. The mRNA transcripts were assembled de novo using Trinity (v2.5.1) software and then to obtain unigenes using CD-HIT software (version 4.7). The obtained sequences were then aligned with the bacterial, fungal, and viral sequences from the NCBI nr database using the BLASTX program. For specie annotation, the MEGAN software (v.6.13.1) was employed. QIIME (v.1.7.0) was used to calculate the relative abundance at each classification level. The unigenesets were also aligned with the KEGG and CAZy databases for functional annotation. The threshold value for gene identification was an e-value < 10⁻⁵.

The homologous reference sequences from five families within Bunyavirales were downloaded from the NCBI database to analyze phylogenetic relationships. Amino acid sequences of RP/NP were aligned using MUSCLE in MEGA-X38 and trimmed. IQ-TREE (v2.3.6)39 was used to construct phylogenetic trees with 1000 bootstrap replicates (-bb 1000) and the ModelFinder function (-m MFP). Phylogenetic trees were visualized and edited using Interactive Tree Of Life (iTOL)40.

Genomic DNA isolation, total RNA isolation, DNase/RNase treatment, and PCR analysis
S. japonicum genomic DNA (adult worms, 28 dpi), murine liver genomic DNA, or planarian D. japonica genomic DNA was isolated using an AxyPrep Mutisource Genomic DNA Miniprep Kit (Axygen, USA) according to the manufacturer’s instructions. The isolated DNA (500 µg) was treated with DNase (1.5 µL) or RNase (1.5 µL) at 37°C for 30 min. Then, ethylenediaminetetraacetic acid (EDTA) was added to a final concentration of 15 mM and heated at 75°C for 10 min to inactivate DNase. Next, PCR was performed in a 20 μL reaction mixture containing 1 μL treated DNA, 10 μL 2× SYBR Primer Ex Taq II (TaKaRa, China), 8.5 μL H2O, and 0.5 μL primers (10 μM) using the following thermal cycling profile: 95°C for 5 min, followed by 40 cycles of amplification (95°C for 5 s, 57°C for 30 s, 72°C for 20 s). The primer sequences that synthesized by Invitrogen were listed in Supplementary Table 1. The PCR products were separated on 2% agarose gels and visualized by staining with GoldView (Yeasen, China).

Total RNAs from S. japonicum, O. hupensis, or murine livers were isolated using TRIzol (Thermo Fisher Scientific) according to the manufacturer’s instructions. Isolated RNA (500 µg) was subject to treatment with DNase (1.5 µL) or RNase (1.5 µL) as described above. Treated RNA was further purified using phenol-chloroform extraction, and reverse transcription reactions were performed in a total volume of 10 μL, containing 6 μL total RNA (100 ng), 2 μL 5× Prime-Script Buffer (TaKaRa), 0.5 μL PrimeScript RT Enzyme Mix I (TaKaRa), 1 μL of RNase inhibitor (TaKaRa), and 0.5 μL random primers (10 μM) or specific strand primer (Supplementary Table 1). PCR and subsequent agarose gel analysis were performed as described above. Similar experiments were performed for planarian D. japonica using the following primers (forward: 5′CGGCACTAGTCAATCAACCCTCT3′, reverse: 5′TGTTTGGTTTACTTGTTGGCTCATTC3′, synthesized by Invitrogen).

5′/3′ RACE
To obtain full-length ngRNA sequences, 5′- and 3′-RACE assays were conducted. Total RNA was extracted from adult S. japonicum worms using TRIzol reagent (Thermo Fisher Scientific) following the manufacturer’s protocol. Gene-specific primers for 5′- and 3′-RACE were designed based on previously assembled ngRNA sequences, respectively (Supplementary Table 1). 5′- and 3′-RACE reactions were performed using the SMARTer RACE cDNA Amplification Kit (Clontech, USA) in accordance with the official instructions. For 5′-RACE outer PCR, the thermal program consisted the following steps: (i) an initial five cycles of 94°C for 5 s, 65°C for 30 s, and 72°C for 3 min; (ii) followed by five cycles of 94°C for 5 s, 60°C for 10 s, and 72°C for 2 min; then (iii) 30 cycles of 94°C for 5 s, 55°C for 10 s, and 72°C for 2 min; and (iv) a final extension at 72°C for 10 min. A diluted aliquot of the outer PCR product was subsequently used as a template for nested PCR, which was run under identical thermal conditions. For 3′-RACE outer PCR, the program was as follows: (i) amplification began with an initial denaturation at 94°C for 3 min; (ii) followed by 20 cycles of 94°C for 30 s (denaturation), 53°C for 20 s (annealing), and iii) 72°C for 3 min (extension), and a final extension at 72°C for 10 min. Similarly, a diluted fraction of the 3′-RACE outer PCR product was subjected to nested PCR using the same conditions. All primers utilized in this study that synthesized by Invitrogen are listed in Supplementary Table 1. Nested PCR products from both 5′- and 3′-RACE reactions were analyzed by 1.2% agarose gel electrophoresis containing ethidium bromide (EtBr). Target bands were excised, purified, and subcloned into the pMD™ 19-T vector (TaKaRa) for Sanger sequencing. Full-length ngRNA sequences were reconstructed by in silico assembly of the 5′- and 3′-RACE-derived sequences.

Analyses of the abundance of ngRNAs in different stages of schistosomes
Total RNA was extracted from schistosomes at various developmental stages (eggs, 7, 14, 21, and 28 dpi) using TRIzol reagent (Thermo Fisher Scientific) following the manufacturer’s protocol. RNA concentrations were quantified with a NanoDrop 1000 spectrophotometer (Thermo Fisher Scientific) for RT-qPCR experiments. Expression levels of each RNA transcript across different stages were determined via RT-qPCR using specific primers (Supplementary Table 1), with the primer targeting S. japonicum nicotinamide adenine dinucleotide dehydrogenase (NADH) gene serving as an internal reference for data normalization. Briefly, RT was carried out in a 10 μL final volume comprising 6 μL of total RNA (400 ng), 2 μL of 5× PrimeScript Buffer (TaKaRa), 0.5 μL of PrimeScript RT Enzyme Mix I (TaKaRa), 1 μL of RNase inhibitor (TaKaRa), and 0.5 μL of 10 μM random primers. RT-qPCR was conducted using SYBR Premix Ex Taq (TaKaRa) in a 15 μL reaction mixture consisting of 1 μL cDNA, 7.5 μL of 2× SYBR Premix Ex Taq II (TaKaRa), 6 μL of nuclease-free H2O, and 0.5 μL of 10 μM primers. Amplification was performed on a Master Cycler Ep Realplex PCR system (Eppendorf, Germany) with the following thermal program: (i) an initial denaturation step at 95°C for 5 min; (ii) followed by 40 cycles of 95°C for 5 s, 60°C for 30 s, and (iii) 72°C for 20 s. Relative mRNA expression levels were calculated using the 2⁻ΔCt method41.

Gonads (testes and ovaries) from adult S. japonicum were dissected following established protocols42,43. Then, total RNA was subsequently extracted from these tissues using TRIzol reagent (Thermo Fisher Scientific) according to the manufacturer’s instructions. The total RNA extracted from testes or ovaries, as well as the total RNAs isolated from AM and AF schistosomes, were reverse transcribed, and qPCR analysis was performed as described above.

Whole mount of in situ hybridization analyses of S. japonicum
Initially, adult males and female schistosomes (28 dpi) were fixed for 15 min in 4% paraformaldehyde (Ted Pella, USA) diluted in PBS containing 0.3% Triton-X100 (hereafter PBSTx). Following fixation, worms were treated with a solution of 50 mM dithiothreitol, 1% NP-40, and 0.1% sodium dodecyl sulfate at 37°C for 5–10 min, then dehydrated through a methanol gradient (100% PBSTx, 50% methanol, 100% methanol) with 5-min incubations per step. Bleaching was performed overnight under light using 6% hydrogen peroxide (diluted from a 30% stock in methanol). Subsequently, worms were rinsed twice with 100% methanol, incubated in 50% methanol for 10 min, and then transferred to PBSTx for an additional 10-min incubation. After permeabilization, worms were rinsed with PBSTx twice. Prehybridization began with incubation for 10 min in hybridization buffer (50% deionized formamide, 1.3× SSC (diluted from 20× SSC from Sigma, pH 5.5 with citric acid), 5 mM EDTA, 50 µg/µL torula yeast RNA, 0.2% Tween, 0.5% CHAPS, and 10% Heparin) mixed with PBSTx at a 1:1 ratio for 10 min. Then, the Worms were subjected to a 10-min wash with hybridization buffer at 54–56°C in a hybridization oven (HBAID). The hybridization buffer was then removed, and the worms were incubated in fresh hybridization buffer for at least 2 h at 54–56°C. Next, the digoxigenin-labeled riboprobe (5 µM) was added to the hybridization buffer, and the worms were incubated with the mixture overnight at 54–56°C. Washes were then performed with solutions preheated to 54–56°C in the following order: prehybridization solution mixed with 2× SSC at a 1:1 ratio, 2× SSC with 0.1% Triton-X100, and 0.2× SSC with 0.1% Triton-X100 in the hybridization oven with agitation. All washes were carried out twice for 30 min each at 54–56°C. After stringency washes, schistosomes were cooled to RT and rinsed twice with maleic acid buffer (100 mM maleic acid, 150 mM NaCl, 0.1% Tween-20, pH 7.5). Blocking was performed with 10% horse serum for 2 h, followed by an overnight incubation at 4°C with anti-DIG-AP-Fab fragments (Merk, catalog#: 11093274910, Lot#: 65752120) diluted 1:2000 in blocking solution (Roche, Switzerland). Next, worms were washed for 2 h in maleic acid buffer (buffer was changed every 20 min). Then, parasites were incubated in alkaline phosphatase (AP) buffer (100 mM Tris [pH 9.5], 100 mM NaCl, 50 mM MgCl2, 0.1% Tween-20 brought up to the final volume with 10% polyvinyl alcohol solution) for 10 min and were then incubated in NBT/BCIP solution for 3 h or overnight. Then, worms were further washed with PBST and post-fixed. Finally, worms were mounted in 80% glycerol under glass coverslips and imaged.

As to probe preparation, cDNA fragments corresponding to each target RNA segment were amplified from cDNAs synthesized from total RNA extracted from adult schistosomes, using the primers listed in Supplementary Table 1. PCR products were purified and then utilized to generate digoxigenin-labeled riboprobes using the MEGAscript Kit (Thermo Fisher Scientific) and digoxigenin-11-UTP (Roche), following the manufacturer’s instructions.

Single-cell RNA sequencing
S. japonicum worms (16, 26 dpi) (100 worms for each gender) were transferred into 15 mL conical tubes and rinsed twice with PBS pre-warmed at 37°C and dissociated into single-cell suspensions in 4 mL of 0.25% trypsin in HBSS for 15 min. Cell suspensions were passed through a nylon mesh (Corning, USA) and immediately doused with 10 mL of PBS supplemented with 1% BSA. The sample was then centrifuged at 500 × g for 5 min at 4°C. Pellets were gently resuspended in PBS with 10% FBS, passed through a 40 μm nylon mesh (Corning, USA), and stained with diamine blue for 30 min, and then were washed and resuspended in PBS with 0.5% BSA and loaded on a SONY SH800S cell sorter. Dead cells were removed using the Dead Cell Removal Kit (Miltenyi Biotec, Germany).

Cell suspension and fluorescent dye AO/PI were mixed at a ratio of 1:1 and incubated for 30 s. The mixture was added to the slide and assessed cell viability using CountStar (Shanghai, China). Viability enrichment was carried out using the Dead Cell Removal Kit (Miltenyi Biotec) according to the manufacturer’s protocol. The cells were resuscitated to a concentration of 700–1200 cells/µl (viability ≥85%) in a final solution of 1× PBS 0.04 BSA before loading on the 10× Genomics Chromium platform. Over 10,000 cells were prepared for scRNA-seq libraries. Chromium Single Cell 3’ Library and Gel Bead Kit V3.1 (10× Genomics, PN1000268) was employed to generate single-cell gel beads in emulsion (GEM). The obtained cells were lysed, and the released RNA was performed for reverse transcription using primers containing polyT, barcode, unique molecular identifiers (UMIs), and read 1 primer sequence in GEMs. Barcoded cDNA was further purified and subsequently amplified by PCR. The adapter was ligated to add the sample index and the read 2 primer sequence.

Upon quality control, the libraries were sequenced using the Illumina Novaseq 6000 platform in a 150 bp pair-ended manner (Berry Genomics Corporation, Beijing, China). Raw reads of fastq files were assembled from the Raw BCL files using Illumina’s bcl2fastq converter. Then, the raw data were firstly processed through primary quality control using the parameters: (i) contain N more than 3; (ii) the proportion of bases with quality value below 5 is more than 20%; (iii) adapter sequence. All the downstream analyses were based on the obtained clean data.

Mapping of single-cell RNA sequencing reads
For reads mapping, a modified reference genome was constructed based on the S. japonicum v2 genome (HuSjv2; from PRJNA520774 WBPS16; https://parasite.wormbase.org), where mitochondrial genes (NCBI Reference Sequence NC_002544.1; naming as “Sj_”) and 4 ngRNA genes (naming as “SjPS”) were manually added. Raw sequencing data were mapped to the constructed reference genome using Cell Ranger v5.0 from 10× Genomics (https://www.10xgenomics.com/support/software/cell-ranger/).

Reference gene information
Orthologues between S. mansoni (v7) and S. japonicum (HuSjv2), protein sequences, and Gene Ontology terms were obtained from WormBase Parasite (https://parasite.wormbase.org/44; using the BioMart tool (access date: 2021-11-25). KEGG mapping was performed using KAAS with the BBH method (https://www.genome.jp/kegg/kaas/; access date: 2021-11-27). Single-cell markers were obtained from published studies for adult S. mansoni (filtered by p_val < 10⁻²⁰)45 and schistosomula46 (filtered by AUC > 0.7). Corresponding S. japonicum gene IDs were generated based on the ortholog information.

Single-cell RNA seq data processing
Cell Ranger filtered matrix files were imported to Seurat (v5.0.1; https://satijalab.org/seurat/) and filtered by nFeature_RNA > 500 & nCount_RNA > 1000 & percent.mt < 5. Each dataset was further processed with doublet detection using scDblFinder (v1.2.0), and only singlets were kept for further analysis. All four datasets were merged using Seurat and integrated using the Harmony Integration method. Briefly, each layer of data was normalized to find the top 2000 variable features, which were scaled (with percentage of mitochondrial expression regressed out) and used for dimensional reduction using the command RunPCA (npcs=30), and IntegrateLayer (method=HarmonyIntegration). Neighborhood graphs were constructed on harmony reduction using dims=1:30, and cell clusters were detected using the default algorithm with resolution 1.3. Afterwards, the layers were joined using the command JoinLayers.

Cell type annotation
To facilitate cell type annotation of Seurat clusters, the above-mentioned S. mansoni single-cell marker genes were used to calculate module scores in the current S. japonicum dataset. This was performed using the Seurat command AddModuleScore using cluster marker as features. Cell type information was further assessed by checking the marker genes, which were calculated using the FindAllMarkers command with test.use=‘roc’. Expression of canonical cell type markers from S. mansoni was also examined. Candidate marker genes were validated using whole amount in situ hybridization, as stated in additional publications.

RNA-seq data processing and differential expression analyses
Parasites treated with rp1 siRNA, np1 siRNA, and control siRNA were subject to total RNA isolation for library preparation. Eight transcriptome libraries (two biological replicates for each treatment or control) were sequenced using the BGISEQ-500 platform. Clean reads were obtained after removing low-quality reads, adapters, as well as any reads containing poly-N from raw reads and mapped to the self-constructed reference genome using STAR_2.6.1d (alignIntronMin 1047). Read counts for genes were summarized using featureCounts v1.6.448). For downstream analysis, the web tool iDEP 1.12 (http://bioinformatics.sdstate.edu/idep/)49 was used using the DESeq250 pipeline. Differentially expressed genes were selected using padj < 0.05 and 1.3-fold change. For functional enrichment of Gene Ontology, only terms with 5–2000 genes in the category “Biological Process” and FDR < 0.05 were selected.

Worm culture and siRNA treatment
S. japonicum were obtained from the infected mice (24–28 dpi), and females were separated and cultured in 12-well flat bottom plates containing 2 mL complete RPMI-1640 medium (Thermo Fisher Scientific) supplemented with 2 g/L glucose, 0.3 g/L L-glutamine, 2.0 g/L NaHCO3, 15% fetal bovine serum (heat inactivated), and 5% penicillin/streptomycin (10,000 U penicillin and 10 mg/streptomycin in 0.9% NaCl). The worms were incubated in a humidified 5% CO2 chamber at 37°C. siRNA duplexes and control siRNAs (3 µg per experiment, chemically synthesized in Shanghai GenePharma, China; Supplementary Table 2) were electroporated (125 V, 20 ms, 1 pulse in 200 µL RPMI 1640 medium) into cultured parasites. Subsequently, worms were then transferred into a 12-well cell culture plate containing 2 mL fresh media. Eggs were collected and counted under an inverted microscope (Olympus, Japan) at the indicated times. Treated worms were collected at 96 h of post electroporation for RT-qPCR analysis and Carmine Red staining. Whole-mount parasites were visualized using confocal microscopy.

Confocal microscopy
At 96 h after electroporation, worms were preserved, stained, and mounted and subjected to morphological observation via a Nikon confocal microscope (Nikon, Japan)51. Briefly, parasites were fixed in a solution containing formalin (10%), alcohol (48%), and glacial acetic acid (2%), stained with Carmine Red, cleared with a 0.5% hydrochloric alcohol solution, and stored as whole mounts. Images were captured with a Nikon CLSI laser confocal microscope (Nikon) under a 488 nm He/Ne laser.

Fast blue BB and EdU staining
Fast Blue BB staining was carried out using the following procedure. Briefly, worms were fixed in 4% formaldehyde in PBSTx for 4 h, then stained for 5 min in freshly prepared, filtered 1% Fast Blue BB solution in PBSTx. Parasites were washed with PBSTx three times (5 min per wash) and mounted on slides with 80% glycerol for microscopic observation.

EdU staining was conducted based on a protocol52. Worms were in vitro labeled with 10 µM EdU (Merck) for 20 h. Subsequently, the parasites were harvested and incubated in 6 M MgCl2 for 1 min to terminate viability. Then, the parasites were fixed in 4% formaldehyde in PBSTx for 4 h, and dehydrated in Methanol. Next, worms were rehydrated in 50% methanol in PBSTx, bleached under bright light for 5–6 h, then treated with 5 µg/ml Proteinase K (Thermo Fisher Scientific) in 1× PBSTx for 30 min. The parasites were treated with 4% formaldehyde in PBSTx for 10 min for post-fixation. Worms were stained with EdU Detection solution (789 µL 1× PBS, 10 µL 100 mM CuSO4, 1 µL 10 mM Azide-fluor 488, 200 µL 500 mM Asorbic Acid) for 30 min in the dark. Parasites were then stained with the DAPI solution (1 µg/ml in PBSTx) for 1 h. The worms were further cleared with 80% glycerol and mounted in Vectashield (Yeasen, China). Images were acquired using a Nikon confocal microscope (Nikon) at 456 and 488 nm.

Egg viability
The viability of schistosome eggs was assessed with the acridine orange staining method, as judged by the intensity of the fluorescence of the stained eggs53. Briefly, the eggs were stained with 0.01% acridine orange solution (Merck) at 37°C and shielded from light for 2 h. The fluorescence of the stained eggs was subsequently recorded and counted under an inverted fluorescence microscope with a CCD camera. Cell viability in eggs was determined using a CellTiter-Lumi Luminescent cell viability assay kit (Beyotime Biotechnology, China) according to the manufacturer’s instructions. Subsequently, the obtained luciferase activity was normalized with protein concentration using the Pierce BCA protein assay kit combined with the Compat-Able protein assay preparation reagent set (Thermo Fisher Scientific).

RT-qPCR analysis of the abundance of ngRNAs in electroporated schistosomes
At 96 h post-electroporation, total RNA was extracted from treated worms using TRIzol (Thermo Fisher Scientific) according to the manufacturer’s instructions. First-strand cDNA was produced from 100 ng RNA using a PrimeScript First Strand cDNA Synthesis Kit (TaKaRa) with random primers. PCR was then conducted using the primers listed in Supplementary Table 1. Gene expression was normalized to S. japonicum NADH as an internal control, and relative expression levels were calculated using 2−ΔCt method41.

Screening of a YTH for analyzing Rp1 RdRp interactors
To identify the interactors of RdRp proteins, we cloned the cDNA encoding the conserved domain of rp1 RdRP (forward: 5′gcatatggccatggaggccgaattcCAGTTCCCCATGATAGAGACTC3′, reverse: 5′gcggccgctgcaggtcgacggatccTCAAGGGTTTTTGCATGCCAACATA3′, synthesized by Invitrogen) and fused it to the binding domain of the pGBKT7 vector to serve as the bait. YTH screening was conducted using a previously prepared S. japonicum cDNA library by the Yeastmaker Yeast Transformation System 2 (TaKaRa). The YTH assay was carried out54. The plasmids were isolated from all of the positive clones for sequencing. For the retransformation validation of the YTH assay, the independently prepared positive clone plasmids and RdRp plasmid were co-transformed into YTH Gold yeast competent cells (TaKaRa), respectively. Yeast transformation was performed using the lithium acetate method at 30°C55. Transformants were first selected on medium lacking Trp and Leu. Positive colonies were then subsequently transferred to medium lacking selection media supplemented with SD/-Trp/-Leu/-His/-Ade with or without X-α-gal.

Inhibition of Rp1 RdRp interactors
Double-stranded RNA (dsRNA) targeting Lin-10 and Lp1 was used to inhibit these Rp1 RdRp interactors. Briefly, PCR was performed using specific primers to generate the templates for in vitro transcription (Supplementary Table 1). The dsRNA was synthesized using the MEGAscript T7 Transcription Kit (Thermo Fisher Scientific). The reaction contained 100 ng of PCR product, 2 µL of 10× reaction buffer, 2 µL of each rNTP, 2 µL of T7 RNA polymerase, and nuclease-free H2O to 20 µL. Ten 26 dpi females were cultured in each well using ABC169 medium and then treated with 30 µg/mL dsRNA. This dsRNA was replenished during medium changes at 1, 3, and 6 days post-treatment. Egg number was monitored daily. Fast Blue BB staining was performed at 6 days post-treatment.

As to RT-qPCR analyses, total RNA was isolated from treated worms using TRIzol reagent (Thermo Fisher Scientific) and reverse-transcribed using the PrimeScript™ RT Reagent Kit (TaKaRa). PCR reactions were performed in 10 µL volumes containing 5 µL ChamQ Universal SYBR qPCR Master Mix (Vazyme, China), 0.2 µL of each 10 µM primer (Supplementary Table 1), 0.5 µL cDNA, and added nuclease-free H2O to 10 µL under the condition at 95°C for 30 s, followed by 40 cycles of 95°C for 10 s and 60°C for 30 s. Data were analyzed as described above.

PCR analysis for determining the complementarity of RNAs
A strand-specific RT-PCR approach was used to detect the complementary strand of ngRNAs56. Total RNA was isolated from adult S. japonicum as described above. RT was performed in a final volume of 10 µL, containing 6 µL total RNA (400 ng), 2 µL 5× Prime-Script Buffer (TaKaRa), 0.5 µL PrimeScript RT Enzyme Mix I (TaKaRa), 1 µL of RNase inhibitor (TaKaRa), and 0.5 µL of specific primers (10 µM) (Supplementary Table 1, the primers were synthesized by Invitrogen) matching the complementary strand of viral RNA. Next, PCR was performed in a 15 µL reaction mixture containing 1 µL cDNA, 7.5 µL 2× SYBR Primer Ex Taq II (TaKaRa), 6 µL H2O, and 0.5 µL primer mixture (10 µM) (Supplementary Table 1) using the following thermal cycling profile: 95°C for 5 min, followed by 40 cycles of amplification (95°C for 5 s, 57°C for 30 s, 72°C for 20 s). The PCR products were analyzed using 2% agarose gels and visualized by staining with GoldView (Yeasen).

Plasmid construction, cell culture, and transfection
The reversed complementary Gaussia luciferase reporter gene embedded into 5′UTR and 3′UTR of Rp1 and reversed complementary np1 flanked by BamHI and HindIII sites were synthesized from GENERAY (GENERAY, China). These fragments were subcloned into the pcDNA3.1(+) vector, respectively. Additionally, the conserved domain of Rp1 RdRP was also subcloned into the pcDNA3.1(+) vector.

293T cells were purchased from the American Type Culture Collection (ATCC, USA; catalog# CRL-3216). 293T cells were maintained in Dulbecco’s Modified Eagle’s Medium (DMEM) (Thermo Fisher Scientific), supplemented with 10% fetal bovine serum (FBS) (Thermo Fisher Scientific) and 1% penicillin/streptomycin (Gibco, MA, USA). Cells were cultured in a humidified incubator at 37°C with 5% CO2. 293T cells were seeded in a 48-well plate after passaging. When cell density reached 70% confluence, transfection was performed. The rp1 RdRp plasmid and either the embedded Gaussia luciferase plasmid or the np1 plasmid at a 1:4 ratio (50 ng: 200 ng) were transfected into the cells and then incubated for 8 h. Meanwhile, the indicated concentrations of Favipiravir were added to the cells. After 24 h of post-transfection and/or treatment, cells were collected for luciferase assay and RT-qPCR analyses as described below.

Preparation of cell lysate, luciferase assay and RT-qPCR analyses
The cells were rinsed with phosphate-buffered saline (PBS, pH 7.4) and lysed using a lysis buffer from Promega Dual-Luciferase Assay (Promega, USA). The lysates were centrifuged at 16,000 g at 4°C for 15 min, and the supernatants were carefully transferred to new tubes. Luciferase activity in cell lysates was measured using the Promega Dual-Luciferase Assay, normalized to protein concentrations in worm lysates using the Pierce BCA protein assay kit and Compat-Able protein assay preparation reagent set (Thermo Fisher Scientific). Partial cells were collected for total RNA isolation. RT-qPCR was performed to determine the expression of np1 using the specific primers (Supplementary Table 1).

Dot blot assay and immunofluorescent analysis
For the dot blot assay, 500 ng dsRNA was dropped on the nylon membrane and then was cross-linked using a UV cross-linking device (1200 mJ/cm2, 3 min). Upon the blocking of the membrane using 5% skimmed milk powder, the membrane was incubated with primary antibody J2 (1:1000 dilution) (Scicons, Hungary; catalog# 10010200, batch# 72-1710) overnight at 4°C. The secondary antibody was used HRP conjugated Goat Anti-Mouse IgG (H+L) (1:5000 dilution) (Servicebio, China; catalog# GB23301, lot# AC231007). Imaging signals were obtained by developing with enhanced chemiluminescence (ECL) reagents (Thermo Fisher Scientific).

Worms were permeabilized in Tris-Triton X-100 buffer with 5% β-mercaptoethanol and incubated for 20 h in a rocker at 37°C. Then the worms were washed and put into collagen IV (1 mg/ml) prepared in PBS for 35 min. Worms were quenched in cold antibody buffer composed of 1× PBS, 0.1% BSA, Triton X-100, and 0.05% sodium azide for 15 min. Worms were washed and incubated with J2 antibody (dilution: 1:1000) overnight at 4°C. After washing, worms were incubated with goat anti-mouse secondary antibody conjugated with Alexa Fluor 568 (dilution: 1:2500; Thermo Fisher Scientific; catalog# A-11004, lot# 3141857) for 2 h. Cell nuclei were counter-stained with 4′,6-diamidino-2-phenylindole (1 µg/ml). Finally, worms were mounted with an antifading mounting medium, and images were acquired with a confocal laser scanning microscope (ZEISS)57.

Structural modeling and drug-protein docking of RdRp
To predict the rp1/2 RdRp structure, we employed the alphafold server (https://alphafoldserver.com/) to predict their structures for the drug-protein docking (DPD) using the default settings. The obtained protein molecule PDB file and the small molecules (Favipiravir and its derivative FTP, in mol2 format) were subjected to the DPD program using AutoDock software (https://autodock.scripps.edu/) following the official protocol. The docking result was further processed by AutoDockTools (https://autodocksuite.scripps.edu/adt/) to obtain a drug-protein complex, which was further visualized by PyMol software (https://pymol.org/). The docking affinity of protein and ligands was quantified by the free binding energy.

RdRp inhibitor treatment
S. japonicum were collected from mice 24–28 days of post-infection, and females were separated and cultured in 12-well flat bottom plates containing 2 mL complete RPMI-1640 medium (Thermo Fisher Scientific) supplemented with 2 g/L glucose, 0.3 g/L L-glutamine, 2.0 g/L NaHCO3, 15% fetal bovine serum (heat inactivated), and 5% penicillin/streptomycin (10,000 U penicillin and 10 mg/streptomycin in 0.9% NaCl) in a humidified 5% CO2 chamber at 37°C. Favipiravir was added to the culture medium at the indicated concentrations. The number of egg production was monitored during the time course. Egg viability was determined as described above.

Expression and purification of np1 in Escherichia coli and identification of recombinant protein
The cloning and expression of np1 were carried out based on the standard molecular procedures. In brief, the ORF of np1 was cloned by PCR-based amplification and the insertion of the fragment into a pET28a(+) expression vector (Novagen, Germany). Restriction recognition sites of EcoRІ and XhoI were introduced by PCR using the specific primers (forward: 5′CCGGAATTCATGTGTAATAAAGCCAGGTC3′; reverse: 5′CCGCTCGAGTTAAACCATTTTCCTCTTAG3′, synthesized by Invitrogen) with template of prepared cDNA. Recombinant np1 was expressed in Escherichia coli BL21 cells transformed with pET28a(+)-np1 plasmids. The recombinant plasmids were induced under 1.0 mM IPTG at 37°C. Recombinant NP1 (rNP1) was purified by a Ni-NTA His*Band Purification Kit according to the manufacturer’s instructions (Novagen). The purified proteins were analyzed by 15.0% sodium dodecyl sulfate–polyacrylamide gel electrophoresis (SDS–PAGE), and the recombinant NP1 was verified by mass spectrometry with ESI-Q-TOF (Applied Biosystems, USA). The results of mass spectrometry were searched against NCBI nr database using MASCOT under default parameters.

Generation of rabbit antiserum
Two rabbits were intramuscularly immunized with 100 µg recombinant NP1 that was emulsified with 206 adjuvants (Sigma, USA), followed by two boosts at a 3-week interval. Sera were collected after 7 days of the last boost and stored at −20°C.

Sucrose ultracentrifugation
Sucrose solutions were prepared (w/w) in PBS. Sucrose gradients were prepared by sequentially adding equal amounts of 30% sucrose solution, 45% sucrose solution and 60% sucrose solution. Schistosomes (28 days old, 0.5 g weight) were lysed in 5 mL PBS, and lysates were then centrifuged at 2000 × g for 5 min. The supernatants were carefully transferred into fresh tubes and then centrifuged at 8000 × g for 5 min. The supernatants were also carefully transferred into fresh tubes and centrifuged at 10,000 × g for 10 min. Supernatants were then filtered using a 0.22-µM filter (EMD Millipore) and further centrifuged at 210,000 × g for 2 h to spin down particles.

Electron microscopy
Schistosome and schistosome cell lysates were prepared as described above. Schistosome cell lysates from sucrose ultracentrifugation were deposited onto 200-mesh formvar-coated grids (Agar Scientific) and left to dry at RT. Grids were rinsed with H2O and subsequently treated with 1% uranyl acetate (System Biosciences) for a 5-min incubation for staining. Then, the grids were washed once in 70% ethanol followed by four times of washes with H2O. Then, grids were loaded onto the sample holder of a transmission electron microscope (Hitachi H-7600; Hitachi, Tokyo, Japan) and exposed to an 80 kV electron beam for morphological observation and image capture.

Ultrathin tissue sections for TEM and Immune TEM
Adult schistosomes (28 dpi) were subjected to overnight fixation at 4°C using 2.5% glutaraldehyde diluted in 0.1 M cacodylate buffer. The worms were washed three times with PBS (3× 10 min). Subsequently, the parasites were rinsed with 0.1 M PBS (pH 7.4) and post-fixed with 1% osmium tetroxide in PBS (pH 7.4) at 4°C for 1 h. Dehydration was performed through a graded ethanol series, followed by transfer to propylene oxide, and the samples were embedded in Eponate 812 resin (Ted Pella, USA). Longitudinal ultra-thin sections (60–70 nm thick) were prepared, contrasted sequentially with uranyl acetate and lead citrate, and imaged under a transmission electron microscope with an accelerating voltage of 80 kV.

Ultrathin sections (90 nm) encompassing the regions of interest were harvested on uncoated nickel grids (M200-NI, Electron Microscope Sciences). Sections were etched with 1% sodium periodate (NaIO₄) for 5 min at room temperature (RT), followed by rinsing with H2O and 0.1 mol/L glycine solution. Subsequently, sections were blocked with 1% BSA for 15 min, probed with primary antibodies (1:500 dilution) in 1% BSA at 4°C for 48 h, and then incubated with Protein A-conjugated 10 nm gold particles (810.111, Aurion) for 2 h at RT. After incubation, the sections were washed 5× 10 min in PBS buffer. Antibody complexes were immobilized via incubation in 2% glutaraldehyde in phosphate buffer for 10 min at RT. The grids were rinsed three times with H2O and stained with 3% uranyl acetate for 10 min and 3% lead citrate (Leica Ultrostain 2) for 6 min. Sections were observed using a JEM-1400 Plus electron microscope, equipped with a Quemsa TEM CCD camera, followed by image capture via iTEM software (TEM Imaging Platform).

In vivo inhibition of rp1 RdRp in mice infected with S. japonicum
BALB/c mice (male, 7–8 week old, mean body weight 25 ± 2 g) were obtained from the Shanghai SLAC Laboratory Animal Co., Ltd and were housed as described above. Mice (n=18) were randomly allocated into three groups. Each mouse was infected with S. japonicum cercariae (n=100 ± 5) via abdominal skin penetration. At 18 dpi, mice in each group received a tail vein injection of 1 µg siRNA or scrambled siRNA diluted in 200 µL of nuclease-free H2O. Concurrently, one group was administered an equal volume of PBS via the tail vein. Four subsequent injections were given on days 20, 22, 24, and 26 dpi. At 28 dpi, worms were collected via liver perfusion, and liver worm burden and egg deposition in the liver were counted. The expression of rp1 in collected worms was analyzed by RT-qPCR. Equal weight of liver tissues was collected from each mouse and then subjected to egg purification as described previously58, and the egg hatch experiment was carried out. Liver tissue was fixed and analyzed for Hematoxylin and Eosin staining as described below.

Histopathology of animal liver infected with S. japonicum
To evaluate the pathological changes in livers following in vivo rp1 inhibition, mouse livers were harvested, chopped into small cubes (1.0× 1.0× 0.3 cm), and fixed in 10% formalin. Subsequently, the fixed tissues were dehydrated through a graded ethanol series (70, 80, 95, and 100%), cleared with xylene, and embedded in paraffin wax. Finally, the paraffin-embedded specimens were sectioned into 6 µm thick slices and stained using the Hematoxylin and Eosin method according to the standard protocols.

WISH and FISH analyses of RdRp location in D. japonica
WISH analysis of RdRp localizations was performed as described above, except the hybridization at 56°C using the probe as in vitro transcribed using PCR products amplified by the specific primers forward 5′CGGCACTAGTCAATCAACCCTCT3′, reverse 5′TAATACGACTCACTATAGGTGTTTGGTTTACTTGTTGGCTCATTC3′, synthesized by Invitrogen). WISH was carried out as described above. Fluorescence in situ hybridization (FISH) was performed as described below. Briefly, worms were washed with ice-cold 2% HCl for 5 min before fixation, and then were fixed for 30 min at RT in 4% formaldehyde in PBS. Then, worms were washed using PBSTx (1× PBS + 0.3% Triton X-100). After fixation, worms were dehydrated and rehydrated in methanol and then treated with 1.2% H2O2 under light for 4 h. The worms were then rehydrated in methanol and then permeabilized with 10 µg/mL Proteinase K (Thermo Fisher Scientific) for 10 min, and fixed for 10 min immediately. The worms were washed with PBSTx twice and then carried out the hybridization at 56°C. Then, the worms were washed with TNTx (100 mM Tris pH 7.5, 150 mM NaCl, 0.3% Triton X-100) solution twice and then were incubated in the solution containing FISH blocking solution (5% horse serum (Sigma) and 0.5% Western Blocking Reagent (Merk) in TNTx (100 mM Tris pH 7.5, 150 mM NaCl, 0.3% Triton X-100) overnight at 4°C. Then, the worms were incubated with anti-DIG-POD (Merk, catalog# 11207733910, lot# 77798000) for 4 h at RT. Upon the wash of several times using TNTx, the Cy5 TSA Fluorescence System Kit (YEASEN) was used to stain the worms for 15 min. Then, the worms were subjected to washing twice using TNTx and stained with DAPI (10 µg/ml). The worms were further washed with PBSTx and then cleared with 80% glycerol, and mounted in Vectashield (Yeasen). The images were capture under a confocal microscope (Nikon).

RdRp silencing in D. japonica
The DNA template for generating RdRp dsRNA was obtained by PCR using the specific primers (Supplementary Table 1). The dsRNA was synthesized using the T7 RiboMAX™ Express RNAi kit (Promega, USA), then it was purified with 2.5 volumes of ethanol and 1/10 volume of NH4OAc (pH 5.2) after DNase I treatment (TaKaRa). dsRNA was delivered in planarian (5 µg/ml) based on the protocol with minor modification59. Control animals were treated with S. japonicum Mettl3 dsRNA. After 6 days of post-treatment, worms were collected for WISH and RT-qPCR analyses. The primers of S. japonicum Mettl3 were used to synthesize the dsRNA using the T7 RiboMAX™ Express RNAi kit (Promega, USA) as an irrelevant control dsRNA (Supplementary Table 1). The effect of RdRp inhibition was also evaluated by qPCR using the specific primers (Supplementary Table 1).

Statistical analysis
The data were analyzed using GraphPad Prism 8. Statistical analyses were conducted using Student’s t-test or one-way ANOVA. P ≤ 0.05 was considered significant.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "In this study, the authors identify ngRNA in the parasitic flatworm Schistosoma japonicum. These ngRNAs appear to be\nrelated to viral RNA-dependent RNA polymerases and nucleocapsid proteins, suggesting a viral origin of these particles.\nUsing ISH and single cell RNAseq, the authors localize the ngRNAs to schistosome reproductive tissues, including the\ntestes, the ovaries, and the vitellaria. The authors then use RNAi and pharmacologic inhibitors to suggest that these\nngRNAs are essential for ovarian health and egg production in S. japonicum. The authors work towards identifying the\nmechanism by which these ngRNAs might be essential for reproduction using yeast-2-hybrid. They also show that the\nngRNAs are likely functional, evidenced by the presence of RDRP activity in the parasites as well as the presence of virus-\nlike particles in parasites. The authors then use in vivo RNAi to show that inhibiting one of these ngRNAs results in a\nreduction of schistosome-associated pathology. Finally, the authors identify a RDRP in the free-living flatworm Dugesia\njaponica, suggesting that endogenous viruses/viral-like elements may be present throughout flatworms.\nOverall, this is a fascinating study that highlights not only the presence of virus-like elements in flatworms, but it also\nsuggests that these virus-like elements may actually play a beneficial role in these animals. By and large, the data support\nthe conclusions and the authors have performed most reasonable experiments. That being said, there are a number of major\nand minor concerns:\nMajor:\nThe authors' Y2H experiments fail to produce a convincing explanation for the mechanism by which rp1 could be affecting\nthe observed phenotype. The STRING analysis is entirely correlative as is the assertations that LIN-10 and RPLP1 regulate\ndevelopment in schistosomes and C. elegans. The authors should perform RNAi of putative interactors to see if any of the\nputative interactors phenocopy rp1 RNAi. Absent that, the authors can really only speculate on the mechanism of rp1\nfunction and certainly are unable to declare that \"these results collectively establish the molecular mechanisms...\" (line 171)\nThere is no explanation in the text, figure legends, or materials and methods that adequately explains how Figure S5 was\ngenerated. This reviewer surmises that Fig S5A shows predicted delta G values from various AlphaFold model predictions\n(RP1_M0, M1, etc.) calculated by AutoDock. There is no indication which structures/models are shown in Fig S5B.\nFurthermore, this figure is largely uninformative. The authors would be better served by purifying rp1 and showing that active\nFavipiravir inhibits its activity. Purified rp1 would also allow the authors to more definitively show that rp1 is responsible for\nRDRP activity in S. japonicum, complimenting Figure 5A-B. Alternatively, expressing rp1 in a cell line could allow the\nauthors to show that Favipiravir is inhibiting rp1 without worrying about the fact that Favipiravir is a prodrug.\nThe authors describe measuring rp1 expression levels in parasites collected from in vivo RNAi experiments in the methods\nsection, but those data are not present to the best of this reviewer's knowledge. The in vivo data are compelling, but without\nactual proof that in vivo RNAi worked, those data can not be relied upon.\nThe authors mapped their Dugesia japonica reads to GCA_022537955.1, which is a Schmidtea mediterranea genome. This\nlikely explains their very high \"Unknown\" reads in Figure 8A. Given this, describing any reads as mapping to D. japonica is\nalmost certainly incorrect. That 0.3% of the reads map to viral genomes is sufficient to suggest the presence of ngRNA in\nDugesia, but it does not preclude the possibility that these ngRNAs are in fact genomic in Dugesia. Figure 8B helps\nprove that the Djap RdRp is probably not genomic, but the authors should either refrain from mapping to the incorrect genome or\nperform a similar analysis with flatworms with full genomes available (i.e., Schmidtea mediterranea).\nThe authors state that coatamer subunit epsilon interacts with rp1 (line 163), but \"Coatamer subunit\" drops out of the Y2H\nhits in the -His condition, suggesting that it is not a bona fide interaction\nThere is no data presented about upregulated genes (Fig 4B-C). These data need to be available to help the reader interpret\nthe consequences of rp1/np1 RNAi.\nMinor:\nThe authors may wish to indicate that RPLP1 is a known host factor in flavivirus infection (PMID: 27974556) (line 162)\nThe authors mapped their Sjap reads to PRJNA520774 whereas the most up to date Sjap genome project is PRJNA724792\n(which is ~35Mb larger). While this probably doesn't matter in the context of most experiments, because the authors are\nattempting to show that the ngRNA sequences don't exist in the Sjap genome, they should probably use the most up to\ndate/comprehensive genome available. That being said, there may be reasons to use an older version in terms of annotation\nquality, etc. If there is a justification for using the older genome, the authors should explain it.\nIn figure 4A, it would be helpful for the authors to label anterior/posterior to help the reader orient themselves with respect to\nthe ovary.\nIt is not clear what the lower panel in Fig 8C represents.\nDugesia has an epidermis, not a tegument (line 239)\nIt is generally difficult to identify the exact gene that the authors are referring to (e.g., Figure S6 lists pGatT7-SJCHGC05978,\nbut SJCHGC05978 is not a current gene accession on Wormbase Parasite). It would be convenient if the authors could\nprovide up-to-date gene accessions for all genes mentioned in the manuscript in a supplemental file.\nOther:\nThe authors stated that they attempted to infect a number of different cell lines to no success. The authors may wish to try\nBGE cells as they are probably the closest thing to schistosomes that's available. This is beyond the scope of this\nmanuscript and the authors do not need to address this.", "reviewer_2": "This paper reports the discovery and biological impact of several viral-like dsRNAs in the pathogenic fluke Schistosoma\njaponicum. While the impact of viral-like elements in protists and fungi is increasingly being recognized, there are fewer\nreports of these in metazoans and especially invertebrates. The data show clearly their presence, properties and functional\nconsequences on worm development, and potentially provide new opportunities for therapies as well as expanding our\nunderstanding of the interrelationships of virus-like elements and pathogenicity. As such, this is potentially an important\npaper. The technical quality and most of the statistical analysis is good.\nHaving said that, there are a great many concerns mostly involving the presentation or analysis of the datasets. These\nshould be addressed not only in the response to the editors, but in the revised text. The order below is arbitrary for the most\npart and not necessarily in priority of their importance.\n1. The many text contains too many experimental details, which could be moved to the supplemental data without loss of\nclarity. The analysis of such dsRNAs is routine in protists and fungi. I think the diagrams of standard analysis included could\nbe mostly omitted as well.\n2. The analysis of the final contig set is inadequate. No information is given about the number of contigs found, especially\nbeyond the 4 rp1/2 np1/2 focus. This is relevant for completeness as well as information relevant to the potentially 'missing'\nglycoprotein. How were the contigs screened? blast or more sophisticated searches? this needs to be spelled out in much\ngreater detail and the key aspects emphasized in the text. It should be noted that absence of the putative glycoprotein\nsegment is not uncommon in many metatranscriptomic datasets, but in few cases has its absence simply because of\ndivergence been ruled out.\n3. The properties of the core 4 contigs are woefully under-described, appearing only as a map in the supplemental data.\nPotential important aspects not discussed are key motifs found within the predicted coding regions, the percent nucleotide\nand/or amino acid similarity between rp1 and 2, or np1 and 2, etc. The map should show the location of the siRNAs\nultimately used for functional studies and some discussion about potential cross-RNAi between rp1 and 2 which seem to be\nquite similar. To me this is a main text figure.\n4. The phylogenetic tree to me deserves more priority than relegation to supplemental data. First the labelling using only\nNCBI numbers is not very helpful - better to have fewer relevant taxa and some names of the more common or type species\nto promote comprehension.\n5. The bunyaviral nomenclature is inconsistent and often incorrect. 'Bunyavirus' was formerly a generic term but after\nmultiple taxonomic revisions the overall group is termed 'Bunyavirales' with various sub-taxa therein. I suggest consulting\nthe most recent ICTV taxonomy to get it correct. Also I believe that Tospoviruses are not Tenuiviruses, in contrast to\nstatements in text or figure, and this should be corrected and made consistent.\n6. The evolutionary relationships amongst RP1/2 and NP1/2 are very superficially discussed. For RP1/2 they appear closely\nrelated to each other but very distant from the other BV taxa provided in the figure - there is no statistical analysis of the\nstrength of these classifications provided, not even a bootstrap. To me they appear more likely to be a unique independent\nclade than a member of existing BV taxa.\nFor NP1/2, they appear very different from each other, but unlike RP1/2 closer to existing taxa. This suggests NP1/2 may\nhave a distinct origin and/or evolutionary path from RP1/2. Some sort of analysis presenting and discussing this more clearly\nis needed.\n7. In contrast to statements throughout out the text, very very few of the protist or fungal viral-like elements are thought to be\ninfectious (Giardiavirus is one of few exceptions). Instead they are thought to have originated deep in evolution by\nspeculative mechanisms (abortive hybridizations exosomes etc) but forever after inherited vertically through the germ line.\nThis should be corrected. More importantly, the origins of the viruses have direct implications to the question of the\nevolution/co-evolution of the viral segments mentioned in #6 above, which should be discussed.\n8. The discussion states that the virus-like elements are 'universal' in S. japonicum but the results suggest it was found in\n89/~100 NCBI datasets. What is the true story? and some analysis supporting this should be provided. If an evolutionary tree\nof the Sj genomes were conducted would the virus element presence/abundance be clustered or not? Given the biological\nconsequences seen in their experimental tests this is an important aspect of the paper. For comparison, the Plasmodium\nvivax virus NLV1 reported by Holmes is absent in many datasets calling into question its universal impact. In this regards, it\nshould be discussed that the future studies might be warranted on those Sj isolates lacking the viral elements.\n9. In protists and fungi these RNA elements are termed viruses or virus-like. Why shouldn't this be done for the Sj elements?\nIt is warranted and would improve the impact and recognition of their work.\n10. Two reasonable mechanisms for how the viral-like RNAs operate are provided, through direct binding of proteins or\nalterations in transcript abundancece. A third mechanism should be considered, namely the impact of dsRNAs themselves\non RNA sensing pathways. This is how Totiviruses and Cryspovirus impact pathogenicity, although there it is through Type I\ninterferon pathways which I believe are absent in worms. What is known about RNA sensing pathways in worms and\nespecially in the germ line, where in vertebrates there are key differences with the somatic cells?\n11. Does the recombinant RP1 protein have enzymatic activity?\n12, Fig S4 is mislabeled - some of figures should be np rather than rp. Fig 7,8 labeling should say rp1 rather than RDRP."}, "author_response": "In response to Reviewer #1 comments\n1. \"In this study, the authors identify ngRNA in the parasitic flatworm Schistosoma japonicum. These\nngRNAs appear to be related to viral RNA-dependent RNA polymerases and nucleocapsid proteins,\nsuggesting a viral origin of these particles. Using ISH and single cell RNAseq, the authors localize\nthe ngRNAs to schistosome reproductive tissues, including the testes, the ovaries, and the vitellaria.\nThe authors then use RNAi and pharmacologic inhibitors to suggest that these ngRNAs are essential\nfor ovarian health and egg production in S. japonicum. The authors work towards identifying the\nmechanism by which these ngRNAs might be essential for reproduction using yeast-2-hybrid. They\nalso show that the ngRNAs are likely functional, evidenced by the presence of RDRP activity in the\nparasites as well as the presence of virus-like particles in parasites. The authors then use in vivo\nRNAi to show that inhibiting one of these ngRNAs results in a reduction of schistosome-associated\npathology. Finally, the authors identify a RDRP in the free-living flatworm Dugesia japonica,\nsuggesting that endogenous viruses/viral-like elements may be present throughout flatworms.\nOverall, this is a fascinating study that highlights not only the presence of virus-like elements in\nflatworms, but it also suggests that these virus-like elements may actually play a beneficial role in\nthese animals. By and large, the data support the conclusions and the authors have performed most\nreasonable experiments. That being said, there are a number of major and minor concerns:\"\nAuthor response: Thank you for your positive comments and insightful evaluation of our\nmanuscript. We deeply appreciate your recognition of our present study.\nMajor:\n1. \"The authors' Y2H experiments fail to produce a convincing explanation for the mechanism by\nwhich rp1 could be affecting the observed phenotype. The STRING analysis is entirely correlative\nas is the assertations that LIN-10 and RPLP1 regulate development in schistosomes and C. elegans.\nThe authors should perform RNAi of putative interactors to see if any of the putative interactors\nphenocopy rp1 RNAi. Absent that, the authors can really only speculate on the mechanism of rp1\nfunction and certainly are unable to declare that \"these results collectively establish the molecular\nmechanisms...\" (line 171)\"\nAuthor response: Thank you very much for the valuable comments. Based on the comments, we\ncarried out additional experiments to suppress two interactors as identified by Yeast two\nhybridization including Lin-10 and LP1 in adult females and then determine effect of inhibition of\nLin-10 and Lp1 on morphological vitellarium and egg production. The results indicated that the\ndecreased transcript levels of either Lin-10 or Lp-1 led to defected vitellocyte morphology and\nreduced number of egg production. The related results were included in the Fig 4G, 4H and 4I.\nThe following statement was included in the results of revised manuscript \" Next, we verified the\nfunctional roles of these identified interactors using RNAi. The results indicated that suppression of\nLin-10 and LP-1 (Fig. 4G) led to abnormal morphology of vitellocytes (Fig. 4H) and reduced\nnumber of egg production (Fig. 4I) (Supplementary data 1).\" (Line:174-176).\nThe related method information was also included in method section of the revised manuscript (Line:\n650-665).\n2. \"There is no explanation in the text, figure legends, or materials and methods that adequately\nexplains how Figure S5 was generated. This reviewer surmises that Fig S5A shows predicted delta\nG values from various AlphaFold model predictions (RP1_M0, M1, etc.) calculated by AutoDock.\nThere is no indication which structures/models are shown in Fig S5B. Furthermore, this figure is\nlargely uninformative.\"\nAuthor response: Thank the reviewer for the comments. To clarify this point, we not only revised\nthe Supplementary Fig. 6 (as shown below) in the revised manuscript but also included the\nfollowing statement in the method and the legend of Supplementary Fig. 6 \"To predict RP1/2 RdRp\nstructure, we employed the alphafold server (https://alphafoldserver.com/) to predict their structures\nfor the drug-protein docking (DPD) using the default settings. The obtained protein molecule PDB\nfile and the small molecules (Favipiravir and its derivative FTP, in mol2 format) were subjected to\nthe DPD program using AutoDock software (https://autodock.scripps.edu/) following the official\nprotocol. The docking result was further processed by AutoDockTools\n(https://autodocksuite.scripps.edu/adt/) to obtain a drug-protein complex, which was further\nvisualized by PyMol software (https://pymol.org/). The docking affinity of protein and ligands was\nquantified by the free binding energy.\" (Line:720-727).\nSupplementary Figure 6. In silico prediction of the chemical inhibitor binding to ngRNA RdRp.\nAdditionally, we added explicit information that the FTPs were docked into the RPs using the\nAlphaFold-predicted 3D structural models RP1_M2 and RP2_M1, respectively. We updated the\nlegend of Supplementary Fig. 6 to specify these models, enhancing the reproducibility of our\ndocking analyses.\n3. \"The authors would be better served by purifying rp1 and showing that active Favipiravir inhibits\nits activity. Purified rp1 would also allow the authors to more definitively show that rp1 is\nresponsible for RDRP activity in S. japonicum, complimenting Figure 5A-B. Alternatively,\nexpressing rp1 in a cell line could allow the authors to show that Favipiravir is inhibiting rp1\nwithout worrying about the fact that Favipiravir is a prodrug.\"\nAuthor response: Thank the reviewer for the constructive comments. To address the reviewer's\ncomments, we designed and carried out additional experiments in the revised manuscript. Briefly,\nwe cloned the conserved Rp1 RdRp domain into pcDNA3.1 vector. The reversed complementary\nGaussia luciferase reporter gene that was fused into 5`UTR and 3`UTR of Rp1 RdRp were also\nconstructed into pcDNA3.1 vector. These recombinant plasmids were transfected into HEK293T\ncells and luciferase assay indicated the activity of RdRp. Moreover, upon the adding of Favipiravir,\nthe luciferase activity was decreased, further supporting the activity of RdRp. Similar results were\nalso observed in the cells co-transfected with the recombinant plasmids containing reversed\ncomplementary np1. These results were included in Fig. 5C and 5D (as shown below) in the revised\nmanuscript.\nThe following statement was included in the results of revised manuscript \" Next, we determined\nwhether that RP1 RdRp is responsible of ngRNA replication and its activity can be inhibited by\nFavipiravir. Upon co-transfection of recombinant plasmids for expressing RP1 RdRp domain and\nreversed complementary Gaussia luciferase reporter gene embedding into 5`UTR and 3`UTR of\nrp1 or reversed complementary np1 into HEK293T cells, luciferase assay and RT-qPCR analyses\nindicated that RP1 RdRp can drive the expressions of either Gaussia reporter gene or np1 gene (Fig.\n5C and 5D) (Supplementary data 1). When applying Favipiravir in transfected cells, their\nexpressions were significantly suppressed (Fig. 5C and 5D) (Supplementary data 1).\" (Line:193-\n199).\nThe related method information was also included in the revised manuscript (Line:678-700).\n4. \"The authors describe measuring rp1 expression levels in parasites collected from in vivo RNAi\nexperiments in the methods section, but those data are not present to the best of this reviewer's\nknowledge. The in vivo data are compelling, but without actual proof that in vivo RNAi worked,\nthose data can not be relied upon. \"\nAuthor response: To address the reviewer's comments. We isolated total RNA from the stored\nparasites from previous in vivo experiment. Then, RT-qPCR was used to determine the transcript\nlevel of rp1 in the parasites. The results indicated the decreased transcript level of rp1 was observed.\nThe results were included in Fig 7C (as shown below) of the revised manuscript.\nThe following statement was included in the results of revised manuscript \" Additionally, total RNA\nwas extracted from the collected worms and the transcript levels of rp1 were analyzed by RT-qPCR.\nThe result indicated that the transcript level of rp1 was significantly decreased in remaining worms\ncollected from the mice administered with rp1 siRNA (Fig. 7C) (Supplementary data 1).\"\n(Line:237-239).\n5. \"The authors mapped their Dugesia japonica reads to GCA_022537955.1, which is a Schmidtea\nmediterranea genome. This likely explains their very high \"Unknown\" reads in Figure 8A. Given\nthis, describing any reads as mapping to D. japonica is almost certainly incorrect. That 0.3% of the\nreads map to viral genomes is sufficient to suggest the presence of ngRNA in Dugesia, but it does\nnot preclude the possibility that these ngRNAs are in fact genomic in Dugesia. Figure 8B helps\nprove that the Djap RdRp is probably not genomic, but the authors should either refrain from\nmapping to the incorrect genome or perform a similar analysis with flatworms with full genomes\navailable (i.e., Schmidtea mediterranea).\"\nAuthor response: Thank the reviewer for the important comments and suggested correction. In the\nrevised manuscript, we used Dugesia japonica genome (GCA_001938525.1) to re-analyze our\nmetragenmics data. The results indicated that 90.77% of reads (58383269 out of 64319570 reads)\ncan be mapped to Dugesia japonica genome. There are 362 reads that successfully map to viral\ndatabase. The following assemble yielded the rp1 RdRp. The updated results were included in Fig.\n8A of the revised manuscript (Fig. 8A).\n6. The authors state that coatamer subunit epsilon interacts with rp1 (line 163), but \"Coatamer\nsubunit\" drops out of the Y2H hits in the -His condition, suggesting that it is not a bona fide\ninteraction\nAuthor response: Thank the reviewer for the important comments and suggested correction. In the\nrevised manuscript, we remove the Coatamer subunit in Fig. 4D, but also keep all of original results\nin Supplementary Fig 7.\n7. There is no data presented about upregulated genes (Fig 4B-C). These data need to be available\nto help the reader interpret the consequences of rp1/np1 RNAi.\nAuthor response: In the revised manuscript, the upregulated genes were included in the\nSupplementary data 3 and data 4.\nMinor:\n1. \"The authors may wish to indicate that RPLP1 is a known host factor in flavivirus infection\n(PMID: 27974556) (line 162) \"\nAuthor response: The related statement was included in the revised manuscript as follow \"In\naddition, ribosomal protein RPLP1 is a known host factor in flavivirus infection\" (Line:168-169).\nAdditionally, t he indicated literature was also cited in revised manuscript.\n2. \"The authors mapped their Sjap reads to PRJNA520774 whereas the most up to date Sjap genome\nproject is PRJNA724792 (which is ~35Mb larger). While this probably doesn't matter in the context\nof most experiments, because the authors are attempting to show that the ngRNA sequences don't\nexist in the Sjap genome, they should probably use the most up to date/comprehensive genome\navailable. That being said, there may be reasons to use an older version in terms of annotation\nquality, etc. If there is a justification for using the older genome, the authors should explain it.\"\nAuthor response: Thank the reviewer for the important comments. In the revised manuscript, we\nused two version of S. japonicum genomes (ASM636876v1 and ASM2146165v1) to perform\nbioinformatic analyses. After removing reads mapped to S. japonicum genomes as well as the rabbit\nhost, assembly of the remaining reads yielded same number of contigs, eg, 8355, 7246, 7069 contigs,\nrespectively (Fig. 1B).\nWe would like to provide the following table to show the mapping results between two versions of\nS. japonicum genomes.\nReads Reads Reads Reads\nHigh\nmapping mapping mapping mapping to\nSamples quality Genome information\nto to to viral viral\nreads\ngenome genome% database database%\nrabbit (mOryCun1.1, GCA_964237555.1) 37037 1.94%\nSRR31784407/ SjR3 1912254 PRJNA520774 (GCA_006368765.1) 547498 28.63% 726359 37.98%\nPRJNA724792(ASM2521551v1) 548442 28.68%\nrabbit (mOryCun1.1, GCA_964237555.1) 12690 1.93%\nSRR31784413/ SjR2 657352 PRJNA520774 (GCA_006368765.1) 172916 26.30% 276889 42.12%\nPRJNA724792(ASM2521551v1) 173057 26.33%\nrabbit (mOryCun1.1, GCA_964237555.1) 7698 1.76%\nSRR31784414/ SjR1 437447 PRJNA520774 (GCA_006368765.1) 112706 25.76% 182422 41.70%\nPRJNA724792(ASM2521551v1) 112848 25.80 %\n3. \"In figure 4A, it would be helpful for the authors to label anterior/posterior to help the reader\norient themselves with respect to the ovary.\"\nAuthor response: In the revised manuscript, we used arrows to indicate the anterior.\n4. \"It is not clear what the lower panel in Fig 8C represents.\"\nAuthor response: It is the WISH results using control probe. We labelled Fig. 8C in the revised\nmanuscript.\n5. \"Dugesia has an epidermis, not a tegument (line 239)\"\nAuthor response: We corrected it accordingly in the revised manuscript.\n6. \"It is generally difficult to identify the exact gene that the authors are referring to (e.g., Figure S6\nlists pGatT7-SJCHGC05978, but SJCHGC05978 is not a current gene accession on Wormbase\nParasite). It would be convenient if the authors could provide up-to-date gene accessions for all\ngenes mentioned in the manuscript in a supplemental file.\"\nAuthor response: In the revised manuscript, we provided the up-to-date gene and protein accession\nID information in Supplementary Fig. 7B (as shown below).\nIn response to Reviewer #2 comments (Remarks to the Author):\nThis paper reports the discovery and biological impact of several viral-like dsRNAs in the\npathogenic fluke Schistosoma japonicum. While the impact of viral-like elements in protists and\nfungi is increasingly being recognized, there are fewer reports of these in metazoans and especially\ninvertebrates. The data show clearly their presence, properties and functional consequences on\nworm development, and potentially provide new opportunities for therapies as well as expanding\nour understanding of the interrelationships of virus-like elements and pathogenicity. As such, this\nis potentially an important paper. The technical quality and most of the statistical analysis is good.\nHaving said that, there are a great many concerns mostly involving the presentation or analysis of\nthe datasets. These should be addressed not only in the response to the editors, but in the revised\ntext. The order below is arbitrary for the most part and not necessarily in priority of their\nimportance.\n1. \"The many text contains too many experimental details, which could be moved to the\nsupplemental data without loss of clarity. The analysis of such dsRNAs is routine in protists and\nfungi. I think the diagrams of standard analysis included could be mostly omitted as well. \"\nAuthor response: Thank the reviewer for the comments. In the revised manuscript, we revised the\nschematic diagram.\n2. \"The analysis of the final contig set is inadequate. No information is given about the number of\ncontigs found, especially beyond the 4 rp1/2 np1/2 focus. This is relevant for completeness as well\nas information relevant to the potentially 'missing' glycoprotein. How were the contigs screened?\nblast or more sophisticated searches? this needs to be spelled out in much greater detail and the\nkey aspects emphasized in the text. It should be noted that absence of the putative glycoprotein\nsegment is not uncommon in many metatranscriptomic datasets, but in few cases has its absence\nsimply because of divergence been ruled out.\"\nAuthor response: Thank you very much for the valuable comments to improve our bioinformatic\nanalyses. After removing reads mapped to S. japonicum genomes as well as the rabbit host, assembly\nof the remaining reads yielded 8355, 7246, 7069 contigs, respectively (Fig. 1B). This following\nstatement was included in the revised manuscript \" After removing reads mapped to S. japonicum\ngenomes (ASM636876v1 and ASM2146165v1) as well as the rabbit host, assembly of the\nremaining reads yielded 8355, 7246, 7069 contigs, respectively (Fig. 1B) (Supplementary data\n1)\"(Line:65-67)\nAs to bioinformatic analyses for contigs, we used DIAMOND (v2.1.12) to perform BLASTx\nalignment on the assembled contigs and raw reads, respectively. Each sequence was aligned against\nthe NCBI Viral Proteome Database (https://ftp.ncbi.nlm.nih.gov/refseq/release/viral/), and an E-\nvalue threshold was set to screen potential viral sequences, with the E-value cutoff set at less than\n10⁻⁵. In addition, we also re-aligned the viral sequences obtained against the NCBI Internal Non-\nProtein Non-Redundant (NVNR) Protein Database to further eliminate potential false-positive viral\nsequences. The following statements were also included the section of methods \" Sequencing reads\nwere first adaptor- and quality-trimmed using the Trim-galore program (v0.6.10) (-phred33 -length\n100 -stringency 3 -paired, https://zenodo.org/records/7598955). The remaining reads were mapped\nto the S. japonicum reference genomes (NCBI accession: GCA_006368765.1 and\nGCA_021461655.1), rabbit (NCBI accession: GCA_964237555.1), and mice (NCBI accession:\nGCF_000001635.27) genomes using the Bowtie2 program (v2.5.0)33. The remaining reads were\nassembled de novo using MEGAHIT (v1.2.9) with default parameters (-min-contig-len 200) to\ngenerate the longest contigs and singlets34. Next, the contigs and singlet reads were analyzed by\nDIAMOND (v2.1.12) with default parameters. Each sequence was aligned against the NCBI virus\nreference proteome (ftp://ftp.ncbi.nih.gov/refseq/release/viral/) in addition to viral protein\nsequences from the NCBI nr Fasta file (based on annotation taxonomy in Virus Kingdom)32, 35, and\nan E-value threshold of <10⁻⁵ was set to screen potential viral sequences. To further eliminate\npotential false-positive viral sequences, the viral sequences obtained through the above steps were\nre-aligned against the NCBI Internal Non-Protein Non-Redundant (NVNR) Protein\nDatabase.\"(Line: 380-391).\nAs to glycoprotein segment, we included the following statements \" Phylogenetic tree and protein\nmotif analyses revealed that these ngRNAs are closely related to Phenuiviridae family (class\nBunyaviricetes). The Phenuiviridae family typically comprises two to eight single-stranded\nmolecules of negative-sense or ambisense encoding RdRp, viral glycoproteins, nucleocapsid protein,\nand non-structural proteins 24. Although we noted that the ngRNAs were significantly enriched,\naccounting for more than 60 percent of total reads, our subsequent bioinformatic analyses detected\nno sequences encoding viral glycoproteins. Additionally, we enriched micro-sized particles by\nsucrose gradient separation and found that the ngRNAs co-fractionated with the most abundant\nparticle component. Furthermore, sequencing and bioinformatic analyses of the RNA from this\ncomponent did not detect any segments encoding glycoproteins. While our study was in progress,\nShi et al. performed deep transcriptome sequencing on various invertebrate species, and reported\nthe presence of similar RNA segments in schistosomes, which also lacked an associated\nglycoprotein RNA25. Collectively, our findings suggest that these ngRNAs may dispense with a\ncanonical matrix protein for the maintenance of their genomic integrity.\" in the discussion\n(Line:283-294).\n3. \"The properties of the core 4 contigs are woefully under-described, appearing only as a map in\nthe supplemental data. Potential important aspects not discussed are key motifs found within the\npredicted coding regions, the percent nucleotide and/or amino acid similarity between rp1 and 2, or\nnp1 and 2, etc. The map should show the location of the siRNAs ultimately used for functional\nstudies and some discussion about potential cross-RNAi between rp1 and 2 which seem to be\nquite similar. To me this is a main text figure.\"\nAuthor response: Thank you very much for the valuable comments to improve clarity of the\nmanuscript. To clarify this point, the following statements were included in the revised manuscript\n\"Protein domain and motif analyses indicated that RP1/2 contained an RNA dependent RNA\npolymerase motif related to Bunyaviricetes and NP1/2 contained a Phenuiviridae nucleocapsid\nprotein motif (Fig. 1K, as shown below). Additionally, the NP1 shows 51.02% identity to NP2\nwhile RP1 displayed 57.84% identify to RP2.\"(Supplementary Fig. 3) (Line: 99-101).\nThe location of the siRNAs ultimately used for functional studies was highlighted in\nSupplementary Fig. 3 (a part of the image is shown below, red box indicating siRNA targeting\nregions), where are not the potential of cross-RNAi.\n4. \"The phylogenetic tree to me deserves more priority than relegation to supplemental data. First\nthe labelling using only NCBI numbers is not very helpful - better to have fewer relevant taxa and\nsome names of the more common or type species to promote comprehension.\"\nAuthor response: Thank you very much for the valuable comments to improve clarity of the\nmanuscript. In the revised manuscript, we constructed the phylogenetic trees for NP1/2 and RP1/2.\nThe results were included in Fig. 1I and 1J (as shown below) in the revised manuscript. The\nfollowing statement was also included the manuscript \"Phylogenetic analyses revealed that either\nRP1 and RP2 or NP1 and NP2 can cluster together belonging to Phenuiviridae family (Fig. 1I and\n1J). Generally, they also showed to be closely related to Trichobi phenuili virus (Fig. 1I and 1J). \"\n(Line: 96-99).\n5. \"The bunyaviral nomenclature is inconsistent and often incorrect. 'Bunyavirus' was formerly a\ngeneric term but after multiple taxonomic revisions the overall group is termed 'Bunyavirales' with\nvarious sub-taxa therein. I suggest consulting the most recent ICTV taxonomy to get it correct. Also\nI believe that Tospoviruses are not Tenuiviruses, in contrast to statements in text or figure, and this\nshould be corrected and made consistent.\"\nAuthor response: Thank you very much for the valuable comments. In the revised manuscript, we\nchecked the most recent ICTV taxonomy and revised them accordingly.\n6. \"The evolutionary relationships amongst RP1/2 and NP1/2 are very superficially discussed. For\nRP1/2 they appear closely related to each other but very distant from the other BV taxa provided in\nthe figure - there is no statistical analysis of the strength of these classifications provided, not\neven a bootstrap. To me they appear more likely to be a unique independent clade than a member\nof existing BV taxa. For NP1/2, they appear very different from each other, but unlike RP1/2 closer\nto existing taxa. This suggests NP1/2 may have a distinct origin and/or evolutionary path from\nRP1/2. Some sort of analysis presenting and discussing this more clearly is needed.\"\nAuthor response: Thank you for the comments. In revised manuscript, we included phylogenetic\ntrees regarding NP1/2 and RP1/2 in Fig. 1I and 1J. As mentioned above, the following statement\nwas also included the manuscript \"Phylogenetic analyses revealed that either RP1 and RP2 or NP1\nand NP2 can cluster together belonging to Phenuiviridae family (Fig. 1I and 1J). Generally, they\nalso showed to be closely related to Trichobi phenuili virus (Fig. 1I and 1J). \" (Line: 96-99).\nWe also included the following in the discussion section \" Phylogenetic analyses of the proteins\nencoding by these ngRNAs showed that they form a distinct clade. Notably, RP1 and RP2 are most\nclosely related to Trichobi phenuili virus 2 and Wufeng rodent phenuivirus 1, whereas NP1 and NP2\nare more closely related to Trichobi phenuili virus 1. The nucleoproteins associated with S.\njaponicum and Trichobilharzia regenti formed a distinct clade, separate from those of\nOpisthorchiidae liver flukes (Metorchis orientalis). Interestingly, RP1/2 is also closely related to\nWufeng rodent phenuivirus 1, which was identified from rodent internal organs and fecal samples\nvia meta-transcriptomic sequencing 27. At present, we are not able to exclude the possibility of\nflatworm infection in the rodent. Nevertheless, these observations suggest that viral-like RNAs may\nbe widespread among the flatworms. A study on the transmission of tapeworm (Schistocephalus\nsolidus) viruses revealed that most viruses are vertically transmitted from parents to offspring 28.\nThis supports our inference that the ngRNAs may be consistently maintained through long-term co-\nevolution with parasites.\" (Line:306-316).\n7. In contrast to statements throughout out the text, very very few of the protist or fungal viral-like\nelements are thought to be infectious (Giardiavirus is one of few exceptions). Instead they are\nthought to have originated deep in evolution by speculative mechanisms (abortive hybridizations\nexosomes etc) but forever after inherited vertically through the germ line. This should be corrected.\nMore importantly, the origins of the viruses have direct implications to the question of the\nevolution/co-evolution of the viral segments mentioned in #6 above, which should be discussed.\nAuthor response: Thank you very much for the comments. We revised and included the following\nstatements in the revised manuscript \" Phylogenetic analyses of the proteins encoding by these\nngRNAs showed that they form a distinct clade. Notably, RP1 and RP2 are most closely related to\nTrichobi phenuili virus 2 and Wufeng rodent phenuivirus 1, whereas NP1 and NP2 are more closely\nrelated to Trichobi phenuili virus 1. The nucleoproteins associated with S. japonicum and\nTrichobilharzia regenti formed a distinct clade, separate from those of Opisthorchiidae liver flukes\n(Metorchis orientalis). Interestingly, RP1/2 is also closely related to Wufeng rodent phenuivirus 1,\nwhich was identified from rodent internal organs and fecal samples via meta-transcriptomic\nsequencing 27. At present, we are not able to exclude the possibility of flatworm infection in the\nrodent. Nevertheless, these observations suggest that viral-like RNAs may be widespread among\nthe flatworms. A study on the transmission of tapeworm (Schistocephalus solidus) viruses revealed\nthat most viruses are vertically transmitted from parents to offspring 28. This supports our inference\nthat the ngRNAs may be consistently maintained through long-term co-evolution with parasites.\"\n(Line:306-316).\n8. \"The discussion states that the virus-like elements are 'universal' in S. japonicum but the results\nsuggest it was found in 89/~100 NCBI datasets. What is the true story? and some analysis\nsupporting this should be provided. If an evolutionary tree of the Sj genomes were conducted would\nthe virus element presence/abundance be clustered or not? Given the biological consequences seen\nin their experimental tests this is an important aspect of the paper. For comparison, the Plasmodium\nvivax virus NLV1 reported by Holmes is absent in many datasets calling into question its universal\nimpact. In this regards, it should be discussed that the future studies might be warranted on those\nSj isolates lacking the viral elements.\"\nAuthor response: Thank you for the comments. Since we are not able to verify whether these\nngRNAs are universal in S. japonicum. We carefully revised the related statements in the revised\nmanuscript. The following statement was included in the discuss of revised manuscript\" these\nngRNAs are probably present across different geographical strains and host isolates, suggesting they\nare broadly relevant to schistosome populations.\" (Line:264-266).\nAs to S. japonicum genome evolutionary tree, the S. japonicum genome dataset mainly resulted\nfrom China mainland stain of S. japonicum, they have same genome sequences although there may\nhave different SNPs in individual isolation. Currently, no enough information regarding S.\njaponicum genomes is able to use to construct an evolutionary tree of the S. japonicum tree.\n9. \"In protists and fungi these RNA elements are termed viruses or virus-like. Why shouldn't this be\ndone for the Sj elements? It is warranted and would improve the impact and recognition of their\nwork.\"\nAuthor response: In the title and main text of revised manuscript, we used \"virus-like\" term.\n10. \"Two reasonable mechanisms for how the viral-like RNAs operate are provided, through direct\nbinding of proteins or alterations in transcript abundancece. A third mechanism should be\nconsidered, namely the impact of dsRNAs themselves on RNA sensing pathways. This is how\nTotiviruses and Cryspovirus impact pathogenicity, although there it is through Type I interferon\npathways which I believe are absent in worms. What is known about RNA sensing pathways in\nworms and especially in the germ line, where in vertebrates there are key differences with the\nsomatic cells?\"\nAuthor response: To clarify these points, we included the following statements in the revised\nmanuscript \" In mammals, cytoplasmic sensors in host cells detect double-stranded RNA (dsRNA),\na common byproduct of viral replication. This recognition then triggers the interferon (IFN) system,\na first line of defense that limit viral replication. However, how dsRNA triggers immune like\nresponse in flatworms and what regulatory mechanisms control these processes in helminths remain\ncritical knowledge gaps. In Caenorhabditis elegans, Orsay virus (OrV) is the known natural virus\naffecting this worm. This nematode employs multiple, orchestrated defenses against OrV infection,\nincluding the RNAi pathway, viral RNA uridylation to direct RNA for degradation, and the\nexpression of antiviral proteins such as those involved in ubiquitination 29, Conversely, OrV\nreplication also depends on host factors like zinc and lipid metabolism, highlighting a complex\ninterplay between proviral and antiviral mechanisms that maintains dsRNA homeostasis. A recent\nstudy demonstrated adult B. malayi worms activate RNA interference as an innate immune response\nagainst viral replication30. Although RNAi is widely used to silence target genes in S. japonicum, it\nremains unknown whether endogenous RNAi is activated during ngRNA replication, particularly\nwithin the germline cells.\" (Line:318-331).\n11. Does the recombinant RP1 protein have enzymatic activity?\nAuthor response: To address the reviewer's comments,\nwe designed and carried out additional experiments in the revised manuscript. Briefly, we cloned\nthe conserved Rp1 RdRp domain into pcDNA3.1 vector. The reversed complementary Gaussia\nluciferase reporter gene that was fused into 5`UTR and 3`UTR of Rp1 RdRp were also\nconstructed into pcDNA3.1 vector. These recombinant plasmids were transfected into HEK293T\ncells and luciferase assay indicated the activity of RdRp. Moreover, upon adding Favipiravir,\nthe luciferase activity was decreased, further supporting the activity of RdRp. Similar results were\nalso observed in the cells co-transfected with the recombinant plasmids containing reversed\ncomplementary np1. These results were included in Fig. 5C and 5D (as shown below) in the revised\nmanuscript.\nThe following statement was included in the results of revised manuscript \" Next, we determined\nwhether that RP1 RdRp is responsible of ngRNA replication and its activity can be inhibited by\nFavipiravir. Upon co-transfection of recombinant plasmids for expressing RP1 RdRp domain and\nreversed complementary Gaussia luciferase reporter gene embedding into 5`UTR and 3`UTR of\nrp1 or reversed complementary np1 into HEK293T cells, luciferase assay and RT-qPCR analyses\nindicated that RP1 RdRp can drive the expressions of either Gaussia reporter gene or np1 gene (Fig.\n5C and 5D) (Supplementary data 1). When applying Favipiravir in transfected cells, their\nexpressions were significantly suppressed (Fig. 5C and 5D) (Supplementary data 1).\" (Line:193-\n199).\nThe related method information was also included in the revised manuscript (Line:680-703).\n12, Fig S4 is mislabeled - some of figures should be np rather than rp. Fig 7,8 labeling should say\nrp1 rather than RDRP.\nAuthor response: we corrected the label in Fig. 8 as well as Supplementary Fig. 5 in the revised\nmanuscript.\nIn response to Reviewer #1 comments\nReviewer #1 (Remarks to the Author):\nThe authors have largely addressed the reviewers concerns. There are a couple of minor points to\naddress:\n1. Line 251: flatworm rather than parasite (planarians are free living).\nAuthor response: we revised it accordingly.\n2. Should include refs to other planarian virus work PMID 30383829 and PMID: 32907972.\nAuthor response: we included these indicated references in the revised manuscript (Line: 249).\n3. Not clear from 8F why the authors conclude there are defects in epidermal integrity? Do the\nworms lyse similar to Neoblast loss?\nAuthor response: Thank you for the comment. Following RdRp dsRNA treatment, we found that\nthe worms exhibited a pronounced morphological phenotype characterized by the appearance of\nnumerous vesicles in the epidermis under microscopic examination. Prolonging the treatment to 8-\n10 days resulted in the progressive deterioration and eventual breakdown of epidermal integrity.\nConsequently, these observations demonstrate that RdRp impairment leads to essential defects in\nepidermal maintenance.\nIn response to Reviewer #2 comments\nReviewer #2 (Remarks to the Author):\nThis is a revised manuscript. Prior critiques were complementary of the signifcance of the findings\nand the overall quality of the work, although a number of specific concerns about presentation or\nexecution of various methods were provided. The authors have responded to these concerns very\nconstructively and comprehensively. The current manuscript will be an important addition to our\nknowledge of viral-like elements in the development of metazoans generally and specifically in\nparasite worms, where some potential utility as targets may be possible.\nAuthor response: Thank you for your positive comments and insightful evaluation of our\nmanuscript. We deeply appreciate your recognition of our present study."}, {"round": 2, "reviewer_comments": {"reviewer_1": "The authors have largely addressed the reviewers concerns. There are a couple of minor points to address:\nLine 251: flatworm rather than parasite (planarians are free living).\nShould include refs to other planarian virus work PMID 30383829 and PMID: 32907972.\nNot clear from 8F why the authors conclude there are defects in epidermal integrity? Do the worms lyse similar to Neoblast loss?", "reviewer_2": "This is a revised manuscript. Prior critiques were complementary of the signifcance of the findings and the overall quality of the work, although a number of specific concerns about presentation or execution of various methods were provided. The authors have responded to these concerns very constructively and comprehensively. The current manuscript will be an important addition to our knowledge of viral-like elements in the development of metazoans generally and specifically in parasite worms, where some potential utility as targets may be possible."}, "author_response": "In response to Reviewer #1 comments\nReviewer #1 (Remarks to the Author):\nThe authors have largely addressed the reviewers concerns. There are a couple of minor points to\naddress:\n1. Line 251: flatworm rather than parasite (planarians are free living).\nAuthor response: we revised it accordingly.\n2. Should include refs to other planarian virus work PMID 30383829 and PMID: 32907972.\nAuthor response: we included these indicated references in the revised manuscript (Line: 249).\n3. Not clear from 8F why the authors conclude there are defects in epidermal integrity? Do the\nworms lyse similar to Neoblast loss?\nAuthor response: Thank you for the comment. Following RdRp dsRNA treatment, we found that\nthe worms exhibited a pronounced morphological phenotype characterized by the appearance of\nnumerous vesicles in the epidermis under microscopic examination. Prolonging the treatment to 8-\n10 days resulted in the progressive deterioration and eventual breakdown of epidermal integrity.\nConsequently, these observations demonstrate that RdRp impairment leads to essential defects in\nepidermal maintenance.\nIn response to Reviewer #2 comments\nReviewer #2 (Remarks to the Author):\nThis is a revised manuscript. Prior critiques were complementary of the signifcance of the findings\nand the overall quality of the work, although a number of specific concerns about presentation or\nexecution of various methods were provided. The authors have responded to these concerns very\nconstructively and comprehensively. The current manuscript will be an important addition to our\nknowledge of viral-like elements in the development of metazoans generally and specifically in\nparasite worms, where some potential utility as targets may be possible.\nAuthor response: Thank you for your positive comments and insightful evaluation of our\nmanuscript. We deeply appreciate your recognition of our present study."}]

Decision letter