Expression profile of microRNAs in porcine alveolar macrophages after Toxoplasma gondii infection

Background Toxoplasma gondii is an apicomplexan protozoan parasite that can cause serious clinical illnesses in both humans and animals. microRNAs (miRNAs) are non-protein-coding RNAs that can regulate the expression of target genes. A previous study found that many miRNAs were differentially expressed after T. gondii infection and exert significant effects and revealed that both host survival and the virulence of different strains can be regulated by different miRNAs. Macrophages play an important role in T. gondii infection, but few studies have investigated the relationship between miRNAs and porcine alveolar macrophages infected with T. gondii. Methods Porcine alveolar macrophages (3D4-21) were infected with the RH (Type I) and Me49 (Type II) strains of T. gondii for 12 h and 24 h and then harvested. miRNA libraries were generated using the NEBNext® Multiplex Small RNA Library Prep Set for Illumina® (NEB, USA), and the miRNA expression levels were estimated based on transcripts per million reads (TPM). Results Our study generated six miRNA expression profiles from macrophages infected with RH and Me49 compared with the control groups. The comparison of the T. gondii-infected and uninfected samples identified 81 differentially expressed miRNAs, including 36 novel miRNAs and 45 mature miRNAs. The target genes of these differentially expressed miRNAs were predicted using miRanda software, and ssc-miR-127 and ssc-miR-143-3p were predicted to regulate nitric oxide synthase 1 (NOS1) and nitric oxide synthase 3 (NOS3), respectively, which play essential roles in synthesizing nitric oxide (NO) by oxidizing L-arginine. These genes were differentially expressed in both the RH- and Me49-infected groups. A KEGG enrichment analysis indicated that the predicted target genes were involved in multiple signaling pathways, including FcγR-mediated phagocytosis, the AMPK signaling pathway, the mTOR signaling pathway, and the FcγRI signaling pathway, all of which are indispensable for the normal functioning of porcine alveolar macrophages. Conclusions Our results provide data on the miRNA profile of porcine alveolar macrophages infected with T. gondii. To our knowledge, this study provides the first demonstration of the relationship between miRNA and macrophages of swine origin. Understanding the functions of these regulated miRNAs will aid the investigation of T. gondii infectious diseases, and the differentially expressed miRNAs might be candidate drug targets for T. gondii infection in pigs. Electronic supplementary material The online version of this article (10.1186/s13071-019-3297-y) contains supplementary material, which is available to authorized users.


Background
Toxoplasma gondii is an apicomplexan protozoan parasite that infects all nucleated cells and a widely distributed zoonotic parasite that can cause serious clinical illnesses in animals and humans. Although T. gondii has a heteroxenous life-cycle and a broad host range, its sexual reproduction only occurs in the small intestine of felids, the only definitive hosts of T. gondii. Toxoplasma gondii can cause severe disease in warm-blooded animals, including humans; immunocompetent individuals may experience fever, lymphadenopathy and muscle aches, and immunosuppressed patients may present with neurological symptoms, encephalitis and toxoplasma encephalitis. Toxoplasma gondii acquired during pregnancy might cause serious damage to newborn babies, including sight and hearing impairments and central nervous system diseases [1][2][3]. In many countries, T. gondii is one of the major causes of abortion and stillbirth in pig and sheep farms. Pork is the most widely consumed meat per capita in most countries, but pigs pose potential zoonotic risks because they might carry bacteria and other pathogens, including T. gondii. In fact, pigs are among the most important sources of material used in human xenoplastic transplantation, and the high rate of T. gondii infection in pigs is a huge challenge for future organ transplantations [4,5]. Therefore, the prevention and control of swine toxoplasmosis is of great significance to animal husbandry and for the prevention of human T. gondii infection.
Macrophages are myeloid immune cells that form part of the mononuclear phagocyte system. These cells, which are found in many tissues and organs, including the spleen, lungs and lymph nodes, consume and digest dead cells, debris, and parasites in the body [6,7]. During inflammation, macrophages are key players in phagocytosis, antigen presentation, and the secretion of various cytokines, chemokines, and growth factors. Toxoplasma gondii infection induces macrophages to secrete various cytokines, including the upregulated cytokines IFN-γ, TNF-α, IL-1, IL-2 and IL-12 and the downregulated cytokines IL-4, IL-6 and IL-10 [8]. IFN-γ plays a leading role in the host's fight against T. gondii, and IFN-γ induces an increase in the insecticidal ability of macrophages. During the processes of phagocytosis, macrophage and pathogen surface receptors assimilate into a membrane-bound digestive vacuole; therefore, macrophages are considered an important component of the innate immune system's first line of defense against pathogen infection [9]. However, many intracellular pathogens exploit this mechanism to ensure their internalization and survival in phagocytes, particularly macrophages [10]. This mechanism is one of the strategies used by T. gondii to evade host immunity.
miRNAs are a class of endogenous non-protein-coding RNAs that serve as negative regulators of the host response through the miRNA-induced silencing complex (miRISC) in apicomplexan parasites [11]. Since the discovery of miRNAs in 1993, these molecules have attracted wide attention due to their unique functions and provide a new dimension of molecular biology. Increasing evidence suggests that miRNAs play multiple roles in diverse biological processes, including organ development, cell proliferation and division, pathological processes [12][13][14] and in T. gondii infection. For example, following T. gondii infection, the miRNA-30c-1, miRNA-125b-2, miRNA-23b-27b-24-1 and miRNA-17~92 cluster genes bind to the STAT3 promoter in human macrophages, which leads to increased apoptosis of these cells [15]. Early during the T. gondii infection process in mice, three miRNAs (mmu-miR-712-3p, mmu-miR-511-5p and mmu-miR-217-5p) are prominently expressed by both the RH and ME49 strains, and the increase in the expression of these miRNAs is specific to Toxoplasma species [16].
In our study, porcine alveolar macrophages were infected with the RH and Me49 strains of T. gondii. The goal of our study was to explore the changes in miRNA expression in porcine alveolar macrophages after infection with T. gondii. A previous study revealed that miR-NAs have significant effects at the posttranscriptional level on physiology, pathogenesis, and immunology. Our results provide full miRNA profiles of porcine alveolar macrophages infected with T. gondii, and differential miRNA expression was further analyzed by target prediction and pathway analysis. Unveiling the roles of differentially expressed miRNAs in macrophages will help researchers obtain a better understanding of the interactions between T. gondii and macrophages. In addition, an understanding of the functions of these regulated miRNAs will aid the investigation of toxoplasmosis, and the differentially expressed miRNAs might be candidate drug targets for T. gondii infection in pigs.

Parasite infection and sample collection
The RH (type I) and Me49 (type II) strains of T. gondii were cryopreserved in liquid nitrogen in our laboratory. Porcine alveolar macrophages and human foreskin fibroblasts (HFF) were cultured in DMEM with 10% FBS. Porcine alveolar macrophages were infected with the RH and Me49 strains at a macrophage-to-tachyzoite ratio of 1:5 (MOI = 5) to ensure an appropriate infective dose. RH-and Me49-infected macrophages were harvested at 12 and 24 h after infection with T. gondii (herein referred to as RH12 and RH24, and M12 and M24), and samples from the control group were collected at the same time points (referred to as D12 and D24). Thus, six samples were generated during the study: RH12, RH24, M12, M24, D12 and D24. RNA was extracted from all samples using TRIzol, and small RNA sequencing was then performed to characterize the global miRNA transcriptional response to T. gondii infection.

RNA quantification and qualification
The degradation and contamination of RNA were monitored using 1% agarose gels. The RNA purity was checked using a NanoPhotometer® spectrophotometer (Implen, GER), and the RNA concentration was measured using a Qubit® RNA Assay Kit with a Qubit® 2.0 Fluorometer (Life Technologies, USA). The RNA integrity was assessed using an RNA Nano 6000 Assay Kit with an Agilent Bioanalyzer 2100 system (Agilent Technologies, USA).

Preparation of libraries for small RNA sequencing and data analysis
The samples were sent to the Beijing Novogene Bioinformatics Institute for Illumina sequencing. Three micrograms of total RNA was used as the input material to generate a small RNA library used for the analysis of miRNAs. Sequencing libraries were generated using the NEBNext® Multiplex Small RNA Library Prep Set for Illumina® (NEB), and first-strand cDNA was then synthesized using Moloney Murine Leukemia Virus Reverse Transcriptase (M-MuLV-RT). PCR amplification was performed using LongAmp Taq 2× Master Mix, and the PCR products were purified on an 8% polyacrylamide gel. The libraries obtained from the different samples were sequenced on an Illumina HiSeq 2500/2000 platform, and 50-bp single-end reads were generated.
The raw data (raw reads) in FASTQ format were first processed using Python and custom Perl scripts. Clean data (clean reads) were obtained after removing various reads, including poly-N reads, reads without a 3' adapter or an insert tag, reads with 5' adapter contamination, reads with poly A, T, G or C and low-quality reads from the raw data. Moreover, the Q20, Q30 and GC contents of the raw data were calculated. We then selected a specific length for the clean reads to perform all subsequent analyses.
To analyze the expression and distribution of the small RNA tags on the reference sequence, we mapped the small RNA tags to the reference sequence using Bowtie [17] without allowing any mismatch; this process was also used to search for known miRNAs. miRBase20.0 was used as a reference, and potential miRNA and secondary structures were obtained using modified mir-deep2 [18] and srna-tools-cli software. Custom scripts were used to obtain the miRNA counts as well as the base bias on the first positions of the identified miRNAs with a specific length and on each position of all identified miRNAs. All small RNA tags were mapped to the Rfam database, Repeat Masker or the various types of data from the specified species to remove tags originating from repeat sequences, protein-coding genes, rRNAs, snRNAs, tRNAs and snoRNAs.
Novel miRNAs can be predicted based on the hairpin structure of the miRNA precursor. To predict novel miRNAs, we explored the Dicer cleavage site, the secondary structure and the minimum free energy of the small RNA tags that were not annotated in the former steps using the available software mirdeep2 [18] and miREvo [19].
The target genes of the miRNAs were predicted using miRanda software. The measure transcripts per million (TPM) was used to estimate miRNA expression levels based on the following criteria: normalized expression = mapped read count/total reads*1,000,000 [20]. Differential expression analysis between two groups was performed using the DESeq R package (v.1.8.3). P-values were adjusted using the Benjamini & Hochberg method, and by default, a corrected P-value of 0.05 was set as the threshold for significant differential expression. The target gene candidates of differentially expressed miRNAs were subjected to gene ontology (GO) enrichment analysis. The statistical enrichment of the target gene candidates in KEGG pathways was tested using KEGG Orthology Based Annotation System (KOBAS) software.

Validation of miRNA expression
The results of the miRNA expression analysis were validated by tailing-reaction quantitative RT-PCR using an ABI StepOne Sequence Detection System with SYBR Green qPCR Super Mix according to the manufacturer's recommended protocol (Invitrogen, USA). All qRT-PCR reactions were performed in triplicate, and snRNA U6 was used as an internal control for the normalization and quantification of miRNA expression.

RNA quality testing
All RNAs were subjected to quality testing to ensure their availability. The detected contents included OD260/280, OD260/230, 28S/18S and RNA integrity (RIN), and all test results were qualified.

miRNA expression profiles of macrophages after T. gondii infection
Our study generated six miRNA expression profiles from macrophages infected with either RH or Me49 or non-infected macrophages (control) at two different time points, and more than ten million raw sequence reads from each group were mapped to the sRNA library. In addition, more than ten million high-quality clean reads were obtained after low-quality reads, junk sequences and adaptor sequences were removed from the raw sequence reads. The sequencing data are summarized in Tables 1 and 2.
Most (> 70%) of the total clean reads in the control group had a length of 22-24 nucleotides, and the length of the clean reads was significantly changed after T. gondii infection. Compared with the control group, the number of sRNAs between 22-24 nucleotides in length decreasing by 50% in experimental groups, and the length of the majority of the sRNAs in the experimental groups was 19-29 nucleotides (Fig. 1a-d). We removed sRNAs of 24-29 nucleotides and 32-35 nucleotides, and the remaining clean reads were used for the identification of mature miRNAs and for the prediction of novel miRNAs. The comparison of the T. gondii-infected and uninfected group revealed 81 differentially expressed miRNAs, including 36 novel miRNAs and 45 mature miRNAs. In the group infected with the RH strain, 19 mature miRNAs were upregulated and nine mature miRNAs were downregulated. In addition, 33 novel miRNAs were identified in this group, all of which were upregulated (Tables 3, 4). Moreover, an analysis of the Me49-infected group revealed that 13 and six mature miRNAs were upregulated and downregulated, respectively, and that all 23 novel miRNAs detected were upregulated (Tables 5, 6). Detailed miRNA data for the different groups (M12, M24, R12 and R24) are summarized in Additional file 1: Table S1.
An analysis of the miRNAs in the different groups is shown in Fig. 2. The comparison of the porcine alveolar macrophages infected with RH and Me49 for 12 h with the corresponding control group identified seven miRNAs (ssc-miR-486, ssc-miR-451, ssc-miR-4332, ssc-miR-1285, ssc-novel-90, ssc-novel-72 and ssc-novel-417) that were differentially expressed in both infected groups. Moreover, an analysis of groups obtained 24 h after infection revealed 26 miRNAs that were differentially expressed in the RH-and Me49-infected groups compared with the control. In total, five miRNAs (ssc-miR-4332, ssc-miR-1285, ssc-novel-90, ssc-novel-72 and ssc-novel-417) were differentially expressed in the M12, M24, R12 and R24 groups compared with the control group, and these detailed miRNA data are summarized in Table 7.

GO and KEGG analyses
In total, 7462 predicted targets were identified (Additional file 2: Table S2). A gene ontology (GO) enrichment analysis was the performed based on the predicted target genes of the differentially expressed miRNAs to identify the enriched biological processes, molecular functions and cellular components. The enriched GO terms of the cellular component, molecular function and biological process categories are shown in Fig. 3. A KEGG enrichment analysis demonstrated that the target genes were involved in multiple signaling pathways, including FcγR-mediated phagocytosis, AMPK signaling, mTOR signaling, phosphatidylinositol signaling, Fcγ epsilon RI signaling, B cell receptor signaling, T cell receptor signaling, pathways in cancer and VEGF signaling. The top 20 KEGG enrichment pathways are shown in Tables 8 and 9.

Verification of differentially expressed miRNAs by realtime fluorescent quantitative PCR
Five mature miRNAs (ssc-miR-335, ssc-miR-199-5p, ssc-miR-125b, ssc-miR-451 and ssc-miR-486) and one novel miRNA (ssc-novel-90) were randomly selected from the set of identified differentially expressed miRNAs, and their expression levels were verified by qRT-PCR using miRNA-specific primers. The expression trends obtained for all the miRNAs by qRT-PCR were similar to those found with the original small RNA-sequencing data (Fig. 4). The results indicated a low false-discovery rate with respect to the RNA-sequencing data and supported the viability of the miRNA expression profiling study. The RNA-sequencing data obtained in our study have been uploaded to the GEO database with the accession number GSE119643.

Discussion
Previous studies have shown that miRNAs play a critical role in host cells, and many miRNAs are up-or downregulated after infection with T. gondii [15,16]. We used small RNA sequencing and qRT-PCR approaches to analyze the miRNA expression profiles of porcine alveolar macrophages after T. gondii infection. Our results showed that miRNAs are differentially expressed between T. gondii-infected and T. gondii-uninfected porcine alveolar macrophages, regardless of whether the RH or Me49 strain is used for the infection. We compared the miRNA profiles of T. gondii-infected and T. gondii-uninfected mice to confirm that the differentially expressed miRNAs might be important in macrophage resistance to T. gondii infection. Compared with the control group, 19 and nine mature miRNAs were upregulated and downregulated, respectively, in the RH-infected group, whereas in the Me49-infected group, 13 and six mature miRNAs were upregulated and downregulated, respectively. The fact that all the differentially expressed novel miRNAs were upregulated and most novel miRNAs were only predicted in the experimental group indicate that the miRNAs expression profile was changed significantly after T. gondii infection. The potential target genes of the identified miRNAs were predicted using miRanda software, and this analysis identified 7462 predicted target genes. Previous research has shown that macrophages play a crucial role in the course of T. gondii infection because nitric oxide (NO) is produced to control the multiplication of T. gondii in activated macrophages both in vitro and in vivo [21,22]. Among the differentially expressed miRNAs, ssc-miR-127 and ssc-miR-143-3p were predicted to regulate nitric oxide synthase 1 (NOS1) and nitric oxide synthase 3 (NOS3), respectively, and these predicted targets play essential roles in NO synthesis by oxidizing L-arginine [23] and were differentially expressed in both the RHand Me49-infected groups. Based on these findings, we can reasonably conclude that ssc-miR-127 and ssc-miR-143-3p might be involved in macrophage resistance to T. gondii. The potential functions and relationships of the identified target genes were predicted through a KEGG pathway enrichment analysis. An analysis of the top 30 KEGG pathways that were found to be enriched in the infected groups revealed many commonalities: the target genes of the differentially expressed miRNAs were mostly enriched in FcγR-mediated phagocytosis, the mTOR signaling pathway and the TNF signaling pathway.
FcγR-mediated phagocytosis plays an important role in defense against parasites through the processes of antigen recognition and phagocytosis in macrophages. SYK kinases are essential in the process of FcγR-mediated phagocytosis and are considered crucial for downstream biochemical processes, including PIP2 breakdown and PLCγ and PI3-kinase activation [24,25]. The activation of the FcγR signaling pathway requires the PI3-kinase isoform PI3Kβ, members of the Rac small GTPase families and Vav guanine nucleotide exchange factors, and Rac2 and Vav3 in particular play notable roles [26][27][28]. Five aberrantly expressed miRNAs were found to be related to these critical genes: ssc-miR-127 was predicted to regulate SYK; ssc-miR-421-3p was predicted to regulate Rac2; and ssc-miR-143-3p, ssc-miR-199a-5p and ssc-miR-1285 were predicted to regulate Vav3. Therefore, we can reasonably conclude that by binding to their target genes, the abnormally expressed miRNAs participate in the elimination of T. gondii by macrophages.
The TNF signaling pathway has a wide range of functions in innate and adaptive immunity. The signal transduction network of TNF comprises many crosslinked pathways, including the MAPK signaling pathway, the NF-κB signaling pathway and the PI3K signaling pathway. TNF, TNF receptor and downstream members form an immense signaling pathway network. Through binding to TNFR1, TNF receptor activation factor-2 (TRAF2) activates the classical NF-κB pathway to produce cytokines such as IFN-γ, TNF-α, IL-1, IL-2 and IL-12, and these cytokines play a vital role in macrophage resistance to T. gondii infection [29][30][31]. Moreover, the activation of TNFR2 by TNF leads to the recruitment of TRAF2 and the downstream activation of the JNK or p38/ MAPK pathways, and the activation of these signal transduction pathways enhances the production of intracellular cytokines. Four aberrantly expressed miRNAs were predicted to be related to the members of these signal transduction pathways: ssc-novel_72 was predicted to regulate TNFR1, and ssc-miR-125b and ssc-miR-143-3p were predicted to regulate TRAF2/5. Thus, these miRNAs might play a crucial role in TNF signaling.