https://orcid.org/0000-0003-2672-4253
ABSTRACT
The LaSota strain of Newcastle disease virus (NDV) is a commonly used vaccine to control Newcastle disease. However, improper immunization is a common reason for vaccine failure. Hence, it is imperative to thoroughly explore innate immunity-related molecular regulatory responses to the LaSota vaccine. In this text, 140 long non-coding RNAs (lncRNAs), 8 microRNAs (miRNAs), and 1514 mRNAs were identified to be differentially expressed by RNA sequencing analysis in the thymic tissues of Chinese Partridge Shank chickens after LaSota vaccine inoculation. Moreover, 70 dysregulated genes related to innate immunity were identified based on GO, Reactome pathway, and InnateDB annotations and differential expression analysis. Additionally, dysregulated lncRNAs and innate immunity-related mRNAs that could interact with dysregulated miRNAs were identified based on bioinformatics prediction analysis via the miRanda software and differential expression analysis. Among these transcripts, expression patterns of five lncRNAs, seven miRNAs, and six mRNAs were further examined by RT-qPCR assay. Both RNA-seq and RT-qPCR outcomes showed that 10 transcripts (MSTRG.22689.1, ENSGALT00000065826, ENSGALT00000059336, ENSGALT00000060887, gga-miR-6575-5p, gga-miR-6631-5p, gga-miR-1727, paraoxonase 2 (PON2), mitogen-activated protein kinase 10, and cystic fibrosis transmembrane conductance regulator (CFTR) were highly expressed, and 4 transcripts (MSTRG.188121.10, gga-miR-6655-5p, gga-miR-6548-3p, and matrix metallopeptidase 9 (MMP9) were low expressed after NDV infection. Additionally, two potential competing endogenous RNA networks (ENSGALT00000060887/gga-miR-6575-5p/PON2 or MSTRG.188121.10/gga-miR-6631-5p/MMP9) and some co-expression axes (ENSGALT00000065826/gga-miR-6631-5p, MSTRG.188121.10/gga-miR-6575-5p, MSTRG.188121.10/CFTR, ENSGALT00000060887/MMP9) were identified based on RT-qPCR and co-expression analyses. In conclusion, we identified multiple dysregulated lncRNAs, miRNAs, and mRNAs after LaSota infection and some potential regulatory networks for these dysregulated transcripts.
Introduction
Newcastle disease (ND), one of the leading devastating diseases in poultry and wild birds worldwide, has undergone several outbreaks and brought enormous economic losses for the poultry industry since its first official report 90 years ago (https://www.oie.int/en/animal-health-in-the-world/animal-diseases/newcastle-disease/) [Citation1,Citation2]. Moreover, ND is a zoonosis (an animal disease that also can infect humans) that can cause mild conjunctivitis and influenza-like symptoms in humans [Citation2,Citation3]. Newcastle Disease virus (NDV), a negative-sense single-stranded RNA virus, belongs to the genus of Avulavirus in the family of Paramyxoviridae [Citation4]. NDV strains can be categorized as lentogenic (mild), mesogenic (moderate), and velogenic (very virulent) pathotypes according to the difference of mean death time in chicken embryos [Citation5]. Lentogenic NDV strains produce mild subclinical signs with negligible mortality, while velogenic strains can cause lethal hemorrhagic lesions, serious nervous and respiratory disorders, rapid transmission, and high mortality in birds and poultry [Citation5]. Once found, mesogenic and virulent NDV viruses need to be immediately notified to the Office of International Epizootics because virulent strains of NDV are defined as the causative agent of ND by the World Health Organization for Animals (OIE) [Citation6]. To date, more than 250 bird species, especially gallinaceous birds (e.g. quail, chickens, and turkeys), have been reported with NDV infection [Citation7,Citation8].
Vaccination is an effective preventive strategy to protect animals from infectious diseases by stimulating the immune system [Citation9,Citation10]. Over the past decades, the NDV LaSota strain has been widely used as one vaccine formulation to protect poultry from virulent NDV infection and prevent ND in many countries due to its favorable immunogenicity [Citation6,Citation8,Citation11]. However, vaccine failure caused by vaccines incompleteness and/or improper immunization is a common issue [Citation11]. To better utilize vaccines and manage the vaccine failure problem, it is requisite to gain a comprehensive knowledge and understanding of molecular responses to NDV vaccines.
Coding RNAs and non-coding RNAs, including long non-coding RNAs (lncRNAs) and microRNAs (miRNAs), have emerged as crucial players in various biological processes and host responses to stress and pathogens in animals including chickens [Citation12–14]. Recently, high-throughput RNA sequencing (RNA-Seq) and microarrays that can simultaneously capture a host of transcripts and decipher the information on transcriptome (e.g. lncRNAs, miRNAs, and mRNAs) have been widely used to identify differentially expressed genes under different conditions and investigate the host-pathogen interactions [Citation15,Citation16]. Compared to the microarray technology, RNA-seq is more accurate and sensitive with low background noise [Citation17,Citation18]. Moreover, RNA-seq is independent of the availability, expression intensity, and prior annotation of probes. Additionally, RNA-seq can identify novel transcripts/low abundance transcripts/biologically critical isoforms and distinguish different isoforms/allelic expression [Citation17,Citation18]. However, RNA-seq also faces several challenges in library construction and bioinformatic analysis [Citation17]. For instance, to examine whether short reads that are identical to each other are genuine RNA, it is requisite to set up multiple biological replicates. To reduce errors in bioinformatic analysis, it is indispensable to remove the low-quality reads [Citation17].
Although previous studies have performed the transcriptome analysis to explore host immune responses against NDV LaSota vaccine infection in embryos [Citation19], spleen [Citation20], trachea [Citation21] from Leghorn and Fayoumi chickens and chick embryo fibroblasts [Citation22], these documents focused on the investigation of NDV LaSota vaccine on differentially expressed genes including some innate immune genes. To our knowledge, the relationships of lncRNAs, miRNAs, and innate immunity-related mRNAs have not been examined in thymic tissue samples of Chinese Partridge Shank chickens after LaSota vaccine injection at the transcriptomic level. In this text, to gain an in-depth insight into the molecular responses to NDV LaSota vaccine in Chinese Partridge Shank chickens, differentially expressed lncRNAs, miRNAs, and mRNAs after NDV LaSota vaccine injection were identified using RNA-seq technology. Given the vital roles of improper immune responses in vaccine failure and relative conservation of innate immune responses among different species, a total of 70 dysregulated genes related to innate immunity were screened out. Moreover, some vital innate immune genes in response to NDV infection were identified by interaction analysis. To thoroughly investigate the relationships of lncRNAs, miRNAs, and mRNAs, lncRNAs or mRNAs that could interact with dysregulated miRNAs were identified by bioinformatics prediction analysis. To further decipher the lncRNA/miRNA/mRNA regulatory networks related to innate immunity, we screened out the differentially expressed lncRNAs and innate immunity-related mRNAs that could interact with dysregulated miRNAs. Among these lncRNAs, miRNAs and mRNAs, five lncRNAs, seven miRNAs, and six innate immune genes were picked out for further RT-qPCR examination in thymic tissues of 10 chickens with or without NDV infection. Based on the co-expression relationships and competing endogenous RNA (ceRNA) hypothesis, two potential ceRNA regulatory networks (MSTRG.188121.10/gga-miR-6631-5p/matrix metallopeptidase 9 (MMP9) and ENSGALT00000060887/gga-miR-6575-5p/paraoxonase 2 (PON2)) were identified.
Materials and methods
Sample collection
Specific pathogen-free (SPF) Chinese Partridge Shank chickens were raised in biosafety level II facilities because the NDV LaSota strain was lentogenic as previously described [Citation23,Citation24]. At the age of 30 days, the challenged chickens were inoculated with 0.2 ml of 10550% egg-infectious dose (EID50) of Lasota suspension through eyes. At 0 or 48 h post-infection, chickens were euthanized, and the thymic tissue samples were collected in RNase free microtubes, immediately frozen in liquid nitrogen, and then stored at −80°C until RNA extraction. Our study was approved by the Animal Care and Use Committee of the Henan University of Animal Husbandry and Economy. All efforts were made to minimize the suffering of birds.
Total RNA extraction, cDNA library construction, and sequencing
RNA was isolated from thymic tissue samples using Trizol Reagent (Thermo Fisher Scientific, Waltham, MA, USA) and then treated with a TURBO DNA-free Kit (Thermo Fisher Scientific) to remove DNA from RNA samples (Thermo Fisher Scientific) as previously described [Citation25]. RNA sequencing (RNA-seq) samples were prepared, and cDNA libraries were constructed and sequenced as previously described [Citation26–28]. Briefly, the concentration and purity of RNA were measured using NanoDrop 2000 spectrophotometer (Thermo Fisher Scientific), and RNA quality and integrity were further tested using an RNA Nano 6000 Assay Kit (Agilent Technologies, Santa Clara, CA, USA) through Agilent Bioanalyzer 2100 (Agilent Technologies). After the removal of ribosomal RNA using the Ribo-Zero rRNA Removal Kit (Illumina), high-quality samples (RNA concentration ≥400 ng/μl, OD260/280: 1.8–2.2, RNA Integrity number ≥8) and the TruSeq Stranded Total RNA Library Prep Kit (Illumina, San Diego, CA, USA) were used for the construction of cDNA libraries of mRNAs and lncRNAs according to the instructions of the manufacturer. For the construction of cDNA libraries of RNAs containing lncRNAs and mRNAs, the remaining RNAs were purified and fragmented. Next, fragmented mRNAs were reversely transcribed into cDNAs, followed by the conversion of residual overhangs into blunt ends. Subsequently, Illumina PE adapter oligonucleotides were ligated with adenylated DNAs (3’ends). Subsequently, cDNA fragments were enriched by PCR amplification and purified through the Agencourt AMPure XP system (Beckman Coulter, Beverly, CA, USA). Next, cDNA libraries were quantified by the Agilent high-sensitivity DNA assay on Agilent Bioanalyzer 2100 (Agilent Technologies). Finally, single-strand cDNA libraries were sequenced on Illumina HiSeq 2500 instrument (Illumina) by Shanghai Personal Biotechnology Co. Ltd. (Shanghai, China). Small RNA libraries were conducted through NEBNext Multiplex Small RNA Library Prep Set for Illumina (New England Biolabs, Ipswich, MA, USA) following the instructions of the manufacturer. Briefly, RNA samples were ligated to RNA 3’ and 5’ adapters and then reverse-transcribed and amplified into cDNA libraries. After being analyzed by Agilent Bioanalyzer 2100 (Agilent Technologies), libraries were sequenced. The sequence data have been uploaded to the NCBI with the accession number PRJNA600558.
Raw data processing
The information of the reference genome (Gallus_gallus.Gallus_gallus-5.0 (assembly GCA_000002315.3) is shown in . Annotation information of Gallus_gallus.Gallus_gallus-5.0 (assembly GCA_000002315.3) genes in different databases was displayed in . Information on the raw data is presented in . After filtering, high-quality clean data without reads containing 3’ adaptors or average base quality < Q20 were aligned to the reference genome (Gallus_gallus.Gallus_gallus-5.0 (assembly GCA_000002315.3)) using Tophat2 as previously described [Citation26]. The alignment results are presented in .
Table 1. Reference genome information
Table 2. Chicken gene annotation information in different databases
Table 3. Raw data statistics
Table 4. RNAseq Map statistics
Differential expression analysis
Expression levels of lncRNAs and mRNAs were measured using the fragments per kilobase of transcript per million mapped reads (FPKM) method as previously described [Citation26]. Difference expression analysis was performed using the DESeq R package. Differentially expressed genes were screened out using the conditions of |log2FoldChange| >1 and P-value <0.05. The MA plot was drawn using the ggplot2 R package.
Functional annotation and enrichment analysis of differentially expressed genes
Genes were annotated using the Gene Ontology (GO) (http://geneontology.org/), Kyoto Encyclopedia of Genes and Genomes (KEGG) (http://www.kegg.jp/), UniProt Knowledgebase (http://www.uniprot.org/help/uniprotkb/), Enzyme Commission (EC) (http://enzyme.expasy.org/), Evolutionary Genealogy of Genes: Non-supervised Orthologous Groups (eggNOG) (http://eggnog.embl.de/version_3.0/) databases. miRNA annotation analysis was performed using the miRbase database. GO terms and KEGG pathways were regarded to be significantly enriched when the P-value was less than 0.05 [Citation29].
Target prediction and venn analysis
The potential interactions between miRNAs and lncRNAs or mRNAs were predicted by the miRanda software. Venn analysis was performed using the jvenn website (http://jvenn.toulouse.inra.fr/app/index.html).
Reverse transcription-quantitative PCR (RT‐qPCR)
RNA was extracted using Trizol Reagent (Thermo Fisher Scientific) and quantified using NanoDrop 2000 spectrophotometer (Thermo Fisher Scientific) as described above. Next, the reverse transcription and quantitative PCR reactions were performed as previously described [Citation30,Citation31]. Briefly, RNA was reversely transcribed into cDNA first strands using an iScript cDNA Synthesis Kit (Bio-Rad, Hercules, CA, USA) and random primer (for mRNAs and lncRNAs) or stem-loop miRNA RT primers in . Next, cDNA was amplified using ITaq Universal SYBR Green Supermix (Bio-Rad) and specific primers in . Expression levels of lncRNAs, mRNAs, and miRNAs were calculated using the formula of 2−ΔΔCT with GAPDH or 5S rRNA as the internal control [Citation32,Citation33].
Table 5. Reverse-transcription primers for miRNAs (stem-loop primers)
Table 6. Quantitative PCR primers
Statistical analysis
Data were analyzed using GraphPad Prism 7 software (GraphPad Software, Inc., San Diego, CA, USA). The student’s T-test was used to compare the difference between groups. The difference was regarded to be statistically significant when the P-value was less than 0.05.
Results
In this project, 140 lncRNAs, 8 miRNAs, and 1514 mRNAs were identified to be differentially expressed in thymic tissues of chickens after NDV vaccine inoculation. Improper immunization is a common cause of vaccine failure. Thus, we believed that the dysregulation of innate immune genes was closely linked with vaccine failure. In this project, 70 differentially expressed innate immune genes were screened out based on differential expression analysis, GO enrichment analysis, and annotation information in the InnateDB and Reactome Pathway databases. To identify key innate immune genes, the interaction networks of these dysregulated innate immune genes were established, and the number of neighboring genes of each node in the network was calculated. We supposed that genes with more neighboring genes might play vital roles in the innate immune responses to NDV infection. Moreover, the major pathways and biological processes that were significantly enriched by these dysregulated innate immune genes were identified by KEGG and GO enrichment analysis. To decipher the regulatory relationships of these dysregulated lncRNAs, miRNAs, and mRNAs, lncRNAs and mRNAs that could interact with differentially expressed miRNAs were predicted by the miRanda software. Moreover, differentially expressed lncRNAs and innate immune genes that could bind with dysregulated miRNAs were identified. Additionally, the expression patterns of five lncRNAs, seven miRNAs, and six innate immune genes were further examined by RT-qPCR assay in thymic tissues of 10 chickens with or without NDV infection. Furthermore, the co-expression relationships of these lncRNAs, miRNAs, and mRNAs were investigated by correlation analysis. The results showed that gga-miR-6631-5p expression was negatively associated with MSTRG.188121.10 or MMP9 expression, and MSTRG.188121.10 expression was positively correlated with MMP9 expression. Also, there was a negative regulatory relationship between gga-miR-6575-5p and ENSGALT00000060887 or PON2, and there was a positive correlation between ENSGALT00000060887 and PON2. A large body of evidence shows that lncRNAs can act as the molecular sponge of miRNAs to alleviate the inhibitory effects of miRNAs on target genes. Hence, we supposed that MSTRG.188121.10 might function as a molecular sponge of gga-miR-6631-5p to regulate MMP9 expression, and ENSGALT00000060887 might regulate PON2 expression by sequestering gga-miR-6575-5p.
Identification of differentially expressed lncRNAs
RNA-seq results showed that 140 lncRNAs were differentially expressed (|log2FoldChange| >1, P-value <0.05) in thymic tissue samples of chickens after NDV vaccine treatment versus the control group (Supplementary Table 1, ).
Figure 1. The MA plots of differentially expressed lncRNAs in the vaccine-treated group vs the control group.
Identification of differentially expressed miRNAs and potential miRNA targets
Eight miRNAs (3 up-regulated, 5 down-regulated) (|log2FoldChange| >1, P-value <0.05) were found to be differentially expressed in thymic tissues of chickens after vaccine treatment compared to the control group (Data have been shown in Supplementary Table 5 (C vs A group) in our recently published article [Citation34]). In detail, expression levels of gga-miR-6631-5p, gga-miR-6575-5p, and gga-miR-1727 were noticeably increased and expression levels of gga-miR-1712-3p, gga-miR-12273-3p, gga-miR-6655-5p, gga-miR-1744-3p, and gga-miR-6548-3p were markedly reduced in thymic tissue samples of chickens in response to NDV vaccine injection. It is well known that miRNAs can exert their functions by regulating specific target mRNAs. Hence, potential targets of differentially expressed miRNAs were identified by bioinformatics analysis using the miRanda software, and the results were presented in Supplementary Table 2.
Identification of differentially expressed mRNAs and related KEGG enrichment analysis
Moreover, 1016 up-regulated mRNAs and 498 down-regulated mRNAs (|log2FoldChange| >1, P-value <0.05) were identified in thymic tissue samples of chickens at 48 h after NDV vaccine infection compared to the control group (Data have been shown in Supplementary Table 7 in our recently published article [Citation34]). KEGG enrichment analysis revealed that these dysregulated genes in response to NDV vaccine infection might be mainly involved in the regulation of immune system-related pathways, such as systemic lupus erythematosus, Staphylococcus aureus infection, complement and coagulation cascades, and asthma (Data have been shown in Supplementary Table 11 (C vs A group) in our recently published article [Citation34]).
Identification of vaccine-influenced innate immune genes
Furthermore, 17 differentially expressed genes were identified to be significantly enriched in biological processes related to innate immunity, which was screened out through searching for the GO terms including the keyword of innate immune in GO enrichment analysis of the dysregulated genes (GO enrichment analytical results have been shown in Supplementary Table 13 sheet 1 in our recently published article [Citation34]). The information of these genes was shown in Supplementary Table 3. Given the relative conservation of innate immune responses among diverse organisms, 712 (Data have been shown in Supplementary Table 17 in our recently published article [Citation34]), and 319 (Data have been shown in Supplementary Table 19 in our recently published article [Citation34]) innate immune genes were identified from the InnateDB (http://innatedb.sahmri.com/annotatedGenes.do?type=innatedb,) and Reactome Pathway databases (http://reactome.ncpsb.org/PathwayBrowser/#/R-GGA-168249&PATH=R-GGA-168256&DTAB=MT), respectively. Combined with the differential expression data, 50 innate immune genes that were annotated by the InnateDB database (Data have been shown in Supplementary Table 18 sheet 3 in our recently published article [Citation34]) and 15 innate immune genes that were annotated by the Reactome Pathway database (Data have been shown in Supplementary Table 20 sheet 3 in our recently published article [Citation34]) were found to be differentially expressed in the post-inoculation group versus the control group. Taken together, a total of 70 differentially expressed genes related to innate immunity were screened out (Supplementary 4).
Establishment of interaction networks of proteins encoded by these filtered innate immune genes
Next, the interaction relationships of proteins that were encoded by the above-dysregulated innate immune genes were analyzed through STRING: functional protein association networks (https://string-db.org/cgi/network?taskId=blvCWmQYCH9w&sessionId=bzHnhukrXlVB). The image of interaction networks is shown in ). The information of these proteins with an interaction score ≥0.4 was presented in Supplementary Table 5. The number of neighboring genes of each node was also examined, and the results were shown in sheet 3 of Supplementary Table 5. Results showed that some genes (e.g. mitogen-activated protein kinase 10 (MAPK10), nuclear factor kappa B subunit 2 (NF-κB2), MMP9) had more neighboring genes, suggesting that these genes might play vital roles in the interaction networks and NDV-related innate immune responses. Moreover, KEGG enrichment analysis revealed that these 70 dysregulated innate immune genes were significantly enriched in retinoic acid‑inducible gene I (RIG‑I)-like receptor (RLR), NOD-like receptor (NLR), C-type lectin receptor (CLR), Toll-like receptor (TLR), and mitogen-activated protein kinase (MAPK) signaling pathways (Supplementary Table 6 sheet 1, )). As shown in Supplementary Table 6 sheet 1, MAPK10 functioned as crucial roles in RLR, NLR, CLR, TLR, and MAPK signaling pathways and NF-κB2 participated in the regulation of CLR and MAPK signaling pathways. In addition, GO enrichment analysis of these 70 dysregulated innate immune genes showed that these genes mainly participated in the regulation of biological processes, such as complement activation, inflammatory response, phagocytosis, and NIK/NF-kappaB signaling, cycle (Supplementary Table 6 sheet 2). The top 20 KEGG pathways and GO_biological_process terms that were significantly enriched by the 70 dysregulated innate immune genes are shown in , respectively.
Figure 2. Interaction and enrichment analysis of filtered dysregulated innate immune genes. (a) Interaction network of proteins encoded by filtered dysregulated innate immune genes. The confidence score was equal to or higher than 0.4. (b) The top 20 KEGG pathways enriched by the dysregulated innate immune genes. (c) The top 20 GO_biological_process terms enriched by the dysregulated innate immune genes.
Potential interaction networks of dysregulated lncRNAs, miRNAs, and innate immune genes
Moreover, lncRNAs and mRNAs that could interact with miRNAs were predicted by the miRanda software. Based on the prediction outcomes and differential expression data, we screened out differentially expressed lncRNAs that could interact with dysregulated miRNAs. The results were shown in Supplementary Table 7. No dysregulated lncRNAs that could bind to gga-miR-12273-3p were identified. Moreover, combined with the data of putative targets of miRNAs in Supplementary Table 2 and differentially expressed genes, potential dysregulated targets of gga-miR-6631-5p, gga-miR-6575-5p, and gga-miR-1727, gga-miR-1712-3p, gga-miR-12273-3p, gga-miR-6655-5p, gga-miR-1744-3p, or gga-miR-6548-3p were screened out by jvenn analysis and shown in Supplementary Table 8. Additionally, the dysregulated innate immunity-related targets of gga-miR-6575-5p, gga-miR-6631-5p, gga-miR-1727, gga-miR-6655-5p, gga-miR-6548-3p, gga-miR-1744-3p, and gga-miR-1712-3p were identified by jvenn analysis of corresponding miRNA targets in Supplementary Table 2 and filtered dysregulated innate immune genes, and the results were shown in Supplementary Table 9 and ). No common genes were found by jvenn analysis between putative gga-miR-12273-3p target group and filtered dysregulated innate immune gene group.
Figure 3. The network of differentially expressed lncRNAs, miRNAs and innate immune genes. The networks were drawn using the Cytoscape software.
RT-qPCR validation of some differentially expressed lncRNAs, miRNAs, and mRNAs
Next, five lncRNAs (MSTRG.22689.1, ENSGALT00000065826, ENSGALT00000059336, ENSGALT00000060887, and MSTRG.188121.10) with higher basemean values were selected for further RT-qPCR detection. Results showed that the expression levels of MSTRG.22689.1, ENSGALT00000065826, ENSGALT00000059336, and ENSGALT00000060887 were notably up-regulated, and MSTRG.188121.10 expression level was markedly down-regulated in thymic tissue samples of chickens at 48 h upon NDV challenge compared to the control group ()). Also, our data revealed that gga-miR-1727, gga-miR-6575-5p, and gga-miR-6631-5p were highly expressed, while gga-miR-6655-5p and gga-miR-6548-3p were low expressed in thymic tissue samples of chickens after NDV vaccine inoculation ()). However, there was no conspicuous difference in the expression levels of gga-miR-1744-3p and gga-miR-1712-3p between the control group and the post-inoculation group ()). Among the filtered genes in Supplementary Table 9, MAPK10 and cystic fibrosis transmembrane conductance regulator (CFTR), E2F transcription factor 1 (E2F1), NF-κB2, and MMP9 were screened out for further RT-qPCR exploration. Also, the expression alteration of an interested gene PON2 in response to NDV infection was examined by RT-qPCR assay. Results showed that PON2, MAPK10, and CFTR were highly expressed, and MMP9 were low expressed in the thymic tissue samples of chickens after NDV vaccine injection relative to the control group ()). However, there was no obvious alteration in the expression levels of E2F1 and NF-κB2 after NDV vaccine inoculation in the thymic tissues of chickens ()).
Figure 4. RT-qPCR validation of some differentially expressed lncRNAs, miRNAs and mRNAs. (a-c) Expression levels of several lncRNAs (MSTRG.22689.1, ENSGALT00000065826, ENSGALT00000059336, ENSGALT00000060887 and MSTRG.188121.10) (a), miRNAs (gga-miR-1744-3p, gga-miR-1712-3p gga-miR-6575-5p, gga-miR-6631-5p, gga-miR-1727, gga-miR-6655-5p and gga-miR-6548-3p) (b) and mRNAs (PON2, MAPK10, CFTR, MMP9, E2F1, and NF-κB2) (c) were further examined through RT-qPCR assay in thymic tissue samples of 10 random chickens at 0 h and 48 h post NDV injection. Note: paraoxonase 2 (PON2), mitogen-activated protein kinase 10 (MAPK10), E2F transcription factor 1 (E2F1), cystic fibrosis transmembrane conductance regulator (CFTR), matrix metallopeptidase 9 (MMP9), nuclear factor kappa B subunit 2 (NF-κB2).
Co-expression analysis of dysregulated lncRNAs, miRNAs, and mRNAs
Next, Pearson correlation analysis revealed that gga-miR-6631-5p expression level was negatively associated with the expression level of ENSGALT00000065826 or MSTRG.188121.10 in thymic tissue samples of chickens after NDV vaccine injection ()). Moreover, there was a negative correlation between gga-miR-6575-5p expression and ENSGALT00000060887 or MSTRG.188121.10 expression in thymic tissue samples of chickens after NDV vaccine challenge ()). In addition, MSTRG.188121.10 expression was positively related to MMP9 expression but was inversely correlated with CFTR expression in thymic tissue samples of chickens following NDV vaccine treatment ()). ENSGALT00000060887 expression was positively associated with PON2 or MMP9 expression in thymic tissue samples of chickens after NDV vaccine administration ()). Moreover, there was a negative correlation between gga-miR-6631-5p and MMP9 expression in thymic tissue samples of chickens after NDV vaccine administration ()). Additionally, gga-miR-6575-5p expression was inversely correlated with PON2 expression in thymic tissue samples of chickens after NDV vaccine challenge ()). Based on the co-expression data of ENSGALT00000060887, gga-miR-6575-5p, and PON2, we supposed that ENSGALT00000060887 might function as a molecular sponge of gga-miR-6575-5p to sequester gga-miR-6575-5p from PON2 in thymic tissues of chickens after NDV vaccine injection. MSTRG.188121.10 might regulate MMP9 expression through acting as a ceRNA of gga-miR-6631-5p in thymic tissues of chickens after NDV vaccine injection.
Figure 5. Pearson correlation analysis of differentially expressed lncRNAs, miRNAs and mRNAs in thymic tissue samples of 10 random chickens at 48 h after NDV inoculation. Note: matrix metallopeptidase 9 (MMP9), cystic fibrosis transmembrane conductance regulator (CFTR), paraoxonase 2 (PON2).
Discussion
In this text, expression levels of 58 lncRNAs, 3 miRNAs, and 1016 mRNAs were identified to be notably up-regulated, and expression levels of 82 lncRNAs, 5 miRNAs, and 498 mRNAs were found to be markedly down-regulated in thymic tissues of Chinese Partridge Shank chickens after LaSota NDV vaccine inoculation compared with the control group without vaccine injection.
Live vaccines are known for their stronger protective effects due to their capacity to efficiently induce a series of robust immune responses [Citation6,Citation11,Citation35]. Given the conservation of innate immune responses among different organisms [Citation36], plenty of innate immune genes were identified through GO, InnateDB and Reactome Pathway analysis. Among these innate immune genes, we found that 70 innate immune genes (31 up-regulated, 39 down-regulated) were differentially expressed in the vaccine treatment group versus the control group. KEGG enrichment analysis of the dysregulated innate immune genes disclosed that these genes were involved in the regulation of RLR, NLR, CLR, TLR, and MAPK signaling pathways, while these pathways have been identified as crucial players in innate immune responses to pathogens [Citation37–40]. For instance, the activation of TLR4 and TLR3 pathways could enhance the immune responses to the NDV vaccine (R2B-mesogenic strain) in chicken [Citation41]. RIG-I overexpression could protect cells from NDV infection [Citation42]. KEGG enrichment analysis also showed that MAPK10 participated in the regulation of RLR, NLR, CLR, TLR4, and MAPK signaling pathways and NF-κB2 were involved in CLR and MAPK signaling pathways. Moreover, the protein–protein interaction analysis of these dysregulated innate immune genes with known names showed that MAPK10, MMP9, and NF-κB2 might be crucial players in the innate immune responses to NDV infection. MAPK10 [Citation43], MMP9 [Citation44,Citation45] and NF-κB2 [Citation46,Citation47] also have been reported to be involved in the regulation of innate immune responses. These data suggested the vital roles of MAPK10, MMP9, and NF-κB2 in the immune responses to NDV infection. In this project, we demonstrated that MAPK10 was highly expressed and MMP9 was low expressed in thymic tissue samples of chickens after NDV vaccine challenge. However, RT-qPCR assay showed that there was no obvious difference in NF-κB2 expression between the NDV treatment group and the uninfected group. Additionally, both RNA-seq and RT-qPCR outcomes demonstrated that PON2 and CFTR expression levels were notably increased in response to NDV vaccine infection in thymic tissues of chickens. PON2, a member of the paraoxonase (PON) family, has been found to be implicated in innate immune responses [Citation48]. For instance, Devarajan et al. demonstrated that PON2 loss hindered bacterial clearance and facilitated macrophage phagocytosis in mice infected with Pseudomonas aeruginosa [Citation49]. CFTR also have been identified to be vital players in innate immune responses [Citation50,Citation51].
Additionally, both RNA-seq and RT-qPCR outcomes showed that expression levels of seven transcripts (i.e. MSTRG.22689.1, ENSGALT00000065826, ENSGALT00000059336, ENSGALT00000060887, gga-miR-6575-5p, gga-miR-6631-5p, and gga-miR-1727) were notably increased, and expression levels of 3 transcripts (i.e. MSTRG.188121.10, gga-miR-6655-5p, and gga-miR-6548-3p) were markedly reduced in thymic tissue samples of chickens after NDV challenge relative to the control group. Furthermore, we demonstrated that there was a negative association between gga-miR-6631-5p and ENSGALT00000065826/MSTRG.188121.10/MMP9, between gga-miR-6575-5p and ENSGALT00000060887/MSTRG.188121.10/PON2, between MSTRG.188121.10 and CFTR in thymic tissue samples of 10 chickens after NDV vaccine injection. ENSGALT00000060887 expression was positively correlated with the expression of PON2 or MMP9 in thymic tissue samples of 10 chickens after NDV vaccine challenge. MSTRG.188121.10 expression was positively associated with MMP9 expression in thymic tissue samples of 10 chickens after NDV vaccine inoculation. Recently, the ceRNA hypothesis has attracted much attention from researchers in elucidating the molecular basis of non-coding RNAs including lncRNAs [Citation52]. The ceRNA hypothesis proposes that lncRNAs can perform as the molecular sponge of miRNAs to attenuate the inhibitory effect of miRNAs on target genes [Citation53–55]. Based on the co-expression analysis, we supposed that MSTRG.188121.10 might act as a ceRNA of gga-miR-6631-5p to regulate MMP9 expression, and ENSGALT00000060887 might regulate PON2 expression by functioning as a molecular sponge of gga-miR-6575-5p in chickens after NDV vaccine inoculation.
Conclusions
In this study, 140 lncRNAs, 8 miRNAs, and 1514 mRNAs were identified to be differentially expressed in thymic tissues of Chinese Partridge Shank chickens after LaSota NDV vaccine inoculation. Moreover, 70 dysregulated innate immune genes after LaSota injection were identified. Furthermore, dysregulated lncRNAs and innate immune genes that could interact with dysregulated miRNAs were filtered out. Additionally, two potential ceRNA networks (ENSGALT00000060887/gga-miR-6575-5p/PON2 and MSTRG.188121.10/gga-miR-6631-5p/MMP9) and multiple possible lncRNA/miRNA and miRNA/mRNA regulatory axes were identified for the first time. These data could deepen our understanding of molecular responses to LaSota infection and contribute to better management of vaccine failure problems.
Highlights
1. NDV infection triggered the dysregulation of 140 lncRNAs, 8 miRNAs, and 1514 mRNAs.
2. A total of 70 innate immune genes were dysregulated after NDV infection.
3. Ten co-expression axes of dysregulated lncRNAs, miRNAs, and mRNAs were identified.
4. MSTRG.188121.10 might regulate MMP9 expression by sponging gga-miR-6631-5p.
5. ENSGALT00000060887 might regulate PON2 expression by sequestering gga-miR-6575-5p.
Supplemental Material
Download Zip (1.2 MB)Acknowledgements
I would like to express my gratitude to all those who have helped me during the writing of this thesis. I gratefully acknowledge the help of the Project of Henan Science and Technology Department that funded our research. Also, I would like to thank Furong Nie, Jingfeng Zhang, MengyunLi, Xuanniu Chang, Haitao Duan, Haoyan Li, Jia Zhou, Yudan Ji, who contributed to the research work.
Disclosure statement
No potential conflict of interest was reported by the author(s).
Data availability statement
The data displayed in this manuscript are available from the corresponding author upon reasonable request.
Supplementary material
Supplemental data for this article can be accessed here.
Additional information
Funding
References
- Rahman A, Habib M, Shabbir M. A daptation of Newcastle disease virus (NDV) in feral birds and their potential role in interspecies transmission. Open Virol J. 2018;12:52–68.
- Ashraf A, Shah M. Newcastle disease: present status and future challenges for developing countries. Afr J Microbiol Res. 2014;8(5):411–416.
- PDQ Integrative, Alternative, and Complementary Therapies Editorial Board. Newcastle disease virus (PDQ®): health professional version. 2018 Aug 22. In: PDQ cancer information summaries [Internet]. Bethesda (MD): National Cancer Institute (US); 2002. PMID: 26389195.
- AfonsoC. Virulence during Newcastle disease viruses cross species adaptation. Viruses. 2021;13:1.
- Dortmans J, Koch G, Rottier P, et al. Virulence of Newcastle disease virus: what is known so far? Vet Res. 2011;42(1):122.
- Bello M, Yusoff K, Ideris A, et al. Diagnostic and vaccination approaches for Newcastle disease virus in poultry: the current and emerging perspectives. BioMed Res Int. 2018;2018:7278459.
- Choi K. Newcastle disease virus vectored vaccines as bivalent or antigen delivery vaccines. Clin Exp Vaccine Res. 2017;6(2):72–82.
- Brown V, Bevins S. A review of virulent Newcastle disease viruses in the United States and the role of wild birds in viral persistence and spread. Vet Res. 2017;48(1):68.
- Vetter V, Denizer G, Friedland L. Understanding modern-day vaccines: what you need to know. Ann Med. 2018;50(2):110–120.
- Burakova Y, Madera R, McVey S, et al. Adjuvants for Animal Vaccines. Viral Immunol. 2018;31(1):11–22.
- Dimitrov K, AfonsoCL, Yu Q, et al. Newcastle disease vaccines—A solved problem or a continuous challenge? Vet Microbio. 2017;206:126–136.
- Kosinska-Selbi B, Mielczarek M. Review: long non-coding RNA in livestock. Animal. 2020;14(10):2003–2013.
- Duan X, WangL, Sun G, et al. Understanding the cross-talk between host and virus in poultry from the perspectives of microRNA. Poult Sci. 2020;99(4):1838–1846.
- Goel A, Ncho C, Choi Y. Regulation of gene expression in chickens by heat stress. J Anim Sci Biotechnol. 2021;12(1):11.
- Stark R, Grzelak M, Hadfield J. RNA sequencing: the teenage years. Nat Rev Genet. 2019;20(11):631–656.
- Westermann A, Barquist L, Vogel J. Resolving host-pathogen interactions by dual RNA-seq. PLoS Pathog. 2017;13(2):e1006033.
- Wang Z, Gerstein M, Snyder M. RNA-Seq: a revolutionary tool for transcriptomics. Nat Rev Genet. 2009;10(1):57.
- Rai M, Tycksen E, Sandell L, et al. Advantages of RNA-seq compared to RNA microarrays for transcriptome profiling of anterior cruciate ligament tears. J Orthop Res. 2018;36:484–497.
- Schilling M, Katani R, Memari S, et al. Transcriptional innate immune response of the developing chicken embryo to Newcastle disease virus infection. Front Genet. 2018;9:61.
- ZhangJ, Kaiser M, Deist M, et al. Transcriptome analysis in spleen reveals differential regulation of response to Newcastle disease virus in two chicken lines. Sci Rep. 2018;8(1):1278.
- Deist M, Gallardo R, Bunn D, et al. Novel mechanisms revealed in the trachea transcriptome of resistant and susceptible chicken lines following infection with Newcastle disease virus. Clin Vaccine Immunol. 2017;24(5):e00027–00017.
- LiuW, Qiu X, Song C, et al. Deep sequencing-based transcriptome profiling reveals avian interferon-stimulated genes and provides comprehensive insight into Newcastle disease virus-induced host responses. Viruses. 2018;10(4):162.
- Liu W, Tian M, Wang Y, et al. The different expression of immune-related cytokine genes in response to velogenic and lentogenic Newcastle disease viruses infection in chicken peripheral blood. Mol Biol Rep. 2012;39:3611–3618.
- Shirvani E, Varghese B, Paldurai A, et al. A recombinant avian paramyxovirus serotype 3 expressing the hemagglutinin protein protects chickens against H5N1 highly pathogenic avian influenza virus challenge. Sci Rep. 2020;10:2221.
- Trofimova I, Chervyakova D, Krasikova A. Transcription of subtelomere tandemly repetitive DNA in chicken embryogenesis. Chromosome Res. 2015;23:495–503.
- Zhang Y, LiD, Han R, et al. Transcriptome analysis of the pectoral muscles of local chickens and commercial broilers using Ribo-Zero ribonucleic acid sequencing. Plos One. 2017;12:e0184115.
- Lim W, KimK, Kim J, et al. Identification of DNA-Methylated CpG Islands associated with gene silencing in the adult body tissues of the ogye chicken using RNA-Seq and reduced representation bisulfite sequencing. Front Genet. 2019;10:346.
- YuZ, GaoX, LiuC, et al. Analysis of microRNA expression profile in specific pathogen-free chickens in response to reticuloendotheliosis virus infection. Appl Microbiol Biotechnol. 2017;101:2767–2777.
- Zhang L, Liu M, Li Q, et al. Identification of differential gene expression in endothelial cells from young and aged mice using RNA-Seq technique. Am J Transl Res. 2019;11(10):6553–6560.
- Yuan L, Wang Y, Ma X, et al. Sunflower seed oil combined with ginseng stem-leaf saponins as an adjuvant to enhance the immune response elicited by Newcastle disease vaccine in chickens. Vaccine. 2020;38:5343–5354.
- Ma X, Chi X, Yuan L, et al. Immunomodulatory effect of ginseng stem-leaf saponins and selenium on Harderian gland in immunization of chickens to Newcastle disease vaccine. Vet Immunol Immunopathol. 2020;225:110061.
- Feng Y, Chen J, Huang P, et al. Expression analysis of differentially expressed miRNAs in male and female chicken embryos. Genet Mol Res. 2014;13:3060–3068.
- Lin J, Xia J, Zhang K, et al. Genome-wide profiling of chicken dendritic cell response to infectious bursal disease. BMC Genomics. 2016;17:878.
- Guo L, Mu Z, Nie F, et al. Thymic transcriptome analysis after Newcastle disease virus inoculation in chickens and the influence of host small RNAs on NDV replication. Sci Rep. 2021;11(1):10270.
- Kapczynski D, Afonso C, Miller P. Immune responses of poultry to Newcastle disease virus. Dev Comp Immunol. 2013;41(3):447–453.
- Shokal U, Eleftherianos I. Evolution and function of thioester-containing proteins and the complement system in the innate immune response. Front Immunol. 2017;8:759.
- Creagh E, O’Neill L. TLRs, NLRs and RLRs: a trinity of pathogen sensors that co-operate in innate immunity. Trends Immunol. 2006;27:352–357.
- Chen S, Cheng A, Wang M. Innate sensing of viruses by pattern recognition receptors in birds. Vet Res. 2013;44:82.
- Bermejo-Jambrina M, Eder J, Helgers L, et al. C-Type Lectin Receptors in Antiviral Immunity and Viral Escape. Front Immunol. 2018;9:590.
- Yu J, Sun X, Goie J, et al. Regulation of Host Immune Responses against Influenza A Virus Infection by Mitogen-Activated Protein Kinases (MAPKs). Microorganisms. 2020 8 7 ;1067.
- Kannaki T, Priyanka E, Reddy M. Co-administration of toll-like receptor (TLR)-3 and 4 ligands augments immune response to Newcastle disease virus (NDV) vaccine in chicken. Vet Res Commun. 2019;43:225–230.
- Fournier P, Wilden H, Schirrmacher V. Importance of retinoic acid-inducible gene I and of receptor for type I interferon for cellular resistance to infection by Newcastle disease virus. Int J Oncol. 2012;40:287–298.
- Kenney MC, Chwa M, Atilano S, et al. Molecular and bioenergetic differences between cells with African versus European inherited mitochondrial DNA haplogroups: implications for population susceptibility to diseases. Biochim Biophys Acta. 2014;1842(2):208–219.
- Bratcher P, Weathington N, Nick H, et al. MMP-9 cleaves SP-D and abrogates its innate immune functions in vitro. PloS One. 2012;7(7):e41881.
- Hong J, Greenlee K, Lee S, et al. Anti-Microbial Function of MMP2 and MMP9 Mediated through Activation of Innate Immunity in Response to Streptococcus pneumoniae. Am J Resp Crit Care. 2009;179:A5753.
- Bonizzi G, Karin M. The two NF-κB activation pathways and their role in innate and adaptive immunity. Trends Immunol. 2004;25(6):280–288.
- Jin J, Hu H, Li H, et al. Noncanonical NF-κB pathway controls the production of type I interferons in antiviral innate immunity. Immunity. 2014;40(3):342–354.
- Shih D, Lusis A. The roles of PON1 and PON2 in cardiovascular disease and innate immunity. Curr Opin Lipidol. 2009;20:288–292.
- Devarajan A, Bourquard N, Grijalva V, et al. Role of PON2 in innate immune response in an acute infection model. Mol Genet Metab. 2013;110:362–370.
- Fiorotto R, Scirpo R, Trauner M, et al. Loss of CFTR affects biliary epithelium innate immunity and causes TLR4–NF-κB—mediated inflammatory response in mice. Gastroenterology. 2011;141(4):1498–1508.
- Vij N, Mazur S, Zeitlin P. CFTR is a negative regulator of NFκB mediated innate immune response. PloS One. 2009;4(2):e4664.
- Salmena L, Poliseno L, Tay Y, et al. A ceRNA hypothesis: the Rosetta Stone of a hidden RNA language? Cell. 2011;146(3):353–358.
- Yang C, Wu D, Gao L, et al. Competing endogenous RNA networks in human cancer: hypothesis, validation, and perspectives. Oncotarget. 2016;7(12):13479.
- Tay Y, Rinn J, Pandolfi P. The multilayered complexity of ceRNA crosstalk and competition. Nature. 2014;505(7483):344.
- Huang Z, Su G, Bi X, et al. Over-expression of long non-coding RNA insulin-like growth factor 2-antisense suppressed hepatocellular carcinoma cell proliferation and metastasis by regulating the microRNA-520h/cyclin-dependent kinase inhibitor 1A signalling pathway. Bioengineered. 2021;12:6952–6966.