Transcriptomic Establishment of Pig Macrophage Polarization Signatures

Macrophages are the foremost controllers of innate and acquired immunity, playing important roles in tissue homeostasis, vasculogenesis, and congenital metabolism. In vitro macrophages are crucial models for understanding the regulatory mechanism of immune responses and the diagnosis or treatment of a variety of diseases. Pigs are the most important agricultural animals and valuable animal models for preclinical studies, but there is no unified method for porcine macrophage isolation and differentiation at present; no systematic study has compared porcine macrophages obtained by different methods. In the current study, we obtained two M1 macrophages (M1_IFNγ + LPS, and M1_GM-CSF) and two M2 macrophages (M2_IL4 + IL10, and M2_M-CSF), and compared the transcriptomic profiles between and within macrophage phenotypes. We observed the transcriptional differences either between or within phenotypes. Porcine M1 and M2 macrophages have consistent gene signatures with human and mouse macrophage phenotypes, respectively. Moreover, we performed GSEA analysis to attribute the prognostic value of our macrophage signatures in discriminating various pathogen infections. Our study provided a framework to guide the interrogation of macrophage phenotypes in the context of health and disease. The approach described here could be used to propose new biomarkers for diagnosis in diverse clinical settings including porcine reproductive and respiratory syndrome virus (PRRSV), African swine fever virus (ASFV), Toxoplasma gondii (T. gondii), porcine circovirus type 2 (PCV2), Haemophilus parasuis serovar 4 (HPS4), Mycoplasma hyopneumoniae (Mhp), Streptococcus suis serotype 2 (SS2), and LPS from Salmonella enterica serotype minnesota Re 595.


RNA-Seq
After the induction, macrophage cells were collected and used for total RNA extraction using Trizol (Invitrogen, San Francisco, CA, USA, 15596026), and the RNA-seq libraries were constructed using the NEBNext ® UltraTM RNA Library Prep Kit (NEB, Ipswich, MA, USA, 7530). The paired-end RNA-seq sequencing libraries were further sequenced by the Illumina Novaseq6000 platform (PE150), yielding 151Gb raw data and an average of 804 million 150-bp paired-end raw reads (Novogene Co. Ltd., Tianjing, China).
Sequenced reads were aligned to the pig reference genome (Sus_scrofa.Sscrofa11.1.104.gtf and Sus_scrofa.Sscrofa11.1.cdna.all.fa.gz). Gene expression was quantified by using Kallis-tov0.44.0 [30] and obtaining read counts for each transcript. We standardize read counts to TPM (per million transcript readings) and differential gene expression analysis was performed using the Edge R package. KEGG and GO annotation were performed using the online tool Metascape and hub genes were identified by using Metascape (https: //metascape.org/gp/index.html, accessed on 4 November 2022), KEGG pathway analysis, and GO biological process analysis was performed for different genes. In order to compare the protein functions of macrophages between different methods, the protein-protein interaction network (PPI) was constructed by String (https://string-db.org, accessed on 17 November 2022) and the key Hub genes were identified.
Raw and processed RNA-seq data were deposited in the NCBI GSE202115.
Since the GSEA method was originally developed for analyzing microarray data, we normalized the raw count for standard GSEA by TPM and transformed TPM into a GCT format gene expression data set. Genes were ranked based on the correlation between their expression and class distinction, by evaluating if an a priori-defined set of genes (M1_IFNγ + LPS, M1_GM-CSF, M2_IL4 + IL10, and M2_M-CSF) were randomly distributed or were primarily associated with a tested class.

GSE34544
The piglets were infected with SS2 and HPS4 28 days later. Transcriptome analysis of infected and uninfected PAM.

Comparing M1 with M2 to Generate Porcine Macrophage Gene Signatures
To fully characterize the specificity of two macrophage phenotypes polarized by different methods, we applied RNA-seq and compared the transcriptomic commonalities and differences across phenotypes and methods within phenotype. The principal components analysis (PCA) plot indicated that the macrophages clearly separated between phenotypes and within phenotype ( Figure 1A). DEGs were evaluated by GO and KEGG after being transformed into gene IDs. The BP (biological process) components of the GO annotations of DEGs were used to examine the functional enrichment of DEGs. To ascertain the connection between DEGs and signaling pathways, KEGG analysis was performed. The DEGs of M1 (M1 IFN + LPS, and M1 GM-CSF) compared to M2 (M2 IL4 + IL10, and M2 M-CSF) were found to be enriched in biological processes related to immune response, such as Cell activation, Regulation of cell activation, Inflammatory response, Innate immune response, and pathways such as Hematopoietic cell lineage, Cytokine-cytokine receptor interaction (Figures 1B and S2 and Additional file S1). Moreover, we observed that the upregulated DEGs mainly enriched into phagosome pathways, pathways in cancer, and disease-related pathways, while the down-regulated genes enriched into tissue development pathways. To verify the reliability of polarized macrophages, 49 classical macrophage marker genes were retrieved from the previous literature, containing 29 and 20 macrophage markers for M1 and M2, respectively (Supplementary Table S2). Consistently, all 29 classical M1 and 20 M2 markers are highly expressed in M1 and M2, respectively ( Figure 1C, Additional file S1). Using a p value < 0.01 and absolute log2FC value > 1 as criteria, 730 differential expressed genes (DEGs) were identified. Comparing M1 to M2 620 and 110 genes were significantly upregulated and downregulated, respectively (Additional file S1). To examine more closely the correlations between the genes of the core response to M1 macrophages (M1_IFNγ + LPS, and M1_GM-CSF) and M2 (M2_IL4 + IL10, and M2_M-CSF) macrophages, we ran a protein-protein interaction analysis using String to identify the hub genes. Ranked by connectivity and betweenness centrality, MPO, S100A12, CTSG, CCR2, CAMP, S100A9, KIT, CXCR4, and CYBB were identified as the major hub genes ( Figure 1D, and Additional file S1). 1). Using a p value < 0.01 and absolute log2FC value > 1 as criteria, 730 differential expressed genes (DEGs) were identified. Comparing M1 to M2 620 and 110 genes were significantly upregulated and downregulated, respectively (Additional file 1). To examine more closely the correlations between the genes of the core response to M1 macrophages (M1_IFNγ + LPS, and M1_GM-CSF) and M2 (M2_IL4 + IL10, and M2_M-CSF) macrophages, we ran a protein-protein interaction analysis using String to identify the hub genes. Ranked by connectivity and betweenness centrality, MPO, S100A12, CTSG, CCR2, CAMP, S100A9, KIT, CXCR4, and CYBB were identified as the major hub genes ( Figure 1D, and Additional file 1).

Revealing Transcriptomic Differences of Macrophages within Phenotypes
We also observed that the specific macrophage phenotypes were separately clustered by induction methods, M1_IFNγ + LPS versus M1_GM-CSF and M2_IL4 + IL10 versus M2_M-CSF ( Figure 2). Comparing M1_IFNγ + LPS to M1_GM-CSF, 207 and 391 genes were significantly upregulated and downregulated, respectively (Additional file S2). DEGs' annotation indicates that they are enriched in biological processes related to immune response, such as Regulation of cell activation, Inflammatory response, Response to the bacterium, and pathways such as Staphylococcus aureus infection, Cytokine-cytokine receptor interaction, Pathways in cancer ( Figures 2C and S2, and Additional file S2). To examine more closely the correlations between the genes of the core response to M1_IFNγ + LPS and M1_GM-CSF, we incorporated a network-based protein-protein interaction analysis approach by String. Ranked by connectivity and betweenness centrality, CD68, MRC1, CD28, HLA-DRA, PPARG, A2M, S100A9, SELP, RETN, and CAMP were identified as the primary hub and bottleneck genes of macrophages ( Figure 2E, and Additional file S2).
Comparing M2_IL4 + IL10 to M2_M-CSF, 391 and 211 genes were significantly upregulated and downregulated, respectively (Additional file S3). Running GO and KEGG annotations, the DEGs enriched in tissue development-related biological processes such as Tube morphogenesis, Reproductive structure development, Regulation of vasculature development, and pathways such as Osteoclast differentiation, Hematopoietic cell lineage, Cell adhesion molecules ( Figures 2D and S3 and Additional file S3). IL6, CD4, IGF1, IL10, PTGS2, FOS, PPARG, NTRK1, ADIPOQ, FLT1, and HMOX1 were identified as the primary hub and bottleneck genes of macrophages induced by colony factors and cytokines' exposure ( Figure 2F and Additional file S3). The interaction between each protein pair is represented by lines, and the size of the circle is directly proportional to the degree of interaction. We adopt a node number greater than 1; (F) The DEGs are studied by using the String online tool, M2_IL4 + IL10, and M2_M-CSF. The interaction between each protein pair is represented by lines, and the size of the circle is directly proportional to the degree of interaction. We adopt a node number greater than 2.

Application of M1_IFNγ + LPS, M1_GM-CSF, M2_IL4 + IL10, and M2_M-CSF Signatures to the Identification of Swine Disease and Its Molecular Mechanism
As macrophages play a key role in determining the activation or resolution of immune responses and can determine the fate of pathogen infection, we evaluated the robustness of our macrophage signatures (M1_IFNγ + LPS, M1_GM-CSF, M2_IL4 + IL10, M2_M-CSF) in discriminating macrophages infected with specific pathogens based on the enrichment analysis of selected genes. The association between the chosen clinical trials where infected pig macrophages were examined in 308 transcriptomes' data obtained from the GEO ( Table 2 and Additional file 4). The DEGs were studied by using the String online tool, M1_IFNγ + LPS, and M1_GM-CSF. The interaction between each protein pair is represented by lines, and the size of the circle is directly proportional to the degree of interaction. We adopt a node number greater than 1; (F) The DEGs are studied by using the String online tool, M2_IL4 + IL10, and M2_M-CSF. The interaction between each protein pair is represented by lines, and the size of the circle is directly proportional to the degree of interaction. We adopt a node number greater than 2.

Application of M1_IFNγ + LPS, M1_GM-CSF, M2_IL4 + IL10, and M2_M-CSF Signatures to the Identification of Swine Disease and Its Molecular Mechanism
As macrophages play a key role in determining the activation or resolution of immune responses and can determine the fate of pathogen infection, we evaluated the robustness of our macrophage signatures (M1_IFNγ + LPS, M1_GM-CSF, M2_IL4 + IL10, M2_M-CSF) in discriminating macrophages infected with specific pathogens based on the enrichment analysis of selected genes. The association between the chosen clinical trials where infected pig macrophages were examined in 308 transcriptomes' data obtained from the GEO ( Table 2 and Additional file S4).
The results indicated that genes from the SS2 infected macrophages were most significantly enriched in the M1_IFNγ + LPS set, genes from the PPRSV, ASFV, Mhp, and LPS infected were most significantly enriched in the M1_GM-CSF, genes from the PCV2 and T. gondii Me49 infections were most significantly enriched in the M2_M-CSF set, and genes from the different species of T. gondii infections were most significantly enriched in the M2_IL4 + IL10 set ( Figure 3).
Together, our results indicate that our macrophage signatures could characterize microarrays and RNA-seq from specific pathological scenarios.
The results indicated that genes from the SS2 infected macrophages were most significantly enriched in the M1_IFNγ + LPS set, genes from the PPRSV, ASFV, Mhp, and LPS infected were most significantly enriched in the M1_GM-CSF, genes from the PCV2 and T. gondii Me49 infections were most significantly enriched in the M2_M-CSF set, and genes from the different species of T. gondii infections were most significantly enriched in the M2_IL4 + IL10 set ( Figure 3).
Together, our results indicate that our macrophage signatures could characterize microarrays and RNA-seq from specific pathological scenarios.

Discussion
Growing evidence shows macrophages have high plasticity and extensive polarization, which hinders the definition of macrophages obtained by different methods. The establishment of the porcine macrophage model is important for pig health and disease research. Transcriptomics is one of the primary tools in this investigation, however, we are aware of its shortcomings. Since the transcriptome is type-specific and changes over time and in response to stimulation, and some tissue-specific disorders

Discussion
Growing evidence shows macrophages have high plasticity and extensive polarization, which hinders the definition of macrophages obtained by different methods. The establishment of the porcine macrophage model is important for pig health and disease research. Transcriptomics is one of the primary tools in this investigation, however, we are aware of its shortcomings. Since the transcriptome is type-specific and changes over time and in response to stimulation, and some tissue-specific disorders cannot be diagnosed by RNA-seq, several studies have demonstrated that not all genes are expressed solely in particular cells [32]. Of course, considering the operability of the laboratory and the better differentiation effect of young piglets, we only used the 7-day-old DLY, which has certain limitations. In the present study, we generated four porcine macrophage phenotypes, namely M1_IFNγ + LPS, M1_GM-CSF, M2_IL4 + IL10, and M2_M-CSF and created four macrophage molecular signatures by transcriptomic analysis. First, we identified 730 DEGs by comparing M1 (merging M1_IFNγ + LPS and M1_GM-CSF) with M2 (mering M2_IL4 + IL10 and M2_M-CSF). In line with previous reports, the porcine macrophages had similar gene profiles to human and mouse classical macrophage phenotypes and alternative phenotypes, respectively. For example, porcine M1 highly expressed CXCR4, S100A9, FCRL3, JAK3, and SLAMF1, while M2 highly expressed CCL11, POSTN, PTX3 MFAP4, PTH1R, and AGTR1. Apart from the well-known macrophage markers, we noticed that porcine M1 highly expressed MPO, S100A12, CTSG, CCR2, CAMP, KIT, and CYBB, and M2 highly expressed ANK2, ALDH1A1, PTHLH, and CACNA1H. MPO is a marker of macrophage aggregation in inflammatory species and promotes the recognition of macrophage scavenger receptors [33]. S100A12, CSTG, CAMP, and KIT are antimicrobial response genes. CSTG, CAMP, and CYBB are involved in macrophage phagocytosis [34][35][36]. Both CCR2 and CXCR4 are pro-inflammatory genes involved in macrophage migration [37]. ANK2 promotes the growth and invasion of pancreatic carcinoma, and the peptide derived from ANK2 is an effective and specific autophagy inhibitor binding to ATG8 [38,39]. Studies have shown that the expression level of ALDH1A1 is reduced in the inflammatory state, which is part of early inflammation [40]. PTHLH is involved in cell and organ growth, development, migration, and survival, and can be used as an independent marker of prognosis [41]. These genes may serve as specific porcine macrophage markers and potential therapeutic targets. It is worth mentioning that four MHC haplotype genes were detected expressing in porcine bone marrow-derived macrophages, such as CIITA, HLA-DOB, HLA-DRA, SLA-DMB, and two of them, HLA-DRA and SLA-DMB, are differentially expressed, comparing M1 versus M2 [42].
After revealing the transcriptomic difference between the two main macrophage phenotypes, we conducted a transcriptomic comparison of M1_IFNγ + LPS with M1_GM-CSF, and M2_IL4 + IL10 with M2_M-CSF. We found that porcine M1_IFNγ + LPS highly expressed S100A9, SELP, RETN, CAMP, and M1_GM-CSF highly expressed CD68, MRC1, CD28, HLA-DRA, PPARG, and A2M. S100A9 is related to the CD68 regulating macrophage function pathway and promoting macrophage migration, which can induce neutrophil chemotaxis and promote macrophage migration and adhesion under inflammatory conditions [43]. RETN has been reported to induce the production of pro-inflammatory cytokines and chemokines in PBMC [44]. CAMP can inhibit the phagocytosis of macrophages through the CAMP-1-activated CAMP-effect-exchange protein [45]. CD68 is highly expressed in macrophages and belongs to the scavenger family. It has the functions of clearing cell debris, promoting phagocytosis, and mediating the recruitment and activation of macrophages. MRC1 and A2M mediate the phagocytosis of macrophages [43]. CD28 can enhance the expression of RANTES and MIP-1α in T cells and MIP-1β, and increase the number and differentiation of macrophages in the wound healing stage [9,46]. HLA-DRA has been proven to inhibit the phagocytosis of macrophages in order to protect the intracellular niche from phagocytosis and killing of host macrophages, which is positively related to the regulation of GM-CSF [47][48][49].
We have investigated the variations among various phenotypes and associated them with disease outcomes because we are interested in investigating the potential application of macrophage phenotypes. We proposed a consensus collection of markers describing major macrophages' activation phenotypes, namely M1_IFNγ + LPS, M1_GM-CSF, M2_IL4 + IL10, and M2_M-CSF, able to characterize robustly controlled in vitro and in vivo scenarios for macrophage induction. Our study confirmed that macrophages cultured in vitro have different enrichment pathways from macrophages infected with different pathogens. At the same time, it is confirmed that macrophages infected with the same pathogen in vivo or in vitro also have different enrichment pathways. Since the description of the polarization state of macrophages has not been unified, our study provides a framework to guide the definition of the phenotype of porcine macrophages in the disease state. Future research into the macrophage model in a disease setting will help us create medications and diagnostic tools for particular disorders. Nevertheless, we can never ignore the bias using public downloaded transcriptomic data.
We conclude here that there are transcriptomic differences between and within two macrophage phenotypes in our porcine model. In general, we found the four types of macrophages (M1_IFNγ + LPS, M1_GM-CSF, M2_IL4 + IL10, and M2_M-CSF) have different functions. Compared with M1_GM-CSF, M1_IFNγ + LPS has a weaker phagocytic capacity, but stronger antibacterial and migratory capacity; M2_IL4 + IL10 has a stronger tissue repair function, while M2_M-CSF has a stronger wound healing ability [66]. At the same time, we used the four established gene characteristics to identify various pig infectious diseases with prognosis and predictive values. Pigs' immunological parameters overlap with either mice or humans' in particular ways more than mice and humans' do with each other [21,[23][24][25][26][27]. Indeed, it would be an interesting research direction to establish similar models of pig, mouse, and human macrophages for future studies.
Supplementary Materials: The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/cimb45030151/s1, Figure S1: The KEGG or GO Analysis of M1 and M2; Figure S2: The KEGG or GO Analysis of M1_IFNγ + LPS and M1_GM-CSF; Figure S3: The KEGG or GO Analysis of M2_IL4 + IL10 and M2_M-CSF; Table S1: Data output quality list. Table S2: The observed marker genes of M1 and M2 macrophages. References  are cited in the supplementary materials.

Institutional Review Board Statement:
All the experiments were conducted following these protocols and were approved by the Experimental Animal Ethics Committee of Sichuan Agricultural University under permit number of 20210167. All experimental steps were performed following relevant guidelines and regulatory requirements.