Key Hub and Bottleneck Genes Differentiate the Macrophage Response to Virulent and Attenuated Mycobacterium bovis

Mycobacterium bovis is an intracellular pathogen that causes tuberculosis in cattle. Following infection, the pathogen resides and persists inside host macrophages by subverting host immune responses via a diverse range of mechanisms. Here, a high-density bovine microarray platform was used to examine the bovine monocyte-derived macrophage transcriptome response to M. bovis infection relative to infection with the attenuated vaccine strain, M. bovis Bacille Calmette–Guérin. Differentially expressed genes were identified (adjusted P-value ≤0.01) and interaction networks generated across an infection time course of 2, 6, and 24 h. The largest number of biological interactions was observed in the 24-h network, which exhibited scale-free network properties. The 24-h network featured a small number of key hub and bottleneck gene nodes, including IKBKE, MYC, NFKB1, and EGR1 that differentiated the macrophage response to virulent and attenuated M. bovis strains, possibly via the modulation of host cell death mechanisms. These hub and bottleneck genes represent possible targets for immuno-modulation of host macrophages by virulent mycobacterial species that enable their survival within a hostile environment.

Introduction

Mycobacterium bovis, the causative agent of bovine tuberculosis (BTB) is a facultative pathogen that has a genome sequence 99.95% identical to that of M. tuberculosis, the etiological agent responsible for human tuberculosis (TB) (1). The pathogenesis of BTB is broadly similar to that of human TB and many of the characteristics of M. tuberculosis infection are thought to apply to M. bovis infection in cattle (2, 3). Both pathogens enter the host principally via inhalation of infected aerosol droplets, where they are phagocytosed by alveolar macrophages and dendritic cells. Pathogenic mycobacteria have evolved mechanisms for subverting the host immune response, including prevention of phagosome maturation and subsequent lysosomal delivery (4), enabling survival and replication within the host macrophage following phagocytosis (5). Infected alveolar macrophages can then remain within the central core of granulomas during latent infection or disseminate to draining lymph nodes and to other host organs during active infection (6, 7).

Mycobacterium bovis Bacille Calmette–Guérin (BCG) is an attenuated strain of M. bovis that has been used as a live vaccine for TB control for almost a century. It is estimated that over 100 million people are immunized with BCG annually, making it the most widely used vaccine in human populations (8). All BCG strains available today are derived from a virulent M. bovis strain that was isolated at the turn of the twentieth century from the milk of an infected cow with tubercular mastitis (9). This strain, named “lait Nocard,” was transferred to the Institut Pasteur in Lille in 1901, where it was cultured on glycerol-soaked potato slices supplemented with ox bile by Albert Calmette and Camille Guérin as part of their studies of TB pathogenesis. Serial passage of the lait Nocard strain on this culture medium resulted in the isolation of a strain, in 1908, that displayed reduced virulence in guinea pigs and calves (9). Over the next 13 years, Calmette and Guérin performed 230 in vitro serial passages of this partially attenuated M. bovis strain, eventually leading to the first successful vaccination of a human infant with a fully attenuated BCG isolate in July 1921 (10, 11).

All currently existing BCG substrains originate from the initial attenuated strain developed by Calmette and Guérin during this 13-year period. It was during this time that BCG underwent the first of two attenuation phases. Early studies identified a region of difference 1 (RD1) that is deleted from all BCG strains (12, 13) and was likely lost during the initial attenuation phase between 1908 and 1921. RD1 encodes, among other genes, the early secretory antigen target (ESAT-6) protein secretion system (ESX-1), which is present in virulent strains of M. bovis and M. tuberculosis (14, 15). A second attenuation phase occurred after 1924 due to dissemination of the BCG vaccine from the Institut Pasteur to other countries (13). During this time further genomic deletions occurred, including the RD2 deletion (8, 9).

Virulent M. bovis and attenuated BCG differ in their ability to generate disease. Infection with M. bovis can lead to subversion of the host immune response and formation of tuberculous granulomas with characteristic BTB pathology (16). In contrast, the BCG vaccine strain does not generally cause progressive infection in immune-competent hosts, but can induce a protective pro-inflammatory T_H1 response via activation of interferon-gamma (IFN-γ) secreting CD4⁺ T-cells (17); however, the precise mechanisms of how protection is afforded against virulent challenge remain unclear.

Despite decades of research, the host innate immune mechanisms elicited following infection with virulent M. bovis relative to attenuated BCG remain largely unknown. To our knowledge, a systems biology approach to understanding innate immune responses to virulent M. bovis and attenuated BCG has not previously been described. Systems analysis of transcriptomics data generated from host macrophages infected with the two strains enables reconstruction of dynamic molecular interaction networks that capture the differential innate immune responses to virulent M. bovis and attenuated BCG (18). Network features such as hub and bottleneck genes – a consequence of scale-free behavior – can also be defined, cataloged, and further explored to identify novel biological attributes relevant to the pathogenesis of mycobacterial infections.

Previous studies have demonstrated that the bovine monocyte-derived macrophage (MDM) is a useful model cell type for understanding the early host immunogenetic response to M. bovis and M. bovis-derived antigens (19–22). In the present study, the Affymetrix^® GeneChip^® Bovine Genome Array was used to examine the transcriptome of MDM from seven age-matched Holstein-Friesian female cattle infected in vitro with M. bovis and BCG across an infection time course of 2, 6, and 24 h. Differentially expressed (DE) genes were identified using established methods and a list of interactions were generated using the InnateDB resource (23) and visualized as interaction networks with the Cytoscape software package (24).

Materials and Methods

Ethics Statement

All animal procedures were carried out according to the provisions of the Cruelty to Animals Act (Irish Department of Health and Children license number B100/3939) and ethical approval for the study was obtained from the UCD Animal Ethics Committee (protocol number AREC-P-07-25-MacHugh).

Animals

Seven age-matched (3-year-old) unrelated Holstein-Friesian females were used in this study. All animals were maintained under uniform housing conditions and nutritional regimens at the UCD Lyons Research Farm (Newcastle, County Kildare, Ireland). The animals were selected from an experimental herd without a recent history of BTB infection and all animals tested negative for the single intradermal tuberculin test (SICTT). These cattle were also negative for infection with Brucella abortus, M. avium subsp. paratuberculosis, Salmonella Typhimurium, bovine herpesvirus 1 (BHV-1), and bovine viral diarrhea (BVD) virus.

Monocyte Extraction and Culture of MDM

Monocyte extraction and culture of bovine MDM was performed as described by us previously (20). Briefly, for monocyte isolation, 300 ml of whole blood was collected in acid citrate dextrose buffer (Sigma-Aldrich Ireland Ltd., Dublin, Ireland) in sterile bottles. Blood was layered onto Accuspin™ tubes containing Histopaque^® 1077 (Sigma-Aldrich Ireland Ltd., Dublin, Ireland), and following density gradient centrifugation (500 g for 20 min) performed at room temperature, peripheral blood mononuclear cells (PBMC) were collected. Contaminating red blood cells (RBC) were removed following resuspension and subsequent incubation of the PBMC in RBC lysis buffer (10 mM KHCO₃, 150 mM NH₄Cl, 0.1 mM EDTA pH 8.0) for 5 min at room temperature. After incubation, PBMC were washed twice with sterile phosphate-buffered saline (PBS; Invitrogen™, Life Technologies Corporation, Paisley, UK) before resuspending cells in PBS containing 1% bovine serum albumin (BSA; Sigma-Aldrich Ireland Ltd., Dublin, Ireland). Monocytes were then isolated using the MACS^® protocol and MACSH MicroBeads conjugated to mouse anti-human CD14 antibodies (Miltenyi Biotec Ltd., Surrey, UK), which has been shown to be cross-reactive with bovine monocytes (25). The MACS^® protocol was performed according to the manufacturer’s instructions.

The identity and purity of monocytes was confirmed by flow cytometry using an anti-CD14 fluorescein-labeled antibody (data not shown). This method has been previously shown by us to yield a purity of CD14⁺ cells ≥99% (19). Purified monocytes were seeded at 1 × 10⁶ per well on 24-well tissue culture plates in RPMI 1640 medium (Invitrogen™, Life Technologies Corporation, Paisley, UK) containing 15% heat inactivated fetal calf serum (FCS; Sigma-Aldrich Ireland Ltd., Dublin, Ireland), 1% non-essential amino acids (NEAA; Sigma-Aldrich Ireland Ltd., Dublin, Ireland), gentamicin (5 mg/ml; Sigma-Aldrich Ireland Ltd., Dublin, Ireland), and incubated at 37°C, 5% CO₂. Following 24 h incubation (day one), the media was replaced with 1 ml fresh antibiotic-containing media to remove any non-adhered cells. On day three, media was replaced with 1 ml antibiotic-free culture media (RPMI 1640 medium containing 15% heat inactivated FCS and 1% NEAA only). To ensure that the same number of MDM were subjected to mycobacterial challenge, cells were dissociated on day five using 1× non-enzymatic cell dissociation solution (Sigma-Aldrich Ireland Ltd., Dublin, Ireland), counted, and then re-seeded at 2 × 10⁵ cells per well in 24-well tissue culture plates (Sarstedt Ltd., County Wexford, Ireland) using antibiotic-free culture media. By day eight, 80–100% confluent monolayers of MDM were generated that displayed the characteristic macrophage morphology as confirmed by Giemsa staining (data not shown). On day eight, MDM were used for the in vitro challenge experiments with M. bovis and M. bovis BCG.

Culture of M. bovis and BCG and Infection of Bovine MDM

The culturing of M. bovis (strain M2137; spoligotype SB0142) has been described by us previously (20). Culturing of BCG (Pasteur strain) was performed in a Biosafety Containment Level 3 (CL3) laboratory and conformed to the national guidelines on the use of Hazard Group 3 infectious organisms. BCG stocks for the in vitro MDM infection experiments were cultured in Middlebrook 7H9 medium (Difco™, Becton, Dickinson Ltd., Oxford, UK) containing 10% (vol/vol) BBL™ Middlebrook OADC Enrichment (Difco™, Becton, Dickinson Ltd., Oxford, UK), 0.05% Tween 80 (Sigma-Aldrich, Dublin, Ireland), and 0.50% (weight/vol) glycerol (Sigma-Aldrich Ltd., Dublin, Ireland) at 37°C. Bacterial cultures were grown to mid-logarithmic phase as determined by spectrophotometric analysis prior to the challenge experiments using the purified bovine MDM. Colony-forming unit (cfu) counting was performed using Middlebrook 7H11 medium (Difco™, Becton-Dickinson Ltd., Dublin, Ireland) containing 10% (vol/vol) Middlebrook (ADC) enrichment (Difco™, Becton-Dickinson Ltd., Dublin, Ireland) and 0.50% (vol/vol) glycerol (Sigma-Aldrich, Dublin, Ireland).

Monocyte-derived macrophage infections with BCG were performed as previously described for MDM infections with M. bovis (20). In brief, MDM (seeded at 2 × 10⁵ cells per well) were challenged with BCG (4 × 10⁵ cells per well) prepared in antibiotic-free culture media (RPMI 1640 medium containing 15% heat inactivated FCS and 1% NEAA only). BCG cell counts were performed using a Petroff Hausser chamber (Fisher Scientific Ltd., Dublin, Ireland) following sterile filtering to prevent clumping using a 5 μm filter (Millipore Ireland Ltd., County Cork, Ireland). BCG cell numbers were then adjusted and 4 × 10⁵ bacilli were added to the appropriate wells giving a multiplicity of infection (MOI) of 2:1. Subsequent cfu for BCG counting also yielded a mean MOI of 2:1. Once challenged, MDM were incubated at 37°C, 5% CO₂ for 2, 6, and 24 h.

RNA Extraction and Microarray Analysis

RNA extraction and microarray analysis has been described by us previously (20). Briefly, global gene expression was analyzed using the pan-genomic high-density Affymetrix^® GeneChip^® Bovine Genome Array (Affymetrix UK Ltd., High Wycombe, UK¹). This array contains 24,072 probe sets representing over 23,000 gene transcripts and includes approximately 19,000 Uni-Gene clusters. cDNA labeling, hybridization, and scanning for the microarray experiments were performed by Almac Diagnostics Ltd. (Craigavon, Co. Armagh, Northern Ireland) using a one-cycle amplification/labeling protocol.

Statistical Analysis of Microarray Data

Affymetrix^® GeneChip^® Bovine Genome Array data were analyzed using BioConductor (26)² contained within the R statistical package³. All raw data were normalized using the factor analysis for robust microarray summarization (FARMS) package. The FARMS package uses only perfect match (PM) probes and a quantile normalization procedure, providing both P-values and signal intensities (27). Normalized data were then further subjected to filtering for informative probes sets using the I/NI-calls package in R (28). This defines a probe set as being informative when many of its probes reflect the same change in mRNA concentration across arrays. DE genes for a paired comparison between M. bovis versus BCG infections at each time point were extracted using a moderated paired t-test within the Linear Models for Microarray Data (LIMMA) R package (29)⁴. Genes displaying differential expression were annotated using the Affymetrix^® bovine gene annotation. The Benjamini–Hochberg multiple-testing correction method (30) was applied to all DE genes to minimize the false discovery rate (FDR) and adjusted P-values for all DE genes were calculated. Only DE genes with an adjusted P-value of ≤0.01 were used for subsequent network reconstruction and analysis. All data are MIAME compliant and have been submitted to the NCBI gene expression omnibus (GEO) database with experiment series accession number GSE59774.

Generation of InnateDB Interaction Lists

Affymetrix^® probe IDs for the DE gene lists at each time point were mapped to human Ensembl IDs, using the BioMart search tool on the Ensembl genome database⁵. Each of the three DE gene lists were uploaded separately to the InnateDB⁶ curated, interaction database (23), where a list of interactions between the uploaded molecules were generated. At the time of this analysis, the InnateDB resource contained over 195,000 experimentally determined interactions, of which more than 18,000 have been manually curated, allowing for a detailed systems-level analysis of the innate immune response. The interactions, for each time point, were then viewed as a reconstructed network using the Cytoscape software (24)⁷. Small clusters of nodes that did not connect to the main network were removed from the analysis.

24 h Interaction Network Analysis

The 24-h interaction network was analyzed as an undirected, Markov, network. Statistical analysis of the network was performed using the NetworkAnalyzer plugin of Cytoscape (31). Hub nodes were identified by calculating the degree of connectivity (DOC) for each node in the network using the Degree Sorted Circle Layout in Cytoscape. DOC is a measure of the number of interactions a node has within the network. Nodes with the highest DOC are classed as hub nodes. Bottleneck nodes were identified by calculating the betweenness centrality index (BCI) for each node in the network, also using the NetworkAnalyzer plugin. BCI is a measure of how often a node appears on the shortest path between nodes in the network and is a reflection of the amount of control a node exerts over the interactions of other nodes in the network (32). A Spearman rank correlation between the DOC and log₂ fold-change in gene expression (M. bovis-infected MDM relative to BCG-infected MDM) was performed in the SPSS statistics version 20 package (IBM, New York, NY, USA).

Combination of interactions and differentially expressed gene lists

Several Perl scripts were implemented to combine gene lists and interactions in various ways. One was used to calculate the numbers of genes that are DE and linked to a set of hub genes by a certain distance of edges/interactions, another one to extract expression values for interactors with bottleneck genes at all time points. All scripts are available through GitHub at: https://github.com/kkillick/Systems-biology-paper.git.

Gene ontology of 24 h network

Ingenuity^® systems pathway analysis (IPA) was used to identify gene ontology (GO) functional categories over-represented within the DE genes found in the 24-h interaction network. The Ingenuity^® knowledge base contains the largest database of manually curated and experimentally validated physical, transcriptional, and enzymatic molecular interactions. Furthermore, each interaction in the Ingenuity^® knowledge base is supported by previously published information. For IPA GO analysis, the Ingenuity^® knowledge base was used as a reference set and only genes within the 24-h interaction network were uploaded. IPA then performed an over-representation analysis that categorized the DE genes within the uploaded list into functional GO categories. Each GO category in IPA is ranked based on the number of DE genes falling into each functional group. Right-tailed Fisher’s exact tests were used to calculate an overlap P-value for each of the biological functions assigned to the list of DE genes.

Expression analysis of the interactors of the top hub and bottleneck genes

The top hub and bottleneck genes (IKBKE, MYC, NFKB1, and EGR1) were specified in InnateDB to generate lists of all human interactions. DE genes for a paired comparison of M. bovis versus BCG infections at each time point, using an adjusted P-value of ≤0.01, were extracted from the complete interaction lists and visualized using the Circos software package⁸(33). Interaction nodes colored in red indicate that the gene was up-regulated within the dataset, while green indicate that the gene was down-regulated.

Identification of scale-free and small-world properties within the 24-h network

If a network has a node degree distribution that can be fitted with a power law distribution, it can be a sign that the network has a scale-free architecture. Node degree distribution of the 24-h network was calculated using the NetworkAnalyzer plugin of cytoscape and a power law distribution was fitted using Microsoft Excel.

To determine whether a network displays small-world properties it must show $L\widetilde{>}{L}_{\mathrm{random}}$ but $C\gg C_{\mathrm{random}}$. Where L is the characteristic path length, defined as the number of edges in the shortest path between two nodes, averaged over all pairs of nodes, and C, the clustering coefficient, is defined as follows: node n has k_n max number of edges. C_n denotes the fraction of these edges that actually exist, while C is defined at the average of C_n over all n (34).

The R package igraph⁹ was used to generate 1,000 random networks, each with 716 nodes and 1,785 edges. The mean characteristic path length and the mean clustering coefficient were calculated across the 1,000 networks also using the average path length and transitivity functions of the igraph R package. All R code used here is available through GitHub at: https://github.com/kkillick/Systems-biology-paper.git.

Real Time Quantitative Reverse Transcription (qRT)-PCR Validation of Microarray Results

A panel of 12 genes that were DE across all three time points based on the microarray data was selected for microarray validation using real time quantitative reverse transcription PCR (qRT-PCR) analysis. The laboratory and statistical methods used to perform this analysis are described in Section “Real Time Quantitative Reverse Transcription (qRT)-PCR Validation of Microarray Results” in Supplementary Material, and Table S1 in Supplementary Material details the qRT-PCR primers used in this study.

Results

Increasing Differential Gene Expression between M. bovis- Versus BCG-Infected MDM Over Time

Comparison of gene expression profiles between paired M. bovis versus BCG in vitro infections revealed 702 DE genes at 2 h and 1,000 and 2,674 DE genes at 6 and 24 h post-infection, respectively (adjusted P-value ≤0.01). Of the 702 DE genes found at the 2-h time point, 334 genes displayed increased expression in M. bovis infections relative to BCG infections (hereafter referred to as up-regulated), while 368 genes showed decreased expression in M. bovis infections relative to BCG infections (hereafter referred to as down-regulated). At the 6-h time point, 569 genes were up-regulated, with 431 genes displaying down-regulation, while at the 24-h time point 1,152 and 1,522 genes were up- and down-regulated, respectively. Tables S2–S4 in Supplementary Material contain DE gene results for the 2, 6, and 24 h post-infection time points, respectively. Figure 1 depicts the number of DE genes and the relative change in expression across the infection time course.

FIGURE 1

Figure 1. Number of DE genes found for a paired comparison of virulent M. bovis versus attenuated BCG infections of bovine MDM (n = 7) across a time course of 2, 6, and 24 h post-infection. DE gene numbers are reported here after filtering using an adjusted P-value of ≤ 0.01. Red bars indicate the number of genes displaying increased relative expression in M. bovis compared to BCG (up-regulated) and green bars indicate the number of genes displaying decreased relative expression in M. bovis compared to BCG (down-regulated).

InnateDB Interaction Analyses

InnateDB is an open-source database containing experimentally verified gene and protein-interactions involved in the human, mouse, and bovine innate immune response (23). As a larger number of curated interactions are available for human data, the DE genes detected here were mapped to corresponding human Ensembl gene IDs prior to analysis with InnateDB. In total, 362 DE genes at 2 h, 554 DE genes at 6 h, and 1,503 DE genes at 24 h mapped to human Ensembl gene IDs. These gene lists were examined separately using InnateDB to generate a catalog of interactions for each infection time point. At 2 h post-infection, 205 interactions were identified, while 413 and 1,785 interactions were detected for 6 and 24 h, respectively. As the 24-h time point generated the largest interaction network, it was chosen for detailed downstream analyses.

24 h Interaction Network Reconstruction and Analysis

Removal of small groups of nodes not linked to the main cluster left 615 nodes for analysis in the 24-h interaction network. These were connected via 1,670 edges. Two hundred eight of these nodes had self-interaction loops. The network diameter was 10 and the average shortest path length was 4.315. The IKBKE, MYC, NFKB1, HDAC5, and TRAF2 genes were identified as the top hub nodes within the network. A degree of connectivity (DOC) of 57 was found for IKBKE, 57 for MYC, 55 for NFKB1, 52 for HDAC5, and 50 for TRAF2. The DOC for these five nodes was markedly higher than the mean number of connections found for nodes within the network of 3.857, suggesting that they are hub genes within the 24-h network. Also, for the 24-h interaction network, IKBKE, MYC, and EGR1 displayed the highest betweenness centrality index (BCI) of 0.1661, 0.1360, and 0.1297, respectively (Figure S1 in Supplementary Material). This finding suggests that IKBKE, MYC, and EGR1 are also key bottleneck nodes within the 24-h interaction network.

No significant correlation between the DOC and the log₂ fold-change in gene expression for nodes was found, indicating that nodes with the greatest changes in expression between groups may not have the biggest influence on other nodes within the network (Spearman rho value 0.066, P = 0.108) (Figure 2). Figures S2 and S3 in Supplementary Material show the biological interaction networks for the 2 and 6-h infection time points, respectively.

FIGURE 2

Figure 2. Top 30 gene nodes in the 24-h interaction network that displayed the highest DOC (red line). The log₂ fold-change values for these genes are also shown (blue line).

Scale-Free and Small-World Properties of the M. bovis Versus BCG 24 h Interaction Network

A power law distribution was fitted to the node degree distribution of the 24-h network, generating an R² value of 0.867 (Figure S4 in Supplementary Material). This result suggests that the 24-h interaction network has a scale-free architecture, whereby many nodes have a small number of connections yet a small number of hub nodes have many connections. Scale-free networks tend to exhibit the small-world network property of high clustering and short average path lengths, $L\widetilde{>}{L}_{\mathrm{random}}$ but $C\gg C_{\mathrm{random}}$, (where L is the characteristic path length and C is the clustering coefficient) (34). This was tested on the 24-h interaction network by comparing it to 1,000 random networks generated with the igraph R package using the same number of nodes and edges. The random networks yielded a mean C = 0.007, with a maximum C = 0.012 and L = 4.251. Conversely, for the observed 24 h interaction network, the C and L values were 0.112 and 4.315, respectively. With a C value 15.9-fold higher than the mean C value for the 1,000 random networks, and the L value only marginally higher, the 24-h interaction network exhibits small-world properties.

Highlighting the gene neighborhoods at increasing distances from a main hub node (MYC) in the 24-h interaction network demonstrates the scale-free architecture of the network. Figure 3 shows that almost all other nodes in the network (587 out of 615) are connected to MYC through less than five edges. The combined effect of the three main hub gene nodes (IKBKE, MYC, and NFKB1) is presented in Figure 4A. It shows that nearly all DE genes in the network (96%) are linked through a maximum of four edges to one of the three key hub genes. This is significantly higher than expected by chance (P-value < 0.01), based on 100,000 permutation tests. In addition, the intersection between the linked genes from all hubs rises quickly from less than 10% at a distance of two edges to nearly 90% at a distance of four edges (Figure 4B), indicating complex intertwinement and a capacity for rapid dissemination of signals throughout the network.

FIGURE 3

Figure 3. The 24-h interaction network demonstrating scale-free network properties. The main hub gene node, MYC, is shown in black. Gene nodes that are one edge/interaction from MYC are colored purple; gene nodes that are two edges/interactions from MYC are colored red; gene nodes that are three edges/interactions from MYC are colored orange; and gene nodes that are four edges/interactions from MYC are colored yellow. Gene nodes colored gray are greater than four edges/interactions from MYC.

FIGURE 4

Figure 4. Numbers of DE genes linked to at least one of the three key hub genes (IKBKE, MYC, and NFKB1) up to a certain distance in terms of edges/interactions. (A) Union of genes linked to the hubs and median average of 100,000 permutations, including 95% of the distribution. (B) Intersection of genes linked to the hubs and median average of 100,000 permutations, including 95% of the distribution.

Gene Ontology of the M. bovis Versus BCG 24 h Interaction Network

The total number of DE genes that could be mapped to molecules in the Ingenuity^® Knowledge Base was 507, out of the 615 genes in the uploaded dataset. GO analysis facilitates discovery of molecular functions that are enriched within a dataset. IPA GO analysis of molecules within the 24-h interaction network revealed Cell death to be a major over-represented functional category. Two hundred thirty-seven of the molecules within the 24-h network were involved in Cell death and an adjusted P-value of 4.18 × 10⁻²⁵ was obtained for this category. Further inspection showed that the top 15 hub nodes and the top 15 bottleneck nodes were also involved in the Cell death GO category.

Expression Analysis of the Interactors of the Top Hub and Bottleneck Genes

Expression changes of the interactors with the top hub and bottleneck genes IKBKE, NFKB1, MYC, and EGR1 were examined at each of the three time points. Using an adjusted P-value of ≤0.01, it was observed that the majority of the interactors were DE at 24 h for a paired comparison of M. bovis versus BCG infections. This is in contrast to the earlier time points where only small numbers of interactors exhibited significant expression changes (adjusted P-value ≤0.01). For example, the number of DE genes that interact with IKBKE increased from 6 to 12 between time points 2 and 6 h, but expanded to 68 at 24 h. These results suggest that differential regulatory activity at the IKBKE hub gene and its interactors increases over time in MDM infected with M. bovis relative to BCG. A similar outcome was observed for the NFKB1 and MYC hub genes and the EGR1 hub/bottleneck gene, where the equivalent numbers of DE genes for the 2, 6, and 24-h time points were 19, 27, and 42; 15, 22, and 72; and 23, 26 and 54, respectively (Figure 5).

FIGURE 5

Figure 5. All interactions within InnateDB for the IKBKE, MYC, and NFKB1 hub genes plus the top bottleneck gene node (EGR1) that show significantly differential expression in at least one of the three time points using an adjusted P-value of ≤ 0.01. Color intensity indicates the degree of up-regulation (red) or down-regulation (green) in the M. bovis- versus the BCG-infected bovine MDM (n = 7) at each time point.

Real Time Quantitative Reverse Transcription PCR (qRT-PCR) Validation of Microarray Results

The results from real time qRT-PCR validation of the microarray results are provided in Table S5 in Supplementary Material.

Discussion

The ability of pathogenic mycobacteria to cause disease has been largely attributed to their capacity to survive and replicate in macrophages via a diverse range of mechanisms that subvert the host innate immune response (35). It is believed that mycobacterial virulence is, in part, governed by several secreted virulence factors, such as members of the ESAT-6 protein family encoded by the RD1 genome region. Notably, the RD1 region is present in all virulent M. tuberculosis and M. bovis strains and is one of several regions of difference deleted from the genome of all strains of attenuated BCG (36). The presence or absence of these virulence factors are thought to result in differential host innate immune responses to virulent and attenuated mycobacteria that can be identified through the analysis of host gene expression data (37–39). In the present study, we have analyzed the global differential bovine macrophage transcriptional response to virulent M. bovis and attenuated BCG across an in vitro infection time course of 2, 6, and 24 h. Furthermore, we have applied a systems biology approach to reconstruct the complex molecular interactions and networks that illustrate the differential macrophage responses to these mycobacterial strains. Further work will be necessary to provide insights into whether key genes identified in this study represent aspects of attenuation regulated by the RD1 region.

Comparison of the MDM response to infection with M. bovis relative to BCG revealed the largest number of DE genes at 24 h post-infection (2,674 DE genes), which was 3.8- and 2.7-fold higher than the number of DE genes at the observed at 2 h (702 DE genes) and 6 h time points (1,000 DE genes), respectively (Figure 1). Notably, comparison between M. bovis- and BCG-infected MDM at the 2-h time point did not reveal significant differences in the expression profiles of key immune-related genes known to be important during mycobacterial infection. For example, mRNA transcripts were detected for genes encoding macrophage pattern recognition receptors (PRRs) such as TLR2 and TLR4, which are involved in the recognition of evolutionarily conserved pathogen-associated molecular patterns (PAMPs), and inflammatory cytokines including IL1B, IL6, and IL10 following infection; however, these were not DE between treatments (adjusted P-value of 0.01).

Systems analysis of the DE genes identified at the 24-h time point generated a biological interaction network that displayed both scale-free and small-world properties (Figure 3). In a scale-free network, nodes are highly clustered and can be connected to each other by relatively few steps due, in part, to the presence of a small number of highly connected hub nodes (34, 40). Hub nodes, identified as those nodes that have the highest DOC within a network, enable scale-free networks to be highly robust. The majority of nodes in the network are non-hub nodes with only a small number of connections; therefore, scale-free networks can withstand the removal of non-hub nodes and still retain the capacity to transmit information across the whole network. As such, these networks are said to be “error tolerant.” There are many examples of both small-world and scale-free networks, including metabolic networks (41–43), neural networks in the brain (44, 45), and the Internet (46). However, scale-free networks are vulnerable: should key hub nodes be removed, the network will be seriously compromised. Therefore, although scale-free and small-world networks are highly robust against random removal of non-hub nodes, they are vulnerable to targeted attacks on key hub nodes (47).

In the current study, IKBKE, MYC, and NFKB1 were the top three ranking hub genes within the 24-h network (Figure 2). The importance of these genes is highlighted by the effect they have on the network based on the number of DE genes that are close neighbors (up to four edges away, Figure 4). There is not enough information available from the interaction databases to infer directionality for all interactions in the network, but the type of products associated with these hubs (inhibitors and activators, see below) suggest that most of their interactors lie downstream within a signaling cascade. The fact that nearly all of the DE genes in the network are linked by a maximum of four edges to one of the three main hub genes emphasizes the rapid and extensive response that can be initiated by change of expression in a very small number of genes. Furthermore, although the first- and second-degree DE neighbors of the hubs are mostly separate sub-groups, these sets of genes merge quickly into one large group after the next two extensions of the network. In other words, the overall number of DE genes (but not necessarily the type of expression change) could potentially be derived from changes at any one of the hubs, which points toward a single response captured by the network, compared to a group of separate responses. This correlates well with the detection of one common category (Cell death) found to be significantly over-represented in the GO analysis; however, observable differences in macrophage viability or morphology between BCG and M. bovis infections were not examined. The identification of three hubs instead of just one could be indicative of (i) redundancy, which has evolved for more robustness in the response to infections; or (ii) fine-tuning derived from interplay between different regulators, which could provide greater control of the response.

The products of the IKBKE and NFKB1 genes are central to the NF-κB signaling pathway that regulates the expression of many immune-related genes following infection. IKBKE encodes a subunit of the inhibitor of κB (IκB) kinase (IKK) complex, which phosphorylates the NF-κB inhibitors (IKBs), resulting in the release of the NF-κB transcription factor complex, a subunit of which is encoded by the NFKB1 gene (48). Notably, at the 24-h time point, both IKBKE and NFKB1 were up-regulated in the M. bovis-infected MDM relative to BCG-infected MDM, suggesting that there is increased activation of NF-κB signaling following infection with the virulent mycobacterial strain. This is further supported by the increased relative expression of several NF-κB-inducible genes known to play a role in the host response to M. bovis, such as IL1A, IL1B, IL6, and TNF (49–51).

The up-regulation of the IKBKE and NFKB1 hub genes indicates that activation of NF-κB and consequent downstream signaling and interactors of these hub genes (Figure 5) are crucial for the differential success of virulent and attenuated M. bovis strains within the host macrophage. Indeed, several pathogens, including Helicobacter pylori and Chlamydia pneumoniae, have been shown to manipulate the host NF-κB signaling pathway to ensure their survival (52). Evidence of the modulation of monocyte NF-κB signaling following mycobacterial infection has also been demonstrated. Dhiman et al. (53) showed that human monocytes infected with virulent M. tuberculosis exhibited increased NF-κB activity that resulted in up-regulation of NF-κB-inducible anti-apoptotic genes and increased survival relative to monocytes infected with avirulent M. tuberculosis. In the current study, we also observed increased expression of several NF-κB-inducible anti-apoptotic genes including BIRC3, CFLAR, and TRAF2, which were among the downstream interactors of the IKBKE and NFKB1 hub genes (Figure 5). Furthermore, GO analysis of all nodes in the 24-h network revealed enrichment for genes involved in apoptosis, which is increasingly regarded as a host innate immune mechanism that contains and limits mycobacterial growth following infection (54–56). However, the induction of anti-apoptotic genes may represent an immuno-modulation strategy employed by the pathogen that inhibits apoptosis, enabling survival, and replication within the macrophage (55). Therefore, it is possible that virulent M. bovis targets key members of the NF-kB pathway causing the induction of anti-apoptotic mechanisms that result in prolonged viability within the bovine macrophage.

MYC, another highly connected hub gene identified within the 24-h network, displayed the same DOC as IKBKE. MYC encodes a multi-functional transcription factor that serves as a key regulator of cell proliferation and apoptosis and has been shown to suppress or activate the expression of its target genes (57, 58). A recent study has shown that MYC expression was induced in primary human blood-derived macrophages following in vitro infection with pathogenic and non-pathogenic mycobacterial species. Moreover, it was demonstrated that the MYC transcription factor mediated the suppression of intra-macrophage mycobacterial growth via the activation of cytokines including TNF and IL-6 (59). In the current study, MYC was down-regulated in M. bovis-infected MDM relative to BCG-infected MDM within the 24-h network. Furthermore, 51 of 73 (69.8%) of the downstream MYC interactors at the 24-h time point were also down-regulated (Figure 5). Previous work by our group using the same M. bovis-infected MDM described here demonstrated that MYC was down-regulated relative to non-infected MDM from the same animals (20, 22). It is also important to note that MYC was not DE in the BCG-infected MDM analyzed here relative to the same control MDM (unpublished results). Collectively, these results suggest that suppression of MYC following M. bovis infection results in the suppression of host pro-inflammatory responses leading to intracellular microbial survival and growth. Hence, MYC may serve as a key target for mycobacterial modulation of innate immune mechanisms, facilitating persistence within host macrophages.

Bottleneck nodes in the 24-h network were identified as nodes with a high BCI (Figure S1 in Supplementary Material). BCI has been shown to be associated with constrained evolutionary rates for gene or protein nodes in biological networks; studies involving eukaryote protein-interaction networks have demonstrated that bottleneck nodes evolve more slowly than non-bottleneck nodes (60, 61). Analyses of protein-interaction networks have also shown that bottleneck nodes are more likely to be essential for organismal survival than non-bottleneck nodes (62, 63).

In the current study, EGR1, together with IKBKE and MYC, was identified as a top bottleneck node within the 24-h interaction network (Figure S1 in Supplementary Material), with the majority of significantly DE genes interacting with EGR1 showing changes at this time point (55 out of 85, Figure 5). EGR1 encodes a zinc-finger transcription factor that regulates the expression of a large number of genes involved in cellular differentiation and mitogenesis. In addition, EGR-1 has also been shown to regulate the transcription of pro-inflammatory cytokines in murine macrophages in response to stimulation with bacterial antigens, including TNF and IL-6 (64, 65). To our knowledge, only one publication exists associating EGR1 with mycobacteria–macrophage interactions (66), although it has recently been identified as part of a conserved macrophage core transcriptional response module that is DE in response to intracellular bacterial pathogens (67). Here, identification of EGR1 as a key bottleneck gene suggests that this gene is a component of a differential bovine macrophage response module to virulent and attenuated mycobacterial strains.

It has been shown that many of the species and subspecies in the M. tuberculosis complex exhibit specific host association (68). It is believed that this host specificity is an evolutionary consequence of reciprocal adaptive genetic changes that have occurred due to strong selection pressures between the interacting pathogenic mycobacteria and their hosts (69, 70). Several examples of reciprocal traits involved in host-mycobacterial pathogen co-evolution have been described, including the ability of the host immune system to clear infection versus the ability of virulent pathogens to subvert and evade host immune responses (71). It is likely that the small-world properties of the 24-h network reflect host–pathogen co-evolution. Small-world networks are characterized by the presence of hub and bottleneck genes that interact with all other nodes within the network via a small number of steps. This model predicts that the expression and products of these hub and bottleneck genes enable the host to regulate an appropriate and efficient immune response to virulent or attenuated mycobacterial strains. Alternatively, key hub and bottleneck genes may serve as targets for the subversion of host immune responses by virulent pathogens that underlie successful infection, such as the suppression of apoptosis-related mechanisms by host transcription factors as demonstrated by Dyer et al. (72).

At the individual gene level, from the results presented here, it is not possible to determine whether host macrophage gene expression changes are due to regulation of host gene networks by the pathogen, the host, or a combination of both. In this regard, previous studies have demonstrated that pathogen-encoded protein virulence factors and, more recently, pathogen-encoded small regulatory RNA molecules manipulate and modulate host macrophage responses to facilitate intracellular survival and dissemination of pathogenic mycobacteria (73–79). Consequently, to fully understand the complex transcriptional regulatory interplay underlying host–mycobacterium interactions, it will be necessary to perform parallel high-throughput transcriptional profiling of mRNA and microRNA in both host and microbial cells (80, 81) and complement these with experimental detection of host–pathogen molecular interactions (82, 83).

Conflict of Interest Statement

The Guest Associate Editor Kieran Meade declares that, despite having collaborated with the authors, Kate E. Killick, David A. Magee, Stephen D. E. Park, John A. Browne, Nicolas C. Nalpas, Eamonn Gormley, Stephen V. Gordon, David E. MacHugh, Karsten Hokamp the review process was handled objectively. The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

This work was supported by Investigator Grants from Science Foundation Ireland (Nos: SFI/01/F.1/B028 and SFI/08/IN.1/B2038), a Research Stimulus Grant from the Department of Agriculture, Fisheries, and the Marine (No: RSF 06 405) and a European Union Framework 7 Project Grant (No: KBBE-211602-MACROSYS). Kate E. Killick was supported by the Irish Research Council (IRC) funded Bioinformatics and Systems Biology thematic PhD program¹⁰. We would like to thank Mr. Eamon Costello and staff at the Irish Department of Agriculture, Food and the Marine Backweston Laboratory Campus (Celbridge Co., Kildare) for providing the isolate of M. bovis. We also would like to thank the staff at the UCD Lyons Research Farm for assistance with blood collection and we thank Prof. Joseph Keane (TCD) for valuable comments.

Supplementary Material

The Supplementary Material for this article can be found online at http://www.frontiersin.org/Journal/10.3389/fimmu.2014.00422/abstract

Footnotes

1. ^www.affymetrix.com
2. ^www.bioconductor.org
3. ^www.r-project.org
4. ^www.bioconductor.org/packages/release/bioc/html/limma.html
5. ^www.ensembl.org/biomart/martview
6. ^www.innatedb.com
7. ^www.cytoscape.org
8. ^http://circos.ca
9. ^http://igraph.sourceforge.net/doc/R/igraph.pdf
10. ^http://bioinfo-casl.ucd.ie/PhD

References

1. Garnier T, Eiglmeier K, Camus JC, Medina N, Mansoor H, Pryor M, et al. The complete genome sequence of Mycobacterium bovis. Proc Natl Acad Sci U S A (2003) 100:7877–82. doi:10.1073/pnas.1130426100

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

2. Pollock JM, Neill SD. Mycobacterium bovis infection and tuberculosis in cattle. Vet J (2002) 163:115–27. doi:10.1053/tvjl.2001.0655

CrossRef Full Text

3. Waters WR, Palmer MV, Thacker TC, Davis WC, Sreevatsan S, Coussens P, et al. Tuberculosis immunity: opportunities from studies with cattle. Clin Dev Immunol (2011) 2011:768542. doi:10.1155/2011/768542

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

4. Clemens DL, Horwitz MA. Characterization of the Mycobacterium tuberculosis phagosome and evidence that phagosomal maturation is inhibited. J Exp Med (1995) 181:257–70. doi:10.1084/jem.181.1.257

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

5. Hestvik AL, Hmama Z, Av-Gay Y. Mycobacterial manipulation of the host cell. FEMS Microbiol Rev (2005) 29:1041–50. doi:10.1016/j.femsre.2005.04.013

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

6. Saunders BM, Britton WJ. Life and death in the granuloma: immunopathology of tuberculosis. Immunol Cell Biol (2007) 85:103–11. doi:10.1038/sj.icb.7100027

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

7. Lin PL, Flynn JL. Understanding latent tuberculosis: a moving target. J Immunol (2010) 185:15–22. doi:10.4049/jimmunol.0903856

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

8. Liu J, Tran V, Leung AS, Alexander DC, Zhu B. BCG vaccines: their mechanisms of attenuation and impact on safety and protective efficacy. Hum Vaccin (2009) 5:70–8. doi:10.4161/hv.5.2.7210

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

9. Oettinger T, Jorgensen M, Ladefoged A, Haslov K, Andersen P. Development of the Mycobacterium bovis BCG vaccine: review of the historical and biochemical evidence for a genealogical tree. Tuber Lung Dis (1999) 79:243–50. doi:10.1054/tuld.1999.0206

CrossRef Full Text

10. Calmette A. Preventive vaccination against tuberculosis with BCG. Proc R Soc Med (1931) 24:1481–90.

11. Bonah C. The ‘experimental stable’ of the BCG vaccine: safety, efficacy, proof, and standards, 1921-1933. Stud Hist Philos Biol Biomed Sci (2005) 36:696–721. doi:10.1016/j.shpsc.2005.09.003

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

12. Mahairas GG, Sabo PJ, Hickey MJ, Singh DC, Stover CK. Molecular analysis of genetic differences between Mycobacterium bovis BCG and virulent M. bovis. J Bacteriol (1996) 178:1274–82.

Pubmed Abstract | Pubmed Full Text

13. Behr MA, Small PM. A historical and molecular phylogeny of BCG strains. Vaccine (1999) 17:915–22. doi:10.1016/S0264-410X(98)00277-1

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

14. Harboe M, Oettinger T, Wiker HG, Rosenkrands I, Andersen P. Evidence for occurrence of the ESAT-6 protein in Mycobacterium tuberculosis and virulent Mycobacterium bovis and for its absence in Mycobacterium bovis BCG. Infect Immun (1996) 64:16–22.

Pubmed Abstract | Pubmed Full Text

15. Fortune SM, Jaeger A, Sarracino DA, Chase MR, Sassetti CM, Sherman DR, et al. Mutually dependent secretion of proteins required for mycobacterial virulence. Proc Natl Acad Sci U S A (2005) 102:10676–81. doi:10.1073/pnas.0504922102

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

16. McNair J, Welsh MD, Pollock JM. The immunology of bovine tuberculosis and progression toward improved disease control strategies. Vaccine (2007) 25:5504–11. doi:10.1016/j.vaccine.2007.02.037

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

17. Flynn JL, Chan J, Triebold KJ, Dalton DK, Stewart TA, Bloom BR. An essential role for interferon gamma in resistance to Mycobacterium tuberculosis infection. J Exp Med (1993) 178:2249–54. doi:10.1084/jem.178.6.2249

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

18. Ideker T, Krogan NJ. Differential network biology. Mol Syst Biol (2012) 8:565. doi:10.1038/msb.2011.99

CrossRef Full Text

19. Taraktsoglou M, Szalabska U, Magee DA, Browne JA, Sweeney T, Gormley E, et al. Transcriptional profiling of immune genes in bovine monocyte-derived macrophages exposed to bacterial antigens. Vet Immunol Immunopathol (2011) 140:130–9. doi:10.1016/j.vetimm.2010.12.002

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

20. Magee DA, Taraktsoglou M, Killick KE, Nalpas NC, Browne JA, Park SD, et al. Global gene expression and systems biology analysis of bovine monocyte-derived macrophages in response to in vitro challenge with Mycobacterium bovis. PLoS One (2012) 7:e32034. doi:10.1371/journal.pone.0032034

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

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

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

22. Nalpas NC, Park SD, Magee DA, Taraktsoglou M, Browne JA, Conlon KM, et al. Whole-transcriptome, high-throughput RNA sequence analysis of the bovine macrophage response to Mycobacterium bovis infection in vitro. BMC Genomics (2013) 14:230. doi:10.1186/1471-2164-14-230

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

23. Lynn DJ, Winsor GL, Chan C, Richard N, Laird MR, Barsky A, et al. InnateDB: facilitating systems-level analyses of the mammalian innate immune response. Mol Syst Biol (2008) 4:218. doi:10.1038/msb.2008.55

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

24. Shannon P, Markiel A, Ozier O, Baliga NS, Wang JT, Ramage D, et al. Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res (2003) 13:2498–504. doi:10.1101/gr.1239303

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

25. Jacobsen CN, Aasted B, Broe MK, Petersen JL. Reactivities of 20 anti-human monoclonal antibodies with leucocytes from ten different animal species. Vet Immunol Immunopathol (1993) 39:461–6.

Pubmed Abstract | Pubmed Full Text

26. Gentleman RC, Carey VJ, Bates DM, Bolstad B, Dettling M, Dudoit S, et al. Bioconductor: open software development for computational biology and bioinformatics. Genome Biol (2004) 5:R80. doi:10.1186/gb-2004-5-10-r80

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

27. Hochreiter S, Clevert DA, Obermayer K. A new summarization method for affymetrix probe level data. Bioinformatics (2006) 22:943–9. doi:10.1093/bioinformatics/btl033

CrossRef Full Text

28. Talloen W, Clevert DA, Hochreiter S, Amaratunga D, Bijnens L, Kass S, et al. I/NI-calls for the exclusion of non-informative genes: a highly effective filtering tool for microarray data. Bioinformatics (2007) 23:2897–902. doi:10.1093/bioinformatics/btm478

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

29. Smyth GK. Linear models and empirical bayes methods for assessing differential expression in microarray experiments. Stat Appl Genet Mol Biol (2004) 3:Article3. doi:10.2202/1544-6115.1027

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

30. Benjamini Y, Hochberg Y. Controlling the false discovery rate – a practical and powerful approach to multiple testing. J R Stat Soc B (1995) 57:289–300.

31. Assenov Y, Ramirez F, Schelhorn SE, Lengauer T, Albrecht M. Computing topological parameters of biological networks. Bioinformatics (2008) 24:282–4. doi:10.1093/bioinformatics/btm554

CrossRef Full Text

32. Brandes U. A faster algorithm for betweenness centrality. J Math Sociol (2001) 25:163–77. doi:10.1080/0022250X.2001.9990249

CrossRef Full Text

33. Krzywinski M, Schein J, Birol I, Connors J, Gascoyne R, Horsman D, et al. Circos: an information aesthetic for comparative genomics. Genome Res (2009) 19:1639–45. doi:10.1101/gr.092759.109

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

34. Watts DJ, Strogatz SH. Collective dynamics of ‘small-world’ networks. Nature (1998) 393:440–2. doi:10.1038/30918

CrossRef Full Text

35. Ottenhoff TH. The knowns and unknowns of the immunopathogenesis of tuberculosis. Int J Tuberc Lung Dis (2012) 16:1424–32. doi:10.5588/ijtld.12.0479

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

36. Forrellad MA, Klepp LI, Gioffre A, Sabio y Garcia J, Morbidoni HR, de la Paz Santangelo M, et al. Virulence factors of the Mycobacterium tuberculosis complex. Virulence (2013) 4:3–66. doi:10.4161/viru.22329

CrossRef Full Text

37. Wedlock DN, Kawakami RP, Koach J, Buddle BM, Collins DM. Differences of gene expression in bovine alveolar macrophages infected with virulent and attenuated isogenic strains of Mycobacterium bovis. Int Immunopharmacol (2006) 6:957–61. doi:10.1016/j.intimp.2006.01.003

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

38. Portevin D, Gagneux S, Comas I, Young D. Human macrophage responses to clinical isolates from the Mycobacterium tuberculosis complex discriminate between ancient and modern lineages. PLoS Pathog (2011) 7:e1001307. doi:10.1371/journal.ppat.1001307

CrossRef Full Text

39. Koo MS, Subbian S, Kaplan G. Strain specific transcriptional response in Mycobacterium tuberculosis infected macrophages. Cell Commun Signal (2012) 10:2. doi:10.1186/1478-811X-10-2

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

40. Barabasi AL, Albert R. Emergence of scaling in random networks. Science (1999) 286:509–12. doi:10.1126/science.286.5439.509

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

41. Jeong H, Tombor B, Albert R, Oltvai ZN, Barabasi AL. The large-scale organization of metabolic networks. Nature (2000) 407:651–4. doi:10.1038/35036627

CrossRef Full Text

42. Krajnc B, Marhl M. The small world in biophysical systems structural properties of glycolysis and the TCA cycle in Escherichia coli. Cell Mol Biol Lett (2002) 7:129–31.

Pubmed Abstract | Pubmed Full Text

43. Ma H, Zeng AP. Reconstruction of metabolic networks from genome data and analysis of their global structure for various organisms. Bioinformatics (2003) 19:270–7. doi:10.1093/bioinformatics/19.2.270

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

44. Gallos LK, Makse HA, Sigman M. A small world of weak ties provides optimal global integration of self-similar modules in functional brain networks. Proc Natl Acad Sci U S A (2012) 109:2825–30. doi:10.1073/pnas.1106612109

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

45. Gallos LK, Sigman M, Makse HA. The conundrum of functional brain networks: small-world efficiency or fractal modularity. Front Physiol (2012) 3:123. doi:10.3389/fphys.2012.00123

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

46. Faloutsos M, Faloutsos P, Faloutsos C. On power-law relationships of the internet topology. Comp Comm Rev (1999) 29:251–62. doi:10.1145/316194.316229

CrossRef Full Text

47. Strogatz SH. Exploring complex networks. Nature (2001) 410:268–76. doi:10.1038/35065725

CrossRef Full Text

48. Yamamoto Y, Gaynor RB. IkappaB kinases: key regulators of the NF-kappaB pathway. Trends Biochem Sci (2004) 29:72–9. doi:10.1016/j.tibs.2003.12.003

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

49. Flynn JL, Goldstein MM, Triebold KJ, Sypek J, Wolf S, Bloom BR. IL-12 increases resistance of BALB/c mice to Mycobacterium tuberculosis infection. J Immunol (1995) 155:2515–24.

Pubmed Abstract | Pubmed Full Text

50. Bergeron A, Bonay M, Kambouchner M, Lecossier D, Riquet M, Soler P, et al. Cytokine patterns in tuberculous and sarcoid granulomas: correlations with histopathologic features of the granulomatous response. J Immunol (1997) 159:3034–43.

Pubmed Abstract | Pubmed Full Text

51. Ladel CH, Blum C, Dreher A, Reifenberg K, Kopf M, Kaufmann SH. Lethal tuberculosis in interleukin-6-deficient mutant mice. Infect Immun (1997) 65:4843–9.

Pubmed Abstract | Pubmed Full Text

52. Yanai A, Hirata Y, Mitsuno Y, Maeda S, Shibata W, Akanuma M, et al. Helicobacter pylori induces antiapoptosis through buclear factor-kappaB activation. J Infect Dis (2003) 188:1741–51. doi:10.1086/379629

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

53. Dhiman R, Raje M, Majumdar S. Differential expression of NF-kappaB in mycobacteria infected THP-1 affects apoptosis. Biochim Biophys Acta (2007) 1770:649–58. doi:10.1016/j.bbagen.2006.11.016

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

54. Keane J, Balcewicz-Sablinska MK, Remold HG, Chupp GL, Meek BB, Fenton MJ, et al. Infection by Mycobacterium tuberculosis promotes human alveolar macrophage apoptosis. Infect Immun (1997) 65:298–304.

Pubmed Abstract | Pubmed Full Text

55. Behar SM, Martin CJ, Booty MG, Nishimura T, Zhao X, Gan HX, et al. Apoptosis is an innate defense function of macrophages against Mycobacterium tuberculosis. Mucosal Immunol (2011) 4:279–87. doi:10.1038/mi.2011.3

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

56. Early J, Fischer K, Bermudez LE. Mycobacterium avium uses apoptotic macrophages as tools for spreading. Microb Pathog (2011) 50:132–9. doi:10.1016/j.micpath.2010.12.004

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

57. Adhikary S, Eilers M. Transcriptional regulation and transformation by Myc proteins. Nat Rev Mol Cell Biol (2005) 6:635–45. doi:10.1038/nrm1703

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

58. Meyer N, Kim SS, Penn LZ. The Oscar-worthy role of Myc in apoptosis. Semin Cancer Biol (2006) 16:275–87. doi:10.1016/j.semcancer.2006.07.011

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

59. Yim HC, Li JC, Pong JC, Lau AS. A role for c-Myc in regulating anti-mycobacterial responses. Proc Natl Acad Sci U S A (2011) 108:17749–54. doi:10.1073/pnas.1104892108

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

60. Hahn MW, Kern AD. Comparative genomics of centrality and essentiality in three eukaryotic protein-interaction networks. Mol Biol Evol (2005) 22:803–6. doi:10.1093/molbev/msi072

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

61. Chae L, Lee I, Shin J, Rhee SY. Towards understanding how molecular networks evolve in plants. Curr Opin Plant Biol (2012) 15:177–84. doi:10.1016/j.pbi.2012.01.006

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

62. Yu H, Greenbaum D, Xin Lu H, Zhu X, Gerstein M. Genomic analysis of essentiality within protein networks. Trends Genet (2004) 20:227–31. doi:10.1016/j.tig.2004.04.008

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

63. Yu H, Kim PM, Sprecher E, Trifonov V, Gerstein M. The importance of bottlenecks in protein networks: correlation with gene essentiality and expression dynamics. PLoS Comput Biol (2007) 3:e59. doi:10.1371/journal.pcbi.0030059

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

64. Min IM, Pietramaggiori G, Kim FS, Passegue E, Stevenson KE, Wagers AJ. The transcription factor EGR1 controls both the proliferation and localization of hematopoietic stem cells. Cell Stem Cell (2008) 2:380–91. doi:10.1016/j.stem.2008.01.015

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

65. Zhang J, Xie S, Ma W, Teng Y, Tian Y, Huang X, et al. A newly identified microRNA, mmu-miR-7578, functions as a negative regulator on inflammatory cytokines tumor necrosis factor-alpha and interleukin-6 via targeting Egr1 in vivo. J Biol Chem (2013) 288:4310–20. doi:10.1074/jbc.M112.351197

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

66. Barthel R, Tsytsykova AV, Barczak AK, Tsai EY, Dascher CC, Brenner MB, et al. Regulation of tumor necrosis factor alpha gene expression by mycobacteria involves the assembly of a unique enhanceosome dependent on the coactivator proteins CBP/p300. Mol Cell Biol (2003) 23:526–33. doi:10.1128/MCB.23.2.526-533.2003

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

67. McDermott JE, Archuleta M, Thrall BD, Adkins JN, Waters KM. Controlling the response: predictive modeling of a highly central, pathogen-targeted core response module in macrophage activation. PLoS One (2011) 6:e14673. doi:10.1371/journal.pone.0014673

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

68. Smith NH, Gordon SV, de la Rua-Domenech R, Clifton-Hadley RS, Hewinson RG. Bottlenecks and broomsticks: the molecular evolution of Mycobacterium bovis. Nat Rev Microbiol (2006) 4:670–81. doi:10.1038/nrmicro1472

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

69. Woolhouse ME, Webster JP, Domingo E, Charlesworth B, Levin BR. Biological and biomedical implications of the co-evolution of pathogens and their hosts. Nat Genet (2002) 32:569–77. doi:10.1038/ng1202-569

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

70. Gagneux S. Host-pathogen coevolution in human tuberculosis. Philos Trans R Soc Lond B Biol Sci (2012) 367:850–9. doi:10.1098/rstb.2011.0316

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

71. Durmus Tekir S, Cakir T, Ulgen KO. Infection strategies of bacterial and viral pathogens through pathogen-human protein-protein interactions. Front Microbiol (2012) 3:46. doi:10.3389/fmicb.2012.00046

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

72. Dyer MD, Murali TM, Sobral BW. The landscape of human proteins interacting with viruses and other pathogens. PLoS Pathog (2008) 4:e32. doi:10.1371/journal.ppat.0040032

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

73. Pieters J. Mycobacterium tuberculosis and the macrophage: maintaining a balance. Cell Host Microbe (2008) 3:399–407. doi:10.1016/j.chom.2008.05.006

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

74. DiChiara JM, Contreras-Martinez LM, Livny J, Smith D, McDonough KA, Belfort M. Multiple small RNAs identified in Mycobacterium bovis BCG are also expressed in Mycobacterium tuberculosis and Mycobacterium smegmatis. Nucleic Acids Res (2010) 38:4067–78. doi:10.1093/nar/gkq101

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

75. Harding CV, Boom WH. Regulation of antigen presentation by Mycobacterium tuberculosis: a role for Toll-like receptors. Nat Rev Microbiol (2010) 8:296–307. doi:10.1038/nrmicro2321

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

76. Huynh KK, Joshi SA, Brown EJ. A delicate dance: host response to mycobacteria. Curr Opin Immunol (2011) 23:464–72. doi:10.1016/j.coi.2011.06.002

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

77. Zhou Y, Xie J. The roles of pathogen small RNAs. J Cell Physiol (2011) 226:968–73. doi:10.1002/jcp.22483

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

78. Miotto P, Forti F, Ambrosi A, Pellin D, Veiga DF, Balazsi G, et al. Genome-wide discovery of small RNAs in Mycobacterium tuberculosis. PLoS One (2012) 7:e51950. doi:10.1371/journal.pone.0051950

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

79. Obregon-Henao A, Duque-Correa MA, Rojas M, Garcia LF, Brennan PJ, Ortiz BL, et al. Stable extracellular RNA fragments of Mycobacterium tuberculosis induce early apoptosis in human monocytes via a caspase-8 dependent mechanism. PLoS One (2012) 7:e29970. doi:10.1371/journal.pone.0029970

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

80. Westermann AJ, Gorski SA, Vogel J. Dual RNA-seq of pathogen and host. Nat Rev Microbiol (2012) 10:618–30. doi:10.1038/nrmicro2852

CrossRef Full Text

81. Bent ZW, Tran-Gyamfi MB, Langevin SA, Brazel DM, Hamblin RY, Branda SS, et al. Enriching pathogen transcripts from infected samples: a capture-based approach to enhanced host-pathogen RNA sequencing. Anal Biochem (2013) 438:90–6. doi:10.1016/j.ab.2013.03.008

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

82. Driscoll T, Dyer MD, Murali TM, Sobral BW. PIG – the pathogen interaction gateway. Nucleic Acids Res (2009) 37:D647–50. doi:10.1093/nar/gkn799

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

83. Gillespie JJ, Wattam AR, Cammer SA, Gabbard JL, Shukla MP, Dalay O, et al. PATRIC: the comprehensive bacterial bioinformatics resource with a focus on human pathogenic species. Infect Immun (2011) 79:4286–98. doi:10.1128/iai.00207-11

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text