High-throughput sequencing of D. discoideum cells infected with M. marinum and L. pneumophila
In order to characterize the early transcriptional regulation of host genes after infection by M. marinum and L. pneumophila respectively, we performed high-throughput sequencing of poly(A) enriched RNA from infected and non-infected (control) cells. To obtain RNA for transcriptional studies of M. marinum infected D. discoideum, we used a high multiplicity of infection (MOI of 200) in order to acquire a strong and synchronized transcriptional signature of infected cells already 2.5 h post infection (hpi). Furthermore, we aimed for similar proportions of infected cells (around 60%) as for the L. pneumophila infection (see below). Flow cytometry analysis revealed that approximately 65% of the D. discoideum cells carried M. marinum at this time point (Fig. 1a). Cell viability could be a concern at this high MOI. Hence, to assay cell death during infection, we challenged D. discoideum cells with M. marinum as described above but with different MOI of M. marinum. The fraction of dead cells were assayed by propidium iodide staining followed by flow cytometry [21]. The results clearly showed that while the proportion of infected cell increased with higher MOI, the cell viability was not affected to great extent since the fraction of dead cells only increased from ~ 2% for uninfected D. discoideum cells up to ~ 8% at MOI 200 (Additional file 1: Figure S1).
The early host response to L. pneumophila infection has previously been investigated in D. discoideum using microarray transcriptome analyses [22, 23]. However, one limitation of these studies is that the microarrays only covered approximately 5400 [22] or 8600 [23, 24] genes out of more than 12,200 protein coding genes in D. discoideum [25]. Therefore, we used high-throughput RNA-seq to further investigate the transcriptional response to L. pneumophila infection. This also allowed us to do a global comparison of regulated genes triggered by M. marinum and L. pneumophila infections respectively. RNA-seq was performed on RNA collected 1 and 6 h after L. pneumophila infection as well as on RNA prepared from non-infected D. discoideum cells. Notably, the RNA used for our RNA-seq study had previously been isolated by Li and coworkers who performed microarray analysis on the same batch of RNA isolated from non-infected cells and cells collected 6 h post L. pneumophila infection [23]. Hence, this also allowed us to perform an evaluation of the two different methods: microarray versus high-throughput RNA-seq analysis (see below).
Each high-throughput sequencing yielded a mean of 18.6 and 18.8 million reads from D. discoideum non-infected cells or cells infected with M. marinum respectively, that mapped to the genome after quality control and filtering steps. The same analyses for L. pneumophila infected and non-infected cells yielded 11.4 and 11.5 million reads.
Principal component (PC) analyses were performed for each type of infection (biological replicates), including their respective non-infected controls (separate for each infection experiment). The result clearly showed that the infected and non-infected samples separated along principal component 1 (PC1) in response to both M. marinum and L. pneumophila (Fig. 1b, c).
Large transcriptional responses early after infection
Differential expression analysis of infected versus non-infected samples was performed for each infection, M. marinum or L. pneumophila, using DESeq2 [26] and genes with a false discovery rate (FDR) < 0.05 were considered to be differentially regulated. For both M. marinum 2.5 hpi and L. pneumophila one hpi, approximately 400 genes were found to be differentially regulated while more than 1300 genes were differentially expressed 6 h post L. pneumophila infection (Additional files 2 and 3). In cells infected with M. marinum, the great majority of the regulated genes showed increased expression while a more even distribution between up- and down-regulated genes was observed for L. pneumophila infected D. discoideum cells (Fig. 2a–c). Separate reverse transcription quantitative PCR (RT-qPCR) was performed on the two RNA-seq replicates to validate the regulation of 12 genes that were up-, down-, and non-regulated in the RNA-seq analysis of M. marinum infected cells and all tested genes showed comparable levels of regulation with both methods (Fig. 2d). The infection was repeated three times and RT-qPCR confirmed the differential expression induced by M. marinum infection, indicating a robust and repeatable gene expression response. This was also apparent when the new RT-qPCR data was compared to the RNA-seq analyses (Additional file 1: Figure S2a-c).
In summary, high-throughput sequencing of RNA from D. discoideum infected by M. marinum and L. pneumophila shows that many genes are differentially expressed already at early time points after uptake of either bacterium. In particular, a dramatic response is set off 6 h after infection with L. pneumophila at which time more than 1300 genes are differentially expressed.
High throughput RNA-seq and microarray analyses yield similar results
Next we compared the gene regulation 6 h after L. pneumophila infection detected by RNA-seq with the previously reported differential gene expression identified by microarray, using the same batch of RNA (Additional file 3) [23]. The RNA-seq analysis, representing all ~ 12,200 genes in D. discoideum, showed differential regulation of 1300 genes (FDR < 0.05, see above), while ~ 900 of the 8600 genes on the microarray were reported as differentially expressed (p-value < 0.05 and log2(FC) > 1 or < − 1) [23]. In order to compare the result from the two methods, we compared the fold changes for the genes identified as significantly differentially regulated by microarray [23] with the changes for the same genes in the RNA-seq data. More than 60% showed similar regulation with a log2(FC) > 1 or < − 1 also in the RNA-seq analysis (Fig. 2e, marked in blue), while approximately 30% showed similar but weaker regulation, including some that appeared unregulated in the RNA-seq analyses (Fig. 2e, marked in black). Less than 9% showed opposite regulation between the two methods (Fig. 2e, marked in green). When we compared the regulation of differentially expressed genes as defined by RNA-seq (FDR < 0.05) to the regulated genes on the microarray (as defined above), more than 99% (446 out of 450) genes showed similar regulation (Additional file 1: Figure S3, Additional file 3).
Notably, of the 1300 genes identified as differentially expressed by RNA-seq, ~ 600 genes had not previously been reported as associated with transcriptional response upon L. pneumophila infection of D. discoideum. In part, this can be explained by the fact that more than 260 of these genes were not included in the microarray design.
Taken together, the RNA-seq and microarray analyses give highly similar results when the host gene expression response to L. pneumophila is compared, which is in line with previously reported comparisons of microarray and RNA-seq transcriptomics data [27].
D. discoideum response to M. marinum is enriched for genes involved in intracellular trafficking, autophagy and phagosome maturation
In order to interpret the transcriptional response triggered by M. marinum, we performed gene ontology (GO) term enrichment analysis for up- and down-regulated genes, respectively. Full list of enriched GO-terms and associated genes are available in Additional file 4. Additional results and key genes involved in the different processes can be found in Additional file 1: Additional results and Table S1.
GTP-binding proteins and actin
Among the up-regulated genes we detected an enrichment of genes coding for GTP binding proteins (Fig. 3, Additional file 4). The majority of these genes are small GTPases belonging to the Ras superfamily, which are known to be important regulators involved in a wide range of biological processes (reviewed in [28]). In our data, several up-regulated genes belong to the Rab family GTPases, whose members are mainly involved in the regulation of intracellular vesicular transport by e.g. enabling vesicle formation and facilitating vesicle fusion [28]. We also detected increased gene expression of several members of Ras and Rho family GTPases, important regulators of e.g. gene expression [28]. Rho family GTPases are also involved in regulating actin reorganization, which is critical for both phagocytosis [29] and subsequent phagosome maturation [30]. The effect on actin dynamics was underscored by the up-regulation of genes for dynamin GTPases and the increased expression of several actin and actin binding protein genes (Additional file 2).
ESCRT and membranes
GO-term enrichment analysis revealed that genes associated with Endosomal Sorting Complexes Required for Transport (ESCRT) were enriched among the up-regulated genes in response to M. marinum infection. In macrophages, M. tuberculosis interfere with the ESCRT machinery, which in turn prevents normal phagosome maturation [31, 32]. In D. discoideum, three of the main complexes, ESCRT-I, −II and -III, are well conserved [33] and the majority of the ESCRT-I and ESCRT-III associated genes were up-regulated in response to M. marinum infection. The genes for ESCRT-II, which is not essential for the function of the ESCRT machinery [34], were unaffected. In addition, we detected up-regulation of the ESCRT-associated genes involved in e.g. recruitment of ESCRT-I components to cytoplasmic membranes.
Autophagy
The ESCRT machinery is also required for macroautophagy, hereafter referred to as autophagy, however its exact role in this process remains to be determined [35]. Although autophagy was not detected as an enriched GO-term in itself, many genes associated with this process were found in several other enriched GO terms, e.g. membrane, vacuole and protein tag (Additional file 4). The autophagic machinery is involved in several steps of the infectious route of M. marinum in D. discoideum, from MCV rupture to the egress of the bacteria through non-lytic ejection [9, 12, 15, 36]. This also applies to human cells where M. tuberculosis manipulates the autophagic machinery to ensure survival within the host [37, 38].
Most of the regulated genes identified with RNA-seq that are associated with autophagy and their products have previously been individually characterized during M. marinum infection in D. discoideum [9, 12, 15, 36]. However, our data revealed increased expression levels of atg5, atg12 and atg18, which previously have not been associated with M. marinum infection, as well as five ubiquitin genes. The Atg5-Atg12 complex is involved in phagophore membrane elongation [39]. Functional autophagy also relies on receptors which bridge the connection between the phagophore and the cargo marked for degradation [39]. Our data showed that two of the three proposed autophagy receptors in D. discoideum [39] were upregulated upon M. marinum infection.
Genes for transmembrane transporters are downregulated during M. marinum infection
Although differential expression analysis showed that the majority of the affected genes were upregulated in D. discoideum cells infected by M. marinum, a fraction (9%) displayed reduced expression. These genes were mainly enriched for GO-terms involved in transmembrane transport (Fig. 3, Additional file 4) and included genes for ATP binding cassette (ABC) G family transporters and iron transporters (orthologues to natural resistance associated to macrophages 1 (nramp1) and mitoferrin (mcfF)).
Transcriptional response to L. pneumophila infection is established already 1 h post infection
In order to characterize the dynamics of the transcriptional response after L. pneumophila infection, we compared the regulation at 1 and 6 h post infection. Of the 380 differentially regulated genes identified at 1 h post infection, 80% was differentially expressed also at 6 h post infection, indicating that the majority of the regulation induced 1 h post infection is maintained at least until 6 h post infection (Additional file 1: Figure S4a and marked in red in Additional file 1: Figure S4b, Additional file 3). However, a considerably larger response was detected at the later time point (1331 vs 380 regulated genes) (Additional file 1: Figure S4a). Interestingly, more than 95% of the genes affected at 6 h post infection (FDR < 0.05) showed similar regulation at the earlier time point when a less stringent cut off was used (cut off = log2(fold change) +/− 0.5) (Additional file 1: red and black marking in Figure S4b). This indicates that almost the entire response detected at 6 h post infection is induced already after 1 h but becomes more pronounced as infection progresses. Some of the differentially regulated genes are discussed below and an extended description, including gene names, can be found in Additional file 1: Additional results and Table S1.
L. pneumophila infection induces expression of genes related to defense responses in D. discoideum
Similarities in the gene regulation at 1 and 6 h post infection was also observed when GO-term enrichment analyses were performed for the up-regulated genes (Fig. 4, Additional file 5). At both 1 and 6 h post infection, an enrichment of genes involved in ubiquitin-dependent protein catabolic processes were detected, which is in line with previous studies characterizing D. discoideum transcriptional response using microarray [22, 23]. Also in line with previous studies in D. discoideum, we detected an up-regulation of tRNA-synthetases at 6 h post infection [22, 23]. In addition to tRNA-synthetases, a wide range of genes predicted to be involved in several aspects of tRNA metabolism, e.g. tRNA splicing and modification, were also up-regulated mainly 6 h post infection, but also 1 h post infection (Additional file 3, Additional file 5). Furthermore, L. pneumophila infection appears to induce the production of reactive oxygen species (ROS) in D. discoideum. For both time points there was an enrichment for the GO-term L-ascorbic acid binding. In human immune cells, ROS are produced in order to kill off any invading pathogen. This process, known as the oxidative burst, leads to the accumulation of L-ascorbic acid within the cell, which is thought to protect the host from oxidative damage [40]. The ROS production in infected D. discoideum cells is further corroborated by up-regulation of genes for the Toll-Interleukin (TIR) receptor domain-containing protein and NADPH oxidase, previously shown to be required for ROS production, as well as a gene for superoxide dismutase [23, 41]. Altogether, the up-regulation of genes involved in both ROS production and scavenging, indicates that D. discoideum induce ROS production in response to L. pneumophila infection.
Reduced ribosome biogenesis and energy production in L. pneumophila infected cells
Similar to previous reports, a down-regulation of many ribosomal protein genes were detected at 1 and 6 h post infection (Additional file 3) [22, 23]. However, our data also indicate a more global inhibitory effect on the translational machinery in D. discoideum after L. pneumophila infection. Ribosome biogenesis factors such as PeBoW and Noc complex genes are down-regulated already 1 h post infection. Both of these complexes are required for ribosome maturation [42, 43]. Also, L. pneumophila infection appears to affect ribosomal RNA transcription as demonstrated by down-regulation of the RNA polymerase I complex. However, the levels of rRNA could not be determined using the RNA-seq data due to poly(A) selection of the RNA. The inhibition of genes associated with snoRNA binding and function and rRNA primary transcript binding indicates that also rRNA processing and maturation are impaired in infected cells. Taken together, this indicates that L. pneumophila actively starts to inhibit the translational machinery in D. discoideum almost immediately after uptake.
L. pneumophila infection also caused inhibition of genes associated with primary energy metabolism pathways (e.g. GO terms pdh complex and mitochondrion in Fig. 4). However, in contrast to the effect on the translational machinery, down-regulation of energy metabolism was not detected until 6 h after infection (Fig. 4). Inhibition of genes coding for pyruvate dehydrogenase complex proteins, as well as genes for ATP citrate synthase and acetyl-CoA carboxylase A, indicates an impairment in acetyl-CoA metabolism affecting both energy production via the citric acid cycle and synthesis of fatty acids. In addition, the down-regulation of many genes associated with the mitochondrial electron transfer chain further support that energy metabolism is reduced in L. pneumophila infected cells.
Common transcriptional responses to M. marinum and L. pneumophila infection
The characterization of the differentially regulated genes in D. discoideum after infection with M. marinum and L. pneumophila, respectively, indicated that the two pathogens induce very different responses. This is not surprising since L. pneumophila and M. marinum follow different routes within the cell after uptake, a fact that makes it difficult to set a defined time point where both pathogens have reached the same stage of infection.
In order to get a better resolution of the pathogen induced expression, we compared the profiles of differentially expressed genes for each bacterium. Indeed, the majority of the regulation is specific for each pathogen and there was a greater overlap of regulated genes between the two time points of L. pneumophila infected D. discoideum cells than when cells infected by L. pneumophila were compared to M. marinum infected cells (Fig. 5a). However, a substantial overlap of 160 genes, that were differentially expressed in response to both M. marinum and L. pneumophila was also identified (Fig. 5a, section A, B and C). Interestingly, the majority of these genes show similar regulation in response to both pathogens (Fig. 5b). Among these, there are several small GTPases and iron metal transporters, e.g. nramp1 (Additional file 6). Notably, there is also an induction of several RNA interferences (RNAi) machinery genes in response to both pathogens.
The common transcriptional response in D. discoideum also covers genes which previously have been associated with only one of the pathogens. For example, ten of these genes are annotated as either Induced or Repressed after Legionella Infection (ili/rli) [25]. However, our data demonstrate that nine out of ten of these genes are regulated in similar ways upon infection by either pathogen (Additional file 6). Also one of the three vacuolin genes, vacA, was up-regulated in response to both pathogens. The vacuolins in D. discoideum are similar to mammalian late endosome associated flotillins [44] and depletion of vacB has been demonstrated to cause decreased proliferation of M. marinum in D. discoideum cells [12]. Taken together, the comparison of the differentially expressed genes identified after infection by L. pneumophila and M. marinum shows that the two pathogens trigger distinct transcriptional responses. However, there is also a substantial overlap including the RNAi machinery, which might be part of a general, rather than pathogen specific, response to infection.
The response to intracellular infection is distinct from the response triggered by food bacteria
In nature, D. discoideum feeds on bacteria and other microorganisms, which are taken up by phagocytosis and subsequently degraded to supply the amoeba with nutrients [1]. Hence, we considered that the common gene expression response, induced by pathogenic M. marinum and L. pneumophila, may be part of a general process used for uptake of any bacteria. To investigate this, the regulated genes in M. marinum and L. pneumophila infected cells were compared to the transcriptional changes 2 h after the addition of Escherichia coli B/r previously studied by microarray transcriptome profiling (cluster 1, 4, 5 and 7 in [45]) (Additional file 6). E. coli B/r is considered non-pathogenic and is commonly used as food source for D. discoideum in the laboratory [25]. Since the microarray design covered only ~ 70% of the genes in D. discoideum, we limited the comparison to the genes represented on the microarray.
The fractions of regulated genes unique and common to the two L. pneumophila infections and the M. marinum infection were almost identical when the genes not represented on the microarray were excluded from the comparison (Figs. 5 and 6), validating the approach to compare RNA-seq data and the microarray analysis. As for L. pneumophila and M. marinum infection, most of the genes (~ 65%) that were differentially regulated upon challenge with the food bacteria E. coli were unique, i.e. not affected by infection with the pathogenic bacteria. However, a common set of 162 genes were differentially regulated in response to L. pneumophila (one and six hpi) and E. coli (Fig. 6: C, E, G, H, I and J). Also, 42 differentially expressed genes were found to overlap between the M. marinum and E. coli challenge (Fig. 6: G, I, J and K) while only 20 regulated genes were shared between the responses to all three bacteria (Fig. 6: G, I and J).
Although there are overlaps between genes differentially regulated when D. discoideum is challenged with E. coli B/r and either of the two pathogens, closer inspection reveals large discrepancies in the responses to pathogenic and food bacteria. Similar to what we previously described (see above), the great majority of the genes, 82–84%, affected in response to both M. marinum and L. pneumophila show similar regulation also when only genes included on the microarray are included in the comparison (Additional file 1: Figure S5). In contrast, only 24–50% of the genes show similar regulation in response to E. coli and either of the pathogens (Additional file 1: Figure S5, Additional file 6). Of the 20 genes that were differentially expressed in response to all three bacteria, nine had common regulatory response where six were up- and three were down-regulated (Additional file 6). In summary, comparison of transcriptional regulation revealed that the responses triggered by the pathogenic bacteria M. marinum and L. pneumophila are overall distinct from that of the food bacteria E. coli B/r.
D. discoideum transcriptional response to pathogens is evolutionary conserved
Next, we investigated if the transcriptional responses in D. discoideum triggered by infection with M. marinum and L. pneumophila are conserved, i.e. if similar responses can be detected in infected human cells. We first searched for human orthologues to the ~ 12,200 protein coding genes in D. discoideum, resulting in 4649 D. discoideum genes that were orthologous to 6123 genes in the human genome (see Materials and methods for detailed information). Next, we cross-examined these orthologues with the list of all genes differentially regulated in D. discoideum when challenged with L. pneumophila and M. marinum (see above) (Additional file 7). This analysis revealed human orthologues for ~ 30% (125 out of 440) of the genes differentially expressed in D. discoideum when infected with M. marinum (Fig. 7a). For D. discoideum infected with L. pneumophila, human orthologues were found for ~ 40% of the differentially regulated genes at both 1 hpi and 6 hpi, i.e. 154 out of 380 genes and 534 out of 1331 genes, respectively (Fig. 7a). Finally, we investigated if these human orthologues also were differentially expressed in human cells when challenged with pathogens. For this we analyzed available data sets for transcriptional responses in human monocyte derived macrophages (HMDM) infected with M. tuberculosis [19] and L. pneumophila [20]. In total, more than 500 of the human genes orthologous to differentially regulated D. discoideum genes were regulated also in macrophages infected with M. tuberculosis, L. pneumophila or both bacteria (Additional file 7). The majority of these orthologues show similar expression pattern for each pathogen upon infection in both D. discoideum and macrophages (Fig. 7b). This trend is less pronounced for the human orthologues to the common set of 55 D. discoideum genes (Fig. 7a: D, E, F) that are differentially expressed in response to both M. marinum and L. pneumophila in D. discoideum (Fig. 7b). Taken together, human orthologues were identified for approximately 40% of the genes involved in host response to mycobacteria and L. pneumophila infection in D. discoideum and in most cases they showed similar regulation in both hosts.
Conserved genes differentially expressed in D. discoideum and human macrophages upon infection
KEGG pathway analyses based on the differentially regulated human orthologues showed several enriched pathways such as “Endocytosis” in response to mycobacteria and “Pyruvate metabolism” in response to L. pneumophila (Fig. 8, Additional file 7). The majority of these pathways are also related to the enriched GO terms we identified in D. discoideum in response to M. marinum and L. pneumophila infections (Figs. 3 and 4). In addition, more similarities were found when we manually inspected the orthologues regulated in both hosts. Taken together, we identified an up-regulation of small GTPases, e.g. RRAS2, RAB8B, RAB13, and ARHGAP24, which regulates actin rearrangement, in HMDM’s in response to mycobacteria infection similar to the regulation of small GTPases in D. discoideum. Also, many genes involved in autophagy were up-regulated in macrophages e.g. GABARAPL1–3 and WIPI1 as well as several E3 ubiquitin-protein ligases. Finally, an induction was observed for ESCRT-I and ESCRT-III gene VPS37A and CHMP5, as well as other ESCRT associated genes, e.g. PDCD6IP and IST1.
Even though less similarities were found among the down-regulated orthologues in response to mycobacteria in D. discoideum and macrophages, some orthologues such as glutathione S-transferases, involved in detoxification of xenobiotic substances [46], showed similar regulation in both hosts.
In D. discoideum, L. pneumophila infection induce an up-regulation of E3 ubiquitin-protein ligases genes, tRNA synthetase genes and genes involved in ROS production and scavenging (see above). Of these, only E3 ubiquitin-protein ligases, e.g. ARIH1 and TRAF6, are represented among the orthologues with similar regulation in response to L. pneumophila infection in both D. discoideum and macrophages. However, a clear repression of genes associated with ribosome biogenesis was detected in both hosts (Additional file 7). As in D. discoideum, L. pneumophila infection triggers a down-regulation of both cytoplasmic and mitochondrial ribosomal protein genes, e.g. RPS27L and MRPL52, as well as genes associated with ribosome assembly, e.g. RRS1 and NOP10. In addition, DNA-directed RNA-polymerase I components POLR2H and ZNRD1 are down-regulated in both hosts, indicating a decrease in ribosomal RNA transcription in response to L. pneumophila infection. Also the down-regulation of e.g. DDX27 and PeBoW complex genes (PES1 and WDR12) indicates that the maturation of ribosomal RNA is impaired in both hosts. In summary, the regulation of many key genes, e.g. autophagy in response to mycobacteria and ribosome biogenesis in response to L. pneumophila, was recapitulated in human macrophages, supporting the relevance of D. discoideum as an infection model.