Figures
Abstract
Discovering genetic biomarkers associated with disease resistance and enhanced immunity is critical to developing advanced strategies for controlling viral and bacterial infections in different species. Macrophages, important cells of innate immunity, are directly involved in cellular interactions with pathogens, the release of cytokines activating other immune cells and antigen presentation to cells of the adaptive immune response. IFNγ is a potent activator of macrophages and increased production has been associated with disease resistance in several species. This study characterizes the molecular basis for dramatically different nitric oxide production and immune function between the B2 and the B19 haplotype chicken macrophages.A large-scale RNA sequencing approach was employed to sequence the RNA of purified macrophages from each haplotype group (B2 vs. B19) during differentiation and after stimulation. Our results demonstrate that a large number of genes exhibit divergent expression between B2 and B19 haplotype cells both prior and after stimulation. These differences in gene expression appear to be regulated by complex epigenetic mechanisms that need further investigation.
Citation: Irizarry KJL, Downs E, Bryden R, Clark J, Griggs L, Kopulos R, et al. (2017) RNA sequencing demonstrates large-scale temporal dysregulation of gene expression in stimulated macrophages derived from MHC-defined chicken haplotypes. PLoS ONE 12(8): e0179391. https://doi.org/10.1371/journal.pone.0179391
Editor: Yun Zheng, Kunming University of Science and Technology, CHINA
Received: January 13, 2016; Accepted: May 29, 2017; Published: August 28, 2017
Copyright: © 2017 Irizarry et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Data Availability: All data are available from the NCBI database repository (accession: PRJNA394901).
Funding: Funded by United States Department of Agriculture 2008-35204-04806 http://www.reeis.usda.gov/web/crisprojectpages/0214956-impact-of-immune-responses-of-chickens-with-defined-b-haplotypes-on-resistance-to-respiratory-coronavirus-infection.html and Merial Veterinary Scholars Program http://www.merialscholars.com/Pages/home.aspx.
Competing interests: The authors have declared that no competing interests exist.
Introduction
Discovering genetic biomarkers associated with disease resistance and enhanced immunity is critical to developing advanced strategies for controlling viral and bacterial infections in various species.
Disease resistance and susceptibility depends on a variety of factors including genetics. In numerous species, disease resistance has been associated with major histocompatibility complex (MHC) haplotype, as well as polymorphisms in several immune genes such as TGFβ and TNFα[1,2]. Cytokine production, specifically secretion of pro-inflammatory molecules, has also been associated with increased resistance against disease [3,4].
Studies have demonstrated association of MHC-B haplotype in chickens and resistance to a variety of viral pathogens, including AIV, Marek’s disease virus (MDV), avian leukosis virus, Newcastle disease virus and Rous sarcoma virus [5–10] as well as other pathogens [11,12]. B2 haplotype chickens are more resistant to avian coronavirus infection than B19 haplotypes and these differences in disease resistance were observed early after infection in our previous studies [10]. This suggests that innate immunity plays a major role with the macrophage being a key player in this enhanced immune response as evidenced by the B2 haplotype birds’ greater capability to produce nitric oxide (NO) in response to IFNγ and Poly I:C [13]. In addition, B2 macrophages activated T cells more efficiently than macrophages derived from B19 haplotypes [14].
Macrophages are directly involved in cellular interactions with pathogens and demonstrate distinct immune responses from more disease resistant animals in response to infection [15–20]. In addition, macrophages release cytokines activating other immune cells and antigen presentation to cells of the adaptive immune response [21–23]. It has become increasingly clear that dysregulation of macrophage function is involved in inflammatory disease processes such as rheumatoid arthritis, inflammatory bowel disease and cancer [24–26]. Involved in these interactions are crucial molecules such as Toll-like receptors (TLRs) that recognize invading microorganisms, resulting in communication with the adaptive immune system such as increased expression of MHC surface molecules, T cell receptors and secreted cytokines [21,23]. Genetic differences in any of those molecules can potentially account for differences in immune competence and thus provide potential immunogenetic markers for disease resistance to various pathogens.
IFNγ is a potent activator of macrophages and increased production has been associated with disease resistance in multiple species [27–31]. These findings indicate that chickens with enhanced IFNγ production are more resistant to certain infections. IFNγ enhances macrophage activation, expression of MHC and nitric oxide release which aides in killing of pathogens and also increases activity of cytotoxic T cells and secretion of Th1 cytokines [31,13], underscoring how crucial this process is for innate immune competence.
Macrophage TLRs appear to be primed by IFNγ, reprogramming cellular responses to other cytokines, such as type I interferons and IL-10 and activating the Jak-STAT pathway (Janus kinase and signal transduction and activator of transcription) [24, 32, 33]. IFNγ, which increases TLR receptor availability for interaction with its ligands, has been shown to induce TLR2, 4, 6 and 9 [34–37].
The response of macrophages to an immune stimulus is not just dependent on cell surface receptor and cytokine expression. Other factors include the differentiation of monocytes into functional macrophages, a tightly regulated process that influences immune competence [38]. Recent studies demonstrated a critical role for molecules such as A2B adenosine receptor for differentiation and proliferation of monocytes and macrophage function in immunity and inflammation [39, 40]. A2B expression is induced by IFNγ and leads to increase of anti-inflammatory signaling counteracting the inflammatory response activated within the IFNγ pathway.
Taken together, these studies emphasize the genetic basis of the activation of macrophages by IFNγ playing an important role in the innate immune response signaling and providing resistance to disease. In addition to inflammatory signaling, a number of transcription factor pathways and epigenetic mechanisms all contribute to immune function. Dysregulation of any of these events will lead to an impaired innate immune response and consequently, increased susceptibility to disease.
Using an ex vivo model, we investigated the gene expression in macrophages from haplotypes B2 and B19 during differentiation and after stimulation with IFNγ. Our experimental design leveraged an initial 6 day window for monocytes to differentiate into macrophages, which was followed by IFNγ stimulation between 1 and 24 h to further characterize subsequent RNA gene expression and the molecular basis for dramatically different nitric oxide production and immune function between the B2 and the B19 haplotype chicken macrophages
Material and methods
Experimental animals
Animal protocols were performed under the approval of the Institutional Animal Care and Use Committee at Western University of Health Sciences, Pomona, California (WesternU). Fertilized eggs, descended from Modified Wisconsin Line 3, were obtained from Dr. W. Elwood Briles, Northern Illinois University, and incubated and hatched under standard conditions at (38°C/50-65% humidity) [10,13] at WesternU. In addition to daily health monitoring, fresh food and water were provided ad libitum. Experimental animals were euthanized by insufflation of isoflurane gas (Butler, Dublin, OH).
Peripheral blood collection
Whole blood samples were collected via jugular venipuncture in EDTA tubes from age matched chicks at 14–18 weeks old. At no time did the amount of blood harvested from each animal exceed 1% of body weight.
Peripheral blood mononuclear cell (PBMC) isolation
Peripheral blood mononuclear cells (PBMCs) from individual birds were isolated using the differential centrifugation as previously described [41, 42, 13] with slight modifications. Briefly, blood was mixed with an equal volume of phosphate buffered saline (PBS) and slowly layered 2:1 on a Ficoll-Hypaque gradient (density 1.083) (Sigma-Aldrich, St. Louis, MO). Samples were centrifuged for 35 min (400 x g; 23°C; brake off) for retrieval of mononuclear cells. Isolated cells were washed 3x in 10 ml PBS at low speed to remove thrombocytes (180 x g; 10 min, 23°C), counted and viability confirmed based on the exclusion of 0.1% trypan blue dye (≥ 90%). PBMCs were re-suspended in PBS to a final concentration of 5 x 107 cells/ml.
Macrophage cell culture
One milliliter of PBMC suspension (5 x 107 cells/ml) was incubated (37°C/5% CO2) for 3 h in each well of a 12-well plate containing RPMI w/o Phenol Red supplemented with, 10% heat inactivated fetal bovine serum (FBS); non-essential amino acids, (0.1mM/ml) (Invitrogen, Carlsbad, CA), L-glutamine (2 mM) (Sigma-Aldrich, St. Louis, MO), 2-mercaptoethanol (55 μM/ml) (Sigma-Aldrich, St. Louis, MO), penicillin (50 U/ml) (Invitrogen, Carlsbad, CA), and streptomycin (50 μg/ml) (Invitrogen, Carlsbad, CA). Following removal of non-adherent cells with warm PBS, medium was replenished and cells were incubated for differentiation with the exception of the -6 day sample which was lysed with 400 μl Trizol (Thermo Scientific, Waltham, MA) and stored at -80°C. Prior to the replacement of medium, adherent cell cultures were washed in warm PBS. Monocytes were cultured for 6 days to allow maturation and differentiation of cells; with medium changes occurring every 3–4 days thus ensuring that optimal nutrient requirements were met. Additionally, -3 day (t-3) samples were lysed with Trizol and stored at -80°C. The morphology of adherent cells was observed daily under bright field microscopy (20x objective).
Purity of monocyte cultures using this culture method was confirmed by IFA and FACS using monoclonal antibody KUL01 as previously described as part of a different aspect of this study [13].
IFNγ stimulation
A 50 ρg/ml ch-IFNγ solution (Invitrogen, Carlsbad, CA) was prepared in RPMI w/o Phenol Red culture medium (Invitrogen, Carlsbad, CA). After washing the cells twice with warm PBS, macrophage cultures were stimulated with 1 ml of RPMI-ch-IFNγ mixture [13].
Nitric oxide assays
Nitric oxide production was measured [10, 43, 44] to confirm macrophage stimulation in assays by interferon (data not shown). Stimulation was evaluated as yes/no based on previously published results from B2 and B19 IFNγ stimulated macrophages (10)
Sample collection and RNA sequencing
A total of 145 gigabytes of RNA sequence data was obtained from B19 and B2 haplotypes of chickens. Two birds from each haplotype were selected for inclusion in the sequencing. Each bird provided blood for extraction and isolation of peripheral blood mononuclear cells. Purified monocytes were cultured for differentiation and cell samples were collected from nine time points for each bird. Samples were collected for sequencing on the day they were cultured (Day t-6), as well as on Day -3 (t-3), Day 0 (called 0 hours), and then six additional times over a 24 hour period corresponding to 1 hour, 2 hours, 4 hours, 8 hours, 16 hours and 24 hours after interferon stimulation. Cells were lysed in wells with RLT buffer containing beta-mercaptoethanol (Qiagen, Valencia, CA) and stored at -80°C. RNA was processed with the Qiashredder and RNAeasy kit from Qiagen (Valencia, CA) according to manufacturer’s instructions and sent on dry ice to Dr. Calvin Keeler at the University of Delaware for generation of libraries and sequencing with an Illumina HiSeq 2000.
An RNA sequence library was constructed from purified RNA. The library was fragmented in order to generate appropriately sized RNA fragments suitable for templates in random primed first-strand cDNA synthesis. Second strand synthesis was completed in accordance with specifications for sequencing with Illumina’s HiSeq2000 platform.
The samples corresponding to each time point from each bird were sequenced and the data was stored in a unique file for each sequenced sample and time point. Forty FASQ files were generated from the data totaling 145 gigabytes. The average file size was 3.65 gigabytes and the standard deviation was 2.25 gigabytes. The sequencing data provided 933,107,885 reads across the biological samples and time points (Table 1). Across all time points for the two B2 samples, one produced 298,903,517 reads and the other produced 165,589,594 reads. Similarly, across all time points, the B19 samples produced 285,392,384 reads and 183,222,390. For each time point (across all four birds) sequencing reads ranged from a low of approximately 78 million reads to a high of just over 171 million reads, with most time points producing over 88.4 million reads each and a few producing over 100 million reads each.
Mapping reads to reference genome and identification of splice junctions / exons
The chicken reference genome WASHUC2, corresponding to Ensembl release 70, was downloaded from Ensembl.org (http://www.ensembl.org/info/data/ftp/index.html). Annotation files included the small RNA annotation files obtained from miBase release 19 (http://www.mirbase.org/). Sequenced reads were filtered to remove low quality sequences from the data. Filtered sequences were aligned to the reference genome using Bowtie and Tophat, available along with the software package Cufflinks, from John Hopkins University Center for Computational Biology (https://ccb.jhu.edu/software.shtml). The aligned reads generated by Bowtie produced gapped alignments on the reference genome which Tophat used to identify splice junctions flanking exons. The resulting aligned reads were analyzed by Cufflinks to construct transcripts corresponding to mRNA sequences. Next, Cufflinks was employed to estimate transcript specific expression levels across the transcripts and genes within the reference genome based on the number of sequence reads for each mRNA. The sequence read data was normalized using the fragments per kilobase of transcript per million mapped reads (FPKM) method to more accurately determine expression levels. The resulting transcriptome data was loaded into the MySQL relational database to more effectively manage, explore, mine and annotate the data.
Hierarchical clustering of genes and production of heat map visualization
Gene expression data was hierarchically clustered using 1-Pearson correlation on the rows and keeping the column order conserved. The resulting clustered data set was visualized as a heat map with black representing lack of gene expression, and darker shades of blue indicating lower expression values. Dark purple represents higher expression values than any shade of blue while bright pink represents the highest expression values. For visualization purposes, the heat maps were generated with maximum heat map color assigned to expression set lower than the absolute maximum expression value contained in the entire data set, subsequently all values of expression greater than or equal to the assigned expression threshold (for example, 1000) shared the same color on the heat map (regardless of whether the actual expression level was 1000, 2000, 20,000 or 90,000). This setting provided the optimal visualization of both high and low expressed genes in the heat maps.
Gene enrichment calculations were performed using the DAVID bioinformatics database tool version 6.8 (https://david.ncifcrf.gov/). The analysis was performed using comparisons of successive time points within the B2 haplotype data to identify sets of genes that were enriched (p-value < 0.05). The B2 haplotype represents the robust macrophage phenotype as characterized by nitric oxide production compared to the B19 haplotype. Subsequently, the gene enrichment was performed on the B2 data. Gene enrichment was determined using three distinct databases: gene ontology biological process, KEGG pathways, and reactome pathways corresponding to S2, S3 and S4 Tables respectively. Because a large number of enrichment annotation terms were produced, a subset of representative highlights from each of these three enrichment analyses was chosen for inclusion in the results. Highlights were selected to provide examples of the biological process annotations, KEGG pathways annotations and reactome annotations.
PCR validation of target genes
Realtime PCR was performed on a selected number of target genes to validate RNA sequencing results. RNA was taken from macrophages stimulated with IFNγ as described above for 2 and 4 hours, unstimulated samples (0h) served as control. For Realtime RT-PCR, cDNA synthesis was performed using SuperScript III First Strand Synthesis kit (Invitrogen, Waltham, MA), according to manufacturer’s instructions. PCR conditions were as follows: 95°C for 10 min-hot start, 40 cycles of 95°C for 15 sec, 60 or 63°C depending on gene (see primers) for 30 sec according to manufacturer instructions for the Biotool 2x Sybr Green qPCR Mix (Biotool, Houston, Tx). Primer sequences were designed using Primer 3 (ATP6VOC, LITAF, IL18R, TLN-1.) Primer sequences previously published were used for TLR2, TLR4, TLR5, TLR6 and TLR7 [45]. ATP6VOC (annealing 60°C) forward TGTTGTCATGGCGGGTATTA, reverse ACAAATAACCTGGGCTGCTG; LITAF(annealing 60°C) forward ATCCTCACCCCTACCCTGTC, reverse GACGTGTCACGATCATCTGG; IL18R (annealing 63°C) forward CTCTTCGTGCCTCCATTGAT, reverse ACCAAGTTCAACTGGCCAAA; TLN-1(annealing 60°C) forward TCAAGCAGAAGTTGCACACC, reverse GGGAGCCATTAAGGATGTCA. PCR analysis was done using the ΔΔ method with 18s serving as housekeeping gene control. Statistics were done using graphpad software (PRISM version 7), paired t-test, two-tailed.
Comparison of IFNγ stimulated vs. cytomegalovirus stimulated macrophage gene expression
A total of 179 gene expression measurements were extracted from a published paper describing the fold change in expression levels of genes induced after 4 h exposure to cytomegalovirus [46]. The data was converted to a tab-delimited text file containing the official gene symbol and the reported expression level. The file was loaded into a MySQL relational database and joined to the expression data produced from the B2 and B19 cells. The data was joined on the gene symbol and a set of 54 genes were identified. The fold change for the B2 and B19 expression data was calculated by taking the log-2 (4 h expression / 0 h expression). B2 and B19 genes having expression = 0 for the initial time point were converted to 0.1 to prevent division by zero. Additionally, the fold-change reported for IL6, 280.8, was changed to 35, in order to preserve the scale of the graphs and legibility of the resulting data represented in the histograms. The fold-change in expression for the B-haplotype birds and the published data was plotted using Microsoft Excel.
Results
Differential gene expression patterns
A set of 13,618 unique genes from among all mapped sequencing reads was generated from the 4 birds across all 9 time points. Next, we analyzed the expression data to determine the number of genes expressed in each haplotype within each time point. Within the minus 6-day (t-6) time point, representing the time point after plating and adherence of monocytes and the start of differentiation into mature macrophages, 11,785 genes were expressed in the B2 birds while 12,089 were observed in the B19 birds, with 11,216 genes expressed in both. Interestingly, 4770 genes were off in both B19 and B2 haplotypes while just 569 genes exhibited expression in only the B2 chickens and 873 genes were expressed only in the B19 birds.
Similar relationships were detected in each of the remaining eight time points. The t-3 day time point, representing 3 days of differentiation in cell culture, exhibited the greatest expression of genes with a total of 11,429 expressed in both B19 and B2 birds while just 4068 genes lacked evidence of expression in both haplotypes. Also, during the t-3 day time point the greatest number of genes (1118) exhibit evidence of expression in the B2 birds while lacking evidence of expression in the B19 birds. At the t0 time point, after 6 days of differentiation and immediately before stimulation with interferon, 10,975 genes were expressed in both haplotypes while 4547 genes were not expressed in macrophages of either haplotype. Likewise, the 1 h and 2 h time points exhibited 11,349 genes and 10,789, on in both haplotypes, respectively. It is worth noting that the time point with the most genes off in both haplotypes is 16 h with 5238 genes.
Overall the data indicates that approximately 10,000 to 11,000 genes are on in both haplotypes at each time point while roughly 4000 to 5300 genes are off in both haplotypes at each time point. The number of genes on in one haplotype, while off in the other haplotype, ranges from about 400 to 1140 depending upon the haplotype and time point (Fig 1)
Visual representation of genes within B2 and B19 haplotypes at each of the time points. Figure includes genes expressed in common, genes expressed only in B2, genes expressed only in B19, all genes expressed in B2, all genes expressed in B19, and total non-redundant genes expressed in either B2 or B19 haplotypes.
Differences in numbers of genes expressed in B19 versus B2 haplotype birds
In order to better understand the cell biology underlying differences in macrophage differentiation and activation between B19 and B2 birds, we searched for genes exhibiting statistically significant differences between different time points within a single B-haplotype haplotype as well within the same time point between haplotypes.
When comparing the expression profiles between the B2 and B19 haplotypes, we identified 210 genes exhibiting differential expression at the t-6 day time point. These genes represent 198 genes with higher expression in the B2 birds and just 12 genes for which expression was greater in the B19. After three days, at the t-3 day time point, thousands of genes exhibited altered expression patterns between the two groups. Surprisingly, 7000 genes showed higher expression in the B19 birds while only 14 genes were expressed at higher levels in the B2 birds.
By t0 hrs, which corresponds to 6 days of monocyte differentiation into macrophages, we observed 955 genes with significant expression patterns between the haplotypes. Of these genes, 544 exhibited greater expression in the B2 haplotype while 411 exhibited higher expression in the B19 haplotype.
Cells were stimulated with IFNγ immediately following RNA collection at the t0 hr time point. At 1 h (t1) post-stimulation 665 genes show evidence of significant patterns of expression between the haplotypes where B19 birds had 109 genes expressed to higher levels while the B19 haplotype was associated with 556 genes having greater expression compared to t0. This pattern of increased expression in the B19 group is reversed by the 2 h time point.
At 2 h after IFNγ treatment, the B2 cells show a global increase in expression for 5989 genes while the B19 cells have just 18 genes on at higher levels than the B2 birds. By 4 hours after stimulation, the B2 birds still exhibit greater expression for 1029 genes while the B19 birds exhibit higher expression for 12 genes. This trend changes by 8 hours after treatment, at which time the slower responding B19 group begin showing increased expression in 797 genes while the B2 cells have greater expression for just 15 genes. By 16 hours after stimulation, only 66 genes are differentially expressed between the two haplotype groups. And, at the 24 hour mark, 406 genes show evidence of statistically significant differences in expression between them with the B2 cells exhibiting greater expression for 339 genes while the B19 cells have higher expression for 67 genes (Table 2).
Different temporal gene expression in B19 versus B2 haplotype birds
The B2 and B19 haplotype birds represent distinct genetic variation within the B-locus on chromosome 16. Subsequently, patterns of gene expression variation of the genes located within this region were investigated. Among the seventeen genes exhibiting statistically significant differences in expression between the B2 and B19 birds, many displayed divergent gene expression patterns prior to IFNγ stimulation. In the B2 cells, gene expression peaks on day t-6 and expression is effectively inhibited by day t-3. This is not the case in the B19 cells. Rather than reach maximum expression levels in a single day, the B19 cells don’t achieve maximum expression until day t-3 (Fig 2).
B-locus haplotypes in chickens provide a mechanism for genetically perturbing the cluster of immunologically important genes on chromosome 16 and producing phenotypic variation affecting infectious disease susceptibility and resistance. The heat map allows visualization of gene expression between the two genetically distinct haplotypes. Each row represents a gene within the B-locus (listed on the right) and each column corresponds to a particular time point when cells were collected for RNA sequencing. Black pixels indicate zero gene expression for a particular gene at a specific point in time, and dark blue corresponds to very low expression, while brighter blue indicates the next higher levels. Dark purple represents higher expression levels than blue colors, and pink represents the highest levels of gene expression. Monocytes were obtained from each haplotype of chicken and allowed to differentiate into macrophages in vitro for seven, days beginning on day minus 6 (t-6). RNA was sampled on day t-6, day t-3, and again three days later which is denoted as 0 hours (t0), when IFNγ was initially added to the cultures. On t0, RNA was sampled immediately before stimulation with IFNγ. Subsequent time points correspond to the time following interferon stimulation, in hours (1 hour, 2 hours, 4 hours, 8 hours, 16 hours and 24 hours). As visible on the heat map, there are distinct differences in gene expression between the B2 and B19 cells. The most dramatic difference occurs on day t-6. B2 cells exhibit a rapid burst of gene expression, indicated as a single column of pink on the left most edge of the heat map. In contrast, the B19 cells appear to undergo a much slower and prolonged gene expression program that was not as rapidly down regulated as in genes in the B2 cells. Additional gene expression data for a number of proteins involved in cell growth and apoptosis, is shown in the bottom half of the figure to highlight a similar pattern in gene expression and kinetics. The green border indicates the B2 haplotype expression pattern and the red border corresponds to the B19 expression pattern.
For example TRIM7, TRIM27.1, BF2, TPN, and TRIM41 exhibit strong expression on day t-6 in the B2 cells while the same genes exhibited prolonged expression over day t-6 and day t-3 in the B19 cells. Members of the TRIM (tripartite motif) family have been implicated in antiviral immune defense and several are ubiquitin ligases [47, 48]. TPN (Tapasin) is a co-factor for MHC I critical for antigen presentation to cytotoxic T-cells and chickens express the single MHCI locus termed BF-2 which is working with TPN in antigen presentation and it has been shown that there are differences in the selection of high affinity peptides in B19 vs B15 haplotypes [49] highlighting their critical role in immune competence. Additional genes within the B-locus display a similar pattern of pre-stimulatory differences in gene expression between the two different haplotypes, including genes involved in differentiation, cell growth and apoptosis such as PTPN2 (tyrosine protein phosphatase non-receptor2) and NFKB. Gene expression decreases to approximately baseline levels by time point t0 hours.
A second distinction in the gene expression patterns between B19 and B2 cells is that B2 cells exhibited a fairly robust expression at 2 and 4 hours after interferon stimulation. Unlike the B2 haplotype, the B19 haplotype appears incapable of generating such a rapid, robust and coherent gene expression profile. In contrast, the B19 cells generate a delayed, weak and uncoordinated lower level of expression that extends up to 8 hours, and in some cases even 16 hours. Overall, this global pattern of temporally dysregulated gene expression represents a re-occurring theme with the B19 monocytes and macrophages.
The divergent timing of gene expression observed in the B-locus genes is mirrored in many other genes as well, including members of the TLR signaling pathway, cellular mediators of apoptosis and cell survival, and components of cytokine signaling.
B2 and B19 display different patterns of gene expression during differentiation
The global dysregulation of gene expression among 700 genes at the t-3 day time point, as well as the expression pattern of 6000 genes exhibiting altered expression led us to explore the pattern of gene expression changes within each haplotype group over all of the time points. At the onset of the study, the B2 cells were actively expressing a diverse set of genes, however by the day t-3, most of those genes displayed reduced expression in the B2 group. Even so, the B19 haplotype cells continue to express these 7000 genes at higher levels than the B2 birds. After stimulation, B2 macrophages again show different patterns of expression compared to B19 cells in regards to timing of peak expression and coherence of expression. Four distinct patterns of divergent gene expression were identified between the B2 haplotype birds and the B19 haplotype birds (Fig 3).
Four distinct patterns were identified as representative of the types of divergent gene expression that re-occur across many genes involved in macrophage differentiation, activation and function in B2 versus B19 macrophages. 1. Day t-6: B2 high vs B19 low. This divergent pattern exhibits strong expression of genes on day -6 in the B2 birds while relatively low levels of expression are observed in the B19 birds at the same time point. Genes of interest include an adenosine receptor (P2RY12) 2. Day t-6: B2 = 1 day vs. B19 = 3 days. This example of divergent patterns is the single peak of day t-6 gene expression in the B2 haplotype cells compared to the prolonged multiple day expression until day t-3 in the B19 haplotype cells. Genes of interest include macrophage differentiation gene GATA, adenosine receptor A2A and macrophage podosome markers VCL and GSN. 3. Maximum IFNγ Stimulation of B2 at 2–4 h versus 4–8 h in B19 macrophages. Another interesting divergent gene expression pattern observed between the two haplotypes occurs after stimulation by IFNγ. There is a four-hour difference in peak expression timing for a large number of induced genes. In the B2 haplotype macrophages, the peak expression occurs between 2 and 4 hours, while in the B19 macrophages, the peak expression occurs between 4 and 8 hours. 4. Maximum IFNγ Stimulation: B2 = Coherent vs. B19 = Non-Coherent Another discernable difference in post-stimulatory induction of genes between the B2 macrophages compared to the B19 macrophages is one of coherence. Specifically, there are a number of genes for which the B2 macrophages are able to rapidly turn on and reach relatively high levels of expression within 2 to 4 hours of IFNγ stimulation. In contrast, these same genes fail to exhibit a coherent peak of expression, even after 4 to 8 hours, in the B19 cells. Instead, they exhibit a dispersed “smear” of gene expression extending from approximately 1 hour after stimulation to 16 hours post-stimulation.
The first interesting divergent pattern shows strong gene expression on day t-6 in the B2 birds while relatively low levels of expression are observed in the B19 birds on the same day. This pattern is of interest because it represents a group of genes that are differentially regulated at the onset of the experimental time course. Specifically, these genes include the macrophage M1 marker PTGS2, as well as the B-locus gene cyp21. Other genes exhibiting this pattern include secreted interleukin ligands IL-1β, IL4I1, and IL6, along with genes associated with inhibition of cellular processes including IRG1 and MIP-3α. Interestingly, the adenosine receptor also displays this pattern of expression. These genes may represent initial modulators of divergent monocyte to macrophage differentiation between the B2 and B19 cells.
The second example of divergent expression patterns is the single peak of day t-6 expression in the B2 haplotype cells compared to the prolonged multiple day expression in the B19 haplotype cells. Some of these genes are macrophage differentiation mediators, like GATA2 [50], and FADD, while others are macrophage podosome (primary matrix structure) markers, including VCL and GSN. Other genes exhibiting this divergent expression pattern include chemokine receptors, like CxCR4, fatty acid transport, such as SLC25A17, and ubiquitin related factors, like DD5, which is associated with proteasomal degradation of gene products.
Additional interesting divergent gene expression patterns were observed between the two haplotypes occurring after stimulation by IFNγ (Fig 3). A notable difference in post-stimulatory induction of gene expression is a four-hour difference in peak expression timing for a large number of induced genes. In the B2 haplotype macrophages, the peak expression occurs between 2 and 4 hours, while in the B19 macrophages, the peak expression occurs between 4 and 8 hours. Some of the most noticeable genes exhibiting this divergent gene expression pattern include LITAF, IL-1β, IL12, and IFIH1, genes involved in macrophage signaling and M1 macrophage polarization [26]. Additionally, a number of genes implicated in invadosome assembly and function also exhibit this temporally displaced pattern of induction such as CD44, RAC1, and SRC.
Another discernable difference in post-stimulatory induction of genes between the B2 macrophages compared to the B19 macrophages is one of coherence (Fig 3). Specifically, there are a number of genes for which the B2 macrophages are able to rapidly turn on and reach relatively high levels of expression within 2 to 4 hours of IFNγ stimulation. In contrast, these same genes fail to exhibit a coherent peak of expression, even after 4 to 8 hours, in the B19 cells. Instead, they exhibit a dispersed “smear” of gene expression extending from approximately 1 hour after stimulation to 16 hours post-stimulation. Some of the most represented genes exhibiting this divergent pattern of expression include molecules involved in lysosome function and phagocytosis. CTTN and ACTR3, genes implicated in FcR mediated phagocytosis, along with lysosomal-associated molecules, like LAPTM5 and LAMP1, as well as the lysosomal transporter molecules ATP6AP1, ATP6V1G1 and ATP6V0C, exhibit this non-coherent pattern of expression in the B19 macrophages.
In contrast, immediately following stimulation, the B2 cells rapidly induce expression of roughly 6000 genes by the 2 h following stimulation; while, at the same time, the cells derived from the B19 birds show no signs of induction among these genes until after 4 h. It is interesting to note that while the B2 birds show a statistically significant increase in expression for 6100 genes between 1 h and 2 h, the B19 cells exhibit increased expression for just 66 genes at this time point. The largest wave of increased gene expression occurs in the B19 cells during the transition from 2 h to 4 h post stimulation, when 1164 genes increase significantly over this time period.
At the transition between 8 h and 16 h, the B2 haplotype group only exhibits differences in expression for 83 genes, with 44 having higher expression at the 16 h time point. Yet, the B19 cells show differences in 386 genes during this same period, but interestingly, 356 of these genes exhibit decreased expression during this same time interval. Taken together, these results suggest that a global disruption of temporal gene expression underlies the observed differences in differentiation, activation and nitric oxide production from macrophages derived from the two different MHC haplotypes.
RT-PCR of B2 and B19 haplotype cells following IFNγ stimulation
Gene expression was measured in separate samples of B2 and B19 cells following stimulation with IFNγ. Change in expression was assessed at 2 hours and 4 hours post stimulation. ATP6V0C exhibited the greatest induction of all genes assayed, showing an increased expression in the B2 cells at 4 hours that was 20 times the initial expression at 0-hours. Expression of ATP6VOC was dramatically less in the B19 birds. Similarly, IL18R exhibited greater than 9 times the initial expression in the B2 cells at 4 hours compared to the B19 cells which exhibited less than 2 times the initial expression at 0-hours. LITAF and TLR2 exhibited more than 7 times the expression at 4 hours in the B2 macrophages, while TLN-1, TLR-5, TLR-6 and TLR-7 exhibited greater than 4 times the initial expression in the B2 macrophages. In contrast, the B19 macrophages failed to exhibit comparable induction of these genes (Fig 4).
Gene expression for ATP6V0C, LITAF, IL18R, TLN-1, TLR2, TLR3, TLR4, TLR5, TLR6, and TLR7 was assessed in B2 and B19 monocytes/macrophages following stimulation with IFNγ. Expression was measured at 0 hours, 2 hours and 4 hours. Expression for transcripts in B2 cells are shown in green and expression for transcripts are shown in red. Standard error is shown for each value. Values were considered statistically significant with p<0.05.
IFNγ stimulated vs. cytomegalovirus stimulated macrophage gene expression
In addition to the RT-PCR validation of gene expression, 54 genes, for which gene expression changes were described following cytomegalovirus stimulation were used as comparisons for the corresponding genes in the B2 and b19 haplotype birds (Fig 5). A total of 25 published genes exhibited decreases in expression following cytomegalovirus stimulation while 29 genes exhibited increased expression following stimulation. Interestingly, all but one gene (FEZ1) in the B2 cells exhibited increased expression following IFNγ stimulation. In contrast, ten genes displayed decreased expression in the B19 cells. Of the ten exhibiting fold-change < 0 in the B19 cells, 70% also exhibited decreased expression in the cytomegalovirus stimulated cells. In total, 28 genes (52%) expressed in the B2 cells matched the direction of the fold change reported in the published data while 33 genes (61%) corresponded between the B19 cells and the published data. Of the ten published genes reported as having greater than 5-fold increased expression, 90% of the B2 genes exhibited fold-change in the same direction. Overall, this data, in conjunction with the RT-PCR data, provides a comprehensive set of validation data providing evidence that the B2 and B19 gene expression data is reproducible and similar to expression patterns observed in cells stimulated towards macrophage activation pathway.
54 genes, for which gene expression changes were previously described following cytomegalovirus stimulation were used as comparisons for the corresponding genes in the B2 and B19 haplotype birds. A total of 25 published genes exhibited decreases in expression following cytomegalovirus stimulation while 29 genes exhibited increased expression following stimulation. All but one gene (FEZ1) in the B2 cells exhibited increased expression following IFNγ stimulation. In contrast, ten genes displayed decreased expression in the B19 cells. Of the ten exhibiting fold-change < 0 in the B19 cells, seven exhibited decreased expression in the cytomegalovirus stimulated cells. Twenty-eight genes (52%) expressed in the B2 cells matched the direction of the fold change reported in the published data while 33 genes (61%) corresponded between the B19 cells and the published data. Of the ten published genes reported as having greater than at least 5-fold increased expression, 90% of the B2 genes exhibited fold-change in the same direction.
Divergent non-coding RNA expression in B2 and B19 macrophages
Visualization of gene expression via heat maps facilitated the identification of distinct expression patterns between the B2 and B19 haplotypes. Because any initial differences in gene expression existing 6 days before IFNγ stimulation represent candidates responsible for the observed phenotypic differences between the two haplotypes. Genes exhibiting divergent gene expression patterns between B2 and B19 birds on day -6 were identified (Fig 6). The genes cluster into four major clades (clade1, clade2, clade3, and clade4 with a singleton labelled clade 5). Among these genes, represented in clade1 and clade2, are a number of miRNAs exhibiting strong expression in B2 cells (mir-147, mir-146b, mir-1618, mir-200a, mir-1649, and mir-1648a) compared to the B19 samples. Likewise, miRNAs contained in clade3 and clade4 exhibit greater expression in B19 cells (mir-1627, mir-222b, mir-1633, and mir-19a).
Visualization of divergent gene expression patterns between the B2 and B19 haplotypes. A subset of genes exhibiting divergent gene expression were identified and visualized in heat following hierarchical clustering of the genes (rows), but not the time points (columns). The genes cluster into four major clades (clade1, clade2, clade3, and clade4) with a singleton gene (labelled clade 5). Among these genes, represented in clade1 and clade2, are a number of miRNAs exhibiting strong expression in B2 cells (mir-147, mir-146b, mir-1618, mir-200a, mir-1649, and mir-1648a) compared to the B19 samples. Likewise, miRNAs contained in clade3 and clade4 exhibit greater expression in B19 cells (mir-1627, mir-222b, mir-1633, and mir-19a). Additionally, a number of small nucleolar RNAs (snoRNAs) exhibit similarly dichotomous gene expression patterns (clade4) such that SNORd24, snoZ40, SNORD74, SNORA17, and SNORD12 exhibit substantially higher levels of expression in B19 cells on day -6 compared to B2 cells while B2 cells express such as snoU2_19 (clade2).
A number of small nucleolar RNAs (snoRNAs) exhibit similarly dichotomous gene expression patterns (clade4). For example, SNORd24, snoZ40, SNORD74, SNORA17, and SNORD12 exhibit substantially higher levels of expression in B19 cells on day -6 compared to B2 cells (Fig 6). However, B2 cells also express snoRNAs exhibiting divergent expression patterns between the two haplotypes, such as snoU2_19 (clade2).
In addition to non-coding RNAs, divergent expression patterns are also observed with protein-coding RNAs (Fig 6). For example, clade1 contains IL6, IL18, IL-1β, CCL1, PTPN2 and MMP10, which exhibit higher initial expression on day -6 in the B2 birds. In contrast, the protein-coding genes LAMP2, UBXN7, UBE4B, PK3CA, UBE2W, and CX3CR1, in clade3, exhibit higher initial expression patterns in B19 birds. Clade5 contains the single gene IFNγ, which exhibits relatively low expression early in both B2 and B19 cells, but following stimulation rises to a higher level at 8 hours in the B19 birds.
Considering the diverse expression patterns discovered and the results indicating the involvement of non-coding RNAs, further results including divergent non-coding RNA expression will be described in more detail in further publications.
Gene enrichment analysis
Enrichment analysis of genes exhibiting statistically significant differences in expression between time points and/or haplotypes (Table 3, S2, S3 and S4 Tables) was performed using gene ontology and both KEGG and reactome pathways. The results of the gene enrichment provided a high-resolution perspective of the functional role of mRNA sequenced within the B2 cells across the experimental time points.
Between the -6 day and -3 day time points, a large number of genes exhibit reduced expression. These genes are enriched for biological processes such as transcription, mRNA splicing, tRNA processing and negative regulation NFκB mediated gene expression. Similarly, KEGG pathways enriched include RNA polymerase, RNA degradation, Metabolic pathways and Basal transcription factors. The reactome enriched pathways mirror these results with annotations of RNA polymerase II initiation, mRNA capping, and some immune functions including antigen presentation and cytosolic sensors of pathogen-associated DNA.
The transition from day -3 to t0 (just prior to IFNγ stimulation) correlates with down regulation of genes associated with chromosome segregation, mitotic nuclear division and DNA repair, cell cycle pathways, acetylation and some immunological functions of interleukin 3 and 5 signaling and interleukin receptor signaling. Genes exhibiting increases in expression during this interval were enriched in processes and pathways related to interleukin 4, biosynthesis of amino acids, and metabolic pathways. Within an hour of IFNγ stimulation genes exhibiting increased expression were associated with inflammatory and defense responses, toll-like receptor signaling pathways, cell chemotaxis, MyD88 signaling, Jak-STAT signaling, cell adhesion, FoxO signaling, and raf activation.
Genes exhibiting an increase in expression within two hours of IFNγ stimulation are enriched for biological processes of intracellular protein transport, ER-to-Golgi vesical mediated transport, endosome to lysosomal transport, chromatin remodeling, histone H3 acetylation, regulation of vesicle fusion and protein import into the nucleus. Among the KEGG pathways that exhibit enrichment for these genes are RNA transport, protein export, lysosome, mRNA surveillance pathway, endocytosis and mTOR signaling. Similarly, the enriched reactome pathways mirror these processes and include cytosolic sensors of pathogen-associated DNA, RNA polymerase II initiation and promoter clearance, mTORC1-mediated signaling, MAP2K and MAPK activation, downstream TCR signaling, FCε Receptor 1 mediated MAPK activation, and Clec7A signaling.
Genes exhibiting decreased expression between 2 hours and 4 hours following IFNγ stimulation correspond to reduced expression of cell-cycle pathways and mediators, as well as genes implicated in G1/S transcription, DNA replication, and separation of sister chromatids. Biological processes identified within these genes include mitotic spindle assembly, chromosome segregation, and mitotic spindle checkpoint assembly. Conversely, genes exhibiting increased expression during this same time are enriched for biological processes of inflammatory response, defense response to virus, positive regulation of interleukin 12 production, negative replication of viral genome replication, bacterial defense processes and positive regulation of IL1α secretion. KEGG pathways associated with these genes include cytokine-cytokine receptor interaction and genes implicated in influenza A signaling. Reactome processes identified included stem cell population maintenance and TOR signaling.
Over the remaining time points, from 4 hours to 8 hours, from 8 hours to 16 hours and from 16 hours to 24 hours the B2 cells exhibit a systemic down regulation of the genes that were initially activated during the IFNγ stimulation. Overall, the gene enrichment analysis of the RNA sequence data provides a cellular-level picture of the specific biological processes that occur over time following activation of monocyte-derived macrophages.
Discussion
Previous work in our laboratories investigated the association between chicken haplotype and disease resistance, specifically the enhanced resistance of B2 haplotypes to avian coronavirus IBV [10] and the influence of innate immunity leading to decreased clinical signs of illness. We showed that macrophages play an important role in this enhanced immunity, demonstrating much better activation in response to stimulation [13]. To analyze the gene expression involved in this process leading to increased macrophage nitric oxide release in B2 haplotypes, we stimulated macrophages from B2 and B19 chicks for RNA sequencing. In addition, we had observed different cell morphology when isolated monocytes from B2 and B19 haplotypes were differentiating into macrophages and therefore time points before stimulation, during differentiation of the macrophages, were included in this study.
The rationale for investigating the gene expression differences between the B2 and B19 haplotype birds was to address underlying questions that were raised at the end of our previous studies: 1. Why do the IFNγ stimulated B19 derived macrophages exhibit decreased nitric oxide production compared to the IFNγ stimulated B2 derived macrophages? 2. How do the two lineages of macrophages differ at the gene expression level? 3. What specific patterns of gene expression correlate with divergent macrophage differentiation, activation and function? 4. What is the underlying cause of the divergent gene expression patterns observed between B2 and B19 macrophages?
The data collection, analysis and interpretation of results described herein provide plausible answers to these questions based on bioinformatics and functional genomics approaches. Although these answers are more realistically new hypotheses for further investigation, they do represent significant advances in the understanding of B2 and B19 monocyte differentiation into macrophages and the resulting divergent patterns of B2 and B19 macrophage activation and function.
Ultimately, the findings and interpretations we report must be functionally and experimentally validated. Even so, the use of computational methods to answer these questions represents a valuable first step in deciphering the cellular phenotypes underlying MHC haplotype variation in macrophage cells.
Our results demonstrate that there are large numbers of genes differentially expressed in the two haplotypes, both during differentiation of peripheral monocytes into mature macrophages, as well as after stimulation of differentiated macrophages with interferon.
The answer to the question of why the IFNγ stimulated B2 haplotype cells produce more nitric oxide than the IFNγ stimulated B19 cells lies in the timing of macrophage differentiation and the phenotypic variation that is set up early prior to IFNγ stimulation, such as divergent expression of genes involved in differentiation and immune competence. At day t-6, after plating of monocytes without interferon stimulation, several genes relating to inflammation, interferon responses and differentiation are upregulated in the B2 haplotype. This is to be expected as adherence of the monocytes is actually an activation signal, but it is notable that this signal is not resulting in the same gene expression pattern in the B19 haplotype. This pattern can be observed for genes such as IL1β, PTSG2, IL6, which are mainly associated with the inflammatory M1 phenotype. On the other hand, Adenosine receptor A2B is also showing increased gene expression at this time in B2 cells, and this receptor plays an important role in differentiation, as well as in the inflammatory response. Expressions of these genes are consequently initially increased at the time of adherence, and then again after stimulation with interferon, in the B2 haplotype. Some, but not all of these genes are expressed after stimulation with interferon in the B19 haplotype, but not to the same extent as the B2 macrophages, which appears to relate to the initial lack of expression at t-6 days, the beginning of differentiation.
Another interesting observation was the differential expression of genes at day t-6 versus day t-3 in the two haplotypes, as a large number of genes is highly expressed in both haplotypes at day t-6, but then is completely shut down in B2 haplotypes while showing delayed expression until day t-3 in the B19 birds. This seems to indicate a lack of appropriate regulation in the B19 birds, consequently leading to less coherent initiation of gene expression when stimulated. Some of the genes showing this pattern are macrophage differentiation associated GATA2 and FADD, as well as macrophage podosome markers VCL and GSN. Taken together, these results appear to suggest that the regulation of B2 differentiation from monocyte to macrophage is very tightly regulated with many genes increasing in expression, quickly followed by shutting down this increased gene expression. In contrast, the regulation of B19 gene expression is not well regulated, appearing to “linger” with either delayed or extended gene expression. Consequently, we observed differences in expression of genes after stimulation with interferon. Specifically, genes that were strongly expressed at day t-6 and not expressed (or only weakly expressed) at day t-3 in the B2 birds, were robustly increased at 2 and 4 hours of stimulation. In contrast, the same genes showed weak and delayed expression in B19 birds after stimulation, emphasizing the importance of the regulation of gene expression during differentiation. This relates very well to the differences we previously reported in morphology of B2 and B19 macrophages during differentiation and after stimulation.
Our results provide insight into the complexity associated with macrophage differentiation, activation and function. Thousands of genes are up-regulated and then down-regulated in a 24-hour period following IFNγ stimulation. The coordinated activity of multiple regulatory and gene expression control mechanisms is required to effectively achieve the dramatic changes in internal cellular programing that occur. Although the B2 and B19 birds’ haplotypes differ within the MHC locus, the functional consequences of this genetic difference extend well beyond the genes encoded within the B-locus, including macrophage differentiation, M1 and M2 macrophage markers, lysosomal factors involved in phagocytosis, podosome development, invadosome capabilities, chemotaxis potential and matrix degradation ability. Additionally, during the process of differentiation and activation, thousands of genes associated with basic cellular biology undergo rapid changes in expression in coordination with the expression of factors associated with cell renewal and proliferation such as cell cycle regulators, mitotic spindle components, factors involved in chromatin remodeling, molecules required for chromosome segregation and nuclear division.
Taken together, our data and interpretations provide a framework of possible mechanisms of B2 and B19 macrophage biology in differentiation and activation. As such, our findings offer a number of hypotheses about macrophage cell biology that can be used for subsequent studies aimed at validating our findings. Although we performed RT-PCR on a set of differentially expressed genes between the B2 and B19 macrophages, it is not feasible, or possible, to systematically verify, via RT-PCR, each and every transcript, observed at each time point, in the experiment. Even so, our PCR validation provides independent evidence that the pattern of gene expression we observed in the RNA sequence data, was consistent and reproducible which is also in line with our previous research detailing differences in macrophage activation and function [13]
Conclusions
We have tried to elucidate possible mechanisms involved in enhanced disease resistance and macrophage functions displayed by B2 haplotype chickens compared to B19. This study highlights the complex gene expression patterns involved in macrophage differentiation and activation.
One of the main conclusions from the large number of differences seen in the gene expression of the two haplotypes is the fact that there are not just a few genes or genetic markers that can be readily identified as being the ultimate cause of enhanced macrophage function in B2 chickens. Rather, it appears that events during differentiation of monocytes into macrophages have a significant impact on the subsequent ability for stimulation of immune genes after IFNγ treatment. The differences in gene expression correlate with the previously observed differences in morphology of the two haplotypes, with B2 macrophages having a more typical macrophage appearance [13].
Considering the global temporal dysregulation of many genes in B19 haplotypes compared to the more resistant B2 chicks, it seems likely that a variety of genomic regulatory mechanisms (such as transcription factors, miRNAs, snoRNAs, ubiquitin mediated proteasomal degradation, and epigenetic regulation) might play a major role in this process which will be further detailed in a future publication. It will be of great interest to further elucidate these mechanisms and their connection to enhanced immunity. Ultimately, our detailed model of macrophage differentiation, activation and function following IFNγ stimulation provides a high resolution molecular map of cellular biology which can be leveraged by other investigators to further explore the role of these genes in immunology.
Supporting information
S1 Table. Significant differences in gene expression with P-Values.
Pairwise gene expression differences between samples (B2 and B19 haplotype) and timepoints (-6 days, -3 days, 0 days, 1 hr, 2 hrs, 4 hrs, 8 hrs, 16 hrs, 24 hrs) are provided with p-values in excel format. The analysis described within this manuscript focused on differences between [1] matched time points between B2 and B19 haplotype chickens (such as B2 1 hr versus B19 1 hr) as well as [2] progressive timepoints within the same haplotype group (such as B2 1hr versus B2 2hr and B19 4 hr versus B19 8 hr). Subsequently the data contained in this supplemental file also focuses predominantly on those comparisons. This file contains a total of 163,043 rows including the header line containing field names (geneName—identifier for each gene either as gene symbol or ensemble geneId; locus—chromosome number and the start-end base pair location of the gene; sample1 and sample2—the “paired” samples compared for significant gene expression, note ‘AB’ corresponds to B2 haplotype and ‘EC’ corresponds to B19 haplotype; testStatus—indication that the analysis method performed by cuffdiff program within the cufflinks package was ‘OK’; fpkm1 and fpkm2—the fpkm values for sample 1 and sample 2 respectively; log2fpkm—the log of the ratio of fpkm1 and fpkm2; testStat—the test statistic generated during the statistical analysis; pValue—the p-value corresponding to the difference in expression between sample1 and sample2; qValue—a multiple testing corrected p-value; signif—‘yes’ indicates that the pairwise difference in expression is in fact statistically significant).
https://doi.org/10.1371/journal.pone.0179391.s001
(XLSX)
S2 Table. Significant gene ontology biological process enrichment analysis results.
Enriched gene ontology biological process terms identified within up and down regulated genes within the B2 haplotype chicken samples across progressive time points (-6 days, -3 days, 0 days, 1 hr, 2 hrs, 4 hrs, 8 hrs, 16 hrs, 24 hrs) are provided in excel file format. Since the B2 chickens exhibited the most robust macrophage phenotype these samples were used for the analysis as a means of characterizing the biological processes that were associated with the altered gene expression across the experimental time points. Note ‘AB’ indicates B2 haplotype. The file contains a total of 362 rows including the header line containing field names (Sample Comparison—indicates the specific pair of time points for which gene expression changes were identified; Gene Set—indicates the specific set of differentially expressed genes, ‘Down’ or ‘Up’; Category—indicates the specific subset of terms that were used for the analysis; Term—provides the gene ontology identifier for the identified gene ontology term/annotation; Description—the specific enriched gene ontology biological process term; Gene Count—the number of genes within the differentially expressed genes that are mapped to the particular enriched gene ontology term; %—the corresponding percent associated with the specific number of genes; P-Value—the p-value associated with the gene ontology term enrichment).
https://doi.org/10.1371/journal.pone.0179391.s002
(XLSX)
S3 Table. Significant KEGG pathway enrichment analysis results.
Since the B2 chickens exhibited the most robust macrophage phenotype these samples were used for the analysis as a means of characterizing the KEGG pathways that were associated with the altered gene expression across the experimental time points. Note ‘AB’ indicates B2 haplotype. The file contains a total of 110 rows including the header line containing field names (Sample Comparison—indicates the specific pair of time points for which gene expression changes were identified; Gene Set—indicates the specific set of differentially expressed genes, ‘Down’ or ‘Up’; Category—indicates the specific subset of terms that were used for the analysis; Term—provides the KEGG pathway identifier for the identified pathway term/annotation; Description—the specific enriched KEGG pathway term; Gene Count—the number of genes within the differentially expressed genes that are mapped to the particular enriched KEGG pathway term; %—the corresponding percent associated with the specific number of genes; P-Value—the p-value associated with the KEGG pathway term enrichment).
https://doi.org/10.1371/journal.pone.0179391.s003
(XLSX)
S4 Table. Significant reactome pathway enrichment analysis results.
Since the B2 chickens exhibited the most robust macrophage phenotype these samples were used for the analysis as a means of characterizing the Reactome pathways that were associated with the altered gene expression across the experimental time points. Note ‘AB’ indicates B2 haplotype. The file contains a total of 110 rows including the header line containing field names (Sample Comparison—indicates the specific pair of time points for which gene expression changes were identified; Gene Set—indicates the specific set of differentially expressed genes, ‘Down’ or ‘Up’; Category—indicates the specific subset of terms that were used for the analysis; Term—provides the Reactome pathway identifier for the identified pathway term/annotation; Description—the specific enriched Reactome pathway term; Gene Count—the number of genes within the differentially expressed genes that are mapped to the particular enriched Reactome pathway term; %—the corresponding percent associated with the specific number of genes; P-Value—the p-value associated with the Reactome pathway term enrichment).
https://doi.org/10.1371/journal.pone.0179391.s004
(XLSX)
Acknowledgments
We would like to thank the vivarium staff of Western University of Health Sciences and the Sequencing facilities staff at the University of Delaware for their work. We would also like to thank Rich Upshaw and Richard Applebee at Western University, and Christopher Sullivan, at Oregon State University, who helped acquire, set up and maintain the computational hardware required to complete this project. Additionally, the authors wish to thank Paul Gettler for his expertise and time in helping with figures for this paper.
Author Contributions
- Conceptualization: KI YD EC.
- Data curation: KI YD LG JC ED RB CK CB TC.
- Formal analysis: KI YD.
- Funding acquisition: EC YD KI.
- Investigation: KI YD LG JC ED RB CK CB TC.
- Methodology: YD KI.
- Project administration: EC YD KI.
- Resources: RK CK.
- Supervision: YD KI.
- Validation: YD KI.
- Visualization: YD KI.
- Writing – original draft: YD KI.
- Writing – review & editing: YD KI.
References
- 1. Malek M.; Lamont S.J. Association of INOS, TRAIL, TGF-beta2, TGF-beta3, and IgL genes with response to Salmonella enteritidis in poultry. Genetics, selection, evolution: GSE. 2003;35 Suppl 1:S99–111.
- 2. Shi K.Q.; Cai X.H.; Xiao D.D.; Wu S.J.; Peng M.M.; Lin X.F.; Liu W.Y.; Fan Y.C.; Chen Y.P.; Zheng M.H. Tumour necrosis factor-alpha-857T allele reduces the risk of hepatitis B virus infection in an Asian population. Journal of viral hepatitis. 2012 Feb;19(2):e66–72. pmid:22239528
- 3. Ferro P.J.; Swaggerty C.L.; Kaiser P.; Pevzner I.Y.; Kogut M.H. Heterophils isolated from chickens resistant to extra-intestinal Salmonella enteritidis infection express higher levels of pro-inflammatory cytokine mRNA following infection than heterophils from susceptible chickens. Epidemiol Infect. 2004 Dec;132(6):1029–1037. pmid:15635959
- 4. Swaggerty C.L.; Pevzner I.Y.; Kaiser P.; Kogut M.H. Profiling pro-inflammatory cytokine and chemokine mRNA expression levels as a novel method for selection of increased innate immune responsiveness. Vet Immunol Immunopathol. 2008 Nov 15;126(1–2):35–42. pmid:18656269
- 5. Yoo B.H.; Sheldon B.L. Association of the major histocompatibility complex with avian leukosis virus infection in chickens. Br Poult Sci. 1992 Jul;33(3):613–620. pmid:1322760
- 6. Heinzelmann E.W.; Clark K.K.; Collins W.M.; Briles W.E. Host age and major histocompatibility genotype influence on Rous sarcoma regression in chickens. Poult Sci. 1981 Oct;60(10):2171–2175. pmid:6276876
- 7. Schat K.A.; Taylor R.L. Jr.; Briles W.E. Resistance to Marek's disease in chickens with recombinant haplotypes to the major histocompatibility (B) complex. Poult Sci. 1994 Apr;73(4):502–508. pmid:8202429
- 8. Kaufman J.; Wallny H.J. Chicken MHC molecules, disease resistance and the evolutionary origin of birds. Curr Top Microbiol Immunol. 1996;212:129–141. pmid:8934816
- 9. Hunt H.D.; Jadhao S.; Swayne D.E. Major histocompatibility complex and background genes in chickens influence susceptibility to high pathogenicity avian influenza virus. Avian Dis. 2010 Mar;54(1 Suppl):572–575. pmid:20521696
- 10. Banat G.R.; Tkalcic S.; Dzielawa J.A.; Jackwood M.W.; Saggese M.D.; Yates L.; Kopulos R.; Briles W.E.; Collisson E.W. Association of the chicken MHC B haplotypes with resistance to avian coronavirus. Dev Comp Immunol. 2013 Apr;39(4):430–437. pmid:23178407
- 11. Kim D.K.; Lillehoj H.S.; Hong Y.H.; Park D.W.; Lamont S.J.; Han J.Y.; Lillehoj E.P. Immune-related gene expression in two B-complex disparate genetically inbred Fayoumi chicken lines following Eimeria maxima infection. Poult Sci. 2008 Mar;87(3):433–443. pmid:18281568
- 12. Lamont S.J. The chicken major histocompatibility complex in disease resistance and poultry breeding. J Dairy Sci. 1989 May;72(5):1328–1333. pmid:2568373
- 13. Dawes M.E.; Griggs L.M.; Collisson E.W.; Briles W.E.; Drechsler Y. Dramatic differences in the response of macrophages from B2 and B19 MHC-defined haplotypes to interferon gamma and polyinosinic:polycytidylic acid stimulation. Poult Sci. 2014 Apr;93(4):830–838. pmid:24706959
- 14. Collisson E, Griggs L, Drechsler Y. Macrophages from disease resistant B2 haplotype chickens activate T lymphocytes more effectively than macrophages from disease susceptible B19 birds. Dev Comp Immunol. 2017 Feb;67:249–256. pmid:27746172
- 15. Qureshi T.; Templeton J.W.; Adams L.G. Intracellular survival of Brucella abortus, Mycobacterium bovis BCG, Salmonella dublin, and Salmonella typhimurium in macrophages from cattle genetically resistant to Brucella abortus. Vet Immunol Immunopathol. 1996 Mar;50(1–2):55–65. pmid:9157686
- 16. Bellamy R. Susceptibility to mycobacterial infections: the importance of host genetics. Genes and immunity. 2003 Jan;4(1):4–11. pmid:12595896
- 17. Van Erp K, Dach K., Koch I., Heesemann J., Hoffmann R. (2006) Role of strain differences on host resistance and the transcriptional response of macrophages to infection with Yersinia enterocolitica. Physiol. Genomics 25, 75–84. pmid:16352694
- 18. Castillo-Velazquez U.; Gomez-Flores R.; Tamez-Guerra R.; Tamez-Guerra P.; Rodriguez-Padilla C. Differential responses of macrophages from bovines naturally resistant or susceptible to Mycobacterium bovis after classical and alternative activation. Vet Immunol Immunopathol. 2013 Jul 15;154(1–2):8–16. pmid:23707003
- 19. Raza S, Barnett MW, Barnett-Itzhaki Z, Amit I, Hume DA, Freeman TC. Analysis of the transcriptional networks underpinning the activation of murine macrophages by inflammatory mediators. J Leukoc Biol. 2014 Aug;96(2):167–83. pmid:24721704
- 20. Turchetti AP, da Costa LF, Romão Ede L, Fujiwara RT, da Paixão TA, Santos RL. Transcription of innate immunity genes and cytokine secretion by canine macrophages resistant or susceptible to intracellular survival of Leishmania infantum. Vet Immunol Immunopathol. 2015 Jan 15;163(1–2):67–76. pmid:25466388
- 21. Medzhitov R.; Janeway C.A. Jr. Innate immunity: the virtues of a nonclonal system of recognition. Cell. 1997a Oct 31;91(3):295–298.
- 22. Medzhitov R.; Janeway C.A. Jr. Innate immunity: impact on the adaptive immune response. Curr Opin Immunol. 1997b Feb;9(1):4–9.
- 23. Medzhitov R.; Janeway C. Jr. Innate immune recognition: mechanisms and pathways. Immunol Rev. 2000 Feb;173:89–97. pmid:10719670
- 24. Sweet M.J.; Hume D.A. CSF-1 as a regulator of macrophage activation and immune responses. Archivum immunologiae et therapiae experimentalis. 2003;51(3):169–177. pmid:12894871
- 25. Schroder K.; Spille M.; Pilz A.; Lattin J.; Bode K.A.; Irvine K.M.; Burrows A.D.; Ravasi T.; Weighardt H.; Stacey K.J.; Decker T.; Hume D.A.; Dalpke A.H.; Sweet M.J. Differential effects of CpG DNA on IFN-beta induction and STAT1 activation in murine macrophages versus dendritic cells: alternatively activated STAT1 negatively regulates TLR signaling in macrophages. J Immunol. 2007b Sep 15;179(6):3495–3503.
- 26. Biswas S. K., Chittezhath M, Shalova I. N., Lim J. Y. (2012) Macrophage polarization and plasticity in health and disease. Immunol. Res. 53, 11–24 pmid:22418728
- 27. Hong C.C.; Sevoian M. Interferon production and host resistance to type II avian (Marek's) leukosis virus (JM strain). Appl Microbiol. 1971 Nov;22(5):818–820. pmid:4332041
- 28. Lowenthal J.W.; York J.J.; O'Neil T.E.; Rhodes S.; Prowse S.J.; Strom D.G.; Digby M.R. In vivo effects of chicken interferon-gamma during infection with Eimeria. J Interferon Cytokine Res. 1997 Sep;17(9):551–558. pmid:9335433
- 29. Bharti D, Kumar A, Mahla RS, Kumar S, Ingle H, Shankar H, Joshi B, Raut AA, Kumar H. The role of TLR9 polymorphism in susceptibility to pulmonary tuberculosis. Immunogenetics. 2014 Dec;66(12):675–81 pmid:25248338
- 30. Nguyen CT, Kim EH, Luong TT, Pyo S, Rhee DK. ATF3 confers resistance to pneumococcal infection through positive regulation of cytokine production. J Infect Dis. 2014 Dec 1;210(11):1745–54 pmid:24951825
- 31. Kaspers B.; Lillehoj H.S.; Jenkins M.C.; Pharr G.T. Chicken interferon-mediated induction of major histocompatibility complex class II antigens on peripheral blood monocytes. Vet Immunol Immunopathol. 1994 Dec;44(1):71–84 pmid:7536986
- 32. Hu X.; Chen J.; Wang L.; Ivashkiv L.B. Crosstalk among Jak-STAT, Toll-like receptor, and ITAM-dependent pathways in macrophage activation. J Leukoc Biol. 2007 Aug;82(2):237–243. pmid:17502339
- 33. Hu X.; Chakravarty S.D.; Ivashkiv L.B. Regulation of interferon and Toll-like receptor signaling during macrophage activation by opposing feedforward and feedback inhibition mechanisms. Immunol Rev. 2008 Dec;226:41–56. pmid:19161415
- 34. Darmani H.; Parton J.; Harwood J.L.; Jackson S.K. Interferon-gamma and polyunsaturated fatty acids increase the binding of lipopolysaccharide to macrophages. Int J Exp Pathol. 1994 Oct;75(5):363–368. pmid:7999637
- 35. Mita Y.; Dobashi K.; Shimizu Y.; Nakazawa T.; Mori M. Toll-like receptor 2 and 4 surface expressions on human monocytes are modulated by interferon-gamma and macrophage colony-stimulating factor. Immunology letters. 2001 Sep 3;78(2):97–101. pmid:11672593
- 36. Dalpke A.H.; Eckerle S.; Frey M.; Heeg K. Triggering of Toll-like receptors modulates IFN-gamma signaling: involvement of serine 727 STAT1 phosphorylation and suppressors of cytokine signaling. Eur J Immunol. 2003 Jul;33(7):1776–1787. pmid:12811837
- 37. Schroder K.; Lichtinger M.; Irvine K.M.; Brion K.; Trieu A.; Ross I.L.; Ravasi T.; Stacey K.J.; Rehli M.; Hume D.A.; Sweet M.J. PU.1 and ICSBP control constitutive and IFN-gamma-regulated Tlr9 gene expression in mouse macrophages. J Leukoc Biol. 2007a Jun;81(6):1577–1590.
- 38. Geissmann F, Manz MG, Jung S, Sieweke MH, Merad M, Ley K. Development of monocytes, macrophages, and dendritic cells. Science. 2010 Feb 5;327(5966):656–61. pmid:20133564
- 39. Haskó G, Pacher P, Deitch EA, Vizi ES. Shaping of monocyte and macrophage function by adenosine receptors. Pharmacol Ther. 2007 Feb;113(2):264–75 pmid:17056121
- 40. Haskó G, Csóka B, Németh ZH, Vizi ES, Pacher P. A(2B) adenosine receptors in immunity and inflammation. Trends Immunol. 2009 Jun;30(6):263–70. pmid:19427267
- 41. Drechsler Y.; Bohls R.L.; Smith R.; Silvy N.; Lillehoj H.; Collisson E.W. An avian, oncogenic retrovirus replicates in vivo in more than 50% of CD4+ and CD8+ T lymphocytes from an endangered grouse. Virology. 2009 Apr 10;386(2):380–386. pmid:19237181
- 42. Drechsler Y, Tkalcic S, Saggese MD, Shivaprasad HL, Ajithdoss DK, Collisson EW. A DNA vaccine expressing ENV and GAG offers partial protection against reticuloendotheliosis virus in the prairie chicken (Tympanicus cupido). J Zoo Wildl Med. 2013 Jun;44(2):251–61. pmid:23805542
- 43. Crippen T.L.; Sheffield C.L.; He H.; Lowry V.K.; Kogut M.H. Differential nitric oxide production by chicken immune cells. Dev Comp Immunol. 2003 Jun-Jul;27(6–7):603–610. pmid:12697316
- 44. Singh S.; Toro H.; Tang D.C.; Briles W.E.; Yates L.M.; Kopulos R.T.; Collisson E.W. Non-replicating adenovirus vectors expressing avian influenza virus hemagglutinin and nucleocapsid proteins induce chicken specific effector, memory and effector memory CD8(+) T lymphocytes. Virology. 2010b Sep 15;405(1):62–69.
- 45. He H, Genovese KJ, Nisbet DJ, Kogut MH. Profile of Toll-like receptor expressions and induction of nitric oxide synthesis by Toll-like receptor agonists in chicken monocytes. Mol Immunol. 2006 Mar;43(7):783–9. pmid:16098593
- 46. Chan G, Bivins-Smith ER, Smith MS, Smith PM, Yurochko AD. Transcriptome analysis reveals human cytomegalovirus reprograms monocyte differentiation toward an M1 macrophage. J Immunol. 2008 Jul 1;181(1):698–711. pmid:18566437
- 47. Gautier EL, Ivanov S, Williams JW, Huang SC, Marcelin G, Fairfax K, Wang PL, Francis JS, Leone P, Wilson DB, Artyomov MN, Pearce EJ, Randolph GJ. Gata6 regulates aspartoacylase expression in resident peritoneal macrophages and controls their survival. J Exp Med. 2014 Jul 28;211(8):1525–31. pmid:25024137
- 48. Nisole S, Stoye JP, Saib A. TRIM family proteins: retroviral restriction and antiviral defence. Nat. Rev. Microbiol. 2005 3: 799–808. pmid:16175175
- 49. Si Z, Vandegraaff N, O’Huigin C, Song B, Yuan W, Xu C, Perron M, Li X, Marasco WA, Engelman A, et al. Evolution of a cytoplasmic tripartite motif (TRIM) protein in cows that restricts retroviral infection. Proc. Natl. Acad. Sci. USA 2006 103: 7454–7459. pmid:16648259
- 50. van Hateren A, Carter R, Bailey A, Kontouli N, Williams AP, Kaufman J, Elliot T. A mechanistic basis for the co-evolution of chicken tapasin and major histocompatibility complex class I (MHC I) proteins. J Biol Chem 2013 288: 32797–32808 pmid:24078633