Rhabdomyosarcoma (RD) cells were cultured at 37°C in 10% fetal bovine serum (FBS; Hyclone, United States) Dulbecco’s modified Eagle medium (DMEM; Hyclone, United States) in complete medium with 5% carbon dioxide (CO₂) until cell confluency reached 80–90%. The EV-D68 (Fermon strain) utilized in this work was maintained and propagated on RD cells. Briefly, fused cells were infected with EV-D68 at various multiplicities of infection (MOI) of 1, and the virus was isolated from the supernatant at the relevant infection times.

RD cells were planted in 6-well plates. When the cell density reached 60–80%, EV-D68 inoculum was supplied and incubated at 37°C with 5% CO₂. The supernatant was changed after 4 h with DMEM containing 2% FBS. Cell samples were obtained at the relevant hours based on the various experimental requirements. Viral titer assays were performed by using the TCID₅₀ method. In detail, RD cells were cultured in 96-well plates. When the cell density reached 60–80%, the previous culture media was removed and the 96-well plate was washed twice with PBS. Simultaneously, the EV-D68 virus stock solution was serially 10⁻¹ diluted from 10⁻¹⁰. At each dilution, the diluted virus was implanted in a single vertical row. The cytopathic effect (CPE) in the 96-well plate was measured every 24 h for 3–5 days, and the TCID₅₀ was estimated using Karber’s method. The formula is lgTCID₅₀ = L–D(S–0.5), where L is the highest dilution logarithm, D is the difference between dilution controls, and S is the sum of CPE-positive well ratios.

Six-well tissue culture plates were inoculated at ~1 × 10⁶ RD cells/well. When cells reach approximately 90% confluence, wash twice with PBS and incubate EV-D68 with cells at MOI of 1. Cells were collected after 0 h of infection and 24 h, and cells treated with EV-D68 for 0 h were designated as the control group. For each time point, total RNA was extracted from three different experiments independently and pooled into one set, following the manufacturer’s instructions (Invitrogen, Carlsbad, CA, United States). Total RNA was then detected and quantified using a nanotitrator and an Agilent 2,100 Bioanalyzer (Thermo Fisher Scientific, MA, United States); all samples containing RNA samples with RIN 8 and 28 s/18 s 1 were considered. cDNA segments containing aptamers were amplified by PCR and then purified using ampoules of XP magnetic beads. For quality control, libraries were analyzed on an Agilent T Technologies 2100 Bioanalyzer. The sequences of splice oligonucleotides were denatured by heating and cycling the previously disclosed double-stranded PCR products. Circular single-stranded DNA was used to create completed libraries (ssCir DNA). More than 300 DNA nanoballs (DNBs) were produced by amplification of the final library with phi29 (Thermo Fisher Scientific, MA, United States). DNBs were included in patterned nanoarrays on the BGISEQ500 platform (BGI, Shenzhen, China) to obtain double-ended 100-base readings (BGI, Shenzhen, China). After filtering the sequencing data with SOAPnuke, clean reads for lncRNAs were acquired and stored in FASTQ format (v1.5.2). HISAT2 (v2.0.4) was used to map clean readings to the reference genome. Using Ericscript (v0.5.5) and RMAT (v3.2.5), differentially spliced genes (DSG) and fusion genes were found. Bowtie2 (v2.2.5) was used to align clean reads to the gene set, a database generated by the Beijing Genome Institute (BGI) in Shenzhen that includes known and new, coding and non-coding transcripts, and then RSEM (v1.2.12) was utilized to calculate gene expression levels. The original sequencing data used in the article, and the original data processing are all consistent with the source of the published articles in our laboratory (Si et al., 2021).

To evaluate transcript expression, Cufflinks software (v2.1.1.) was used to determine transcripts per kilobase and per million mapped read fragments (FPKM) based on fragment length and read counts mapped to these fragments. FPKM 0.1 generally implies that the transcript is expressed. Following the experimental design, we screened for differentially expressed genes (DEG) between EV-D68-infected and control groups using the cufflinks software. | Log₂ (fold change) | ≥ 1 and Q. Value < 0.05 serve as indications for distinguishing between the infected and uninfected groups.

The lncRNAs regulate cellular biology by binding competitively to miRNAs as competing endogenous RNAs (ceRNAs), consequently modulating the amounts of encoded proteins (Salmena et al., 2011). The ceRNA network regulatory map was used to predict all possible interactions of the differential lncRNAs, miRNAs, and mRNAs, calculate the correlation of targeting relationships, and visualize them using Cytoscape (v3.8.2). The Coexpression Network (CN) assumes that co-expressed genomes are involved in the same pathway, regulated by the same stimuli, and have similar biological activity (Cen et al., 2022). Based on the Pearson correlation values algorithm, lncRNAs and coding genes were identified. Using modified signal intensities for different mRNA and lncRNA expression levels, co-expression networks were created and displayed with Cytoscape (v3.8.2). LncRNAs can also interact with mRNAs through base complementary pairing, and the functions of lncRNAs were investigated by targeting mRNAs of lncRNAs (Li M. et al., 2022). The software lncTar is employed to predict the targeting relationships of lncRNAs and all the differential mRNAs. All the reciprocal pairs predicted to have targeting relationships are drawn and visualized using Cytoscape (v3.8.2). In addition to interacting with some RNAs by a complementary pairing of base sequences, lncRNAs can also interact with RNA-binding protein (RBP) via specific spatial structures (Wang et al., 2020). The database RBPDB was used, all possible target proteins of the predicted top five lncRNA of up-and down-regulated differential ploidy (without differential screening). All interaction pairs with target relationships are plotted and the network diagram is constructed via Cytoscape (v3.8.2).

This study used Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways based on the Database for Annotation, Visualization, and Integrated Discovery (DAVID) to investigate the potential functions of target gene mRNAs derived from the prediction of differentially expressed lncRNAs. Gene ontology keywords give information regarding biological processes (BP), cellular components (CC), and molecular functions (MF), whereas KEGG analysis is utilized to investigate the pathways and functional categories associated with target genes.

The miR-150-5p inhibitor mimics as well as the corresponding negative controls (NC) were obtained from RIBOBIO Co., Ltd. Small interfering RNAs for SNHG9 (si-SNHG9) were purchased from RIBOBIO Co., Ltd. (Guangzhou, China; Table 1). All the plasmids, miRNA inhibitors, mimics and siRNAs were transfected into cells using Lipofectamine 2000 reagent (Invitrogen, United States) following the manufacturer’s instructions.

The Mir-X miRNA First-Strand synthesis kit (TaKaRa, China) was used for the reverse transcription of total RNA. PrimeScript^™ RT kit (TaKaRa, China) containing gDNA scavenger was used for reverse transcription of total RNA. The mRQ 3′ Primer supplied with the kit is the 3′primer for miRNA RT-qPCR. Tsingke Biotechnology Co., Ltd. manufactured qualified primers (Table 2). The expression levels of lncRNAs and mRNAs were normalized using GAPDH, whereas the expression levels of miRNAs were adjusted to U6. The 2^–△△Ct equation was used to determine the relative expression levels of various RNAs.

The total protein from treated cells was extracted by the Total Protein Extraction Kit (Thermo Fisher Scientific, Waltham, MA, United States) according to the manufacturer’s instructions. Protein concentration in the lysate was measured using the Protein BCA Assay Kit (Beyotime, Shanghai, China). On 10% sodium dodecyl sulfate-polyacrylamide gel electrophoresis, proteins were separated and then transferred to PVDF membranes for protein blotting detection. The membranes were sealed for 1 h at room temperature with 5% skim milk in Tris-buffered saline with Tween 20. The primary antibodies were incubated with EV-D68 VP1 (1:5000, GeneTex, Shanghai, China) and GAPDH (1:8000, Proteintech Group, Wuhan, China) for 1 h at 4°C, followed by horseradish peroxidase-conjugated secondary antibodies. The target protein was monitored and quantified using the Bio hypersensitive Rad ECL chemiluminescence kit (Suzhou New Cell & Molecular Biotechnology Co., Ltd., Suzhou, China). ImageJ software (National Institutes of Health) was utilized for the quantitative analysis of proteins.

A luciferase reporter plasmid containing wild-type (WT) or mutant (MUT) SNHG9/c-Fos with miR-150-5p binding site was constructed. Different promoter sequences containing WT or MUT BS were cloned into the pmiRGLO vector (Promega, Madison, WI, United States) and transfected into target cells using Lipofectamine 3,000 (Invitrogen). Samples were assayed for Luciferase activity using Promega’s Dual-Luciferase Reporter Assay System (E1910). At 48 h post-transfection, the values were read using a GloMax bioluminescence detector for 2 s.

GraphPad Prism version 8.0 (GraphPad Software, San Diego, CA, United States) was used to conduct all statistical analyses. Student’s t-test was used to compare two groups, and a one-way analysis of variance was followed by Tukey’s post hoc test to compare multiple groups. p < 0.05 was determined to be statistically significant.

As depicted in Figures 1A,B, qRT-PCR and western blot were utilized to determine the expression level of VP1 at MOI = 1 after 0, 6, 12, 18, and 24 h of RD cell infection. Both RNA and protein levels of VP1 expression increased considerably with increasing infection time, demonstrating a time dependence. In a conclusion, RD cells are an excellent model for understanding the EV-D68 infection process.

Since differentially expressed genes will be utilized as the basis for further analysis in high-throughput whole transcriptome sequencing, screening for differentially expressed genes is critical. The entire dataset (Supplementary Data Sheet 1) was filtered using the following criteria. (|Log₂(Fold change) | ≥ 1, Q. Value < 0.05). Differential expression analysis plotted in the volcano plot (Figure 2A) shows that 375 differentially expressed lncRNAs were found in the infected group compared to the control group, of which 154 were upregulated and 221 were downregulated. Mapping the chromosomal localization of these lncRNAs indicated that the majority of the highly differently expressed lncRNAs originated from the autosomal chromosome (Figure 2B). In addition, these lncRNAs were utilized for clustering analysis, which revealed the expression pattern of lncRNAs between various samples (Figure 2C). These findings imply that EV-D68 infection generates alterations in several lncRNAs and that lncRNAs co-expressed differentially during EV-D68 infection may be regulatory factors in the pathogenesis of EV-D68.

As the mode of action of lncRNAs is very complex, so the target gene prediction must be based on the various modes of action of lncRNAs.

As shown in Figure 3A, we screened a total of 25 lncRNAs, 22 miRNAs, and 143 mRNAs from differentially expressed genes (Supplementary Data Sheet 2) to develop a regulatory map for the ceRNA network. And the building of a regulatory map for the ceRNA network provides a new perspective for transcriptome research and a more comprehensive and in-depth explanation of certain biological occurrences in RD cells infected with EV-D68. As shown in Figure 3B, 8 mRNAs co-expressed with differential lncRNAs were obtained using expression profiling data (Supplementary Data Sheet 3), and coregulated lncRNAs and mRNAs of the experimental samples were mapped. These co-expressed differential mRNAs and lncRNAs may participate in similar EV-D68 infection pathways.

According to the lncRNA regulatory mechanism, the five most substantially differently expressed lncRNAs were screened, and the differential mRNAs that were complementarily associated with them were predicted (Supplementary Data Sheet 4). Many divergent mRNAs addressed lnc00235, as seen in Figure 3C, indicating that lnc00235 may be a crucial hub lncRNA in the control of EV-D68 infection. LncRNAs perform various functions by forming RNA-protein complexes with proteins such as chromosomal regulatory complexes, transcription factors, and RNP complexes. The lncRNAs with the most significant up- and down-regulation were chosen for RBP prediction (Supplementary Data Sheet 5), as shown in Figure 3D. The SR protein family (SFRS) had the closest targeting to lncRNAs. The most extensively investigated function of this RBP family is its participation in the shearing of RNA and the metabolic processes of mRNAs, including stability and exonucleation (Wan et al., 2022). This also shows that via binding to lncRNAs, SFRS may be involved in controlling EV-D68 infection.

To determine the key functions and pathways of action of differentially expressed mRNAs in the EV-D68 infection model, GO and KEGG enrichment analyses were performed on differential mRNAs screened by lncRNA regulatory mechanisms (Supplementary Data Sheet 6). First, as shown in Figure 4A, the most important pathways enriched by BP, CC, and MF in GO enrichment analysis were “Hematopoietic regulation,” “Monocyte differentiation,” “Rhodopsin-hemoglobin complex,” “Hemoglobin complex” and “DNA-binding transcriptional activator activity, RNA polymerase II specificity,” respectively.

The KEGG pathway enrichment analysis (Figure 4B) revealed that the signaling pathways were predominantly enriched (Supplementary Data Sheet 7) in the “MAPK signaling route,” “Measles,” “Hippo signaling pathway,” and “Herpes simplex virus 1 infection.” Based on the results of KEGG and GO enrichment analyses, we will also study the functions of the genes that are enriched in these pathways.

To ensure the accuracy of the sequencing results and the reliability of the subsequent data analysis, 8 lncRNAs were randomly selected from 375 differentially expressed lncRNAs for qRT-PCR validation, including 4 up-regulated lncRNAs (TALAM1, SNHG9, SLMO2-ATP5E, and LINC01970) and 4 down-regulated lncRNAs (BORCS7-ASMT, MEG8, MIR133A1HG, and ITGA6-AS1), as shown in Figures 5A,B, at MOI = 1, infected RD cells for 24 h. Moreover, as depicted in Figures 5C,D, it was discovered that the expression levels of these up-regulated genes rose dramatically with time and concentration. Figures 5E,F demonstrate that with the prolongation of infection time and the increase of MOI, the expression level of down-regulated genes decreased gradually. All these results indicate the dependability of high-throughput sequencing data and give prediction data and experimental viability for future mechanistic research.

In a preliminary study in our laboratory, c-Fos was demonstrated to have an antiviral effect on RD cells infected with EV-D68 (Si et al., 2021). Overexpression of c-Fos in EV-D68 infection resulted in a significant decrease in viral VP1 expression at both the RNA and protein levels, according to this study. Viral titer experiments also revealed that c-Fos overexpression reduced viral replication capacity. Based on the ceRNA hypothesis, we reasonably hypothesized that lncRNA SNHG9 competes as a ceRNA for c-Fos to bind miR-150-5p, resulting in reduced degradation of c-Fos and thus indirectly regulating the expression level of c-Fos to regulate viral replication. Using qRT-PCR, the expression levels of SNHG9, miR-150-5p, and c-Fos were evaluated in RD cells infected with EV-68. Figure 6A reveals that EV-D68 infection of RD cells increased the expression of SNHG9 and c-Fos, whereas it decreased the expression of miR-150-5p. Therefore, when the miR-150-5p expression is down-regulated in this infection paradigm, SNHG9 and c-Fos expression will exhibit the inverse trend. In addition, considering the reliability of lncRNA SNHG9/miR-150-5p/c-Fos axis regulation, miR-150-5p mimic and inhibitor were transfected into RD cells by miRNA transfection technique. As shown in Figure 6B, qRT-PCR results demonstrated miR-150-5p expression levels. The RNA expression of lncRNA SNHG9 and c-Fos were further examined and the results in Figure 6C showed that when miR-150-5p was overexpressed, the expression of SNHG9 and c-Fos was significantly reduced compared with the control group. In contrast, SNHG9 and c-Fos expression increased when miR-150-5p was knocked down. It is reasonable to speculate whether these phenomena influence EV-D68 replication. Figures 6D,E shows that after EV-D68 infected RD cells with MOI = 1 for 24 h, the expression level of EV-D68 VP1 RNA was significantly increased in the miR-150-5p mimic group compared to the control group. Viral titer assays similarly showed a significant increase when miR-150-5p was overexpressed. This indicates that EV-D68 virus infectivity was enhanced, and replication ability was promoted at this time. When miR-150-5p was knocked down, the opposite change was verified by qRT-PCR and virus titer experiments. The above results suggest that miR-150-5p can regulate the replication of EV-D68 by affecting the expression of c-Fos, i.e., the antiviral effect of Fos.

SNHG9 plasmid-mediated overexpression or siRNA-mediated knockdown was then utilized to confirm the involvement of SNHG9 in RD cells infected with EV-D68. Like Figure 7A, EV-D68 was infected with RD cells for 24 h at MOI = 1, and SNHG9 expression levels were validated by qRT-PCR. As shown in Figures 7B–E, when SNHG9 is overexpressed, TCID₅₀/mL decreases significantly, western blot and qRT-PCR results similarly show a decrease in viral VP1 expression. The elimination of SNHG9 has a profound effect on these variables. These data imply that SNHG9 is involved in the regulation of EV-D68 infection propagation.

To further investigate the mechanism by which SNHG9 regulates EV-D68 replication, the miR-150-5p/c-Fos axis was focused on. As demonstrated by qRT-PCR and western blot in Figures 8A–D, upon EV-D68 infection of RD cells, miR-150-5p expression was decreased when SNHG9 was overexpressed and increased when SNHG9 was knocked down. In contrast, c-Fos expression was increased in RD cells when SNHG9 was overexpressed and decreased when SNHG9 was silenced. The connection between miR-150-5p and SNHG9/c-Fos in RD cells was then studied. Figure 8E depicts miR-150-5p binding sites to SNHG9 and c-Fos. Figure 8F demonstrates, using a dual luciferase reporter, that miR-150-5p mimics dramatically reduced the relative luciferase activity of the WT-SNHG9/c-Fos group. It did not affect the MUT-SNHG9/c-Fos group. Thus, miR-150-5p was capable of binding SNHG9/c-Fos directly in RD cells. As demonstrated in Figure 8G, the miR-150-5p inhibitor greatly decreased the stimulation of viral replication in EV-D68-infected RD cells when SNHG9 was knocked down, and this impact was confirmed at the protein level (Figures 8I,J). In addition, miR-150-5p inhibitors were found to prevent SNHG9 knockdown-mediated c-Fos downregulation (Figures 8H–J). It was inferred that SNHG9 regulates viral replication through the miR-150-5p/c-Fos axis to affect EV-D68 infection of RD cells.