Transcriptomic Analysis of Chicken Lungs Infected With Avian and Bovine Pasteurella multocida Serotype A

Pasteurella multocida (P. multocida) is a common animal pathogen responsible for many animal diseases. Strains from different hosts exhibit disparate degrees of effect in other species. Here, we characterize an avian P. multocida serogroup A strain (PmQ) showing high lethality to chickens and a bovine P. multocida serogroup A strain (PmCQ2) with no lethality to chickens. We used RNA-seq to profile the transcriptomes of chicken lungs infected with PmQ and PmCQ2. A total of 1,649 differentially expressed genes (DEGs) due to PmQ infection (831 upregulated genes and 818 downregulated genes) and 1427 DEGs (633 upregulated genes and 794 downregulated genes) due to PmCQ2 infection were identified. Functional analysis of these DEGs demonstrated that the TNF signaling pathway, the toll-like receptor signaling pathway, complement and coagulation cascades, and cytokine–cytokine receptor interaction were both enriched in PmQ and PmCQ2 infection. STAT and apoptosis signaling pathways were uniquely enriched by PmQ infection, and the NOD-like receptor signaling pathway was enriched only by PmCQ2 infection. Cell-type enrichment analysis of the transcriptomes showed that immune cells, including macrophages and granulocytes, were enriched in both infection groups. Collectively, our study profiled the transcriptomic response of chicken lungs infected with P. multocida and provided valuable information to understand the chicken responses to P. multocida infection.


INTRODUCTION
Pasteurella multocida (P. multocida) is a Gram-negative bacteria, first characterized by Louis Pasteur as the causative agent of fowl cholera in 1881 (1). P. multocida causes many animal diseases, such as fowl cholera, swine atrophic rhinitis, rabbit septicemia, and bovine pneumonia (2)(3)(4). Fowl cholera, caused by avian P. multocida, is a highly contagious disease of various domestic and wild bird species, resulting in great economic losses worldwide (5). In the case of peracute or acute infection, the lung lesions dominated by hemorrhages can be characteristic of fowl cholera (6). Antibiotics are still the main treatment for pasteurellosis due to the lack of effective multi-serotype vaccines (7). In some cases, P. multocida can cause human infection via animal bites and scratches (8,9). Thus, further study on the mechanism of its pathogenesis is still needed.
Pasteurella multocida can be classified into five serotypes, A, B, D, E, and F, according to the specificity of capsular antigens (10). Different types of P. multocida strains often associated with different diseases depend on the infected host; for example, avian cholera is usually caused by P. multocida of serotype A and F, and swine atrophic rhinitis usually results from serotype D (2). Moreover, P. multocida isolated from certain hosts shows differing pathogenicity to other animals (11,12). A P. multocida B:2 strain can kill cattle and buffaloes at a low dose but has no affect on chickens even at very high doses (13).
In a previous study, we isolated two P. multocida strains of serotype A, PmQ from tissues of a dead duck and PmCQ2 from pneumonic lungs of a beef cow (14,15). Both strains showed high virulence to mice, but only PmQ was virulent in chickens and ducks. Here, we have used RNA-seq to investigate the chicken lung transcriptional response to these two P. multocida strains. This study provides comprehensive insights into the chicken immune response to P. multocida and demonstrates that different host-isolated strains show different virulence and may induce diverse immune responses.

Bacteria Strains and Culture Conditions
PmQ (GenBank accession No. CP033597) used in this study was isolated from tissues of a dead duck, and PmCQ2 (GenBank accession No. CP033599) was isolated from pneumonic lungs of a beef cow (16). Colonies were picked from Martin's broth agar plate, cultured overnight in Martin's broth aerobically at 37 • C with shaking at 220 rpm, and counted by plating on Martin agar plates.

Experimental Animals and Ethics Statement
A total of 150 commercial chickens (Ross 308, 7 days old), obtained from a commercial supplier (JinGang Village Farm, Chongqing), were used in our experiments. All animal experiments were carried out with the approval of the Animal Ethics and Research Committee of Southwest University (permit No. 11-1025), Chongqing, China.

Chicken Challenge
Chickens were kept in wire cages separately in groups (n = 10 per group) with free access to food and water. Chickens of experimental groups were intraperitoneally infected with PmQ (PmQ-P) at a dosage of 1 × 10 3 colony-forming units (CFU) in 200 µL saline and PmCQ2 (PmCQ2-P) at a dosage of 2 × 10 9 CFU in 200 µL saline and then monitored for 7 days. The control group (P-C) was treated in parallel with sterile saline of equal volume.
Chickens infected with PmQ were sacrificed at 16 h postinfection (hpi) and lungs were collected. Chickens infected with PmCQ2 were sacrificed at 4, 8, 16, 24, and 48 hpi. Bacterial burdens were determined by homogenizing tissues and diluting ten-fold serially in saline. The different dilutions were plated

Pathological Examinations
For histopathological examination, the lungs of chickens were collected and immediately fixed in 4% paraformaldehyde for 24 h, dehydrated in graded ethanol, and then embedded in paraffin wax. The tissues were sliced at 3 µm thickness and then stained with hematoxylin and eosin (H&E).

RNA Isolation, Library Construction, and Illumina Deep Sequencing
About 50-100 mg of each lung were homogenized in liquid nitrogen. Total RNA was extracted using Trizol reagent (Invitrogen Life Technologies, USA) following the manufacturer's protocol. RNA integrity was confirmed by agarose gel electrophoresis and quantified by NanoDrop (NanoDrop 2000, Thermo Scientific). Following that, the Ribo-Zero rRNA Removal Kit (Illumina, San Diego, CA, USA) was used to remove

Qualification and Identification of Differently Expression Genes (DEGs)
HTSeq 0.6.1p2 (http://www-huber.embl.de/users/anders/ HTSeq) were used to obtain the read count of genes. To eliminate the effects of gene length and sample differences, reads per kilo bases per million (RPKM) was used to qualify the gene expression. Differential expression analysis was performed using estimateSizeFactors in DESeq (version 1.18.0) R package, and p-value and fold-change value were identified using the nbinomTest in DESeq with the R package. P < 0.05 and |log 2 foldchange| >1 were set as the threshold for DEG identification. Hierarchical cluster analysis was used to identify DEGs with certain patterns of expression under three different challenge groups using Heatmap based on R.

Functional Annotation of DEGs
All DEGs were mapped to the gene ontology (GO) and KEGG ontology (KO) databases to obtain respective GO terms and KO classification and then enrichment analysis of the DEGs was performed using R based on hypergeometric distribution in the context of the entire reference genome. Gene clusters were subject to network analysis by using STRING (version 10.0) (17). Target genes were mapped to the STRING database to gain the relationship between each other and then drawn by Cytoscape (www.cytoscape.org).
To determine lung cell type enrichment, the transcript identifiers of the DEGs for each of the three major groups were converted into their associated gene names using BioMart software (http://uswest.ensembl.org/biomart/martview/). The three DEG lists were then input into the CTen database (18) to obtain the enriched cell types based on a highly expressed, cell-specific gene database. One-sided Fisher's exact test was used for enrichment analysis, and the CTen enrichment score was returned from the -log 10 of the Benjamini-Hochberg (BH) adjusted P-values. The CTen database is originally based on human and mouse data and has previously been used for chickens (19,20).  All sequence data discussed in this study have been deposited to NCBI's Sequence Read Archive (SRA) database, and the BioProject number is PRJNA597560.

Quantitative Real-Time PCR (qRT-PCR)
Total RNA of the lungs was extracted, and then cDNA was synthesized using FastKing gDNA Dispelling RT SuperMix (Tiangen, Beijing, China). Thirteen immune response-related genes, including Il6, Il1β, Il22, Il4i1, S100A9, Socs1, Socs3, AVD, AvBD1, AvBD2, AvBD4, TLR2A, and TLR4, were selected to confirm RNAseq results. Criteria for gene selection were based on immune response function and significance in the RNAseq study. Beta-actin was used as a reference gene to normalize transcript levels of the target gene. Specific primers were designed according to the reference sequences in NCBI with Primer-BLAST under the criteria: (a) PCR product size 70-200 bp, (b) Melting temperatures 60 ± 2 • C, and (c) primer must span an exon-exon junction (21). Primer sequences are listed in Table 1. The relative expression levels of genes are shown as a ratio of the target gene to the reference gene using the formula 2 −( Ct) (22).

Statistical Analysis
All data are displayed as means ± standard error of mean (SEM). GraphPad Prism 6.0 software was applied for all statistical analysis. Survival rates of chicks were evaluated using Kaplan-Meier analysis. All the other data were evaluated using unpaired, two-tailed Student's t-test. Significant differences were considered at P < 0.05. * P < 0.05, * * P < 0.01, and * * * P <0.001.

Mortality and Histopathology of Infected Chickens
Chickens infected with PmCQ2 and saline all survived, and PmQ-infected chickens all died rapidly within 24 h (Figure 1A). The bacterial burden of PmCQ2 in chicken lungs was approximately 10 5 CFU/mg at 4 and 8 hpi, maintained a high density at 16 hpi, then decreased to 10 1 CFU/mg at 24 hpi ( Figure 1B). The bacterial load of PmQ in chicken lungs was about 10 4 CFU/mg, on average, at 16 hpi ( Figure 1C). Although chickens infected with PmQ had many inflammatory cell infiltrates, PmCQ2 infected chickens had few inflammatory cells detected (Figure 1D). PmCQ2-infected chickens all recovered with 7 days with no clinical signs of infection.

Evaluation of Transcriptomic Sequencing Quality
A total of nine samples were generated from three biological replicates under three treatments; 44,374,891, 95,465,238, and 99,426,423 average clean reads were obtained from control and two infection group samples, respectively. The percentages mapped to the reference genome were more than 80% in all samples. A total of 19,119 genes were annotated in the Ensembl database for chicken genome. The average number of detected genes for individuals are 14,363, accounting for more than 75% of the reference genes. The percentage of uniquely mapped reads in all these samples was above 96%. Sequence data statistics are shown in Table 2.

Identification of DEGs in the Lungs After P. multocida Challenge
Based on the criteria: |log 2 foldchange| >1 and P < 0.05, there were 1,649 DEGs in P-C vs. PmQ-P, of which 831 were upregulated and 818 were downregulated; 1,427 DEGs (633 upregulated genes and 794 downregulated genes) in P-C vs. PmCQ2-P; and 318 DEGs with 232 upregulated genes and 86 downregulated genes in PmCQ2-P vs. PmQ-P (Figure 2A). The description of total DEGs is displayed in Supplementary Table 1. A total of 899 DEGs were common in P-C vs. PmQ-P and P-C vs. PmCQ2-P (Figure 2B), of which 31 DEGs were also found in PmQ infection vs. PmCQ2. The full gene list is displayed in Supplementary Table 2. Unsupervised hierarchical clustering analysis of global gene expression profiles revealed that each of the three biological triplicates clustered together ( Figure 2C).

Functional Analysis of DEGs
GO enrichment revealed that there were 34 GO terms enriched (of 3,864 annotated) in P-C vs. PmQ-P; 21 GO terms enriched (of 3,448 annotated) in P-C vs. PmCQ2-P ( Table 4). The detailed information is listed in Supplementary Table 3. In both P-C vs. PmQ-P and P-C vs. PmCQ2-P, "defense response" and "immune system process" categories were predominant. Although categories including "inflammatory response, " "response to cytokine, " and "response to chemical" are only enriched in P-C vs. PmQ-P, the enriched GO terms in P-C vs. PmCQ2-P are more often associated with leukocytes and granulocytes, such as "leukocyte migration, " "granulocyte chemotaxis, " "leukocyte chemotaxis, " and "granulocyte migration." KEGG pathway enrichment showed that 33 pathways were enriched in P-C vs. PmQ-P, of which the top 20 enriched pathways are displayed (Figure 3). Similarly, 30 pathways were enriched in P-C vs. PmCQ2-P (Figure 4). Many of these pathways are related to the immune system; signal transduction; the digestive system; and metabolism of lipids, amino acids, and xenobiotics. Both of the two contrasts clustered several conserved pathways, including "TNF signaling pathway, " "Tolllike receptor signaling pathway, " "Complement and coagulation cascades, " "Cytokine-cytokine receptor interaction, " "NF-kappa B signaling pathway, " "Fc gamma R-mediated phagocytosis, " and "Chemokine signaling pathway, " indicating that these pathways are related to P. multocida infection (Supplementary Table 4).
Furthermore, cell type enrichment analysis based on upregulated DEGs showed that different types of macrophages are the most enriched immune cell type in both PmQ and PmCQ2 infected chicken vs. uninfected chicken (Figures 5A,B). Other immune related cell types, such as microglia, osteoclasts, granulocytes, lymph nodes, CD8a+ Dend. cells myeloid, bone marrow, and mast cells, were all significantly enriched (Supplementary Table 3). Specifically, comparing the PmCQ2 infected group with PmQ infection, the most enriched cell type was also focused on macrophages (Figure 5C), which suggests that macrophage may be the key immune cell upon P. multocida infection in chicken. Moreover, the expression level of macrophage related genes was displayed with a heat map (Figure 6). The clustered gene group of each gene showed that TLR4, TLR2A, Macro, Ccr2, Ccr5, MRC1L-B, and MRC1L-D were upregulated in both PmCQ2 and PmQ infected groups. Moreover, genes such as Ifr9, Ifr8, Ifr7, CSF2RA, STAT1, and STAT2 were expressed to a greater extent in the PmQ infection group than PmCQ2, and Il22, Il6, Nos2, Ccl4, Ccl5, Ifr1, Socs1, and Socs3 were increased only in PmQ infected chickens rather than PmCQ2.

DISCUSSION
Pasteurella multocida is a very common animal pathogen and is classified into five serotypes (A, B, D, E, and F) depending on the capsular antigen (10). Among them, serotype A tends to have a wide host range (2,3). However strains isolated from nonavian hosts generally do not infect poultry, and the underling mechanisms are poorly understood (23). An avian P. multocida of serotype A, PmQ, isolated from tissues of a dead duck shows high virulence to poultry as well as to mammals. A bovine strain of serotype A, PmCQ2, from pneumonic lungs of a beef cow can cause severe pneumonia in mice (14,16) but is non-lethal to chickens. It has been suggested that some of the differences in host susceptibility to P. multocida infection may be due to differences in the host response expressed in the lung during the early phase of infection (24). In the present study, we infected chickens with PmCQ2 at an extremely high dosage (2 × 10 9 CFU), via the abdominal cavity to artificially create a chicken "infection model." However, the PmCQ2 is quickly eliminated by the hosts, and a dramatic drop between 16 and 24 hpi of the bacterial load occurs. Whereas, the PmQ infected chickens died rapidly within 24 hpi. Therefore, transcripts of the chicken lungs at 16 hpi were used for RNA-seq. The innate immune system is the first line in host defense and protects the host from bacterial infection (25). Previously, we found TLRs and inflammatory cytokines such as IL17, IL6, TNF-α, IFN-γ were significantly upregulated in mouse lungs infected with PmCQ2 (14). Another study also reports that IL6, IL17, and IL22 expression levels were significantly increased in the lungs of chicken challenged with avian P. multocida strain X73, which succumbed to acute infection (26). Consistent with these studies, the transcriptome of chicken lungs infected with PmQ also demonstrate that cytokines, such as IL6 and IL22, were significantly upregulated in PmQ-P vs. P-C, which implies PmQ infection induces a strong inflammatory reaction in the lungs. IL6 is a pleiotropic cytokine and can be both pro-and anti-inflammatory depending on the context (27). Impaired IL6 function causes enhanced susceptibility of mice to infection of Listeria monocytogenes (28), but excessive and sustained production of IL6 can be life-threatening, including systemic inflammatory response syndrome (SIRS) and chronic inflammatory diseases (29). IL22 is a member of the IL10 family and targets epithelial cells to induce various innate immune responses, which is essential for host defense against pathogen invasion (30). It can not only protect the tissue integrity and maintain mucosal homeostasis, but also act as a proinflammatory cytokine capable of amplifying the inflammatory reaction, and overproduction of IL22 may result in tissue damage (31). IL22 can mediate heme-associated iron scavenging away from pathogens during systematic infection (32,33).
In PmCQ2 infected chickens, the expression level of IL22 and IL6 were unaffected. Similarly, chickens challenged with capsular-deficient strain 775 with low virulence also fail to increase IL6 and IL22 expression (26). Further verification using q-PCR indicated that IL6 and IL22 expression levels were unchanged during PmCQ2 infection (Supplementary Figure 1), which implies that PmCQ2 causes a lower inflammatory response. CCL4, CCL5, and CCL28 were upregulated in PmQ and PmCQ2 infected chickens, but the expression level was increased to a high level in PmQ-P. CCL4, also known as macrophage inflammatory protein-1β (MIP-1β), is a chemoattractant for natural killer cells, monocytes, and a variety of other immune cells (34); CCL5 is chemotactic for T cells, eosinophils, and basophils and plays an active role in recruiting leukocytes into inflammatory sites (35); and CCL28 shows both antimicrobial activity and functions as a chemoattractant for T and B lymphocytes and migration of eosinophils. Together, these chemokines are ascribed to more immune cell infiltration, such as lymphocytes and eosinophils, in the PmQ infected chicken lungs. Whereas, the expression level of CCL19 and CCL26 were only increased in PmQ-P vs. P-C. CCL19 binding to its receptor CCR7 attracts dendritic cells to lymph nodes (36), and CCL26, also called macrophage inflammatory protein 4-alpha (MIP-4alpha), can drive eosinophils and basophil migration by binding to CCR3 (37). In sum, these data suggest a different immune response of the host infected with the two P. multocida strains.
In this study, a number of enrichment GO function categories were identified based on the DEGs of PmQ-P vs. P-C and PmCQ2-P vs. P-C. The protein-protein interactions based on DEGs annotated in immune system response of two infection groups (Figures 7A,B) indicated a different immune response interaction network. In PmQ infected chicken, IL6 and PTPRC (CD45) are in the center of the net. Upregulation of IL6 combined with STAT1 and STAT2/3, attributed to the key points of the Jak-STAT pathway, demonstrated that the STAT pathway was activated during P. multocida infection, which is consistent with our previous study (14). The activation of the STAT pathway can then trigger a signaling cascade required for immune response and inflammatory reaction (38), resulting in production of IFNγ, contributing to bacterial clearance (39,40). CD45, expressed on leukocytes, is a crucial player in the activation of T-cell receptor signaling by controlling the activation of the Src family protein-tyrosine kinases Lck and Fyn. CD45 deficiency can result in T-and B-lymphocyte dysfunction in the form of severe combined immune deficiency (41). Some immunomodulatory proteins, such as UL11 of human cytomegalovirus, can interact with CD45, resulting in function disruption of T-cells (42). Different from P. multocida toxin (PMT) (43,44), lipopolysaccharides isolated from P. multocida serotype A and B stains can induce apoptosis of bovine leukocytes, including lymphocytes, neutrophils, monocytes, and macrophages, resulting in tissue damage (45,46). In the current study, genes (FAS, APAF1, TNFRSF1A, CASP7, and CASP8) related to apoptosis were significantly upregulated in PmQ infected chickens but not in chickens infected with PmCQ2. In contrast to our previous study, the apoptosis pathway is enriched in PmCQ2 infected mice with a deadly outcome (14). Together, these results hint that the PmQ infected chicken attempt to cause cell death of the host, leading to tissue lesions of the lungs.
The NOD-like receptors (NLRs) are a specialized group of intracellular innate immune receptors that mediate IL1β production against bacterial infection (47). In the present study, the NOD-like receptor signaling pathway was enriched alone in the PmCQ2 infected chicken, which is consistent with the previous study of PmCQ2 infected mice (14). NLRC5, which can act as a negative modulator of inflammatory pathways, is critical for LPS-induced IL-10 production in macrophages (48), which upregulated in PmCQ2 infection, which suggests it might contribute to the lower inflammatory reaction in PmCQ2 infection. Our previous study revealed the NLRP3 inflammasome of mouse peritoneal macrophages infected with PmCQ2 played a role in mediating IL1β against P. multocida infection (49). In line with the study, IL1β was also significantly increased in PmCQ2 infected chicken, which might responsible for PmCQ2 clearance.
Macrophages play a crucial role in sensing P. multocida and secretion of inflammatory factors (49). We found that Lserine supplementation can decrease the macrophage-mediated inflammatory reaction in PmCQ2 infected mice and reduce the lung damage (50). In the current study, cell type enrichment analysis based on DEGs showed that macrophages were the most significantly enriched cell types, suggesting that macrophages were also activated in chickens infected with two different host-isolated P. multocida strains (Figure 5). Moreover, when analyzing the DEGs between PmCQ2 and PmQ infected chicken, macrophage was still the most enriched cell type. However, the different outcomes of chicken challenged with PmCQ2 and PmQ hint that macrophage may play different roles. For example, NOS2, a key marker, and PTX3 and IRF5, two important regulators of M1 macrophage activation (51)(52)(53), were significantly overexpressed in PmQ infected chicken, tending to activate to more pro-inflammatory activation. Whereas, in the PmCQ2 infected chicken lung, M2 markers (MRC1LD, MRC1LB, and MACRO) and M1 markers (CD86 and CD80) were upregulated, implying a mixed activation of both M1 and M2 macrophage occurred. In Marek's disease virus (MDV) infection, macrophages derived from differentially sensitive chicken lines undergo different immune strategies, and the macrophages from resistant chicken lines activate different biological pathways and suppression of oncogenic potential (54). Macrophages can be activated into different states depending upon the microenvironment and stimulus encountered (55,56), and the activation styles are often related to its morphology (57).
In addition to the immune response, metabolism was also affected upon P. multocida infection. For example, GAD1, a key enzyme that converts L-glutamate to GABA, used for succinic semialdehyde synthesis and eventually succinate, which is the main route of M1 macrophage and fundamental for stabilization of HIFα (58,59), was significantly upregulated in PmQ infected chicken while GLUL (also known as glutamine synthetase, GS), highly expressed in M2 macrophages and important to M2 phenotype (60), is increased in PmCQ2 infected chicken. However, PmCQ2 infection also caused reduction of all NADH dehydrogenase subunits (ND1∼ND6), implying the oxidative phosphorylation of chicken was impaired, which was unfavorable to M2 macrophage activation (61). A recent study reported that arachidonic acid can bind to myeloid differentiation factor-2 (MD2) and prevent MD2/TLR4 signaling activation, thus, it inhibits inflammatory response (62). Ten genes associated with arachidonic acid metabolism are upregulated in PmQ infected chickens, resulting in overproduction of a number of metabolites, for example, PGD2 (Prostaglandin D 2 ), which can bind to its receptor DP2 and amplify inflammatory response in asthma (63). Overall, these results suggest macrophages might undergo different polarization and contribute to opposite outcomes of chickens.
In conclusion, different host-isolated serotype A P. multocida strains, bovine PmCQ2 and avian PmQ, exhibited opposite pathogenicity to birds. Additionally, this study also revealed several crucial pathways related to anti-infection and metabolisms were enriched and that macrophages might play a key role in P. multocida infection.

DATA AVAILABILITY STATEMENT
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found at: https://www.ncbi.nlm.nih. gov/, PRJNA597560.

ETHICS STATEMENT
This animal study was reviewed and approved by Animal Ethics and Research Committee of Southwest University.

AUTHOR CONTRIBUTIONS
YP and NL conceived the experiments. PL conducted the experiments and wrote the paper. FH, CW, and GZ contributed to the experiments and data analysis. PH, YP, and NL supervised the project and revised the paper. All authors discussed the results and approved the final manuscript.

FUNDING
This work was supported by the earmarked fund of the China Agriculture Research System (Beef/Yak Cattle, CARS-37), Chongqing Science & Technology Commission (cstc2017shms-zdyfx0036 and cstc2017jcyjAX0288). The funding organizations had no role in study design, data collection or analysis, decision to publish, or preparation.