Skip to main content

Temporal transcriptomic changes in long non-coding RNAs and messenger RNAs involved in the host immune and metabolic response during Toxoplasma gondii lytic cycle

Abstract

Background

Long non-coding RNAs (lncRNAs) are important regulators of various biological and pathological processes, in particular the inflammatory response by modulating the transcriptional control of inflammatory genes. However, the role of lncRNAs in regulating the immune and inflammatory responses during infection with the protozoan parasite Toxoplasma gondii remains largely unknown.

Methods

We performed a longitudinal RNA sequencing analysis of human foreskin fibroblast (HFF) cells infected by T. gondii to identify differentially expressed long non-coding RNAs (lncRNAs) and messenger RNAs (mRNAs), and dysregulated pathways over the course of T. gondii lytic cycle. The transcriptome data were validated by qRT-PCR.

Results

RNA sequencing revealed significant transcriptional changes in the infected HFFs. A total of 697, 1234, 1499, 873, 1466, 561, 676 and 716 differentially expressed lncRNAs (DElncRNAs), and 636, 1266, 1843, 2303, 3022, 1757, 3088 and 2531 differentially expressed mRNAs (DEmRNAs) were identified at 1.5, 3, 6, 9, 12, 24, 36 and 48 h post-infection, respectively. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis of DElncRNAs and DEmRNAs revealed that T. gondii infection altered the expression of genes involved in the regulation of host immune response (e.g., cytokine–cytokine receptor interaction), receptor signaling (e.g., NOD-like receptor signaling pathway), disease (e.g., Alzheimer's disease), and metabolism (e.g., fatty acid degradation).

Conclusions

These results provide novel information for further research on the role of lncRNAs in immune regulation of T. gondii infection.

Graphical Abstract

Background

Toxoplasma gondii is an obligate intracellular parasite which can infect many warm-blooded animals and is widespread in the human population [1, 2]. Human infection occurs mainly via ingestion of food or water contaminated with T. gondii oocysts or cysts, respectively [3]. Infected people with a compromised immune system can experience severe symptoms, such as headache, seizures or even death [4]. Unfortunately, there is no effective human vaccine to prevent infection, and drugs used to treat infected individuals have limitations [1, 2].

Toxoplasma gondii divides by endodyogeny inside an intracellular parasitophorous vacuole (PV) [5]. When the number of tachyzoites within the PV reaches a certain level, the tachyzoites egress from infected cells and infect new host cells in a process known as the lytic cycle [6]. Following host cell invasion, the proliferating tachyzoites of T. gondii establish an acute infection. In response to the pressure conferred by the host immune response, tachyzoites differentiate to slowly replicate bradyzoites and establish a lifelong latent infection, which can be reactivated in immunodeficient hosts [7]. The interaction of T. gondii with the host cells occurs at the transcriptional and post-transcriptional levels, which is fundamental for establishing a lytic or latent infection, but the molecular mechanisms that mediate these interactions remain poorly understood.

Several approaches, including transcriptomics [8,9,10], proteomics [11,12,13], and metabolomics [14], have been used to reveal the mechanisms that mediate the interaction between T. gondii and the host cell. Long non-coding RNAs (lncRNAs) have recently received increased attention because they play vital roles in cell development [15], chromatin modification [16], and immune regulation [17]. Changes in lncRNA expression have been associated with developmental defects [18], tumorigenesis [19], and autoimmune diseases [20]. There are a few studies on the expression of lncRNAs and their role in infectious diseases, such as those caused by viruses [21, 22]. Interestingly, deregulation of the expression of 996 lncRNAs and 109 messenger RNAs (mRNAs) has been demonstrated in human foreskin fibroblast (HFF) cells infected by T. gondii by using microarray analysis [23].

In the present study, we used RNA sequencing to obtain greater insight into the temporal changes in the expression of lncRNAs and mRNAs during the lytic cycle of T. gondii infection in HFF cells.

Methods

Parasites and cell culture

The HFF cells, purchased from the American Type Culture Collection (ATCC), were cultured in Dulbecco’s modified Eagle medium (DMEM) supplemented with 15% (vol/vol) fetal bovine serum (FBS), 100 U/ml penicillin, and 100 µg/ml streptomycin. The cell culture was maintained at 37 °C and 5% CO2. The RH (Type I) strain used in the present study was maintained in HFF cells as described previously [24].

Sample collection, RNA extraction, and RNA sequencing (RNA-seq)

The HFF cells were infected with 10 × 106 freshly egressed tachyzoites using a multiplicity of infection (MOI) of 1 (1 tachyzoite to 1 HFF cell). After 1.5 h, the culture medium was removed and replaced with fresh DMEM. The cultured cells were further incubated at 37 °C with 5% CO2. Then, infected cell samples were collected at 1.5, 3, 6, 9, 12, 24, 36 and 48 h post-infection (hpi) with T. gondii. Mock-infected cells were collected at 0 h. The experiment was carried out three times. The collected cell samples were stored at −80 °C until use for RNA extraction. Total RNA was extracted using TRIzol (Invitrogen, Carlsbad, CA, USA) according to the manufacturer’s instructions. The quantity and quality of the extracted RNA were determined by using the NanoDrop spectrophotometer and Agilent 2100 Bioanalyzer (Thermo Fisher Scientific, MA, USA), respectively. Samples with an RNA integrity number (RIN) ≥ 8.0 were used for RNA sequencing. The construction of sequencing libraries and sequencing were performed at Beijing Genomics Institute in Shenzhen using the BGISE500 platform.

Data processing and differential expression analysis

The raw sequencing data were analyzed using SOAPnuke (v1.5.2) [25] to obtain high-quality reads. HISAT2 was then used to map the clean reads to the reference genome [26], followed by the application of Bowtie 2 (V2.2.5) [27] to align the clean reads to the gene set. A novel, coding and non-coding transcript database built by BGI (Beijing Genomic Institute in Shenzhen). RSEM (v1.2.12) was used to calculate gene expression levels [28], and differential expression analysis was performed using DESeq2 (v1.4.5) with Q value ≤ 0.05 [29].

Function prediction of the differentially expressed lnRNAs

Gene Ontology (GO) (http://www.geneontology.org/) and the Kyoto Encyclopedia of Genes and Genomes (KEGG) (http://www.genome.jp/kegg/) enrichment analysis of the annotated host genes corresponding to differentially expressed lnRNAs (DElncRNAs) and DEmRNAs were performed by Phyper (https://en.wikipedia.org/wike/Hypergeometric_distribution) based on the hypergeometric test. The significance of the enriched GO terms and KEGG pathways was corrected by Q value with a threshold (Q value ≤ 0.05) using the Bonferroni method [30].

Quantitative real-time PCR analysis

Validation of the differential expression of lnRNAs and mRNAs identified by RNA-seq was performed by using quantitative real-time PCR (qRT-PCR) analysis of eight randomly selected RNA samples. qRT-PCR was performed using the SYBR assay according to the manufacturer’s instructions in a total of 20 µl reaction volume, containing 1 µl of each primer, 2 µl cDNA, 6 µl nuclease-free H2O, and 10 µl master mix (Takara, Japan). Three replicates were performed for each gene, and the glyceraldehyde-3-phosphate dehydrogenase (GAPDH) gene was used as an internal control. The level of gene expression was performed by using the 2−∆∆CT method. The primers used in the study (Table 1) were designed according to the sequences available in the NCBI database and synthesized at Sangon Biotech (Shanghai, China).

Table 1 The qRT-PCR primers used in the study

Statistical analysis

All statistical analyses were performed using GraphPad Prism 7 software (GraphPad Software, Inc., San Diego, CA, USA) and a value of P < 0.05 was considered statistically significant.

Results

Identification of DElncRNAs and DEmRNAs

To characterize the dynamic expression and function of lncRNA and mRNAs in the lytic cycle of T. gondii, HFFs were collected for transcriptome analysis at 0, 1.5, 3, 6, 9, 12, 24, 36 and 48 hpi. Q value < 0.05 and |log2 fold change|> 1 were used as thresholds to identify differentially expressed genes (DEGs). The number and hierarchical clustering of DElncRNAs (Fig. 1a, c) and DEmRNAs (Fig. 1b, d) in infected cells at different time points compared with the sample at 0 h (control) are shown in Fig. 1. The DE-lncRNAs were clustered into three distinct groups as shown in the dendrogram on the top of the heatmap (Fig. 1c), indicating differences in the transcriptional profiles between the early (0 h, 1.5 h, and 3 h), middle (6 h, 9 h, and 12 h) and late (24 h, 36 h, and 48 h) stages of infection. Principal component analysis (PCA) of all identified lncRNAs (Additional file 1: Fig. S1a) and mRNAs (Additional file 1: Fig. S1b) further verified the reproducibility of these results. The volcano plots of DElncRNAs (Additional file 2: Fig. S2) and DEmRNAs (Additional file 3: Fig. S3) in infected cells at different time points compared with the sample at 0 h (control) are shown in Additional file 2: Fig. S2 and Additional file 3: Fig. S3. Our analysis identified 7,722 differently expressed lncRNAs (DElncRNAs), of which 697, 1,234, 1,499, 873, 1,466, 561, 676, and 716 lncRNAs were differently expressed at 1.5, 3, 6, 9, 12, 24, 36 and 48 hpi, respectively (Additional file 4: Table S1). A total of 16,446 differently expressed mRNAs (DEmRNAs) were identified; of which 636, 1,266, 1,843, 2,303, 3,022, 1,757, 3,088, and 2,531 were differentially expressed at 1.5, 3, 6, 9, 12, 24, 36 and 48 hpi, respectively (Additional file 5: Table S2). In addition, a total of 262 DElncRNAs (Additional file 6: Table S3) and 554 DEmRNAs (Additional file 7: Table S4) were co-expressed at all time points. The results of the RNA-seq were validated by analyzing the expression of four mRNAs (IL-11, IL-32, CCL2, and HOXA10) and four lncRNAs (LINC00941, LUCAT1, Loc107985080, and Loc101927226) at five time points (6, 12, 24, 36 and 48 hpi) using qRT-PCR. The qRT-PCR results showed that the expression profiles obtained by qRT-PCR analysis are consistent with the results of RNA-seq (Fig. 2), demonstrating the correctness of the RNA-seq data.

Fig. 1
figure 1

Differentially expressed (DE) lncRNAs (a) and mRNAs (b), and hierarchical clustering heatmaps of the DElncRNAs (c) and DEmRNAs (d) detected in HFF cells at different time points after T. gondii infection. The x-axis shows the time points after T. gondii infection and the y-axis shows the number of DElncRNAs and DEmRNAs. Red and green represent the upregulated and downregulated DElncRNAs and DEmRNAs, respectively

Fig. 2
figure 2

qRT-PCR-based validation of the expression of representative DElncRNAs (LINC00941, LUCAT1, Loc107985080 and Loc101927226) and DEmRNAs (IL-11, IL-32, CCL2 and HOXA107). The x-axis shows the names of the selected RNAs and the time points after infection. The y-axis shows the relative expression levels. The GADPH gene was used for mRNA and lncRNA normalization. qRT–PCR data are mean of three replicates ± standard deviation. Results obtained by RNA-seq analysis and qRT–PCR produced similar gene expression patterns

The gene targets of DElncRNAs

To identify the biological functions of DElncRNAs, GO enrichment analysis was performed using the gene target of DElncRNAs. The top 20 most-enriched GO terms in the biological processes are shown in Fig. 3. The biological processes of DElncRNAs at different groups were mainly involved in DNA recombination (GO:0006310), reverse transcription involved in RNA-mediated transposition (GO:0032199), transcription, DNA-templated (GO:0006351), and nucleic acid phosphodiester bond hydrolysis (GO:0090305) (Fig. 3). The top 20 most-enriched KEGG pathways of DElncRNAs are shown in Fig. 4. Most of the DElncRNAs, including CLMAT3, TARID, and CARMN, were involved in autophagy and mitophagy at different times. In addition, DElncRNAs were involved in immune-related pathways, such as the TNF signaling pathway (e.g., CFLAR-AS1, CEBPB-AS1, and LOC100130476), and the nuclear factor kappa B (NF-κB) signaling pathway (e.g., C10orf55, CFLAR-AS1, and LOC646626), and were significantly enriched at 1.5, 3 and 6 hpi. However, the metabolism-related pathways, such as fatty acid degradation (e.g., CHKB-DT) and tryptophan metabolism (e.g., TMEM161B-AS1, LOC101929352 and LOC105376957) were enriched at other time points after infection. The upregulated and downregulated lncRNAs were mainly involved in autophagy and mitophagy (Additional file 8: Fig. S4 and Additional file 9: Fig. S5).

Fig. 3
figure 3

GO enrichment analysis of the target genes of the DElncRNAs in HFF cells at different time points after T. gondii infection. ah The top 20 most-enriched GO terms in the biological process at 1.5, 3, 6, 9, 12, 24, 36 and 48 hpi, respectively. The y-axis represents the significantly enriched GO terms and the x-axis denotes the pathway enrichment and the number of DEmRNAs

Fig. 4
figure 4

KEGG pathway analysis of the target genes of the DElncRNAs in HFF cells at different time points after T. gondii infection. ah Scatterplots show the top 20 pathways at 1.5, 3, 6, 9, 12, 24, 36 and 48 hpi, respectively. The x-axis denotes the pathway enrichment. The y-axis represents the names of the significantly enriched pathways. The P-values are indicated by variations from blue to red, with darker blue denoting more significant difference

Functional analysis of the DEmRNAs

The top 20 most-enriched biological processes are shown in Fig. 5. The biological processes of DEmRNAs at different time points were mainly involved in immune- or inflammation-related GO terms, such as the immune response (GO:0002376), inflammatory response (GO:0006954), defense response to virus (GO:0051607) and cytokine-mediated signaling pathway (GO:0019221). The top 20 most-enriched pathways of DEmRNAs are shown in Fig. 6. Most of the DEmRNAs were involved in immune-related pathways, such as cytokine-cytokine receptor interaction (e.g., IL6, CXCL10 and CSF2), IL-17 signaling pathway (e.g., CCL2, NFKB1 and FOSL1), TNF signaling pathway (e.g., MAP3K8, IRF1 and ICAM1), and NOD-like receptor signaling pathway (e.g., GBP1, GBP2 and IFNAR2). The upregulated mRNAs were mainly involved in immune-related pathways (Additional file 10: Fig. S6), and the downregulated mRNAs were mainly involved in metabolism-related pathways and neurodegenerative diseases-related pathways, such as Alzheimer's disease and Huntington's disease (Additional file 11: Fig. S7).

Fig. 5
figure 5

GO enrichment analysis of the DEmRNAs in HFF cells at different time points after T. gondii infection. ah Scatterplots show the top 20 most-enriched GO terms in the biological process at 1.5, 3, 6, 9, 12, 24, 36 and 48 hpi, respectively. The y-axis represents the significantly enriched GO terms and the x-axis denotes the pathway enrichment and the number of DEmRNAs

Fig. 6
figure 6

KEGG pathway analysis of the DEmRNAs in HFF cells at different time points after T. gondii infection. ah Scatterplots show the top 20 pathways at 1.5, 3, 6, 9, 12, 24, 36 and 48 hpi, respectively. The x-axis denotes the pathway enrichment. The y-axis represents the names of the significantly enriched pathways. The P-values are indicated by variations from blue to red, with darker blue denoting more significant difference

Functional analysis of the gene targets of the co-expressed DElncRNAs

GO analysis and KEGG pathway enrichment were performed at different time points to identify the function of host genes corresponding to the DElncRNAs. The top 20 most-enriched biological processes are shown in Fig. 7a. They were mainly involved in reverse transcription linked to RNA-mediated transposition (e.g., EGOT, LYRM4-AS1 and SLCO4A1-AS1), nucleic acid phosphodiester bond hydrolysis (e.g., LYRM4-AS1, MSC-AS1 and GNG12-AS1), and DNA recombination (e.g., PLCE1-AS1, ZNF667-AS1 and LYRM4-AS1). The top 20 most-enriched pathways of DElncRNAs are shown in Fig. 7b. Most of DElncRNAs were involved in autophagy and mitophagy (e.g., LINC00173, LUCAT1 and ITGB1-DT).

Fig. 7
figure 7

GO enrichment and KEGG pathway analyses of the target genes of co-expressed DElncRNAs in HFF cells after T. gondii infection. Scatterplots show a the 20 most-enriched GO terms in the biological process of target genes of DElncRNAs and b the top 20 pathways of the target genes of DElncRNAs. The x-axis denotes the pathway enrichment. The y-axis represents the names of the significantly enriched GO terms and pathways. The P-values are indicated by variations from blue to red, with darker blue denoting more significant difference

Functional analysis of the co-expressed DEmRNAs

The top 20 most-enriched biological processes are shown in Fig. 8a. The biological processes of DEmRNAs were mainly involved in immune-related GO terms, such as the immune system process (e.g., IFI30, TRIM13, and DHX36), defense response to virus (e.g., TRIM22, IFIT5, and IFI6), and inflammatory response (e.g., NLRP3, IFI16, and IL-6). The top 20 most-enriched pathways of DEmRNAs are shown in Fig. 8b. Most DEmRNAs were involved in immune-related pathways such as cytokine-cytokine receptor interaction (e.g., CSF1, CXCL1, and IL1A), TNF signaling (e.g., FAS, IRF1, and CXCL3), and NOD-like receptor signaling (e.g., GBP4, NAMPT and BIRC2).

Fig. 8
figure 8

GO enrichment and KEGG pathway analysis of the co-expressed DEmRNAs in HFF cells after T. gondii infection. Scatterplots show a the 20 most-enriched GO terms in the biological process of DEmRNAs and b the top 20 pathways of DElncRNAs. The x-axis denotes the pathway enrichment. The y-axis represents the names of the significantly enriched GO terms and pathways. The P-values are indicated by variations from blue to red, with darker blue denoting more significant difference

Competing endogenous RNA (ceRNA) networks of lncRNAs-mRNAs

To reveal the correlations between DElncRNAs and DEmRNAs, the ceRNA network was constructed using Cytoscape v3.6.0 (2015), based on the potential target relationship between DElncRNAs and DEmRNAs. One lncRNA was associated with one or more mRNAs (Fig. 9).

Fig. 9
figure 9

The competing endogenous RNA (ceRNA) networks of lncRNAs and mRNAs. Red and green represent upregulated and downregulated genes, respectively. Gray represents both upregulated and downregulated genes at different times post-infection

As shown in Fig. 9, lncRNAs can interact with innate immune-related mRNAs and apoptotic-related mRNAs; for example, HOXA10-AS can interact with MAVS, STAT1, IL10 and BCL2L1. However, the underlying regulatory mechanisms of these DEmRNAs and DElncRNAs need to be elucidated in further studies.

Discussion

The progressive proliferation of T. gondii tachyzoites in the host tissue causes acute infection, disrupting physiological homeostasis and resulting in clinical manifestations such as fever and encephalitis [31]. Understanding how the expression of cellular lncRNAs and mRNAs is altered during T. gondii infection can provide new insight into the mechanisms underlying the early stages of the interaction between T. gondii and its host. In this study, a total of 7722 DElncRNAs and 16,446 DEmRNAs were identified, and DElncRNAs were mainly involved in immuno-inflammation and metabolism. These results suggest that lnRNAs play a significant role in the pathogenesis of T. gondii infection.

Previous studies using transcriptome sequencing detected many DEGs in T. gondii-infected animal models [31,32,33]. These DEGs were mainly involved in immune regulation. The host implements a variety of measures to counter T. gondii infection. For example, IFN-inducible GTPase families, including GBP1, GBP2, GBP3, GBP5, and GBP7, play roles in controlling T. gondii infection by regulating IFN-γ-mediated Irgb6-dependent cell-mediated immunity [34]. In the present study, we found that several genes associated with IFN-γ signaling were overexpressed during T. gondii infection, including GBP1, GBP2, GBP3, and GBP5. We also found that several genes (e.g., TLR3, STAT1, DHX58, IRF7, and Myd88) involved in the Toll-like receptor signaling pathway, NOD-like receptor signaling pathway, and RIG-I-like receptor signaling pathway were upregulated during T. gondii infection (Fig. 6). In addition, inflammation-related genes including NLRP3, IFI16, IL6, IL1A, CXCL2, and CCL6 were upregulated, further confirming that infection by T. gondii induces a robust immune response to limit infection.

Recent studies have revealed that pathogen invasion can alter the host lnRNA expression, which seems to help pathogens evade the host immune response [35, 36] or enhance the host innate immune response [37]. In the present study, GO and KEGG analysis revealed that many of the DElncRNAs play a role in immune regulation (Fig. 4). Toxoplasma gondii can inactivate human plasmacytoid dendritic cells by functional mimicry of IL-10 [38], and lncRNA TMC3-AS1 can negatively regulate the expression of IL-10 [39]. In the present study, lncRNA TMC3-AS1 was downregulated. The deletion of lncRNA LUCAT1 in myeloid cells increases the expression of type I interferon-stimulated genes in response to LPS, and overexpression of LUCAT1 reduces inducible ISG response via their interaction with STAT1 in the nucleus [40]. In the present study, the expression of LUCAT1 was upregulated. Previous studies showed that MAVS plays a role in inducing an innate immune response and many viruses employ various strategies to suppress its activity [41,42,43]. In the present study, T. gondii downregulated the expression of MAVS in infected cells. These findings indicate that T. gondii infection can manipulate the host’s immune response by altering the expression of host lncRNAs (e.g., LUCAT1 and MAVS); however, the exact functions of these lncRNAs during T. gondii infection remain to be elucidated.

Latent T. gondii infection has been associated with various psychiatric disorders in humans [44, 45]. However, the mechanism underlying this connection remains a subject of further investigation [12]. In the present study, T. gondii infection increased the expression of genes associated with Alzheimer's disease, Parkinson's disease, and Toxoplasmosis, and genes that play a role in the neurotransmitter release cycle and transmission across chemical synapses. Alterations in solute carrier family genes, such as SLC2A6, SLC1A5, and SLC22A4, have been implicated in the impairment of the synaptic transmission, which seems relevant to the pathology of schizophrenia [46]. Other studies have revealed that altered expression of DRD1, which was downregulated in the present study, may contribute to the pathophysiology of schizophrenia and affective disorders [47].

Recent studies have enriched our understanding of the metabolism of both T. gondii and its host cells [14, 48, 49]. In the present study, we have shown that T. gondii infection changes the host cell’s metabolism by altering the expression of mRNAs and lncRNAs involved in carbon metabolism, metabolic pathways, fatty acid degradation, biosynthesis of secondary metabolites, and pyruvate metabolism (Figs. 4, 6). The tricarboxylic acid (TCA) cycle is essential for T. gondii growth [48]. In the absence of glucose, T. gondii can utilize the gluconeogenic enzyme fructose bisphosphatase 2 to provide carbon for gluconeogenesis [49, 50]. In the present study, the downregulated lnRNAs were mainly involved in metabolism-related pathways (Additional file 3: Fig. S3), indicating that lncRNAs play a role in regulating host metabolism during T. gondii infection.

Conclusions

Our data revealed significant changes in lncRNA and mRNA expression in HFF cells following T. gondii infection. The identified DElncRNAs and DEmRNAs were mainly involved in metabolism, signal transduction, and immune responses. Further investigation of these DElncRNAs and DEmRNAs is warranted to reveal the exact role of these transcripts in the pathophysiology of T. gondii infection.

Availability of data and materials

The datasets supporting the findings of this article are included within the article and its additional files. The original data were deposited in the NCBI repository under accession number is PRJNA721229.

Abbreviations

lncRNAs:

Long non-coding RNAs

HFF:

Human foreskin fibroblast

DElncRNAs:

Differentially expressed lncRNAs

DEmRNAs:

Differentially expressed mRNAs

hpi:

Hours post-infection

GO:

Gene Ontology

KEGG:

Kyoto Encyclopedia of Genes and Genomes

PV:

Parasitophorous vacuole

ATCC:

American Type Culture Collection

DMEM:

Dulbecco’s modified Eagle medium

FBS:

Fetal bovine serum

MOI:

Multiplicity of infection

qRT-PCR:

Quantitative real-time PCR

GAPDH :

Glyceraldehyde-3-phosphate dehydrogenase

IL-11 :

Interleukin-11

IL-32 :

Interleukin-32

CCL2 :

C–C motif chemokine ligand 2

HOXA10:

Homeobox A10

LUCAT1 :

Lung cancer-associated transcript 1

CLMAT3 :

Colorectal liver metastasis-associated transcript 3

TARID :

Transcription factor 21 antisense RNA-inducing promoter demethylation

CARMN :

Cardiac mesoderm enhancer-associated non-coding RNA

CFLAR-AS1 :

Caspase 8 and Fas associated via death domain like apoptosis regulator antisense RNA 1

CEBPB-AS1 :

CCAAT enhancer binding protein beta antisense RNA 1

C10orf55 :

Chromosome 10 putative open reading frame 55

CHKB-DT :

Choline kinase beta divergent transcript

IL6 :

Interleukin-6

TMEM161B-AS1:

transmembrane protein 161B antisense RNA 1

CXCL10 :

C–X–C motif chemokine ligand 10

MAP3K8:

mitogen-activated protein kinase kinase kinase 8

CSF2 :

Colony-stimulating factor 2

IL-17:

Interleukin-17

NFKB1 :

Nuclear factor kappa B subunit 1

FOSL1 :

FOS-like 1

GBP1 :

Guanylate binding protein 1

GBP2 :

Guanylate binding protein 2

IFNAR2 :

Interferon alpha and beta receptor subunit 2

EGOT :

Eosinophil granule ontogeny transcript

LYRM4-AS1 :

LYR motif-containing 4 antisense RNA 1

SLCO4A1-AS1 :

Solute carrier organic anion transporter family member 4A1 antisense RNA 1

MSC-AS1 :

Musculin antisense RNA 1

GNG12-AS1 :

G protein subunit gamma 12

PLCE1-AS1 :

Phospholipase C epsilon 1 antisense RNA 1

ZNF667-AS1 :

Zinc finger protein 667 antisense RNA 1

ITGB1-DT :

Integrin subunit beta 1 divergent transcript

IFI30 :

Interferon gamma lysosomal thiol reductase

TRIM13 :

Tripartite motif-containing 13

DHX36 :

DEAH-box helicase 36

TRIM22 :

Tripartite motif-containing 22

IFIT5 :

Interferon-induced protein with tetratricopeptide repeats 5

IFI6 :

Interferon-induced protein with tetratricopeptide repeats 6

NLRP3 :

NOD-like receptor family pyrin domain-containing 3

IFI16 :

Interferon-induced protein with tetratricopeptide repeats 16

CSF1 :

Colony-stimulating factor 1

CXCL1 :

C–X–C motif chemokine ligand 1

IL1A :

Interleukin 1 alpha

FAS :

Fas cell surface death receptor

IRF1 :

Interferon regulatory factor 1

CXCL3 :

C–X–C motif chemokine ligand 3

GBP4 :

Guanylate-binding protein 4

NAMPT :

Nicotinamide phosphoribosyltransferase

BIRC2 :

Baculoviral IAP repeat-containing 2

IFN-γ:

Interferon gamma

GBP3 :

Guanylate-binding protein 3

GBP5 :

Guanylate-binding protein 5

GBP7 :

Guanylate-binding protein 7

TLR3 :

Toll-like receptor 3

STAT1 :

Signal transducer and activator of transcription 1

DHX58 :

DExH-box helicase 58

IRF7 :

Interferon regulatory factor 7

Myd88 :

Myeloid differentiation primary response gene 88

CXCL2 :

C–X–C motif chemokine ligand 2

CCL6 :

Chemokine (C–C motif) ligand 6

TMC3-AS1 :

Transmembrane channel-like 3 antisense RNA 1

LPS:

Lipopolysaccharides

ISG:

Interferon-stimulated gene

MAVS :

Mitochondrial antiviral signaling protein

SLC2A6:

Solute carrier family 2 member 6

SLC1A5 :

Solute carrier family 2 member 5

SLC22A4 :

Solute carrier family 22 member 4

DRD1:

Dopamine receptor D1

TCA:

Tricarboxylic acid

ceRNA:

Competing endogenous RNA

References

  1. Elsheikha HM, Marra CM, Zhu XQ. Epidemiology, pathophysiology, diagnosis, and management of cerebral toxoplasmosis. Clin Microbiol Rev. 2021;34:e00115-e119.

    CAS  PubMed  Google Scholar 

  2. Smith NC, Goulart C, Hayward JA, Kupz A, Miller CM, van Dooren GG. Control of human toxoplasmosis. Int J Parasitol. 2021;51:95–121.

    CAS  PubMed  Google Scholar 

  3. Hill D, Dubey JP. Toxoplasma gondii: transmission, diagnosis and prevention. Clin Microbiol Infect. 2002;8:634–40.

    CAS  PubMed  Google Scholar 

  4. Jones-Brando L, Torrey EF, Yolken R. Drugs used in the treatment of schizophrenia and bipolar disorder inhibit the replication of Toxoplasma gondii. Schizophr Res. 2003;62:237–44.

    PubMed  Google Scholar 

  5. Mordue DG, Håkansson S, Niesman I, Sibley LD. Toxoplasma gondii resides in a vacuole that avoids fusion with host cell endocytic and exocytic vesicular trafficking pathways. Exp Parasitol. 1999;92:87–99.

    CAS  PubMed  Google Scholar 

  6. Black MW, Boothroyd JC. Lytic cycle of Toxoplasma gondii. Microbiol Mol Biol Rev. 2000;64:607–23.

    CAS  PubMed  PubMed Central  Google Scholar 

  7. Tu V, Yakubu R, Weiss LM. Observations on bradyzoite biology. Microbes Infect. 2018;20:466–76.

    PubMed  Google Scholar 

  8. Jia B, Lu H, Liu Q, Yin J, Jiang N, Chen Q. Genome-wide comparative analysis revealed significant transcriptome changes in mice after Toxoplasma gondii infection. Parasit Vectors. 2013;6:161.

    CAS  PubMed  PubMed Central  Google Scholar 

  9. Tanaka S, Nishimura M, Ihara F, Yamagishi J, Suzuki Y, Nishikawa Y. Transcriptome analysis of mouse brain infected with Toxoplasma gondii. Infect Immun. 2013;81:3609–19.

    CAS  PubMed  PubMed Central  Google Scholar 

  10. Radke JB, Worth D, Hong D, Huang S, Sullivan WJ Jr, Wilson EH, et al. Transcriptional repression by ApiAP2 factors is central to chronic toxoplasmosis. PLoS Pathog. 2018;14:e1007035.

    PubMed  PubMed Central  Google Scholar 

  11. Yang J, Du F, Zhou X, Wang L, Li S, Fang R, et al. Brain proteomic differences between wild-type and CD44- mice induced by chronic Toxoplasma gondii infection. Parasitol Res. 2018;117:2623–33.

    PubMed  Google Scholar 

  12. Sahu A, Kumar S, Sreenivasamurthy SK, Selvan LD, Madugundu AK, Yelamanchi SD, et al. Host response profile of human brain proteome in toxoplasma encephalitis co-infected with HIV. Clin Proteomics. 2014;11:39.

    PubMed  PubMed Central  Google Scholar 

  13. Lv L, Wang Y, Feng W, Hernandez JA, Huang W, Zheng Y, et al. iTRAQ-based differential proteomic analysis in Mongolian gerbil brains chronically infected with Toxoplasma gondii. J Proteomics. 2017;160:74–83.

    CAS  PubMed  Google Scholar 

  14. Olson WJ, Genova B, Gallego-Lopez G, Dawson AR, Stevenson D, et al. Dual metabolomic profiling uncovers Toxoplasma manipulation of the host metabolome and the discovery of a novel parasite metabolic capability. PLoS Pathog. 2020;16:e1008432.

    CAS  PubMed  PubMed Central  Google Scholar 

  15. Peng L, Paulson A, Li H, Piekos S, He X, Li L, et al. Developmental programming of long non-coding RNAs during postnatal liver maturation in mice. PLoS ONE. 2014;9:e114917.

    PubMed  PubMed Central  Google Scholar 

  16. Angrand PO, Vennin C, Le Bourhis X, Adriaenssens E. The role of long noncoding RNAs in genome formatting and expression. Front Genet. 2015;6:165.

    PubMed  PubMed Central  Google Scholar 

  17. Imamura K, Akimitsu N. Long non-coding RNAs involved in immune responses. Front Immunol. 2014;5:573.

    PubMed  PubMed Central  Google Scholar 

  18. Grote P, Wittler L, Hendrix D, Koch F, Wahrisch S, Beisaw A, et al. The tissue-specific lncRNA Fendrr is an essential regulator of heart and body wall development in the mouse. Dev Cell. 2013;24:206–14.

    CAS  PubMed  PubMed Central  Google Scholar 

  19. Lin C, Wang Y, Wang Y, Zhang S, Yu L, Guo C, et al. Transcriptional and posttranscriptional regulation of HOXA13 by lncRNA HOTTIP facilitates tumorigenesis and metastasis in esophageal squamous carcinoma cells. Oncogene. 2017;36:5392–406.

    CAS  PubMed  Google Scholar 

  20. Satpathy Ansuman T, Chang HY. Long noncoding RNA in hematopoiesis and immunity. Immunity. 2015;42:792–804.

    CAS  PubMed  Google Scholar 

  21. Ma Y, Ouyang J, Wei J, Maarouf M, Chen JL. Involvement of host noncoding RNAs in the pathogenesis of the influenza virus. Int J Mol Sci. 2016;18:e39.

    PubMed  Google Scholar 

  22. Peng X, Gralinski L, Armour CD, Ferris MT, Thomas MJ, Proll S, et al. Unique signatures of long noncoding RNA expression in response to virus infection and altered innate immune signaling. MBio. 2010;1:00206–10.

    Google Scholar 

  23. Liu W, Huang L, Wei Q, Zhang Y, Zhang S, Zhang W, et al. Microarray analysis of long non-coding RNA expression profiles uncovers a Toxoplasma-induced negative regulation of host immune signaling. Parasit Vectors. 2018;11:174.

    PubMed  PubMed Central  Google Scholar 

  24. Zhang ZW, Li TT, Wang JL, Liang QL, Zhang HS, Sun LX, et al. Functional characterization of two thioredoxin proteins of Toxoplasma gondii using the CRISPR-Cas9 system. Front Vet Sci. 2021;7:614759.

    PubMed  PubMed Central  Google Scholar 

  25. Li R, Yu C, Li Y, Lam TW, Yiu SM, Kristiansen K, et al. SOAP2: an improved ultrafast tool for short read alignment. Bioinformatics. 2009;25:1966–7.

    CAS  PubMed  Google Scholar 

  26. Kim D, Langmead B, Salzberg SL. HISAT: a fast spliced aligner with low memory requirements. Nat Methods. 2015;12:357–60.

    CAS  PubMed  PubMed Central  Google Scholar 

  27. Langmead B, Salzberg SL. Fast gapped-read alignment with Bowtie 2. Nat Methods. 2012;9:357–9.

    CAS  PubMed  PubMed Central  Google Scholar 

  28. Li B, Dewey CN. RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome. BMC Bioinform. 2011;12:323.

    CAS  Google Scholar 

  29. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:550.

    PubMed  PubMed Central  Google Scholar 

  30. Abdi H. The Bonferonni and Sˇidák Corrections for multiple comparisons. Encycl Meas Stat. 2007;1:1–9.

    Google Scholar 

  31. Zhou CX, Ai K, Huang CQ, Guo JJ, Cong H, He SY, et al. miRNA and circRNA expression patterns in mouse brain during toxoplasmosis development. BMC Genomics. 2020;21:46.

    PubMed  PubMed Central  Google Scholar 

  32. Hermes G, Ajioka JW, Kelly KA, Mui E, Roberts F, Kasza K, et al. Neurological and behavioral abnormalities, ventricular dilatation, altered cellular functions, inflammation, and neuronal injury in brains of mice due to common, persistent, parasitic infection. J Neuroinflammation. 2008;5:48.

    PubMed  PubMed Central  Google Scholar 

  33. Wang M, Zhang FK, Elsheikha HM, Zhang NZ, He JJ, Luo JX, et al. Transcriptomic insights into the early host-pathogen interaction of cat intestine with Toxoplasma gondii. Parasit Vectors. 2018;11:592.

    CAS  PubMed  PubMed Central  Google Scholar 

  34. Yamamoto M, Okuyama M, Ma JS, Kimura T, Kamiyama N, Saiga H, et al. A cluster of interferon-γ-inducible p65 GTPases plays a critical role in host defense against Toxoplasma gondii. Immunity. 2012;37:302–13.

    CAS  PubMed  Google Scholar 

  35. Wang Y, Wang P, Zhang Y, Xu J, Li Z, Li Z, et al. Decreased expression of the host long-noncoding RNA-GM facilitates viral escape by inhibiting the kinase activity TBK1 via S-glutathionylation. Immunity. 2020;53:1168-1181.e7.

    CAS  PubMed  Google Scholar 

  36. Heaton NS, Cullen BR. Viruses hijack a long non-coding RNA. Nature. 2017;552:184–5.

    CAS  PubMed  Google Scholar 

  37. Zhou Y, Li M, Xue Y, Li Z, Wen W, Liu X, et al. Interferon-inducible cytoplasmic lncLrrc55-AS promotes antiviral innate responses by strengthening IRF3 phosphorylation. Cell Res. 2019;29:641–54.

    CAS  PubMed  PubMed Central  Google Scholar 

  38. Pierog PL, Zhao Y, Singh S, Dai J, Yap GS, Fitzgerald-Bocarsly P. Toxoplasma gondii inactivates human plasmacytoid dendritic cells by functional mimicry of IL-10. J Immunol. 2018;200:186–95.

    CAS  PubMed  Google Scholar 

  39. Ye M, Xie M, Zhu J, Wang C, Zhou R, Li X. LPS-inducible lncRNA TMC3-AS1 negatively regulates the expression of IL-10. Front Immunol. 2020;11:1418.

    CAS  PubMed  PubMed Central  Google Scholar 

  40. Agarwal S, Vierbuchen T, Ghosh S, Chan J, Jiang Z, Kandasamy RK, et al. The long non-coding RNA LUCAT1 is a negative feedback regulator of interferon responses in humans. Nat Commu. 2020;11:6348.

    CAS  Google Scholar 

  41. Zhang W, Wang Q, Yang F, Zhu Z, Duan Y, Yang Y, et al. JMJD6 negatively regulates cytosolic RNA induced antiviral signaling by recruiting RNF5 to promote activated IRF3 K48 ubiquitination. PLoS Pathog. 2021;17:e1009366.

    CAS  PubMed  PubMed Central  Google Scholar 

  42. Jiang X, Tang T, Guo J, Wang Y, Li P, Chen X, et al. Human herpesvirus 6B U26 inhibits the activation of the RLR/MAVS signaling pathway. MBio. 2021;12:e03505-e3520.

    CAS  PubMed  PubMed Central  Google Scholar 

  43. Zeng Y, Xu S, Wei Y, Zhang X, Wang Q, Jia Y, et al. The PB1 protein of influenza A virus inhibits the innate immune response by targeting MAVS for NBR1-mediated selective autophagic degradation. PLoS Pathog. 2021;17:e1009300.

    CAS  PubMed  PubMed Central  Google Scholar 

  44. Elsheikha HM, Büsselberg D, Zhu XQ. The known and missing links between Toxoplasma gondii and schizophrenia. Metab Brain Dis. 2016;31:749–59.

    PubMed  Google Scholar 

  45. Elsheikha HM, Zhu XQ. Toxoplasma gondii infection and schizophrenia: an inter-kingdom communication perspective. Curr Opin Infect Dis. 2016;29:311–8.

    CAS  PubMed  Google Scholar 

  46. Eastwood SL, Harrison PJ. Decreased expression of vesicular glutamate transporter 1 and complexin II mRNAs in schizophrenia: further evidence for a synaptic pathology affecting glutamate neurons. Schizophr Res. 2005;73:159–72.

    CAS  PubMed  Google Scholar 

  47. Kaalund SS, Newburn EN, Ye T, Tao R, Li C, Deep-Soboslay A, et al. Contrasting changes in DRD1 and DRD2 splice variant expression in schizophrenia and affective disorders, and associations with SNPs in postmortem brain. Mol Psychiatry. 2014;19:1258–66.

    CAS  PubMed  Google Scholar 

  48. MacRae JI, Sheiner L, Nahid A, Tonkin C, Striepen B, McConville MJ. Mitochondrial metabolism of glucose and glutamine is required for intracellular growth of Toxoplasma gondii. Cell Host Microbe. 2012;12:682–92.

    CAS  PubMed  PubMed Central  Google Scholar 

  49. Blume M, Nitzsche R, Sternberg U, Gerlic M, Masters SL, Gupta N, et al. A Toxoplasma gondii gluconeogenic enzyme contributes to robust central carbon metabolism and is essential for replication and virulence. Cell Host Microbe. 2015;18:210–20.

    CAS  PubMed  Google Scholar 

  50. Blume M, Rodriguez-Contreras D, Landfear S, Fleige T, Soldati-Favre D, Lucius R, et al. Host-derived glucose and its transporter in the obligate intracellular pathogen Toxoplasma gondii are dispensable by glutaminolysis. Proc Natl Acad Sci U S A. 2009;106:12998–3003.

    CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

We would like to thank BGI-Shenzhen for technical assistance.

Funding

Project support was provided by the Agricultural Science and Technology Innovation Program (ASTIP) of China (Grant No. CAAS-ASTIP-2016-LVRI-03), the Yunnan Expert Workstation (Grant No. 202005AF150041), Veterinary Public Health Innovation Team of Yunnan Province (Grant No. 202105AE160014), the Fund for Shanxi “1331 Project” (Grant No. 20211331-13), and the Special Research Fund of Shanxi Agricultural University for High-level Talents (Grant No. 2021XG001).

Author information

Authors and Affiliations

Authors

Contributions

GHZ, CXZ, HME, and XQZ conceived and designed the study. SSW performed the experiments, analyzed the data, and wrote the manuscript. CXZ participated in the implementation of the study. JJH, WBZ, and FCZ contributed reagents/materials/analysis tools. HME, JJH, WBZ, GHZ, and XQZ critically revised the manuscript. All authors read and approved the final manuscript.

Corresponding authors

Correspondence to Xing-Quan Zhu or Guang-Hui Zhao.

Ethics declarations

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests. The co-corresponding author Prof Xing-Quan Zhu serves as the Subject Editor for the section “Parasite genetics, genomics and proteomics” of Parasites & Vectors.

Additional information

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Figure S1.

Principal component analysis (PCA) of all identified lncRNAs (a) and mRNAs (b).

Additional file 2: Figure S2.

Volcano plots showing the differentially expressed lncRNAs at 1.5, 3, 6, 9, 12, 24, 36, and 48 hpi, respectively (a-h). The negative log10-transformed P-values (y-axis) are plotted against the average log2 fold changes in expression (x-axis). Data points representing lncRNAs that were not differentially expressed are shown in black. Transcripts that are differentially expressed with an absolute |log2 fold change (FC)| more than or less than 1 are shown as red (upregulated) and green (downregulated) dots.

Additional file 3: Figure S3.

Volcano plots showing the differentially expressed mRNAs at 1.5, 3, 6, 9, 12, 24, 36, and 48 hpi, respectively (a–h). The negative log10-transformed P-values (y-axis) are plotted against the average log2 fold changes in expression (x-axis). Data points representing mRNAs that were not differentially expressed are shown in black. Transcripts that are differentially expressed with an absolute |log2 fold change (FC)| more than or less than 1 are shown as red (upregulated) and green (downregulated) dots.

Additional file 4: Table S1.

List of all differentially expressed messenger RNAs (DEmRNAs).

Additional file 5: Table S2.

List of all differentially expressed long non-coding RNAs (DElncRNAs).

Additional file 6: Table S3.

List of co-expressed differentially expressed messenger RNAs (DEmRNAs).

Additional file 7: Table S4.

List of co-expressed differentially expressed long non-coding RNA (DElncRNAs).

Additional file 8: Figure S4.

KEGG pathway analysis of the target genes of the upregulated DElncRNAs at different time points after infection of HFF cells by T. gondii. (ah) Scatterplots show the top 20 pathways at 1.5, 3, 6, 9, 12, 24, 36, and 48 hpi, respectively. The x-axis denotes the pathway enrichment. The y-axis shows the names of the significantly enriched pathways. The P-values are indicated by variations from blue to red, with darker blue indicating more significant difference.

Additional file 9: Figure S5.

KEGG pathway analysis of the target genes of the downregulated DElncRNAs in HFF cells at different time points after T. gondii infection. (ah) Scatterplots show the top 20 pathways at 1.5, 3, 6, 9, 12, 24, 36, and 48 hpi, respectively. The x-axis denotes the pathway enrichment. The y-axis shows the names of the significantly enriched pathways. The P-values are indicated by variations from blue to red, with darker blue indicating more significant difference.

Additional file 10: Figure S6.

KEGG pathway analysis of the upregulated mRNAs in HFF cells at different time points after T. gondii infection. (ah) Scatterplots show the top 20 pathways at 1.5, 3, 6, 9, 12, 24, 36, and 48 hpi, respectively. The x-axis denotes the pathway enrichment. The y-axis shows the names of the significantly enriched pathways. The P-values are indicated by variations from blue to red, with darker blue indicating more significant difference.

Additional file 11: Figure S7.

KEGG pathway analysis of the downregulated DEmRNAs in HFF cells at different time points after T. gondii infection. (ah) Scatterplots show the top 20 pathways at 1.5, 3, 6, 9, 12, 24, 36, and 48 hpi, respectively. The x-axis denotes the pathway enrichment. The y-axis shows the names of the significantly enriched pathways. The P-values are indicated by variations from blue to red, with darker blue indicating more significant difference.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Wang, SS., Zhou, CX., Elsheikha, H.M. et al. Temporal transcriptomic changes in long non-coding RNAs and messenger RNAs involved in the host immune and metabolic response during Toxoplasma gondii lytic cycle. Parasites Vectors 15, 22 (2022). https://doi.org/10.1186/s13071-021-05140-3

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s13071-021-05140-3

Keywords

  • Toxoplasma gondii
  • RNA-seq
  • lncRNA
  • mRNA
  • Immune response