Integrating Analysis of Publicly Available Microarray Data to Study the Immune Response of Cattle to Infection with Mycobacterium bovis

Iranian Journal of Applied Animal Science

دوره 13، شماره 1، خرداد 2023، صفحه 11-19   

نوع مقاله: Research Articles

نویسندگان

Nemat Hedayat-Evrigh^* ؛ R. Khalkhali-Evrigh؛ Y. Ramezani؛ R. Seyed Sharifi؛ M.D. Shakouri

Department of Animal Science, Faculty of Agriculture and Natural Resources, University of Mohaghegh Ardabili, Ardabil, Iran

چکیده

Bovine tuberculosis is one of the serious public health challenges which also causes economic damage in the livestock industry. Understanding the interaction mechanism between the host immune system and the causative pathogen of tuberculosis is one of the essential areas of study for success in designing effective drugs to treat tuberculosis. Here, we used four publicly available microarray data to light up the response of the cattle immune systems to Mycobacterium bovis at the gene level. Our integrating analysis results on microarray data led to identifying 189 (160 up- and 29 down-regulated) differentially expressed genes for infected samples with Mycobacterium bovis against uninfected samples. Gene ontology and pathway analysis indicated that most of the differentially expressed genes are related to the immune system's (especially innate immune system) response to the pathogen. Finally, 122 proteins (108 up-regulated and 14 down-regulated) were included in the constructed protein-protein interaction network among the proteins from differentially expressed genes. We identified 11 genes as hub genes based on three methods using the cyto Hubba plug-in in the Cytoscape. Based on our analysis, most differentially expressed genes are related to the innate immune system. However, considering the impact of time on the microarray data analysis indicated that associated gene expression with the adaptive immune system increased by time.

کلیدواژه‌ها

cattle؛ cytokine؛ immune system؛ macrophage؛ microarray؛ Mycobacterium tuberculosis

اصل مقاله

INTRODUCTION

Bovine tuberculosis is a zoonotic disease, caused by Mycobacterium bovis, is one of the most critical diseases in terms of economic losses and threats to public health (Ejeh et al. 2014). Annually, more than 50 million cattle become infected with Mycobacterium bovis globally, and the economic cost of this disease is about $3 billion (Waters and Palmer, 2015). Although cattle are the primary hosts of Mycobacterium bovis, however, this bacterium is still responsible for about 5% of tuberculosis in humans (Michel et al. 2010). Tuberculosis is a deadly disease that killed about 1.7 million people in 2016 (Alam et al. 2019). Macrophages are essential components of the innate immune system against tuberculosis-causing agents such as Mycobacterium bovis. These intracellular pathogens have acquired some abilities to overcome their host immune systems during evolution. The success of these pathogens depends on factors such as their ability to survive and proliferate inside infected macrophages (Philips and Ernst, 2012). The most important strategy of these bacilli to survive inside macrophages is to manipulate the phagocytosis process and prevent the normal maturation of phagosome into the acidic and hydrolytic components. Through this, pathogens are protected from the direct access of the immune system (Russell, 2001). One of the most effective ways to help researchers design effective drugs to treat tuberculosis is to study gene interaction between pathogens and the host immune system. Several studies (Killick et al. 2011; Magee et al. 2012; Caimi et al. 2013; Killick et al. 2014) have investigated gene expression differences in the macrophages and other immune cells infected with Mycobacterium bovis in cattle. In this study, four publicly available data of the gene expression omnibus (GEO) database have been combined to recognize differentially expressed genes (DEGs) of cattle immune systems in response to Mycobacterium bovis. Also, pathway enrichment analysis, gene ontology (GO), protein-protein interaction (PPI), and detection of hub genes construction were carried out for detected DEGs.

 

MATERIALS AND METHODS

Obtaining and pre-processing microarray data

To identify studies profiling the response of cattle immune system to Mycobacterium bovis infection, we searched the National Center for Biotechnology Information (NCBI) GEO database (http://www.ncbi.nlm.nih.gov/geo/). Mycobacterium bovis was used as a search keyword, and filtration was done based on the organism (Bos taurus), study type (expression profiling by array), entry type (Dataset/Series), and platform (Affymetrix Bovine Genome Array, GPL2112). Based on the mentioned criteria, 4 microarray data sets (GSE33309, GSE33359, GSE39819, GSE59774) were selected, and their CEL files were downloaded for downstream analysis (Table 1). After downloading the mentioned datasets, we found that 21 samples are common between GSE33309 and GSE59774 datasets. A series of duplicate data were removed to avoid bias in the analysis. Data reading and normalization of each experiment were done using a robust multi-array average (RMA) algorithm in the R affy package (Gautier et al. 2004). It should be noted that due to the corruption in 4 files (GSM979839, GSM979840, GSM979841, GSM979842) related to the GSE39819 dataset, we used the deposited gene expression matrix of the GSE39819 dataset (which contained all 10-samples) instead CEL files. In the next step, the quality test was performed using the array quality metrics package (Kauffmann et al. 2009) for all samples based on normalized datasets. Finally, 95 samples were used in the final data set. Affymetrix bovine annotation data (Bovine.db) Bioconductor package was used to annotate (adding official gene symbol and Entrez ID) each normalized dataset.

 

Identification of differentially expressed genes

After data integration of 4 sets, batch effect removal was done using ComBat in sva R package. The batch effect represents the non-biological and systematic technical variations that arise during sample processing and measurement at different times, places and, conditions (Lazar et al. 2013). Principal component analysis (PCA) was performed before and after batch effect removal to visualize clustering patterns of studied sets and samples (Figure 1). LIMMA package (Phipson et al. 2016) was performed to identify differentially expressed genes. Cut-off criteria for final selection of DEGs defined as false discovery rate (FDR) adjusted p-value < 0.05 and |log2 Fold change| ≥ 0.58 (Fold change≥1.5) (Figure 2). Also, based on PCA results, we divided samples into three groups, including G1 (healthy samples), G2 (infected samples 2-6 hours after the infection, as well as infected samples which the time after infection was not specified for them), and G3 (infected samples 24 hours after infection). We compared G2 (43 samples) vs. G1 (38 samples) and G3 (14 samples) vs. G1 to identify DEGs. GO analysis was performed on DEGs to understand the time effect on gene expression profile for the classified samples.

 

Functional analysis

We subjected the list of down- and up-regulated genes for the KEGG pathway and GO analyses using DAVID (Sherman and Lempicki, 2009). The results of GO analysis were categorized into three classes, including Biological Process, Cellular Component, and Molecular Function. Benjamini-Hochberg correction was used for calculated p-values of the multiple testing, and p-adjusted < 0.05 was considered a significant level.

 

Construction PPI network and identification of hub genes

The STRING app from Cytoscape (Shannon et al. 2003) software was used for PPI network construction of DEGs. In the constructed network, nodes edges represent the proteins and the interaction between them, respectively, and the labels of nodes indicate the official symbols of the identified DEGs. PPI networks facilitates the understanding and interpreting of the relationship between proteins in a complex biological process. cytoHubba (Chin et al. 2014) plug-in in the cytoscape was used to identify hub genes, commonly identified genes using maximal Clique Centrality (MCC), maximum neighborhood component (MNC), and degree (based on the number of interactions per gene) methods, were considered as hub genes in this study. Hub genes are defined as genes with a highly connected to other genes in a PPI network. Since each extensive network may have contained several sub-networks, we used the MCODE plug-in in the Cytoscape to identify denser sub-network.

 

Table 1 Summary of datasets presented in this study

 

Figure 1 PCA results for studied samples before (a) and after (b) batch effect removal. Control, infected and infected_T24 are classified as G1, G2, and G3, respectively

 

Figure 2 Volcano plot showing the DEGs between infected and healthy samples. Green (up-regulated) and red (down-regulated) dots indicate DEGs with adjusted P-value < 0.05 and |log2 Fold change| ≥ 0.58 (Fold change≥1.5). The black dots fail to meet our criteria. Genes with the highest fold change and lowest adjusted P-value were marked for up- and down-regulated genes

 

Network clustering was performed based on topology to find densely connected regions. Then, a Cytoscape plug-in named CluGO (Bindea et al. 2009) was carried out for GO analysis and to visualize functionally grouped terms of the obtained sub-network. GO analysis allows us to understand better and express the biological function of genes and it classifies them into three categories, including biological process (BP), molecular function (MF), and cellular component (CC) (Khalkhali-Evrigh et al. 2022).

 

RESULTS AND DISCUSSION

After data collection and running the pre-processing on them, finally, we reached 95 samples which comprised 57 infected samples and 38 healthy individuals as control samples. A total of 189 DEGs were identified after analyzing all mentioned samples (infected vs. control), of which 160 were up- and 29 were down-regulated genes. The fold change of up-regulated genes ranged from 1.50 for the SPX gene to 7.57 for the PTX3 gene, while this range for down-regulated genes was from 1.50 for the CD101 gene to 2.09 for the FOS gene. Also, based on the statistic parameter (adjusted P-value), the most significant up- and down-regulated genes were the PVR gene (adjusted P-value equal to 1.27E-09) and the GADD45A gene (adjusted P-value equal to 2E-06), respectively. Significantly enriched terms from GO analysis on up-regulated genes in BP, CC, and MF categories were 19, 2, and 5 terms, respectively (Table 2). With the mentioned criteria in the material and methods section, we found no significantly enriched terms for down-regulated genes. The number of identified DEGs for G2 vs. G1 and G3 vs. G1 that passed our criteria were 262 (213 up- and 49 down-regulated genes) and 627 (394 up- and 233 down-regulated genes), respectively. Also, the number of shared DEGs between the G2 vs. G1 and G3 vs. G1 was 42 (36 up- and 6 down-regulated genes).

 

Pathway analysis

Obtained results from kyoto encyclopedia of genes and genomes (KEGG) analysis revealed that 34 KEGG pathways significantly enriched for up-regulated genes list (Top ten pathways in Table 3. The KEGG is a collection of manual pathway maps demonstrating our knowledge of the molecular interaction, reaction, and communication networks in different areas of biological science like metabolism, genetic information processing, cellular processes, and organism systems. In the present study, the most enriched pathways for mentioned genes were associated with the immune system (such as Cytokine-cytokine receptor interaction, TNF signaling pathway, Toll-like receptor signaling pathway, etc.) and disease (such as Rheumatoid arthritis, Salmonella infection, Tuberculosis, etc.). For enriched pathways, the p-value adjusted based on Benjamin-Hochberg less than 0.05 was considered the significant cut-off value. No significantly enriched pathway was found for down-regulated genes.

 

PPI network

We used all DEGs (189) for PPI network construction, and 184 genes from this list were included in the network. After removing singleton nodes (nodes without edges) from the constructed network, finally, 122 proteins (belonging to the 108 up- and 14 down-regulated genes) with 804 edges (interaction) were included in the network (Figure 3a). The number of interactions for ten genes with the most interaction is represented in Table 4. Extracted sub-network from the leading network that was the densest subset consisting of 21 proteins and 204 edges (Figure 3b). The GO analysis of 21 genes present in the sub-networks showed that more than77 percent of enriched terms were classified in four gropes including “cellular response to lipopolysaccharide” (40.78%), “cytokine receptor binding” (20.39%), “granulocyte chemotaxis” (8.74%) and “positive regulation of response to external stimulus” (7.38%; Figure 3c). The defense response of the immune system against different pathogens like bacteria and viruses is a routine process in organisms. However, studying how the immune system responds at the level of genes and proteins can help us to control diseases by developing effective diagnostic and treatment methods. Inflammation can be defined as a multi-component response of the immune system to different external threats like pathogens (Medzhitov and Horng, 2009). Macrophages, a prominent part of the innate immune system, are among the first defensive lines against tuberculosis (Queval et al. 2017). Along with the phagocytic activity of macrophages, secretion of chemotactic and inflammatory cytokines, and the beginning of a process called chemotaxis bring together the components of the innate immune system against pathogens (Sokol and Luster, 2015). Here, the results of GO analysis for DEGs were consistent with the behavior of the immune system against pathogens, so most up-regulated DEGs are enriched in “inflammatory response”, “immune response”, “cellular response to interferon-gamma”, “monocytes chemotaxis” and “defense response to bacterium” in BP category; also “chemokine activity” and “cytokine activity” in MF category. Cellular component analysis revealed that the most productive up-regulated DEGs were localized in the extracellular space and the external side of the plasma membrane of cells.

 

Table 2 Significantly gene ontology (GO) terms for up-regulated genes

BP: biological process; MF: molecular function and CC: cellular component.

 

Table 3 Top 10 enriched KEGG pathways for up-regulated genes

 

In the present study, of the 10 most up-regulated genes, the PTX3 gene was in the first rank (Fold change equal to 7.56). Based on an RNA-Seq study on humans, this gene is mainly expressed in the bone marrow tissue (Fagerberg et al. 2014). The product of this gene increases in response to inflammation and stress, and interacting with pathogens and damaged host cells help the innate immune system remove them (Balhara et al. 2013). Two (CCL20 and CXCL2) and four genes (IL1-a, IL1-b, IL6, and TNF) from the mentioned list belonged to chemokines and cytokines, respectively. CCL20, and CXCL2 were produced by macrophage and mast cells and induced the recruitment and migration of innate immune cells to infected tissues (Sokol and Luster, 2015). Production of IL1-b by macrophages due to its antimicrobial properties is considered an effective strategy to prevent disease development (Lee et al. 2019). 

 

Figure 3 DEGs PPI network construction. (a) Constructed PPI network for DEGs including 108 up- (green) and 14 down-regulated (red) genes; the size of nodes is based on the number of connections; larger nodes have more edges. (b) Extracted sub-network from the main PPI network including 21 nodes. (c) Grouped results of GO analysis for 21 genes included in the sub-network

 

Table 4 Ten top genes with the most interaction in the constructed protein-protein interaction (PPI) network

 

However, there is evidence to show that its long-term production can lead to host tissue damage (Zhang et al. 2014). The PVR gene, along with interaction with other genes like DNAM-1 promotes monocyte transendothelial migration (Reymond et al. 2004; Gerhardt and Ley, 2015), as a prerequisite for innate immune system response to pathogens. The encoded protein by the TFPI2 gene can play an essential role in inhibiting serine proteases. Studies on Mycobacterium tuberculosis have shown that this pathogen employs serine protease as a virulence factor to damage and disable the host immune system (Dave et al. 2002; Naffin-Olivos et al. 2017). Therefore, it is possible that increasing the expression of the TFPI2 gene suppresses this ability in the bacteria. Pathway analysis for up-regulated genes revealed mostly enriched terms related to disease (like Tuberculosis) and immune system response. In the present study, some crucial enriched pathway terms were “TNF signaling pathway”, “Cytokine-Cytokine receptor interaction”, “NF-kappa B signaling pathway” and “Toll-like receptor signaling pathway”. In response to mycobacterial pathogens, macrophages produce toll-like receptors (TLRs) to interact with pathogen-associated molecular patterns (PAMPs) and induce TNF secretion (Killick et al. 2013). One of the essential roles of TNF as a pro-inflammatory cytokine is inducing apoptosis and necroptosis in infected macrophages as a preventive strategy against the development of infection (Xu et al. 2014). The NF-kappa B signaling pathway is considered a key and central to the immune system response to many pathogens like mycobacteria. This pathway is activated through pattern-recognition receptors (PRRs) expressed in response to PAMPs and responsible for stimulating the production of cytokines, chemokines, adhesion molecules, and many other molecules required by the innate and adaptive immune system (Liu et al. 2017). We observed that the number of immune-related genes with a significantly different expression between infected and control samples increased over time so that DEGs in G2 vs. G1 and G3 vs. G1 were 262 and 627, respectively. Comparison of GO and pathway analysis results showed that one of the exciting differences between the two groups was the enriched term “antigen processing and presentation of endogenous peptide antigen via MHC class I” in the biological process category and “Antigen processing and presentation” in the KEGG pathway for G3 vs. G1. During antigen processing and presentation as an immunological process, foreign proteins are digested into peptide fragments and transported to the surface of antigen-presenting cells (APCs) for recognition by T cells. Major histocompatibility complex (MHC) class I molecules form a complex with antigen peptides in the endoplasmic reticulum to export to the plasma membrane (Abele and Tampe, 2009). It seems MHC class I is a critical factor against tuberculosis, so disruption of MHC class I antigen processing made the studied mice more susceptible to Mycobacterium tuberculosis infection (Urdahl et al. 2003). Also, evidence suggests that MHC class I serve as a target for suppression by Mycobacterium tuberculosis (Meng et al. 2017). Since antigen presentation by MHC proteins plays a vital role in adaptive immunity (Wieczorek et al. 2017), it seems, over time, genes associated with the innate immune system are replaced by genes that light up the adaptive immune system for a more efficient defeat of Mycobacterium bovis. While, results of GO analysis revealed that in the early hours of infection (G2), genes associated with GO terms such as “cell chemotaxis”, “cytokine activity”, “chemokine activity” and “CXCR chemokine receptor binding” have significantly changed in expression level. Chemokines are chemotactic cytokines that play various vital roles in immune system function (Griffith et al. 2014), especially in the first hours of infection. All 11 identified hub genes (CCL20, CD40, CSF2, CXCL2, GRO1, ICAM1, IFNG, IL1b, IL6, IL10, and TNF) belonged to up-regulated genes. The product of the CD40 gene is a transmembrane protein that acts as a receptor on APCs. Its interaction with CD40L (ligand on T cells) activates macrophage and other APCs antimicrobial properties (Méndez-Samperio et al. 2003) it has a vital role in immunity against intracellular pathogens (Gurunathan et al. 1998) like Mycobacterium bovis. CSF2, also known as a granulocyte-macrophage colony-stimulating factor (GM-CSF) is expressed by nonhematopoietic and hematopoietic cells and restricts the development of tuberculosis (Robinson, 2017). Following our study, several studies have shown that ICAM1 expression increases after infection with Mycobacterium bovis and tuberculosis (DesJardin et al. 2002; Li et al. 2017).One of the crucial roles reported for ICAM1 expressed on APCs is facilitating leukocyte transmigration across the endothelium (Müller, 2019) as an essential part of the immune system's response to pathogens. IL6 which is produced by T cells and macrophages is probably involved in the secretion of IFNG (Goovaerts et al. 2013). IFNG is a cytokine that mediates the activation of macrophages to kill intracellular pathogens and the production of other cytokines (Flynn et al. 1993).

 

CONCLUSION

In the present study, we combined four publicly available microarray data to integrate separately conducted studies. Analysis of the mentioned data led to the detection of 160 up- and 29 down-regulated genes. These genes were significantly enriched in immune-related GO terms and pathways. We highlighted 11 genes defined as hub genes based on three methods, including MCC, MNC, and the number of interactions (degree). Most identified hub genes as cytokines and chemokines play a role in innate immune system response to pathogens. Also, the investigation of the effect of time on the pattern of gene expression in infected samples showed that associated gene expression with the adaptive immune system increase over time. Some of the obtained results in this study can be the basis for designing experiments to understand better the host (cattle) immune response to Mycobacterium bovis. The genes identified in this study can be used as genetic markers to identify the mentioned disease and design effective drugs. However, more studies are needed to confirm these results.

 

ACKNOWLEDGEMENT

The authors thank all the teams who worked on the experiments and provided results during this study.

مراجع

Abele R. and Tampe R. (2009). Peptide trafficking and translocation across membranes in cellular signaling and self-defense strategies. Curr. Opin. Cell Biol. 21(4), 508-515.

Alam A., Imam N., Ahmed M.M., Tazyeen S., Tamkeen N., Farooqui A., Malik M. and Ishrat R. (2019). Identification and classification of differentially expressed genes and network meta-analysis reveals potential molecular signatures associated with tuberculosis. Front. Genet. 10, 932-941.

Balhara J., Koussih L., Zhang J. and Gounni A.S. (2013). Pentraxin 3: an immuno-regulator in the lungs. Front. Immunol. 4, 127-135.

Bindea G., Mlecnik B., Hackl H., Charoentong P., Tosolini M., Kirilovsky A., Fridman W.H., Pagès F., Trajanoski Z. and Galon J. (2009). ClueGO: A Cytoscape plug-in to decipher functionally grouped gene ontology and pathway annotation networks. Bioinformatics. 25(8), 1091-1093.

Caimi K., Blanco F., Soria M. and Bigi F. (2013). Transcriptional response of bovine monocyte-derived macrophages after the infection with different Argentinean Mycobacterium bovis isolates. Biomed. Res. Int. 2013, 458278.

Chin C.H., Chen S.H., Wu H.H., Ho C.W., Ko M.T. and Lin C.Y. (2014). CytoHubba: Identifying hub objects and sub-networks from complex interactome. BMC Syst. Biol. 4, 1-7.

Dave J.A., Van Pittius N.C.G., Beyers A.D., Ehlers M.R. and Brown G.D. (2002). Mycosin-1, a subtilisin-like serine protease of Mycobacterium tuberculosis, is cell wall-associated and expressed during infection of macrophages. BMC Microbiol. 2(1), 1-8.

DesJardin L.E., Kaufman T.M., Potts B., Kutzbach B., Yi H. and Schlesinger L.S. (2002). Mycobacterium tuberculosis-infected human macrophages exhibit enhanced cellular adhesion with increased expression of LFA-1 and ICAM-1 and reduced expression and/or function of complement receptors, FcγRII and the mannose receptor. Microbiology. 148(10), 3161-3171.

Ejeh E.F., Raji M.A., Bello M., Lawan F.A., Francis M.I., Kudi A.C. and Cadmus S.I.B. (2014). Prevalence and direct economic losses from bovine tuberculosis in Makurdi, Nigeria. Vet. Med. Int. 63(1), 41-47.

Fagerberg L., Hallström B.M., Oksvold P., Kampf C., Djureinovic D., Odeberg J., Habuka M., Tahmasebpoor S., Danielsson A., Edlund K. and Asplund A. (2014). Analysis of the human tissue-specific expression by genome-wide integration of transcriptomics and antibody-based proteomics. Mol. Cell. Proteomics. 13(2), 397-406.

Flynn J.L., Chan J., Triebold K.J., Dalton D.K., Stewart T.A. and Bloom B.R. (1993). An essential role for interferon gamma in resistance to Mycobacterium tuberculosis infection. Exp. Med. 178(6), 2249-2254.

Gautier L., Cope L., Bolstad B.M. and Irizarry R.A. (2004). Affy-analysis of Affymetrix GeneChip data at the probe level. J. Bioinform. 20(3), 307-315.

Gerhardt T. and Ley K. (2015). Monocyte trafficking across the vessel wall. Cardiovasc. Res. 107(3), 321-330.

Goovaerts O., Jennes W., Massinga-Loembé M., Ceulemans A., Worodria W., Mayanja-Kizza H., Colebunders R. and Kestens L. (2013). TB-IRIS study group. LPS-binding protein and IL-6 mark paradoxical tuberculosis immune reconstitution inflammatory syndrome in HIV patients. PloS One. 8(11), e81856.

Griffith J.W., Sokol C.L. and Luster A.D. (2014). Chemokines and chemokine receptors: positioning cells for host defense and immunity. Annu. Rev. Immunol. 32, 659-702.

Gurunathan S., Irvine K.R., Wu C.Y., Cohen J.I., Thomas E., Prussin C., Restifo N.P. and Seder R.A. (1998). CD40 ligand/trimer DNA enhances both humoral and cellular immune responses and induces protective immunity to infectious and tumor challenge. J. Immunol. Res. 161(9), 4563-4571.

Kauffmann A., Gentleman R. and Huber W. (2009). ArrayQualityMetrics-a bioconductor package for quality assessment of microarray data. J. Bioinform. 25(3), 415-416.

Khalkhali-Evrigh R., Hedayat N., Ming L. and Jirimutu J. (2022). Identification of selection signatures in Iranian dromedary and Bactrian camels using whole genome sequencing data. Sci. Rep. 12(1), 1-10.

Killick K.E., Browne J.A., Park S.D., Magee D.A., Martin I., Meade K.G., Gordon S.V., Gormley E., O'Farrelly C., Hokamp K. and MacHugh D.E. (2011). Genome-wide transcriptional profiling of peripheral blood leukocytes from cattle infected with Mycobacterium bovis reveals suppression of host immune genes. BMC Genom. 12(1), 1-18.

Killick K.E., Magee D.A., Park S.D., Taraktsoglou M., Browne J.A., Conlon K.M., Nalpas N.C., Gormley E., Gordon S.V., MacHugh D.E. and Hokamp K. (2014). Key hub and bottleneck genes differentiate the macrophage response to virulent and attenuated Mycobacterium bovis. Front. Immunol. 5, 422-432.

Killick K.E., NíCheallaigh C., O'Farrelly C., Hokamp K., Mac Hugh D.E. and Harris J. (2013). Receptor-mediated recognition of mycobacterial pathogens. Cell. Microbiol. 15(9), 1484-1495.

Lazar C., Meganck S., Taminau J., Steenhoff D., Coletta A., Molter C., Weiss-Solís D.Y., Duque R., Bersini H. and Nowé A. (2013). Batch effect removal methods for microarray gene expression data integration: A survey. Brief. Bioinform. 14(4), 469-490.

Lee M.R., Chang L.Y., Chang C.H., Yan B.S., Wang J.Y. and Lin W.H. (2019). Differed il-1 beta response between active tb and ltbi cases by ex vivo stimulation of human monocyte-derived macrophage with tb-specific antigen. Dis. Markers. 2019, 7869576.

Li P., Wang R., Dong W., Hu L., Zong B., Zhang Y., Wang X., Guo A., Zhang A., Xiang Y. and Chen H. (2017). Comparative proteomics analysis of human macrophages infected with virulent Mycobacterium bovis. Front. Cell. Infect. Microbiol. 7, 65-73.

Liu T., Zhang L., Joo D. and Sun S.C. (2017). NF-κB signaling in inflammation. Signal Transduct. Target. Ther. 2(1), 1-9.

Magee D.A., Taraktsoglou M., Killick K.E., Nalpas N.C., Browne J.A., Park S.D., Conlon K.M., Lynn D.J., Hokamp K., Gordon S.V. and Gormley E. (2012). Global gene expression and systems biology analysis of bovine monocyte-derived macrophages in response to in vitro challenge with Mycobacterium bovis. PloS One. 7(2), e32034.

Medzhitov R. and Horng T. (2009). Transcriptional control of the inflammatory response. Nat. Rev. Immunol. 9(10), 692-703.

Méndez-Samperio P., Ayala H. and Vázquez A. (2003). NF-κB is involved in regulation of CD40 ligand expression on Mycobacterium bovis bacillus Calmette-Guérin-activated human T cells. Clin. Diagn. Lab. Immunol. 10(3), 376-382.

Meng L., Tong J., Wang H., Tao C., Wang Q., Niu C., Zhang X. AND Gao Q. (2017). PPE38 protein of Mycobacterium tuberculosis inhibits macrophage MHC class I expression and dampens CD8+ T cell responses. Front. Cell. Infect. Microbiol. 7, 68-75.

Michel A.L., Müller B. AND Van Helden P.D. (2010). Mycobacterium bovis at the animal–human interface: A problem, or not? Vet. Microbiol. 140(3), 371-381.

Müller N. (2019). The role of intercellular adhesion molecule-1 in the pathogenesis of psychiatric disorders. Front. Pharmacol. 10, 1251-1259.

Naffin-Olivos J.L., Daab A., White A., Goldfarb N.E., Milne A.C., Liu D., Baikovitz J., Dunn B.M., Rengarajan J., Petsko G.A. AND Ringe D. (2017). Structure determination of Mycobacterium tuberculosis serine protease Hip1 (Rv2224c). Biochemistry. 56(17), 2304-2314.

Philips J.A. and Ernst J.D. (2012). Tuberculosis pathogenesis and immunity. Annu. Rev. Pathol: Mech. Dis. 7, 353-384.

Phipson B., Lee S., Majewski I.J., Alexander W.S. and Smyth G.K. (2016). Robust hyperparameter estimation protects against hypervariable genes and improves power to detect differential expression. Ann. Appl. Stat. 10(2), 946-954.

Queval C.J., Brosch R. and Simeone R. (2017). The macrophage: a disputed fortress in the battle against Mycobacterium tuberculosis. Front. Microbiol. 8, 2284.

Reymond N., Imbert A.M., Devilard E., Fabre S., Chabannon C., Xerri L., Farnarier C., Cantoni C., Bottino C., Moretta A. and Dubreuil P. (2004). DNAM-1 and PVR regulate monocyte migration through endothelial junctions. Exp. Med. 199(10), 1331-1341.

Robinson R.T. (2017). T Cell Production of GM-CSF protects the host during experimental tuberculosis. mBio. 8, e02087.

Russell D.G. (2001). Mycobacterium tuberculosis: Here today, and here tomorrow. Nat. Rev. Mol. Cell Biol. 2(8), 569-578.

Shannon P., Markiel A., Ozier O., Baliga N.S., Wang J.T., Ramage D., Amin N., Schwikowski B. and Ideker T. (2003). Cytoscape: A software environment for integrated models of biomolecular interaction networks. Genome Res. 13(11), 2498-2504.

Sherman B.T. and Lempicki R.A. (2009). Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nat. Protoc. 4(1), 44-52.

Sokol C.L. and Luster A.D. (2015). The chemokine system in innate immunity. Cold Spring Harb. Perspect. Biol. 7(5), a016303.

Urdahl K.B., Liggitt D. and Bevan M.J. (2003). CD8+ T cells accumulate in the lungs of Mycobacterium tuberculosis-infected Kb−/− Db−/− mice, but provide minimal protection. J. Immunol. Res. 170(4), 1987-1994.

Waters W.R. and Palmer M.V. (2015). Mycobacterium bovis infection of cattle and white-tailed deer: Translational research of relevance to human tuberculosis. Inst. Lab. Anim. Res. 56(1), 26-43.

Wieczorek M., Abualrous E.T., Sticht J., Álvaro-Benito M., Stolzenberg S., Noé F. and Freund C. (2017). Major histocompatibility complex (MHC) class I and MHC class II proteins: Conformational plasticity in antigen presentation. Front. Immunol. 8, 292-305.

Xu G., Wang J., Gao G.F. and Liu C.H. (2014). Insights into battles between Mycobacterium tuberculosis and macrophages. Protein Cell. 5(10), 728-736.

Zhang G., Zhou B., Li S., Yue J., Yang H., Wen Y., Zhan S., Wang W., Liao M., Zhang M. and Zeng G. (2014). Allele-specific induction of IL-1β expression by C/EBPβ and PU.1 contributes to increased tuberculosis susceptibility. PLoS Pathog. 10, 1004426.