- Open Access
Differential expression analysis of mRNAs, lncRNAs, and miRNAs expression profiles and construction of ceRNA networks in PEDV infection
BMC Genomics volume 23, Article number: 586 (2022)
Porcine Epidemic Diarrhea Virus (PEDV) is a coronavirus that seriously affects the swine industry. MicroRNAs and long noncoding RNAs are two relevant non-coding RNAs (ncRNAs) class and play crucial roles in a variety of physiological processes. Increased evidence indicates a complex interaction between mRNA and ncRNA. However, our understanding of the function of ncRNA involved in host-PEDV interaction is limited.
A total of 1,197 mRNA transcripts, 539 lncRNA transcripts, and 208 miRNA transcripts were differentially regulated at 24 h and 48 h post-infection. Gene ontology (GO) and KEGG pathway enrichment analysis showed that DE mRNAs and DE lncRNAs were mainly involved in biosynthesis, innate immunity, and lipid metabolism. Moreover, we constructed a miRNA-mRNA-pathway network using bioinformatics, including 12 DE mRNAs, 120 DE miRNAs, and 11 pathways. Finally, the target genes of DE miRNAs were screened by bioinformatics, and we constructed immune-related lncRNA-miRNA-mRNA ceRNA networks. Then, the selected DE genes were validated by qRT-PCR, which were consistent with the results from RNA-Seq data.
This study provides the comprehensive analysis of the expression profiles of mRNAs, lncRNAs, and miRNAs during PEDV infection. We characterize the ceRNA networks which can provide new insights into the pathogenesis of PEDV.
Porcine epidemic diarrhea (PED) is one of the infectious intestinal diseases of swine, manifested by severe diarrhea, vomiting, and dehydration . Porcine epidemic diarrhea virus (PEDV), the etiological agent of PED, was originally identified in the United Kingdom in the 1970s . PEDV is a single-stranded positive-stranded RNA virus that is a member of the Coronavirus family. The genome of PEDV is roughly 28 kb and contains seven open reading frames (ORF), including the ORF1a, ORF1b, S, ORF3, E, M, and N genes . In 2010, the outbreak of PED caused by highly virulent PEDV strains spread in China, resulting in a mortality rate of nearly 100% in suckling piglets and huge economic losses [4,5,6]. Due to the high frequency of mutations in the PEDV genome, existing vaccines are ineffective in protecting pigs infected with highly dangerous strains . Although the pathogenesis of PEDV has been studied, the mechanisms of host-PEDV interactions are still limited.
Non-coding RNAs (ncRNAs), such as small ncRNAs (sncRNAs) and long ncRNAs (lncRNAs), are transcribed from a large percentage of host genomes . More and more studies have indicated that ncRNAs play crucial roles in many biological processes, including immune responses, post-transcriptional regulation, apoptosis, and autophagy [9,10,11]. MicroRNAs (miRNAs) are a type of abundant sncRNAs (~ 22nt long) that negatively regulate gene expression by binding to the 3’UTR of the mRNA . It is well known that miRNAs are closely related to the infection and proliferation of various viruses [13, 14]. For example, miR-30c-5p inhibits SOCS1 expression by targeting the SOCS1 3'UTR region to facilitate PEDV replication . MiR-129a-3p activates the NF-κB pathway by targeting the 3'UTR of EDA, thereby inhibiting the replication of PEDV .
LncRNAs are a class of RNAs that are more than 200 nucleotides (nt) in length and have no protein-coding capability . Compared with mRNA, lncRNA is more tissue-specific and less conservative . Several studies have shown that viral infection induces lncRNAs to regulate antiviral innate immunity . LncRNA NEAT1 exerts an antiviral effect by promoting IFN response during Hantavirus infection . On the other hand, another study has shown that lncRNA-CMPK2 acts as a negative regulator of IFN response to promote HCV replication . In recent studies, the competing endogenous RNA (ceRNA) theory has been reported as a new regulatory mechanism between lncRNA and mRNA . LncRNAs act as miRNA sponges to indirectly regulate the degradation or translational inhibition of target mRNA . So far, ceRNA theory has become the mainstream method for studying lncRNA-miRNA-mRNA interaction and has been verified by many physiological processes. Lnc-ISG20 binds to miR-326 to enhance the expression of ISG20 and inhibit IAV replication . However, there are no studies on the lncRNA-miRNA-mRNA network during PEDV infection.
Recent studies have analyzed the changes of ncRNA expression profiles during PEDV infection [1, 24, 25]. However, the regulatory mechanism of ceRNA in PEDV-infected cells remains unclear. In this study, we performed the systematic analysis of mRNA, lncRNA, and miRNA expression profiles at two different time points during PEDV infection. We analyzed differentially expressed genes, identified the target genes of DE lncRNAs and DE miRNAs, and characterized the mRNA-miRNA-pathway network and immune-related ceRNA network. Our study aims to explore the potential role of ncRNAs during PEDV infection and expand the understanding of the host-virus interactions.
Overview of RNA-sequencing
To identify the ncRNAs (lncRNA and miRNA) and protein-coding transcripts associated with viral infections, 18 RNA-seq libraries were sequenced. A total of 995 M raw reads were obtained and an average of 983.63 M clean reads per library was produced (Table S1). In addition, at least 89% of transcripts in each library were covered by reads (Table S2). Meanwhile, approximately 71.33 M clean reads were obtained from 9 small RNA sequencing libraries (Table S3). In this study, we identified 17,500 mRNAs, 5,748 lncRNAs, and 659 miRNAs. The RNA analysis suggested that, except for miRNA, the distribution of mRNA and lncRNA in each treatment group was essentially the same (Fig. 1a-c, Fig. S1). The proportion of miRNA expression < 1 was greater than 41% in the 48hpi-vs-Mock group. Following that, the properties of mRNA and lncRNA transcripts were assessed. The number of mRNA exons was noticeably greater than that of lncRNA exons. The findings revealed that mRNAs with more than 7 exons exceeded 60%, while lncRNAs mostly contained only two or three exons (Fig. 1d). The average length of these mRNAs was longer than the average length of lncRNAs, and the majority of the lncRNAs were longer than 2000 bp (Fig. 1e). Furthermore, miRNA length analysis revealed that it was commonly distributed at 22nt (Fig. 1f).
Identification and analysis of differentially expressed mRNAs
In this study, we identified 660 and 669 DE mRNAs in two groups, including 813 up-regulated (61.17%) and 516 down-regulated (38.83%) (Fig. 2a and b, Table S4). Next, the hierarchical cluster analysis revealed that the expression profiles of mRNAs changed significantly during PEDV infection, and 132 common DE mRNAs were found between the two groups through Venny 2.1 (Fig. 2c-e). GO enrichment analysis suggested that these DE mRNAs were significantly enriched in the biosynthetic process, metabolic process, cell–cell signaling, and binding (Fig. 3a and b). Meanwhile, KEGG analysis indicated that several DE mRNAs were enriched in signal transduction pathways, including the Rap1 signal, Wnt signal, cGMP-PKG signal, TGF-beta signal, and NOD-like receptor signaling pathway (Fig. 3c and d).
Identification and analysis of miRNAs
To analyze the characteristics of miRNAs during PEDV infection, miRNA expression levels were compared in different groups. We detected 43.1% mature miRNAs after comparing the clean reads to the GenBank, Rfam database, and reference genome. (Fig. 4a). A total of 616 and 641 miRNAs were detected from 24hpi-vs-Mock and 48hpi-vs-Mock, respectively. Among them, 598 miRNAs were co-expressed, while 18 and 43 miRNAs were specifically expressed in the two groups. According to fold change > 2 and P-value < 0.05, there were 3 and 206 DE miRNAs, including 100 up-regulated (31.9%) and 109 down-regulated (68.1%) were identified (Fig. 4b-d, Fig. S2a and b, Table S5). A total of 7,352 target mRNAs for DE miRNAs were identified, and we identify the potential physiological functions of DE miRNAs through GO and KEGG enrichment analyses (Fig. S3a-d).
Construction of miRNA-mRNA-pathway network
Among intersecting genes, 12 genes related to cellular immunity and cell growth along with their corresponding 120 DE miRNAs were selected (Table 1). CSF2 (Colony Stimulating Factor 2), BMP2 (Bone Morphogenetic Protein 2), WNT11 (Wnt Family Member 11), SFN (Stratifin), and TNFAIP3 (TNF Alpha Induced Protein 3) genes are involved in a variety of pathways. To visualize the interactions and further explore the function of the DE mRNAs and their corresponding DE miRNAs, the miRNA-mRNA-pathway network was constructed using the data in Table 1 (Fig. 5).
Identification and analysis of lncRNAs
In this study, CNCI, CPC, Pfam-scan, and PLEK were used to remove potential coding transcripts, and 5,769 lncRNAs were obtained. Through further filtering, we finally identified 5,748 lncRNAs, which included 37.2% intergenic lncRNAs and 40.4% intronic lncRNAs (Fig. 6a and b). Most of the lncRNA transcripts were found on chromosome 20, however, few lncRNAs were found on chromosome Y (Fig. 6c). A total of 404 and 209 DE lncRNAs were detected in 24hpi-vs-Mock and 48hpi-vs-Mock, respectively, including 269 up-regulated (31.9%) and 344 down-regulated (68.1%) (Fig. 6d-f, Fig. S2c and d, Table S6). In addition, 1,197 protein-coding genes were found to exhibit significant correlations with DE lncRNAs in expression, of which 1,169 were detected as trans-regulated target genes. As shown in Fig. S4, GO terms analysis and KEGG analysis was performed.
Construction of ceRNA network
We used miRanda (v3.3a) to identify the targets of each DE miRNA from DE mRNA and DE lncRNA transcripts. We selected 3 DE miRNAs important in the innate immunity and built a ceRNA network by incorporating miRNA-target interactions. The ceRNA network consisted of 3 DE miRNAs, 82 DE mRNAs, 32 DE lncRNAs, and 171 edges (Fig. 7). The DE mRNAs in the ceRNA network included TNFAIP3, MEFV (MEFV innate immunity regulator, pyrin), TRBV13 (T Cell Receptor Beta Variable 13), DUSP5 (Dual Specificity Phosphatase 5), ID3 (Inhibitor of DNA Binding 3, HLH Protein), TRIM54 (Tripartite Motif Containing 54), and IL2RB (Interleukin 2 Receptor Subunit Beta), etc. Functional enrichment analysis revealed that these DE mRNAs were involved in several immune-related pathways, including apoptosis, NOD-like receptor signaling, PI3K-Akt, NF-kappa B signaling, and RIG-I-like receptor signaling pathway. In addition, we also constructed an extensive ceRNA network to reveal the relationships among all DE genes (Fig. S5).
qRT-PCR validation of gene expression
To further verify the RNA-seq data, we randomly identified 6 DE mRNAs (ADSS1, GAS1, DUSP5, MEFV, IL2RB, and C1R), 6 DE lncRNAs (TCONS_00009568, TCONS_00007263, TCONS_00021035, TCONS_00024125, TCONS_00027963, and TCONS_00040950), and 6 DE miRNAs (mml-miR-30b-5p, mml-miR-1-3p, mml-miR-30c-5p, mml-miR-654-3p, mml-miR-7180-3p, and mml-miR-485-5p) by qRT-PCR. The results showed the expression pattern of these genes was consistent with those obtained by the RNA-Seq (Fig. 8).
PEDV caused huge economic losses to the livestock industry, and commercial vaccines cannot provide complete protection for pigs . Although research on PEDV has made great progress in recent decades, the pathogenic mechanism of PEDV has not been fully studied. Recent research has demonstrated that lncRNAs can act as miRNA sponges, regulating the expression of host genes and hence influencing virus replication. However, the characterization of the ceRNA mechanism in host cells during PEDV infection has not been elucidated. In this study, we detected changes in the expression profiles of mRNAs, lncRNAs, and miRNAs in PEDV-infected MARC-145 cells. The highest number of DE lncRNAs was discovered at 48 hpi, and there were more up-regulated genes were identified. These characteristics reveal that the host response to PEDV infection involves activating a large number of genes. Comparative analysis of the genomic characteristics of mRNAs and lncRNAs reveals that lncRNAs have shorter transcripts, fewer exons, and lower expression levels, which is consistent with the results of transcriptome analysis in different species [26,27,28]. In addition, the expression changes of randomly selected differentially expressed genes were confirmed by qRT-PCR analysis, indicating the accuracy of the sequencing results.
Recent studies suggest a close relationship between cholesterol metabolism and innate immunity [29,30,31]. In addition, the ATP binding cassette (ABC) transporters play a key role in many metabolic processes and contribute to resisting the penetration of xenobiotic . By analyzing these DE mRNAs, we found that the expression levels of ABCG1, ABCC3, LPIN1, DHCR7, HMGCR, and INSIG1 genes were significantly up-regulated during PEDV infection. These genes related to cholesterol metabolism were significantly up-regulated, suggesting that the host may resist PEDV infection by regulating lipid metabolism. We also studied the changes in the expression of genes associated with immunity. The sequencing results revealed that toll-like receptors (TLR2 and TLR10), some inflammatory factors (IL24, IL1RN, and IL37), and several antiviral genes such as IFIT1, ISG15, and TNK1 were differentially expressed. These up-regulated DE mRNAs were mainly focused on the immune-related pathways, such as Toll-like receptor signaling, p38 MAPK signaling, IFN signaling, and TNF signaling pathway. The data presented here may indicate that the host resists viral infection by up-regulating the genes described above in 48 h. As we all know, PEDV can use a variety of strategies to antagonize innate immune responses to ensure self-proliferation. However, no significant changes in mRNA levels of typical immune signaling factors such as MAVS, TRAF3, and IRF3 were detected during the early stages of PEDV infection. Therefore, we hypothesize that this may be one of the ways PEDV evades the host innate immune response. Further in vivo and in vitro studies are required to confirm these findings.
At present, only two studies have focused on lncRNA profile variation in PEDV-infected cells. Yu et al. analyzed PEDV-infected Vero cells at different infection time points and DE lncRNAs were screened . Chen et al. compared the lncRNA profiles in PEDV infected IPEC-J2 cells at different time points and obtained the DE lncRNAs related to immunity . There are several differences between our study and other studies. Firstly, we conducted a synthetic analysis of the mRNA, lncRNA, and miRNA expression profiles of PEDV-infected cells and performed KEGG and GO enrichment analysis. Secondly, we predicted the DE mRNAs and DE lncRNAs interacting with DE miRNAs through the database, and the ceRNA network was built. Subsequently, the enrichment analysis of DE mRNAs and DE lncRNAs showed the enrichment of cell growth, steroid synthesis, fatty acid synthesis, and other lipid metabolism pathways, which is different from previous reports.
In our study, we found that the position of lncRNA in the genome is non-random, which indicates a connection between the lncRNA function and genomic position. Ørom et al. found that a group of human lncRNAs can perform enhancer-like functions and increase the expression of neighboring protein-coding genes through cis-acting . Furthermore, another report clarified that lncRNA-mediated cis-regulation of adjacent gene transcription is common in mammalian gene regulation . Therefore, knowledge of the cis-regulation of lncRNA is essential to reveal its biological effects. In our study, several lncRNAs (TCONS_00031360, TCONS_00030751, TCONS_00005777, TCONS_00021191, etc.) were involved in the cis-regulation of autophagy-related genes during PEDV infection. Our data also indicated that some lncRNAs targeted immunity-related genes such as DDIT3, TXNIP, and FOSL1. Although lncRNA-mediated cis-regulation is ubiquitous, some lncRNAs can also function in trans. For example, the expression of several lncRNAs (TCONS_00039521, TCONS_00019044, TCONS_00024478, and TCONS_00021764) was significantly increased, and their target genes were ADA2, FADS1, MEA1, ADGRE1, MAN1C1, PNPLA3, and RND2, which were located in the trans-position.
Several studies have demonstrated the extensive involvement of miRNAs in viral pathogenesis as well as in multiple immune responses. It is universally known that IRF3 is a member of the IRF family that can induce the expression of type I IFN to regulate the antiviral innate immune response. Interestingly, our sequencing data showed that the expression levels of miR-296-3p, miR-486-3p, miR-940, and miR-1306-5p, whose target gene is IRF3, were significantly down-regulated. Therefore, we hypothesize that the host innate immune system was effectively activated to down-regulate miRNA expression in response to viral infection. Furthermore, several DE miRNAs targeting immune-related genes were identified, but with different expression trends. This finding also indicates that these miRNAs are differentially regulated and perform diverse functions during viral infection. In addition, we found that the target genes of several DE miRNAs were enriched in lipid metabolism and synthesis pathways, which demonstrated that miRNAs may be involved in the regulation of lipid.
More and more research has shown that lncRNAs can act as miRNA sponges, and the construction of ceRNA networks is a new strategy to investigate gene function. Therefore, we selected 3 DE miRNAs to build a ceRNA network to explore the exact roles and mechanisms of DE lncRNAs. Integration analysis of the lncRNA-associated ceRNA network revealed the functions of lncRNAs and their roles in cell proliferation and innate immune response. We predicted many lncRNA-miRNA-mRNA regulatory axes that may be involved in the regulation of PEDV replication, for example lnc_00037724-miR-455-5p-MEFV, lnc_00036554-miR-27b-3p-ID3, and lnc_00010748-miR-145-5p-TNFAIP3. It has been reported that MEFV regulates innate immunity by degrading multiple inflammasome components, including CASP1, NLRP1, and NLRP3 . Other studies have documented that ID3 is involved in several biological processes, including cell proliferation, senescence, and apoptosis . TNFAIP3 is a ubiquitin editing enzyme with immunosuppressive properties, preventing NF-κB activation and TNF-mediated apoptosis [38, 39]. Furthermore, several lncRNA-miRNA-mRNA axes with similar expression patterns were identified, which were enriched in different biological processes. Although the results predicted by the software need to be further investigated, our findings provide new evidence for the involvement of lncRNAs in host antiviral immunity.
In summary, our study provides the systematic description of the expression profiling of mRNAs, lncRNAs, and miRNAs. A total of 17,500 mRNAs, 5,748 lncRNAs, and 659 miRNAs were identified. In addition, we constructed a ceRNA regulatory network which contained 3 DE miRNAs, 82 DE mRNAs, and 32 DE lncRNAs based on the co-expressed transcripts. Therefore, the data acquired in this study provides a beneficial reference for understanding the interaction between the host and PEDV.
Cells culture and virus infection
MARC-145 cells were maintained in Dulbecco’s modified Eagle’s essential medium (DMEM; Hyclone, USA) supplemented with 10% fetal bovine serum (FBS; Gibco, USA) and cultured at 37 °C with 5% CO2. PEDV strain CH/SXYL/2016 (GenBank: MF462814.1) used in this study was isolated as previously described . The cells were infected with PEDV at an MOI of 1 and incubated in DMEM containing 2 μg/mL Trypsin for 1 h at 37 °C with 5% CO2. After incubation, the cells were washed with phosphate-buffered saline (PBS) and cultured in DMEM supplemented with 2% FBS.
The kinetics of PEDV replication in Marc-145 cells showed that the viral RNA levels peaked at 48 hpi and decreased thereafter .Therefore, the samples from 24 and 48 hpi were selected for RNA sequencing. MARC-145 cells were divided into three groups: cells infected with PEDV for 24 h (S24hpi_1, 2, 3), cells infected with PEDV for 48 h (S48hpi_1, 2, 3), and mock-infected cells as controls (Mock_1, 2, 3). Three biological replicates were performed. According to the manufacturer's instructions, the total RNA of PEDV-infected cells and uninfected cells was extracted using the TRIzol® reagent (Invitrogen, USA). RNA quality and concentration were tested by NanoDrop 2000 (Thermo Scientific, USA) and Agilent 2100 Bioanalyzer (Agilent Technologies, USA). The samples with RNA Integrity Numbers (RIN) ≥ 7 were subjected to the subsequent analysis.
Library construction and sequencing
A total of 1 μg RNA from each sample was used as input for RNA-seq library preparation. According to the manufacturer’s instructions, ribosomal RNA (rRNA) was removed using an Epicenter Ribo-Zero rRNA Removal Kit (Illumina, USA), and small RNA libraries were constructed using TruSeq Small RNA Sample Prep Kit (Illumina, USA). Meanwhile, sequencing libraries of mRNA and lncRNA were generated using NEBNext® Ultra™ II Directional RNA Library Prep Kit for Illumina® (NEB, USA). Thereafter, the library products were assessed for purity using an Agilent 2100 Bioanalyzer and then sequenced using the Illumina HiSeq X Ten platform by OE Biotech Co.
Analysis of differentially expressed mRNAs
First, low-quality reads were removed from the raw data using Trimomatic (v0.36)  to obtain clean reads. Then, the clean reads were mapped to the Chlorocebus sabaeus genome using HISAT2 (v2.2.1) , and the mapped reads of each sample were assembled using StringTie (v2.1.1) . FPKM  of each gene was calculated using Cufflinks (v2.2.1) , and the read counts of each gene were obtained by HTSeq-count (v0.9.1) . Analysis of DE mRNAs was performed using the DESeq R package . Fold change > 2 and P-value < 0.05 were set as the thresholds for significantly differential expression.
miRNA identification and analysis
Cutadapt (v1.14) was used to remove adapter dimers, repetitive sequences, and sequences shorter than 18 nucleotides in the raw data. The clean data was aligned by ACGT101-miR (v4.2) to eliminate RNA families (rRNA, tRNA, snRNA, snoRNA) and repeat sequences. Then, the clean reads were mapped to the miRBase database (v22.0) to identify known miRNAs and novel miRNAs. The miRNA expression was calculated and normalized as transcripts per million (TPM), and DE miRNAs were selected according to fold change > 2 and P-value < 0.05.
lncRNA identification and analysis
The raw data were first filtered by Trimomatic (v0.36) to remove low-quality reads and adapter sequences, which were then mapped to the Chlorocebus sabaeus genome using HISAT2 (v2.2.1). Afterward, transcripts that overlapped with known mRNAs and shorter than 200 bp were discarded. Next, we utilized CNCI (v1.0), CPC (0.9-r2), Pfam-scan (v30), and PLEK (v1.2) to predict the coding ability of transcripts. Transcripts with coding potential predicted by the four tools listed above were discarded, and those without coding potential were classified as lncRNAs. The target genes of lncRNAs were predicted using cis- and trans-regulation analyses . Cufflinks (v2.2.1) was used to calculate the expression level of each transcript by calculating FPKM. DE lncRNAs were selected with fold change > 2 and P-value < 0.05 by DESeq R package .
GO and KEGG Enrichment analysis
The DAVID online tool was used to perform Gene Ontology (GO) enrichment analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis of differentially expressed genes and target genes [50, 51]. GO terms and KEGG pathways with P-value < 0.05 were considered to be significantly enriched.
Construction of lncRNA-miRNA-mRNA network
Based on the expression of different mRNAs, miRNAs, and lncRNAs, Pearson’s correlation coefficient and P-value were calculated for miRNA-target pairs. The Pearson correlation coefficient ≥ 0.8 or ≤ -0.8, with P < 0.05 was considered statistically significant. MiRanda (v 3.3a) was employed to predict the relationship between miRNAs and mRNAs as well as between lncRNAs and miRNAs. The threshold parameter was set as described previously: S ≥ 150, ΔG ≤ − 30 kcal/mol and strict 5’ seed pairing . Cytoscape (v3.5.1) was used to create and visualize the lncRNA-miRNA-mRNA networks.
Quantitative real time-PCR (qRT-PCR) verification
To validate the data obtained from RNA-seq, qRT-PCR was performed for the randomly selected immune-related RNAs. Total RNA from PEDV-infected and mock-infected cells was extracted using TRIzol reagent (Invitrogen, USA), and reversely transcribed using the FastKing RT Kit (Tiangen, China). Then, qRT-PCR assays were performed with SYBR Green qPCR Master Mix (Bio-Rad, USA). The β-actin gene and U6 snRNA were selected as the internal control. The relative quantification was determined using the 2−ΔΔct method. The sequence of the primers used in this study is listed in Table S7.
GraphPad Prism 6.0 software was used for statistical analyses. Differences were determined to be statistically significant at values of P < 0.05 by Student’s t-test for two groups and one-way ANOVA analysis for more than two groups. (Data are shown as means ± standard deviation (SD). n = 3).
Availability of data and materials
All raw data has been deposited into the National Center for Biotechnology Information (https://www.ncbi.nlm.nih.gov) under the accession number PRJNA758783.
Porcine epidemic diarrhea virus
Long non-coding RNA
Small nuclear RNA
Small nucleolar RNA
Competing endogenous RNA
Kyoto encyclopedia of genes and genomes
Chen J, Wang H, Jin L, Wang L, Huang X, Chen W, et al. Profile analysis of circRNAs induced by porcine endemic diarrhea virus infection in porcine intestinal epithelial cells. Virology. 2019;527:169–79.
Wood EN. An apparently new syndrome of porcine epidemic diarrhoea. Vet Rec. 1977;100(12):243–4.
Cheun-Arom T, Temeeyasen G, Tripipat T, Kaewprommal P, Piriyapongsa J, Sukrong S, et al. Full-length genome analysis of two genetically distinct variants of porcine epidemic diarrhea virus in Thailand. Infect Genet Evol. 2016;44:114–21.
Sun RQ, Cai RJ, Chen YQ, Liang PS, Chen DK, Song CX. Outbreak of porcine epidemic diarrhea in suckling piglets. China Emerg Infect Dis. 2012;18(1):161–3.
Wang XM, Niu BB, Yan H, Gao DS, Yang X, Chen L, et al. Genetic properties of endemic Chinese porcine epidemic diarrhea virus strains isolated since 2010. Arch Virol. 2013;158(12):2487–94.
Li W, Li H, Liu Y, Pan Y, Deng F, Song Y, et al. New variants of porcine epidemic diarrhea virus, China, 2011. Emerg Infect Dis. 2012;18(8):1350–3.
Li S, Yang J, Zhu Z, Zheng H. Porcine Epidemic Diarrhea Virus and the Host Innate Immune Response. Pathogens. 2020;9(5).
Esteller M. Non-coding RNAs in human disease. Nat Rev Genet. 2011;12(12):861–74.
Hu G, Niu F, Humburg BA, Liao K, Bendi S, Callen S, et al. Molecular mechanisms of long noncoding RNAs and their role in disease pathogenesis. Oncotarget. 2018;9(26):18648–63.
Dhanoa JK, Sethi RS, Verma R, Arora JS, Mukhopadhyay CS. Long non-coding RNA: its evolutionary relics and biological implications in mammals: a review. J Anim Sci Technol. 2018;60:25.
Gong Q, Su G. Roles of miRNAs and long noncoding RNAs in the progression of diabetic retinopathy. Biosci Rep. 2017;37(6).
Kumar A, Wong AK, Tizard ML, Moore RJ, Lefevre C. miRNA_Targets: a database for miRNA target predictions in coding and non-coding regions of mRNAs. Genomics. 2012;100(6):352–6.
Lai FW, Stephenson KB, Mahony J, Lichty BD. Human coronavirus OC43 nucleocapsid protein binds microRNA 9 and potentiates NF-kappaB activation. J Virol. 2014;88(1):54–65.
Zhao X, Song X, Bai X, Fei N, Huang Y, Zhao Z, et al. miR-27b attenuates apoptosis induced by transmissible gastroenteritis virus (TGEV) infection via targeting runt-related transcription factor 1 (RUNX1). PeerJ. 2016;4: e1635.
Wang C, Shan L, Qu S, Xue M, Wang K, Fu F, et al. The Coronavirus PEDV Evades Type III Interferon Response Through the miR-30c-5p/SOCS1 Axis. Front Microbiol. 2020;11:1180.
Qi X, Cao Y, Wu S, Wu Z, Bao W. miR-129a-3p Inhibits PEDV Replication by Targeting the EDA-Mediated NF-kappaB Pathway in IPEC-J2 Cells. Int J Mol Sci. 2021;22(15).
Wang KC, Chang HY. Molecular mechanisms of long noncoding RNAs. Mol Cell. 2011;43(6):904–14.
Djebali S, Davis CA, Merkel A, Dobin A, Lassmann T, Mortazavi A, et al. Landscape of transcription in human cells. Nature. 2012;489(7414):101–8.
Fortes P, Morris KV. Long noncoding RNAs in viral infections. Virus Res. 2016;212:1–11.
Ma H, Han P, Ye W, Chen H, Zheng X, Cheng L, et al. The Long Noncoding RNA NEAT1 Exerts Antihantaviral Effects by Acting as Positive Feedback for RIG-I Signaling. J Virol. 2017;91(9).
Kambara H, Niazi F, Kostadinova L, Moonka DK, Siegel CT, Post AB, et al. Negative regulation of the interferon response by an interferon-induced long non-coding RNA. Nucleic Acids Res. 2014;42(16):10668–80.
Salmena L, Poliseno L, Tay Y, Kats L, Pandolfi PP. A ceRNA hypothesis: the Rosetta Stone of a hidden RNA language? Cell. 2011;146(3):353–8.
Chai W, Li J, Shangguan Q, Liu Q, Li X, Qi D, et al. Lnc-ISG20 Inhibits Influenza A Virus Replication by Enhancing ISG20 Expression. J Virol. 2018;92(16).
Yin L, Shen X, Yin D, Wang J, Zhao R, Dai Y, et al. Characteristics of the MicroRNA Expression Profile of Exosomes Released by Vero Cells Infected with Porcine Epidemic Diarrhea Virus. Viruses. 2022;14(4).
Chen J, Zhang C, Zhang N, Liu G. Porcine endemic diarrhea virus infection regulates long noncoding RNA expression. Virology. 2019;527:89–97.
Ning X, Sun L. Identification and characterization of immune-related lncRNAs and lncRNA-miRNA-mRNA networks of Paralichthys olivaceus involved in Vibrio anguillarum infection. BMC Genomics. 2021;22(1):447.
Li H, Cui P, Fu X, Zhang L, Yan W, Zhai Y, et al. Identification and analysis of long non-coding RNAs and mRNAs in chicken macrophages infected with avian infectious bronchitis coronavirus. BMC Genomics. 2021;22(1):67.
Li B, Yang J, He J, Peng X, Zeng Q, Song Y, et al. Characterization of the whole transcriptome of spleens from Chinese indigenous breed Ningxiang pig reveals diverse coding and non-coding RNAs for immunity regulation. Genomics. 2021;113(4):2468–82.
Mackenzie JM, Khromykh AA, Parton RG. Cholesterol manipulation by West Nile virus perturbs the cellular immune response. Cell Host Microbe. 2007;2(4):229–39.
Petersen J, Drake MJ, Bruce EA, Riblett AM, Didigu CA, Wilen CB, et al. The major cellular sterol regulatory pathway is required for Andes virus infection. PLoS Pathog. 2014;10(2): e1003911.
Blanc M, Hsieh WY, Robertson KA, Watterson S, Shui G, Lacaze P, et al. Host defense against viral infection involves interferon mediated down-regulation of sterol biosynthesis. PLoS Biol. 2011;9(3): e1000598.
Bruckmueller H, Cascorbi I. ABCB1, ABCG2, ABCC1, ABCC2, and ABCC3 drug transporter polymorphisms and their impact on drug bioavailability: what is our current understanding? Expert Opin Drug Metab Toxicol. 2021;17(4):369–96.
Yu L, Dong J, Liu Y, Zhang L, Liang P, Wang L, et al. Genome-wide analysis of long noncoding RNA profiles in Vero cells infected with porcine epidemic diarrhea virus. Arch Virol. 2020;165(9):1969–77.
Orom UA, Shiekhattar R. Long non-coding RNAs and enhancers. Curr Opin Genet Dev. 2011;21(2):194–8.
Yan P, Luo S, Lu JY, Shen X. Cis- and trans-acting lncRNAs in pluripotency and reprogramming. Curr Opin Genet Dev. 2017;46:170–8.
Kimura T, Jain A, Choi SW, Mandell MA, Schroder K, Johansen T, et al. TRIM-mediated precision autophagy targets cytoplasmic regulators of innate immunity. J Cell Biol. 2015;210(6):973–89.
Lasorella A, Benezra R, Iavarone A. The ID proteins: master regulators of cancer stem cells and tumour aggressiveness. Nat Rev Cancer. 2014;14(2):77–91.
Coornaert B, Carpentier I, Beyaert R. A20: central gatekeeper in inflammation and immunity. J Biol Chem. 2009;284(13):8217–21.
Catrysse L, Vereecke L, Beyaert R, van Loo G. A20 in inflammation and autoimmunity. Trends Immunol. 2014;35(1):22–31.
Sun P, Wu H, Huang J, Xu Y, Yang F, Zhang Q, et al. Porcine epidemic diarrhea virus through p53-dependent pathway causes cell cycle arrest in the G0/G1 phase. Virus Res. 2018;253:1–11.
Zhou H, Zhang Y, Wang J, Yan Y, Liu Y, Shi X, et al. The CREB and AP-1-Dependent Cell Communication Network Factor 1 Regulates Porcine Epidemic Diarrhea Virus-Induced Cell Apoptosis Inhibiting Virus Replication Through the p53 Pathway. Front Microbiol. 2022;13: 831852.
Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30(15):2114–20.
Kim D, Langmead B, Salzberg SL. HISAT: a fast spliced aligner with low memory requirements. Nat Methods. 2015;12(4):357–60.
Pertea M, Pertea GM, Antonescu CM, Chang TC, Mendell JT, Salzberg SL. StringTie enables improved reconstruction of a transcriptome from RNA-seq reads. Nat Biotechnol. 2015;33(3):290–5.
Roberts A, Trapnell C, Donaghey J, Rinn JL, Pachter L. Improving RNA-Seq expression estimates by correcting for fragment bias. Genome Biol. 2011;12(3):R22.
Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, et al. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol. 2010;28(5):511–5.
Anders S, Pyl PT, Huber W. HTSeq–a Python framework to work with high-throughput sequencing data. Bioinformatics. 2015;31(2):166–9.
Anders S, Huber W. Differential expression analysis for sequence count data. Genome Biol. 2010;11(10):R106.
Han D, Zhang Y, Chen J, Hua G, Li J, Deng X, et al. Transcriptome analyses of differential gene expression in the bursa of Fabricius between Silky Fowl and White Leghorn. Sci Rep. 2017;7:45959.
Kanehisa M, Goto S. KEGG: kyoto encyclopedia of genes and genomes. Nucleic Acids Res. 2000;28(1):27–30.
da Huang W, Sherman BT, Lempicki RA. Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nat Protoc. 2009;4(1):44–57.
Zhao Y, Lin Q, Li N, Babu VS, Fu X, Liu L, et al. MicroRNAs profiles of Chinese Perch Brain (CPB) cells infected with Siniperca chuatsi rhabdovirus (SCRV). Fish Shellfish Immunol. 2019;84:1075–82.
This work was supported by grants from the agricultural special fund project of Shaanxi Province (No. 2019SJNYZX45), and the innovation project for agrotechnology of Shaanxi Province (No. 2019NY-081), China.
Ethics approval and consent to participate
Consent for publication
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Table S1. Overview of the RNA sequencing data.
Table S2. Alignment information of clean reads with thedatabase.
Table S3. Characteristics of the reads from small RNAsequencing libraries.
Table S4. Differentially expressed mRNA.
Table S5. Differentially expressed miRNA.
Table S6. Differentially expressed lncRNA.
Table S7. Primers used in this study.
Figure S1. Identification of mRNA, miRNA, and lncRNAexpression profiles.
Figure S2. Characteristics of DE miRNAs and DE lncRNAsexpression levels.
Figure S3. Functional enrichment analysis of DE miRNAs.
Figure S4. Functional enrichment analysis of DE lncRNAs.
Figure S5. An overview of the top 1000 ceRNA network.
About this article
Cite this article
Shi, X., Zhang, Q., Wang, J. et al. Differential expression analysis of mRNAs, lncRNAs, and miRNAs expression profiles and construction of ceRNA networks in PEDV infection. BMC Genomics 23, 586 (2022). https://doi.org/10.1186/s12864-022-08805-0
- Functional enrichment
- Signaling pathway
- ceRNA network