Jil Sander1, Susanne V Schmidt2, Branko Cirovic3, Naomi McGovern4, Olympia Papantonopoulou1, Anna-Lena Hardt1, Anna C Aschenbrenner1, Christoph Kreer5, Thomas Quast6, Alexander M Xu7, Lisa M Schmidleithner1, Heidi Theis1, Lan Do Thi Huong1, Hermi Rizal Bin Sumatoh8, Mario A R Lauterbach2, Jonas Schulte-Schrepping1, Patrick Günther1, Jia Xue1, Kevin Baßler1, Thomas Ulas1, Kathrin Klee1, Natalie Katzmarski3, Stefanie Herresthal1, Wolfgang Krebs1, Bianca Martin2, Eicke Latz9, Kristian Händler1, Michael Kraut1, Waldemar Kolanus6, Marc Beyer10, Christine S Falk11, Bettina Wiegmann12, Sven Burgdorf5, Nicholas A Melosh7, Evan W Newell8, Florent Ginhoux8, Andreas Schlitzer13, Joachim L Schultze14. 1. Genomics and Immunoregulation, LIMES-Institute, University of Bonn, 53115 Bonn, Germany. 2. Institute of Innate Immunity, University Hospital Bonn, University of Bonn, Sigmund-Freud-Strasse 25, 53127 Bonn, Germany. 3. Myeloid Cell Biology, LIMES-Institute, University of Bonn, 53115 Bonn, Germany. 4. Agency for Science, Technology and Research (A(∗)STAR), Singapore Immunology Network (SIgN), 138648 Singapore, Singapore; Department of Pathology and Center for Trophoblast Research, University of Cambridge, CB2 1QP Cambridge, UK. 5. Cellular Immunology, LIMES-Institute, University of Bonn, 53115 Bonn, Germany. 6. Molecular Immunology & Cell Biology, LIMES-Institute, University of Bonn, 53115 Bonn, Germany. 7. Department of Materials Science and Engineering, Stanford University, Stanford, CA 94305, USA. 8. Agency for Science, Technology and Research (A(∗)STAR), Singapore Immunology Network (SIgN), 138648 Singapore, Singapore. 9. Institute of Innate Immunity, University Hospital Bonn, University of Bonn, Sigmund-Freud-Strasse 25, 53127 Bonn, Germany; Department of Medicine, University of Massachusetts Medical School, Worcester, MA 01605, USA; German Center for Neurodegenerative Diseases, 53127 Bonn, Germany. 10. Genomics and Immunoregulation, LIMES-Institute, University of Bonn, 53115 Bonn, Germany; Molecular Immunology, German Center for Neurodegenerative Diseases (DZNE), Sigmund-Freud-Str. 27, 53127 Bonn, Germany. 11. Institute of Transplant Immunology, Integrated Research and Treatment Center Transplantation, Hannover Medical School, 30625 Hannover, Germany. 12. Department of Cardiothoracic, Transplantation and Vascular Surgery, Hannover Medical School, 30625 Hannover, Germany. 13. Myeloid Cell Biology, LIMES-Institute, University of Bonn, 53115 Bonn, Germany; Agency for Science, Technology and Research (A(∗)STAR), Singapore Immunology Network (SIgN), 138648 Singapore, Singapore. Electronic address: andreas.schlitzer@uni-bonn.de. 14. Genomics and Immunoregulation, LIMES-Institute, University of Bonn, 53115 Bonn, Germany; Platform for Single Cell Genomics and Epigenomics (PRECISE) at the German Center for Neurodegenerative Diseases and the University of Bonn, 53127 Bonn, Germany.
Abstract
Human in vitro generated monocyte-derived dendritic cells (moDCs) and macrophages are used clinically, e.g., to induce immunity against cancer. However, their physiological counterparts, ontogeny, transcriptional regulation, and heterogeneity remains largely unknown, hampering their clinical use. High-dimensional techniques were used to elucidate transcriptional, phenotypic, and functional differences between human in vivo and in vitro generated mononuclear phagocytes to facilitate their full potential in the clinic. We demonstrate that monocytes differentiated by macrophage colony-stimulating factor (M-CSF) or granulocyte macrophage colony-stimulating factor (GM-CSF) resembled in vivo inflammatory macrophages, while moDCs resembled in vivo inflammatory DCs. Moreover, differentiated monocytes presented with profound transcriptomic, phenotypic, and functional differences. Monocytes integrated GM-CSF and IL-4 stimulation combinatorically and temporally, resulting in a mode- and time-dependent differentiation relying on NCOR2. Finally, moDCs are phenotypically heterogeneous and therefore necessitate the use of high-dimensional phenotyping to open new possibilities for better clinical tailoring of these cellular therapies.
Human in vitro generated monocyte-derived dendritic cells (moDCs) and macrophages are used clinically, e.g., to induce immunity against cancer. However, their physiological counterparts, ontogeny, transcriptional regulation, and heterogeneity remains largely unknown, hampering their clinical use. High-dimensional techniques were used to elucidate transcriptional, phenotypic, and functional differences between human in vivo and in vitro generated mononuclear phagocytes to facilitate their full potential in the clinic. We demonstrate that monocytes differentiated by macrophage colony-stimulating factor (M-CSF) or granulocyte macrophage colony-stimulating factor (GM-CSF) resembled in vivo inflammatory macrophages, while moDCs resembled in vivo inflammatory DCs. Moreover, differentiated monocytes presented with profound transcriptomic, phenotypic, and functional differences. Monocytes integrated GM-CSF and IL-4 stimulation combinatorically and temporally, resulting in a mode- and time-dependent differentiation relying on NCOR2. Finally, moDCs are phenotypically heterogeneous and therefore necessitate the use of high-dimensional phenotyping to open new possibilities for better clinical tailoring of these cellular therapies.
Recently, transcriptomic, epigenetic, functional, and fate-mapping studies established the identity of three mononuclear cell lineages within the myeloid cell network: macrophages (Macs), dendritic cells (DCs), and monocytes (MOs) (Schlitzer et al., 2015a). During murine embryogenesis, progenitors colonize developing tissues differentiating into tissue-resident Macs, which are long lived and self-maintained in most tissues (Ginhoux and Guilliams, 2016). DCs can be separated into plasmacytoid DCs, conventional DC 1 (cDC1s), and cDC2s (Merad et al., 2013). cDC1s and cDC2s arise from specialized precursors within the bone marrow forming their functional specialization during development (Grajales-Reyes et al., 2015, Schlitzer et al., 2015b, See et al., 2017).The third component of the mononuclear phagocyte (MP) network are MOs. In mice, MOs are divided into Ly6clo and Ly6chi MOs (Varol et al., 2015), while in human blood, three different monocyte subsets are identified by expression of CD14 and CD16, in which CD14+CD16− MOs correspond to murine Ly6chi MOs and CD14−CD16+ MOs to murine Ly6clo MOs (Auffray et al., 2009). Yet little is known whether this classification relates to functional specialization of distinct subsets.Prior work suggests functional overlap between DCs, MOs, and Macs, including basic mechanisms of phagocytosis (Fraser et al., 2009), anti-bacterial activity (Fehlings et al., 2012), antigen uptake, processing, and presentation (Cella et al., 1997), the capacity to activate adaptive immune cells (Toujas et al., 1997), cellular motility (Rossi et al., 2011), and activation programs (Frucht et al., 2000). This overlap has proven difficult to parse, but new knowledge concerning the distinct ontogeny of these cells provided the opportunity to reassess the division of labor between DCs, MOs, and tissue Macs.Currently, assigning functions to subsets of the MP system requires use of simplified in vitro systems to focus on basic molecular aspects. Murinegranulocyte macrophage colony-stimulating factor (GM-CSF)- or macrophage colony-stimulating factor (M-CSF)-driven bone marrow-derived DCs and Mac cultures are frequently used to elucidate and assign molecular mechanisms of functions to subsets of MPs. However, these in vitro cultures create heterogeneous cultures, making attribution of distinct cellular functions difficult (Helft et al., 2015). This conundrum highlights the need for a detailed investigation of cellular identity and the regulation thereof in such in vitro cultures (Xue et al., 2014).Sallusto and Lanzavecchia (1994) described the in vitro generation of human MO-derived DCs (moDCs) by culturing peripheral blood MOs with GM-CSF and IL-4. Here, the term moDCs has been attributed to an activated MO population with DC-like function based on morphological and functional criteria. Similar functionally converging phenotypes are observed in human in vitro systems of MO-derived M-CSF-driven Macs (moMacs) (Akagawa et al., 2006) or GM-CSF (Xue et al., 2014). Systems biology-based definitions of MP function and nomenclature have been established, yielding insights about identity, regulation, and developmental origin of those cells (Xue et al., 2014). However, studies directly addressing their relationships to MPs observed in vivo remain limited (Ohradanova-Repic et al., 2016). Understanding such relationships and linking this knowledge to cellular ontogeny is crucial considering the interest in using in vitro generated MPs for immunotherapy (Brower, 2005, Garg et al., 2017). Therefore, the functional convergence, plasticity, and heterogeneity of MO-derived MPs paired with the clinical interest raises several important questions. What are the in vivo counterparts of in vitro MO derivatives? Do MOs integrate cytokine signaling in a temporal fashion and how is it regulated molecularly? Lastly, how heterogeneous are human in vitro MO cultures?Computational analysis of MP transcriptomes and analysis of cellular phenotype, function, and perturbation experiments elucidated the relationship of human moDCs and moMacs to the in vivo MP system. The differentiation of MO in vitro culture systems is multifaceted, integrating time-dependent signals delivered by GM-CSF and IL-4 and orchestrated by nuclear receptor corepressor 2 (NCOR2). Finally, mass cytometry (MC) revealed cellular heterogeneity of moDCs with several subsets being identified. These results uncover the in vivo counterparts of MO derivatives and identify a novel regulator of MO differentiation and plasticity.
Results
In Vitro Differentiated Human MO-Derived MPs Are Transcriptionally Similar to MO-Derived Inflammatory MPs
Human MOs differentiated with M-CSF are used as models for human Macs (Akagawa et al., 2006), whereas MOs differentiated with GM-CSF and IL-4 are models for human DCs (Sallusto and Lanzavecchia, 1994). For clarity and in light of recent findings concerning DC, MO, and Mac ontogeny (Guilliams and van de Laar, 2015), we term differentiated MOs according to their activation, e.g., MOs activated with M-CSF are named MOs-M-CSF and MOs differentiated for a specified duration (0–72 hr; 0–144 hr) with GM-CSF and IL-4 are MOs-GM-CSFIL-4. To establish transcriptional similarity between ex vivo isolated cells and differentiated MOs, we compared blood CD14+ MOs, CD1c+ DCs, CD141+ DCs (Haniffa et al., 2012), and T, B, and NK cells alongside CD45+lin−HLA-DRhi lung derived cells, to MOs-M-CSF, MOs-GM-CSF, MOs-GM-CSFIL-4(0-72h), and MOs-GM-CSFIL-4(0-144h) (Figure S1A).Principle component analysis (PCA) revealed that T, B, and NK cells formed one of three clusters (Figures 1A and 1B, green), all ex vivo isolated MPs formed a second cluster (red and yellow), and both these clusters were most distinct from a third cluster containing ex vivo polarized MO-derived MPs (blue, purple, cyan, lilac). This three-cluster structure was validated by hierarchical clustering (HC) of the 1,000 most variable genes (Figure 1C; Table S1) and Pearson correlation coefficient matrix (PCCM) analysis (Figure 1D), showing that in vitro differentiated MOs are transcriptionally unique compared to ex vivo human peripheral blood cells.
Figure 1
Relationship of In Vitro Activated Monocyte-Derived Cells
(A and B) PCA (21,250 present probes); displayed principal components (PCs): (A) 1 versus 2 and (B) 1 versus 3.
(C) Heatmap of 1,000 most variable genes in dataset. Log2-expression values, z-transformed, scaled (−2 [blue] to 2 [red]).
(D) Heatmap, Pearson correlation values (PCV) calculated pairwise between all cell types (top 1,000 most variable genes).
(E) PCA (23,952 present probes).
(F) Relative fractions of MO, BDCA1+ DC, infM, and infDC gene signatures in CD14+ MOs, different MO-derived cells, DCs.
(G and H) Heatmaps of genes specifically expressed in (G) infDCs compared to infMs, BDCA1+ DCs, MOs (dataset 1), and MOs-GM-CSFIL-4(0-72h) versus MO-derived cells, CD14+ MOs, and DCs (dataset 2), or in (H) infMs compared to infDCs, BDCA1+ DCs, and MOs (dataset 1), and in MOs-GM-CSF and MOs-M-CSF versus MOs-GM-CSFIL-4(0-72h), CD14+ MOs, and DCs (dataset 2). PCVs between indicated group patterns of dataset 1 versus dataset 2, barplot next to heatmaps (correlation cutoff > 0.4). Genes analyzed in (J) and (K) highlighted in red. Log2-expression values, z-transformed, scaled (−2 [blue] to 2 [red]).
(I) Histograms (flow cytometry analysis), relative expression CD226, MARCO, VSIG4, and CCR7 (representative data, n = 4).
(J) Analysis of cell culture supernatants of MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h) for CCL22 and CCL2 using ELISA (n = 3, 2 technical replicates each, mean + SEM, one-way RM [repeated-measures] ANOVA, Tukey’s method for multiple test correction, with ∗p < 0.05, ∗∗p < 0.01, and ∗∗∗p < 0.001; n.d., not detected).
(K) MMP12 quantification (relative) in CD14+ MOs, MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h) by immunoblot (n = 3, mean + SEM, one-way ANOVA and Tukey’s method for multiple test correction, with ∗p < 0.05).
Please also see Figure S1.
Relationship of In Vitro Activated Monocyte-Derived Cells(A and B) PCA (21,250 present probes); displayed principal components (PCs): (A) 1 versus 2 and (B) 1 versus 3.(C) Heatmap of 1,000 most variable genes in dataset. Log2-expression values, z-transformed, scaled (−2 [blue] to 2 [red]).(D) Heatmap, Pearson correlation values (PCV) calculated pairwise between all cell types (top 1,000 most variable genes).(E) PCA (23,952 present probes).(F) Relative fractions of MO, BDCA1+ DC, infM, and infDC gene signatures in CD14+ MOs, different MO-derived cells, DCs.(G and H) Heatmaps of genes specifically expressed in (G) infDCs compared to infMs, BDCA1+ DCs, MOs (dataset 1), and MOs-GM-CSFIL-4(0-72h) versus MO-derived cells, CD14+ MOs, and DCs (dataset 2), or in (H) infMs compared to infDCs, BDCA1+ DCs, and MOs (dataset 1), and in MOs-GM-CSF and MOs-M-CSF versus MOs-GM-CSFIL-4(0-72h), CD14+ MOs, and DCs (dataset 2). PCVs between indicated group patterns of dataset 1 versus dataset 2, barplot next to heatmaps (correlation cutoff > 0.4). Genes analyzed in (J) and (K) highlighted in red. Log2-expression values, z-transformed, scaled (−2 [blue] to 2 [red]).(I) Histograms (flow cytometry analysis), relative expression CD226, MARCO, VSIG4, and CCR7 (representative data, n = 4).(J) Analysis of cell culture supernatants of MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h) for CCL22 and CCL2 using ELISA (n = 3, 2 technical replicates each, mean + SEM, one-way RM [repeated-measures] ANOVA, Tukey’s method for multiple test correction, with ∗p < 0.05, ∗∗p < 0.01, and ∗∗∗p < 0.001; n.d., not detected).(K) MMP12 quantification (relative) in CD14+ MOs, MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h) by immunoblot (n = 3, mean + SEM, one-way ANOVA and Tukey’s method for multiple test correction, with ∗p < 0.05).Please also see Figure S1.To analyze MP subsets, T, B, and NK cells were excluded from the previous dataset and remaining cell types were analyzed by PCA (Figures 1E, S1B, and S1C). In vitro generated cells clustered separately compared to in vivo MPs, with lung-derived CD45+lin−HLA-DRhi and CD14+ MOs forming a separate cluster away from a DC cluster and a cluster for MOs-M-CSF, MOs-GM-CSF, MOs-GM-CSFIL-4(0-144h), and MOs-GM-CSFIL-4(0-72h). Additionally, in vitro MO-derived MPs shared a common set of differentially expressed genes (DEGs) in comparison to ex vivo CD14+ MOs (Figures S1D and S1E; Table S1), supporting transcriptional difference to homeostatic cells. These findings raised the question of which cells found in vivo represent the MO model systems.We assessed whether the in vitro MO-derived MPs reflected human MO-derived inflammatory DCs (infDCs) and inflammatory Mac (infMs) in vivo utilizing a published dataset (Segura et al., 2013). We generated a signature matrix representing infDCs, infMs, BDCA1+ (CD1c) DCs, CD14+CD16− MOs, and CD14dimCD16+ MOs and assessed the relative signature enrichment in our dataset using linear support vector regression (Figure 1F; Table S1; Newman et al., 2015). MOs-M-CSF showed highest enrichment of infM genes, while infDC genes were highly enriched in MOs-GM-CSFIL-4, suggesting that in vitro polarization conditions reflected in vivo infDCs and infMs as shown for the mouse (Helft et al., 2015). Control gene sets derived from CD14+CD16− MOs were most highly enriched in ex vivo CD14+ MOs and lung-derived CD45+lin−HLA-DRhi cells, while the BDCA1+ DC signature enriched in both ex vivo myeloid DC subsets. Gene set enrichment analyses (GSEA) confirmed transcriptional similarities between infMs and MOs-M-CSF but also MOs-GM-CSF and between infDCs and MOs-GM-CSFIL-4 (Figure S1F). Collectively, we defined four groups in both datasets, describing comparable cell subsets (see STAR Methods). Then, we performed PCCM by comparing expression patterns in both datasets based on the four groups and visualized genes with highest correlation scores between infDCs and MOs-GM-CSFIL-4(0-72h/144h) (Figure 1G; Table S1) and between infMs and both MOs-M-CSF and MOs-GM-CSF (Figure 1H; Table S1), which included several surface markers and secreted molecules (Figures S1G–S1J; Table S1). Many genes associated with activated DCs (CCL22, MMP12, CD226, CCR7) were highly elevated in both infDCs and MOs-GM-CSFIL-4 (Figure 1G), while typical Mac genes (MARCO, CCL2, VSIG4) were most highly expressed in infMs, MOs-M-CSF, and MOs-GM-CSF (Figure 1H). Furthermore, differential regulation of the proteins CD226, MARCO, VSIG4, CCR7, CCL2, CCL22, and MMP12 was validated (Figures 1I–1K). Collectively, polarization of MOs in vivo (infDCs, infMs) (Segura et al., 2013) and in vitro leads to similar transcriptomic identities, including cell surface and effector molecules, allowing us to use these models to understand the role of M-CSF, GM-CSF, and IL-4 for inflammation-associated MO differentiation.
GM-CSF + IL-4 but Not GM-CSF or M-CSF Alone Enforce a Unique Transcriptional Signature in Human CD14+ MOs
Next, we wanted to understand the similarities and differences in MO activation induced by M-CSF, GM-CSF, and IL-4. Previous work suggested important differences between MOs-M-CSF and MOs-GM-CSF (Lacey et al., 2012). However, these studies did not answer the overall relationship between all three activation conditions (Figure 2A). Using the well-established surface markers CD14, CD11b, and CD209, we assessed the differences between MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h) (Figure S2A). This revealed that CD14 marked MOs, MOs-M-CSF, and MOs-GM-CSF, but not MOs-GM-CSFIL-4(0-72h) (Sallusto and Lanzavecchia, 1994). CD209 was exclusively expressed by MOs-GM-CSFIL-4(0-72h). CD11b did not discriminate between the cell populations. However, when reassessing the overall relationship between ex vivo isolated CD14+ MOs, MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h) on transcriptome level (Figure S2B), MOs-M-CSF and MOs-GM-CSF clustered together, while MOs and MOs-GM-CSFIL-4(0-72h) clustered separately as demonstrated by PCA (Figure 2B), PCCM (Figure S2C), and HC (Figure 2C; Table S2). MOs-M-CSF and MOs-GM-CSF were marked by high expression of Mac genes (CD81, VSIG4, SIGLEC1, MARCO, FPR3) while MOs and MOs-GM-CSFIL-4(0-72h) cell populations formed separated gene clusters marked by expression of MO-associated (AHR, SELL, CLEC4D) or DC-associated (CD1C, ZBTB46) genes, respectively. Gene-level analysis of the present surfaceome (see STAR Methods) of ex vivo isolated CD14+ MOs, MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h) revealed only a small number of DEGs for MOs-M-CSF and MOs-GM-CSF but a markedly different expression profile of surface markers related to MOs-GM-CSFIL-4(0-72h) (Figure 2D). Additionally, profiling the expression of pattern recognition receptors (PRRs) on these cell populations revealed clear differences (Figure S2D). MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h) decreased expression of components of the inflammasome signaling complex (CASP1, NLRP1, 2, and 3) but increased expression of the intracellular PRR NOD1. NOD2 expression was maintained only by MOs-GM-CSF. MOs-GM-CSFIL-4 (0-72h) displayed a unique set of PRRs (high CD209 and CLEC10A, loss of Toll-like receptor [TLR] 7 and 5). To determine transcriptional differences, we performed co-expression network analysis based on the union of DEGs between MOs and the three differentiation conditions and mapped differential gene expression onto the network (Figure 2E; Table S2). Within the network topology, a central large MO-related gene cluster was surrounded by separate clusters for each of the three differentiation conditions, indicating that despite an overall close relationship, MOs-M-CSF and MOs-GM-CSF are characterized by signal-specific subclusters of regulated genes (Figure 2E). Identification of DEGs between MOs and MOs-M-CSF and MOs-GM-CSF further supported a close overall relationship, but also indicated differently regulated genes in only one or the other condition (Figures S2E and S2F; Table S2). Gene ontology enrichment analysis (GOEA) revealed that enriched terms in MOs-GM-CSF are associated with immune response and regulation of protein metabolism, whereas most GO terms enriched in MOs-M-CSF relate to metabolism and G-protein-coupled receptor signaling (Figure S2G; Table S2).
Figure 2
MOs-GM-CSFIL-4 Are Most Distinct from MOs-M-CSF and MOs-GM-CSF
(A) Schema describing the questions addressed here and in Figure S2.
(B) PCA (18,318 present probes).
(C) Heatmap of 1,000 most variable genes in dataset. Log2-expression values, z-transformed, scaled (−2 [blue] to 2 [red]).
(D) Heatmap of specifically expressed genes in a single out of the three MO-derived cells versus CD14+ MOs and the other two MO-derived cell types. Log2-expression values, z-transformed, scaled (−1.5 [blue] to 1.5 [red]).
(E) Co-expression networks (union of 2,086 DEG [fold-change > 2 or < −2 and FDR-adjusted p value < 0.05]) between each of three MO-derived cells types versus CD14+ MOs. Fold-change of respective cell type versus overall mean mapped onto networks and displayed blue (negative fold-change) over white to red (positive fold-change). Based on fold-change patterns, networks were divided into four clusters, each cluster representing one of four cell types.
(F) Co-expression network (411 TRs expressed in dataset). For each cell type, fold-change of respective cell type versus overall mean mapped onto network. Cell type-specific clusters of upregulated regulators were generated, indicated by color-coded shadings behind network. TRs highlighted in red were predicted as unique master regulators of corresponding cell type. Prediction performed on all genes highlighted in red (fold-change > 1.5 over overall mean) in corresponding cell type-specific cluster in (E).
(G) t-SNE display of CD14+ MOs, MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-144h) analyzed by MC (n = 3).
(H) Heatmap and HC of mean surface marker expression analyzed using MC. Normalized intensity values, z-transformed, scaled (−6 [blue] to 6 [red]). Color code depicts cluster assignment according to culture condition. Color code as in (G).
Please also see Figure S2.
MOs-GM-CSFIL-4 Are Most Distinct from MOs-M-CSF and MOs-GM-CSF(A) Schema describing the questions addressed here and in Figure S2.(B) PCA (18,318 present probes).(C) Heatmap of 1,000 most variable genes in dataset. Log2-expression values, z-transformed, scaled (−2 [blue] to 2 [red]).(D) Heatmap of specifically expressed genes in a single out of the three MO-derived cells versus CD14+ MOs and the other two MO-derived cell types. Log2-expression values, z-transformed, scaled (−1.5 [blue] to 1.5 [red]).(E) Co-expression networks (union of 2,086 DEG [fold-change > 2 or < −2 and FDR-adjusted p value < 0.05]) between each of three MO-derived cells types versus CD14+ MOs. Fold-change of respective cell type versus overall mean mapped onto networks and displayed blue (negative fold-change) over white to red (positive fold-change). Based on fold-change patterns, networks were divided into four clusters, each cluster representing one of four cell types.(F) Co-expression network (411 TRs expressed in dataset). For each cell type, fold-change of respective cell type versus overall mean mapped onto network. Cell type-specific clusters of upregulated regulators were generated, indicated by color-coded shadings behind network. TRs highlighted in red were predicted as unique master regulators of corresponding cell type. Prediction performed on all genes highlighted in red (fold-change > 1.5 over overall mean) in corresponding cell type-specific cluster in (E).(G) t-SNE display of CD14+ MOs, MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-144h) analyzed by MC (n = 3).(H) Heatmap and HC of mean surface marker expression analyzed using MC. Normalized intensity values, z-transformed, scaled (−6 [blue] to 6 [red]). Color code depicts cluster assignment according to culture condition. Color code as in (G).Please also see Figure S2.To identify the transcriptional regulators (TRs) involved in MO differentiation, we predicted upstream transcription factors for each of the four condition-specific clusters identified in Figure 2E and generated a co-expression network of TRs expressed in the dataset to identify specific clusters of upregulated TRs for CD14+ MOs (yellow), MOs-M-CSF (turquoise), MOs-GM-CSF (light blue), and MOs-GM-CSFIL-4(0-72h) (dark blue) (Figure 2F; Table S2). Finally, we mapped the predicted master transcription factors onto the co-expression network and identified NFIL3, ATF4, and ETS2 among others to be putative regulators of CD14+ MOs. TCF12, MEF2C, and ARID3A were predicted to specifically regulate MOs-M-CSF, whereas ESR1, MTF1, and SREBF1 were anticipated to regulate MOs-GM-CSF identity. RELB, implicated as important for mouse DC differentiation (Wu et al., 1998), was predicted as a central regulator of the transcriptional identity of MOs-GM-CSFIL-4(0-72h), highlighting the uniqueness of the transcriptional identity induced by GM-CSF and IL-4.Since traditional surface markers (CD14, CD11b, CD209) are not informative to discriminate in vitro polarized MO subsets, we designed a comprehensive MC panel incorporating markers defined by our transcriptomic approach and compared ex vivo isolated blood CD14+ MOs to in vitro M-CSF-, GM-CSF-, and GM-CSF + IL-4-polarized CD14+ MOs. Dimensionality reduction using t-distributed neighbor embedding (t-SNE) (Maaten and Hinton, 2008) of the CD45+lin−HLA-DR+ cell fraction comparing CD14+ MOs, MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-144h) revealed donor-independent separation into four cellular clusters (Figures 2G, 2H, and S2H; Table S7). Overlaying their differentiation history on the t-SNE topology shows that the identified clusters corresponded to the four differentiation cues used, validating their transcriptomic differences (Figure 2B). Widely used markers for the delineation of MOs, Macs, and DCs, such as CD11b, CD68, CD11c, and HLA-DR were expressed uniformly across all four clusters, showing that only a high-dimensional phenotyping approach enables robust detection of polarized subsets across all four differentiation conditions (Figure 2H). CD14+ MOs showed a high expression of CLA and CD64, whereas MOs-GM-CSF displayed a high expression of MARCO. VSIG4 was expressed by MOs-GM-CSF and MOs-M-CSF, whereas MOs-M-CSF expressed high amounts of CD163, CD169, and MERTK. The MOs-GM-CSFIL-4(0-144h) cluster specifically expressed SEPP1, FcεR1, CD1c, and CD48. Taken together, MC enabled us to identify transcriptionally validated markers, facilitating separation of different transcriptomic entities on the protein level.
In Vitro Differentiated MO Subsets Are Functionally and Metabolically Different
To understand how the transcriptomic and phenotypic differences between the Mo-derived cells influence their ability to phagocytize and to secrete cytokines in response to PRR stimulation, their motility, and their metabolic profile, we first measured receptor-mediated uptake of GFP-labeled yeast or YG beads. MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h) phagocytosed GFP+ yeast buds after 1 hr of incubation, similarly indicating no differential induction of receptors and signaling pathways involved in yeast uptake (Figures 3A and S3A). MOs-M-CSF displayed an up to 12 times enhanced uptake of YG beads in comparison to MOs-GM-CSFIL-4(0-72h) (Figures 3B and S3B), strongly suggesting that M-CSF but not GM-CSF drives phagocytic capacity in MOs. When assessing cell motility (Figures 3C, S3C, and S3D), MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h) showed very little, intermediate, and high motility, respectively, as assessed by distance and velocity analysis suggesting that migratory capacity was linked to GM-CSF and further potentiated by IL-4 signaling (Figures 3C, S3C, and S3D). Metabolically, MOs-M-CSF and MOs-GM-CSF showed a similar rate of oxidative phosphorylation (OXPHOS), extracellular acidification (ECAR), ATP production, and glycolytic capacity (Figures 3D, 3E, and S3E–S3I). MOs-GM-CSFIL-4(0-72h) displayed an increase of OXPHOS, ATP production, and glycolytic capacity alongside an elevated maximal respiration capacity, indicative of increased energetic fitness paired with higher metabolic flexibility induced by IL-4-specific signaling (Figures 3D, 3E, and S3E–S3I).
Figure 3
Prediction of Differentiated MO Functionality
(A and B) Flow cytometry analysis of MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h) after incubation with (A) GFP-expressing yeast (1 hr, histogram: representative result [n = 3], bar plot [n = 3], mean + SEM, one-way RM ANOVA, with p > 0.05) or (B) YG beads (4 hr, n = 5–6, mean + SEM, one-way RM ANOVA and Tukey’s method for multiple test correction, with ∗∗∗p < 0.001).
(C) Migration tracks (3 hr) of MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h) (representative result, n = 3).
(D and E) OCR (D) and ECAR (E) of MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h).
(F) Heatmap of mean secreted cytokine concentrations (n = 4). Expression values, z-transformed, scaled (–3 [blue] to 3 [red], raw data, Table S3).
Please also see Figure S3.
Prediction of Differentiated MO Functionality(A and B) Flow cytometry analysis of MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h) after incubation with (A) GFP-expressing yeast (1 hr, histogram: representative result [n = 3], bar plot [n = 3], mean + SEM, one-way RM ANOVA, with p > 0.05) or (B) YG beads (4 hr, n = 5–6, mean + SEM, one-way RM ANOVA and Tukey’s method for multiple test correction, with ∗∗∗p < 0.001).(C) Migration tracks (3 hr) of MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h) (representative result, n = 3).(D and E) OCR (D) and ECAR (E) of MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-72h).(F) Heatmap of mean secreted cytokine concentrations (n = 4). Expression values, z-transformed, scaled (–3 [blue] to 3 [red], raw data, Table S3).Please also see Figure S3.Linking these phagocytosis, migratory, and metabolic data back to their transcriptomes, we identified key genes involved in the regulation of these processes. RAB10 (Cardoso et al., 2010), MSR1 (Bonilla et al., 2013), and DAB2 (Tumbarello et al., 2013) were implicated as crucial regulators of phagocytosis in immune cells. These genes alongside other regulators of phagocytosis (RAPH1, RILPL2, TNS1, SCARB2) were markedly upregulated in MOs-M-CSF (Figure S3J; Table S3). Essential migration and cell motility molecules (Lymphotoxin β [LTB] [Yu et al., 2002], CCL13 [Stellato et al., 1997], CCL22 [Godiska et al., 1997], ASAP1 [Curtis et al., 2015]) were upregulated in MOs-GM-CSFIL-4(0-72h) (Figure S3K; Table S3). Genes regulating glycolysis (PFKL, PFKP) were highly upregulated in MOs-GM-CSFIL-4(0-72h) corresponding to the higher glycolytic capacity shown before (Figure S3L; Table S3). Furthermore, UQCRC1, SDHA, ATP5D, COX10, and ATP5I, genes of the respiratory chain and IDH3G, a molecule involved in the TCA cycle were highly upregulated in MOs-GM-CSFIL-4(0-72h), further linking the transcriptomic and functional level.Finally, we stimulated Mo-derived cells with LPS, LPS + interferon-γ (IFN-γ), CL097, or Flagellin and measured cytokine release in response to PRR ligation (Figure 3F; Table S3). IL-10, IL-1B, and MCP1 were secreted only by MOs-M-CSF upon activation with either of the four stimuli, demonstrating their similarities to in vivo infMs. Conversely, IL-12p70 was secreted only by MOs-GM-CSFIL-4(0-72h) upon LPS, LPS + IFN-γ, and CL097 stimulation, indicating a functional overlap with DCs regarding the induction of T helper 1 (Th1) cell responses. IL-23, a major inflammatory driver essential for Th17 T cell induction was produced only by cells differentiated by GM-CSF. Both MOs-GM-CSF and MOs-GM-CSFIL-4(0-72h) produced IL-23 following LPS + IFN-γ and CL097 stimulation, respectively. This is in line with their similarity to infDCs (Segura et al., 2013).
IL-4 Regulates Transcriptomic and Functional Polarization of moDCs and MO-Derived “M2-like” Macs
Presuming that IL-4 induced a functional and phenotypic convergence of MOs and DCs, we next asked whether MOs-GM-CSFIL-4(0-144h) were distinct from what was previously described as M2 Macs, better described as MO-derived Macs further activated by IL-4, termed here MOs-GM-CSFIL-4(72-144h). To reduce variables to a minimum, we generated MOs-GM-CSFIL-4(72-144h) with GM-CSF, so that the only differences to MOs-GM-CSFIL-4(0-144h) cells were the onset and duration of IL-4 exposure. As controls, MOs polarized for only 3 days with either GM-CSF (MOs-GM-CSFIL-4(0h)) or GM-CSF+IL-4 (MOs-GM-CSFIL-4(0-72h)) were included (Figure 4A). PCA determined that MOs-GM-CSFIL-4(72-144h) were surprisingly distinct from MOs-GM-CSFIL-4(0-72h), (0-144h) (Figures 4B, S4A, and S4B). Co-expression network analysis (Figure 4C), HC using the most variable genes (Figure 4D; Table S4), and PCCM analysis (Figure S4C) supported these results. Surprisingly, CD23 (Figure S4D)—a marker formerly associated with MOs-GM-CSFIL-4(72-144h) (Mantovani et al., 2002)—and CD209 (Figure S4E)—which has been linked to MOs with DC functionality (Geijtenbeek et al., 2000)—were expressed similarly between MOs-GM-CSFIL-4(72-144h) and MOs-GM-CSFIL-4(0-72h), (0-144h). No difference in the expression of MMP-12 (Figures S4F and S4G) nor in the release of CCL22 (Figure S4H) were detected. In contrast, MOs-GM-CSFIL-4(0-144h) took up more yeast within 60 min of exposure time (Figure 4E) but phagocytosed similar amounts of YG bead (Figure S4I) when compared to MOs-GM-CSFIL-4(72-144h). MOs-GM-CSFIL-4 (0-144h) were more motile than MOs-GM-CSFIL-4(72-144h) (Figure 4F), with significantly higher accumulated distance (Figure S4J) and velocity (Figure S4K). Investigation of metabolic parameters of MOs-GM-CSFIL-4(0-144h) and MOs-GM-CSFIL-4(72-144h) revealed no differences in their rate of OXPHOS, ATP production, and glycolysis (data not shown). Collectively, this strongly suggest that the IL-4 signal integrates in a time-dependent manner representing a critical checkpoint for MO differentiation.
Figure 4
MOs-GM-CSFIL-4(0-144h) Differ from MOs-GM-CSFIL-4(72-144h) MO-Derived Cells
(A) Schema describing questions addressed herein and in Figure S4.
(B) PCA (18,857 present probes).
(C) Co-expression network (13,691 present genes) describing relationships between CD14+ MOs and four types of MO-derived cells.
(D) Heatmap of 1,000 most variable genes in dataset. Log2-expression values, z-transformed, scaled (−2 [blue] to 2 [red]). Highly expressed genes grouped together (black boxes) according to cell type. Corresponding group-related cell types highlighted, left side of heatmap. Important genes of each cluster depicted, right side of heatmap.
(E) Flow cytometry analysis of MOs-GM-CSFIL-4(72-144h) and MOs-GM-CSFIL-4(0-144h) after incubation with GFP-expressing yeast (1 hr, n = 4–6, mean + SEM, Student’s t test with ∗p < 0.05).
(F) Migration tracks (3 hr) of MOs-GM-CSFIL-4(72-144h) and MOs-GM-CSFIL-4(0-144h) (representative result, n = 3).
Please also see Figure S4.
MOs-GM-CSFIL-4(0-144h) Differ from MOs-GM-CSFIL-4(72-144h) MO-Derived Cells(A) Schema describing questions addressed herein and in Figure S4.(B) PCA (18,857 present probes).(C) Co-expression network (13,691 present genes) describing relationships between CD14+ MOs and four types of MO-derived cells.(D) Heatmap of 1,000 most variable genes in dataset. Log2-expression values, z-transformed, scaled (−2 [blue] to 2 [red]). Highly expressed genes grouped together (black boxes) according to cell type. Corresponding group-related cell types highlighted, left side of heatmap. Important genes of each cluster depicted, right side of heatmap.(E) Flow cytometry analysis of MOs-GM-CSFIL-4(72-144h) and MOs-GM-CSFIL-4(0-144h) after incubation with GFP-expressing yeast (1 hr, n = 4–6, mean + SEM, Student’s t test with ∗p < 0.05).(F) Migration tracks (3 hr) of MOs-GM-CSFIL-4(72-144h) and MOs-GM-CSFIL-4(0-144h) (representative result, n = 3).Please also see Figure S4.
Timing of IL-4 Stimulation Determines Transcriptional Regulation MOs-GM-CSFIL-4
Differences between MOs-GM-CSFIL-4(0-144h) and MOs-GM-CSFIL-4(72-144h) could be explained either by a dichotomous model with MOs differentiating into MOs with DC or Mac functionality, or a continuum model that integrates time of exposure suggesting plasticity of MO-derived cells. To determine the best-fitting model, we performed a time kinetics experiment, by adding IL-4 at the start of the culture (Mo-GM-CSFIL-4(0-144h)) or at 12 (MOs-GM-CSFIL-4(12-144h)), 24 (MOs-GM-CSFIL-4(24-144h)), 48 (MOs-GM-CSFIL-4(48-144h)), or 72 (MOs-GM-CSFIL-4(72-144h)) hr after initiation of differentiation with GM-CSF (Figure 5A). Transcriptomes were assessed after 144 hr and MOs-GM-CSFIL-4(0h) and MOs-GM-CSFIL-4(0-72h) were used as controls. Assessing CD14 (Figures 5B and S5A) and CD209 (Figures 5C and S5B) expression, we observed a dichotomous distribution for cells differentiated with IL-4 being CD14loCD209hi while cells differentiated only by GM-CSF were CD14+CD209−/lo. In contrast, transcriptomic analysis revealed a different model (Figure S5C). PCA (Figure 5D), HC (Figure 5E; Table S5), and self-organizing map (SOM) clustering (Figure 5E) revealed a gradual ordering of samples corresponding to the exposure time to IL-4 indicating a differentiation continuum. Mapping gene expression information for each time point onto a co-expression network (Figure 5F; Table S5) revealed a dense network with two major clusters, one characterized by genes elevated in MOs-GM-CSFIL-4(0) (0 hr, red: upregulated; blue: downregulated), the other one defined by IL-4 exposure (Figure 5G). Adjacent time points showed partially overlapping gene sets suggesting a plastic continuum integrating IL-4 signaling over time, arguing against the dichotomous model of polarization. Collectively, these data suggest that IL-4 signaling differentiates MOs along a transcriptomic continuum with MOs-GM-CSFIL-4(0h) and MOs-GM-CSFIL-4(72-144h) being at the extreme ends.
Figure 5
Timing of IL-4 Determines Transcriptional Regulation in Differentiated MOs
(A) Schema describing the IL-4 time kinetic experiment.
(B and C) Histograms, relative expression of CD14 and CD209 analyzed by flow cytometry.
(D) PCA (12,794 present genes).
(E) Heatmap of 1,000 most variable genes across dataset. Log2-expression values, z-transformed, scaled (−2 [blue] to 2 [red]). Below, SOM-clustering (12,794 present genes across cell types).
(F) Co-expression networks (union of 2,775 DEGs, fold-change > 1.5 or < −1.5 and FDR-corrected p value < 0.05) between MOs-GM-CSFIL-4 and MOs-GM-CSF. For each cell type, fold-change of respective cell type versus overall mean mapped onto networks, displayed in blue (fold-change ≤ 1.5) or red (fold-change ≥ 1.5).
(G) Example of genes located in condition-related clusters depicted in (F) and in first column (fold-change ≥ 1.5). First column: condition-specific genes; following columns: genes shared between clusters of two consecutive time points.
Please also see Figure S5.
Timing of IL-4 Determines Transcriptional Regulation in Differentiated MOs(A) Schema describing the IL-4 time kinetic experiment.(B and C) Histograms, relative expression of CD14 and CD209 analyzed by flow cytometry.(D) PCA (12,794 present genes).(E) Heatmap of 1,000 most variable genes across dataset. Log2-expression values, z-transformed, scaled (−2 [blue] to 2 [red]). Below, SOM-clustering (12,794 present genes across cell types).(F) Co-expression networks (union of 2,775 DEGs, fold-change > 1.5 or < −1.5 and FDR-corrected p value < 0.05) between MOs-GM-CSFIL-4 and MOs-GM-CSF. For each cell type, fold-change of respective cell type versus overall mean mapped onto networks, displayed in blue (fold-change ≤ 1.5) or red (fold-change ≥ 1.5).(G) Example of genes located in condition-related clusters depicted in (F) and in first column (fold-change ≥ 1.5). First column: condition-specific genes; following columns: genes shared between clusters of two consecutive time points.Please also see Figure S5.
NCOR2 Is a Transcriptional Regulator of MOs-GM-CSFIL-4
To understand how IL-4 enforces the unique transcriptional program in MOs-GM-CSFIL-4(0-72h), a co-expression TR network was generated first using TRs differentially expressed between MOs, MOs-GM-CSF, and MOs-GM-CSF differentiated with different durations of IL-4 (Figures 6A, 6B, and S6A). TRs were then filtered according to differential expression in MOs-GM-CSFIL-4(0-72h) and ranked based on absolute expression (Figures 6A and 6B, red: upregulated, blue: downregulated; Table S6). We identified seven eligible TRs with NCOR2 showing the highest expression (Figure 6C). We confirmed NCOR2 expression in MOs-GM-CSF and MOs-GM-CSFIL-4(0-72h) (n = 2) using real-time PCR showing that MOs-GM-CSFIL-4(0-72h) expressed significantly more NCOR2 and CD209 mRNA (Figures S6B and S6C). Analysis of NCOR2 protein expression by immunoblot analysis (Figures 6D and 6E) or confocal microscopy (Figures 6F and 6G) revealed a significant enrichment of NCOR2 in the nucleus of MOs-GM-CSFIL-4(0-144h) but not in MOs-GM-CSF. Using nanostraw technology (Figure S6D), we introduced anti-NCOR2 siRNAs for the last 24 hr of the MOs-GM-CSFIL-4(0-72h) differentiation. After silencing, mean NCOR2 mRNA expression in MOs-GM-CSFIL-4(0-72h) was reduced to 65% relative to the control (Figures S6E and S6F), reflecting effective downregulation of NCOR2 transcription considering its long half-life of more than 24 hr (Raghavan et al., 2002). siRNA silencing of NCOR2 in MOs-GM-CSFIL-4(0-72h) also reduced CD209 mRNA (Figures S6E and S6F) and protein (Figure S6G). To understand the impact of NCOR2 on IL-4-mediated transcriptional regulation in MOs-GM-CSFIL-4(0-72h), we performed a global transcriptome analysis of anti-NCOR2 siRNA-treated MOs-GM-CSFIL-4(0-72h) versus scrambled siRNA-treated MOs-GM-CSFIL-4(0-72h) (Figures S6H and S6I). NCOR2 silencing resulted in 1,834 variable genes (Figure S6J; Table S6). Classification of NCOR2-regulated genes, based on a literature-derived IL-4 signature (GEO: GSE13762, GSE35304, GSE32164; 457 induced, 498 repressed genes) revealed that the large majority of genes regulated by NCOR2 are IL-4 signature genes (Figures 6H and S6J; Table S6). Additionally, GSEA was performed separately for up- and downregulated genes within the IL-4 signature to identify signature enrichments in the transcriptomes of MOs-GM-CSFIL-4(0-72h) treated with scrambled siRNAs or with anti-NCOR2 siRNAs (Figures 6I and 6J). We found a statistically significant enrichment of genes upregulated in the IL-4 signature in MOs-GM-CSFIL-4(0-72h) treated with scrambled siRNAs (Figure 6I) whereas genes downregulated in the IL-4 signature were enriched in MOs-GM-CSFIL-4(0-72h) treated with anti-NCOR2 siRNAs (Figure 6J), establishing NCOR2 as a key regulator for IL-4-induced MO differentiation.
Figure 6
NCOR2 Is a Transcriptional Regulator of MOs-GM-CSFIL-4(0-72/144h)
(A and B) Co-expression networks (267 TRs) for (A) MOs-GM-CSF, (B) MOs-GM-CSFIL-4(0-72h), fold-change versus CD14+ MOs mapped onto network. MOs-GM-CSFIL-4(0-72h)-specific cluster of elevated regulators (dark blue).
(C) TR Heatmap, specifically upregulated in MOs-GM-CSFIL-4(0-72h), MOs-GM-CSFIL-4(0-144h) versus CD14+ MOs, MOs-GM-CSF. Log2-expression values, z-transformed, scaled (−1.15 [blue] to 1.15 [red]).
(D) Immunoblot of NCOR2, Lamin A/C, and β-tubulin in cytoplasm and nucleus of MOs-GM-CSF and MOs-GM-CSFIL-4(0-144h) (representative result, n = 4).
(E) Quantification of relative enrichment of NCOR2 in nuclear fractions isolated from MOs-GM-CSF and MOs-GM-CSFIL-4(0-144h) (n = 4, ∗p < 0.05, mean ± SEM).
(F) Representative confocal microscopy images of medial nuclear region (MNR) of MOs-GM-CSF and MOs-GM-CSFIL-4(0-144h) (n = 3, green: NCOR2; blue: DAPI; red: tubulin).
(G) Quantification of mean fluorescence intensity of NCOR2 in MNR identified by confocal microscopy of GM-CSF and MOs-GM-CSFIL-4(0-144h) (n = 3, green: NCOR2; blue: DAPI; red: tubulin, ∗p < 0.05, mean ± SEM).
(H) Scatterplot (1,834 variable genes in dataset) containing cells treated with αNCOR2 (y axis), scrambled siRNA (x axis, log2-mean expression values). Highlighted genes determined (see Figure S6F) induced (red) or repressed (blue) by IL-4.
(I and J) GSEA of genes upregulated (I), downregulated (J) in IL-4 signature in MOs-GM-CSFIL-4(0-72h) treated with scrambled or αNCOR2 siRNA (nominal p value, empirical phenotype-based permutation test [p < 0.05, FDR < 0.25], 1,000 samples permutations).
Please also see Figure S6.
NCOR2 Is a Transcriptional Regulator of MOs-GM-CSFIL-4(0-72/144h)(A and B) Co-expression networks (267 TRs) for (A) MOs-GM-CSF, (B) MOs-GM-CSFIL-4(0-72h), fold-change versus CD14+ MOs mapped onto network. MOs-GM-CSFIL-4(0-72h)-specific cluster of elevated regulators (dark blue).(C) TR Heatmap, specifically upregulated in MOs-GM-CSFIL-4(0-72h), MOs-GM-CSFIL-4(0-144h) versus CD14+ MOs, MOs-GM-CSF. Log2-expression values, z-transformed, scaled (−1.15 [blue] to 1.15 [red]).(D) Immunoblot of NCOR2, Lamin A/C, and β-tubulin in cytoplasm and nucleus of MOs-GM-CSF and MOs-GM-CSFIL-4(0-144h) (representative result, n = 4).(E) Quantification of relative enrichment of NCOR2 in nuclear fractions isolated from MOs-GM-CSF and MOs-GM-CSFIL-4(0-144h) (n = 4, ∗p < 0.05, mean ± SEM).(F) Representative confocal microscopy images of medial nuclear region (MNR) of MOs-GM-CSF and MOs-GM-CSFIL-4(0-144h) (n = 3, green: NCOR2; blue: DAPI; red: tubulin).(G) Quantification of mean fluorescence intensity of NCOR2 in MNR identified by confocal microscopy of GM-CSF and MOs-GM-CSFIL-4(0-144h) (n = 3, green: NCOR2; blue: DAPI; red: tubulin, ∗p < 0.05, mean ± SEM).(H) Scatterplot (1,834 variable genes in dataset) containing cells treated with αNCOR2 (y axis), scrambled siRNA (x axis, log2-mean expression values). Highlighted genes determined (see Figure S6F) induced (red) or repressed (blue) by IL-4.(I and J) GSEA of genes upregulated (I), downregulated (J) in IL-4 signature in MOs-GM-CSFIL-4(0-72h) treated with scrambled or αNCOR2 siRNA (nominal p value, empirical phenotype-based permutation test [p < 0.05, FDR < 0.25], 1,000 samples permutations).Please also see Figure S6.
MC Identifies Heterogeneity in Surface Molecules Expressed by MOs-GM-CSFIL-4(0-144h)
Within the MC-defined populations (MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-144h)) (Figure 2G), we recognized a sub-cluster structure indicating further heterogeneity within each in vitro subset (Figures 7A and 7B). MOs-M-CSF revealed four subpopulations (clusters 2, 5, 10, 11) characterized by co-expression of tissue Mac markers MERTK, CD64, CD169, and CD163 and by variable expression of L-selectin (CD62L, low in cluster 2, 11) and CD26 (low in cluster 2), indicating different migration and maturation profiles (Figure 7B). MOs-GM-CSF were defined by two subclusters (clusters 7 and 8). Both clusters expressed CD14, CD64, CD68, and CD206, but not Mac markers such as MERTK, CD169, and CD163, further corroborating their difference from MOs-M-CSF. Cluster 7 additionally expressed the molecules FcεR1 and CD1a, supporting an activated status of this subpopulation. MOs-GM-CSFIL-4(0-144h) showed a similar heterogeneity as MOs-M-CSF, but with more pronounced phenotypic differences (cluster 1, 3, 4, 9). Therefore, we isolated the MOs-GM-CSFIL-4(0-144h) dataset and analyzed it using Rphenograph (Figures 7C–7E). This analysis revealed 11 phenotypically different clusters within MOs-GM-CSFIL-4(0-144h) (Figures 7C and 7D). All 11 subpopulations expressed the moDC markers CD1c, CD226, CD48, and CD11c. The biggest differences across all different subpopulations were seen in the expression of activation and antigen presentation-associated molecules, such as CD1a, HLA-DR, FcεR1, CD62L, and CD86. Marker expression differed between the clusters 6, 1, 10, 11 showing high expression of CD1a and CD62L, whereas only cluster 5 and 2 expressed high amounts of the co-stimulatory molecule CD86. Expression of CD11b, a marker used to analyze the purity of moDCs in the clinic, varied across the examined moDC population, exemplifying the urgent need for a higher-dimensional phenotyping to improve purity and therapeutic outcome when using these differentiated MOs therapeutically (Figure 7E).
Figure 7
Mass Cytometry Analysis Identifies Unappreciated Phenotypic Heterogeneity in Clinically Relevant Mo-GM-CSFIL-4(0-144h) Cultures
(A) Phenograph of CD14+ MOs, MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-144h) based on MC expression data (n = 3, 36 myeloid-related surface markers). Affiliation of cells to the 11 identified clusters indicated by color coding and visualized in t-SNE plot.
(B) Heatmap and HC of mean surface marker expression of 11 individual clusters. Right side: Differentiating conditions according to (A) and Figure 2G.
(C) Phenograph of MOs-GM-CSFIL4(0-144h) and visualized in t-SNE plot (representative donor; 3,500 cells).
(D) Heatmap and HC of mean surface marker expression of 11 individual clusters.
(E) Expression feature plot of the depicted surface markers in MOs-GM-CSFIL-4(0-144h).
Mass Cytometry Analysis Identifies Unappreciated Phenotypic Heterogeneity in Clinically Relevant Mo-GM-CSFIL-4(0-144h) Cultures(A) Phenograph of CD14+ MOs, MOs-M-CSF, MOs-GM-CSF, and MOs-GM-CSFIL-4(0-144h) based on MC expression data (n = 3, 36 myeloid-related surface markers). Affiliation of cells to the 11 identified clusters indicated by color coding and visualized in t-SNE plot.(B) Heatmap and HC of mean surface marker expression of 11 individual clusters. Right side: Differentiating conditions according to (A) and Figure 2G.(C) Phenograph of MOs-GM-CSFIL4(0-144h) and visualized in t-SNE plot (representative donor; 3,500 cells).(D) Heatmap and HC of mean surface marker expression of 11 individual clusters.(E) Expression feature plot of the depicted surface markers in MOs-GM-CSFIL-4(0-144h).
Discussion
HumanCD14+ MOs differentiated by either M-CSF, GM-CSF, or GM-CSF + IL-4 (formerly moDCs [Sallusto and Lanzavecchia, 1994]) have been extensively studied as in vitro models for Macs, “M2-like” Macs, or DCs (Akagawa et al., 2006, Sallusto and Lanzavecchia, 1994). In the past, in vitro generated Mac and DC models were justified by morphological, phenotypical, and functional similarities to cells identified in vivo. Recently, it became clear that Macs, MOs, and DCs present with very different transcriptomes in vivo, reflecting different ontogeny (Guilliams et al., 2014). Therefore, reassessment of the relationships between M-CSF-, GM-CSF-, and GM-CSF + IL-4-differentiated MOs and their alleged in vivo counterparts on a transcriptomic, phenotypic, and functional level is imperative.Using data of ascites-associated infMs and infDCs (Segura et al., 2013) revealed that MOs-M-CSF aligned closely to infMs on the transcriptional level. This was phenotypically supported by low expression of CD1a and FcεR1 while MOs-GM-CSFIL-4 aligned closely to infDCs with high expression of CD1a, CD206, and FcεR1. Taken together, the currently available in vitro models best resemble inflammatory rather than homeostatic DC and Mac populations, therefore, serving best as reductionist models to study the role of these cells in inflammation. For the future, we encourage researchers to identify culture conditions for MO-derived cells that resemble homeostatic DC or Mac phenotypes guided by transcriptomic analysis, as shown for in vitro cultures of microglia (Butovsky et al., 2014), and pursue the identification of dedicated progenitors of DCs in the human blood and bone marrow (Breton et al., 2015).Transcriptome analysis defined the cellular relationships between these model systems on a much more precise level and revealed a close association of MOs differentiated by M-CSF and GM-CSF, while IL-4 was the major driver for moDC identity. In addition, we found that MOs integrate the GM-CSF and IL-4 signals over time, which necessitates a reassessment of our dichotomous definition of MOs differentiating toward a Mac or DC-like phenotype. The varying time of onset and the variance in overall exposure to IL-4 resulted in gradually changed transcriptional and functional identities of MO-derived cells. The heterogeneity defined by MCs added an additional layer of individual cell reactivity. These observations challenge a dichotomous differentiation model of MOs toward Macs or DCs induced by a single cytokine, and instead support a dynamic model in which cell identity is a function of the duration of signal input. These findings might also help to better understand IL-4-mediated MO-derived cell expansion during Th2 cell-driven inflammatory conditions, such as parasitic infections and tissue repair (Rückerl and Allen, 2014).We identified NCOR2 as a key transcriptional hub linked to IL-4-dependent differentiation of MOs. NCOR2 plays an important role during development and homeostatic processes in muscle, adipose, and liver tissues (Mottis et al., 2013). NCOR2 has been shown to be important for cell fate, differentiation, and lineage progression in the nervous system (Jepsen et al., 2007). Elevated NCOR2 expression was observed in tissues with high OXPHOS activity, similar to our observations of elevated expression of NCOR2 in MOs-GM-CSFIL-4 (Reilly et al., 2010). Additionally, signaling by nuclear receptors, such as peroxisome proliferator-activated receptor-γ (PPAR-γ) or liver X receptor, NCOR2 was linked to the repression of NF-κB target genes in response to LPS stimulation of Macs (Pascual et al., 2005). It has been speculated that IL-4 activation of human MOs leads to endogenous PPAR-γ ligand production (Czimmerer et al., 2012), but further work is necessary to establish a link to NCOR2-mediated gene repression.Unexpectedly, we identified significant heterogeneity within MO cultures, which was revealed only by a set of surface markers as applied by MCs. In particular, within MOs-GM-CSFIL-4 we identified subsets that either expressed HLA-DR and CD86 or CD1a and FcεR1, the former representing a subpopulation with elevated antigen presenting capacity. While induction of FcεR1 on CD1a+ DCs derived from CD34+ stem cells has been reported (Allam et al., 2004), this has not been studied during MO to MO-GM-CSFIL-4 differentiation. Considering findings that CD1a+ and CD1a− MOs-GM-CSFIL-4 differ in their capacity to direct Th cell differentiation (Cernadas et al., 2009), monitoring of these cultures in a clinical setting might be beneficial for optimizing efficiency of cellular vaccines. Furthermore, in studies using bulk MOs instead of CD14+ MOs, the problem of heterogeneity might even be more pronounced and should be addressed when analyzing in vitro MO cultures. Consequently, high-dimensional characterization should be used to optimize culture conditions, generate more homogeneous cell populations, and thereby open avenues for optimizing cellular products derived from human MOs for vaccination strategies.
STAR★Methods
Key Resources Table
Contact for Reagent and Resource Sharing
Further information and requests for resources and reagents should be directed and will be fulfilled by the Lead Contact, Andreas Schlitzer (andreas.schlitzer@uni-bonn.de).
Experimental Model and Subject Details
Human primary cell isolation
Buffy coats were obtained from healthy donors (University hospital Bonn, local ethics vote 203/09) after written consent was given according to the Declaration of Helsinki. Peripheral blood mononuclear cells (PBMC) were isolated by Pancoll (PAN-Biotech) density centrifugation from buffy coats. CD14-, CD56-, CD4- and CD19-specific MACS beads (Miltenyi Biotec) were used for the enrichment of CD14+ MO, CD56+ NK cells, CD4+ T cells and CD19+ B cells, respectively. Lung-derived myeloid cells were isolated from human perfusates of lung transplant recipients with informed consent and immediately sorted for CD45+Lin-HLA-DRhi cells using the FACS Fusion cell sorter (BD, USA).
Human primary cell culture
CD14+ MO were cultured in RPMI1640 medium supplemented with 10% FCS, 1% Penicillin-Streptomycin, 1% Sodium pyruvate and 1% Glutamax (GIBCO) for 3 days. CD14+ MO were differentiated into Mo-M-CSF or Mo-GM-CSF in the presence of 100 IU/ml rhM-CSF or 800 IU/ml rhGM-CSF, respectively. Mo-GM-CSFIL-4 were generated by the addition of 800 IU/ml rhGM-CSF and 500 IU/ml rhIL4 and were incubated for up to 6 days. All cytokines were purchased from Immunotools.
Method Details
Flow cytometry
Cells were washed with ice cold PBS. After FcR blockage (Miltenyi, Germany), cells were stained with the respective antibodies in PBS supplemented with 0.5% FCS, 2.5 mM EDTA for 20min at 4°C. Following antibodies were purchased from Biolegend (USA): CD19 (HIB19), CD11b (CBRM1/5), CD23 (EBVCS-5), CD14 (M5E2); R&D: VSIG4 (polyclonal), MARCO (polyclonal), CCR7 (150503); Becton Dickinson (BD, USA): CD56 (B159), HLA-DR (L243), CD209 (DCN46), CD226 (DX11); Data acquisition was performed with a LSR II (BD). Analyses were performed with FlowJo software (Tree Star).
Mass cytometry
Following culture, cells were washed with PBS (GIBCO, Life Technologies, Carlsbad, CA) and stained with cisplatin (Sigma-Aldrich, St Louis, MO). Then, cells were washed with PBS containing 4% FBS and 0.05% Sodium azide and fixed with 2% paraformaldehyde (PFA; Electron Microscopy Sciences, Hatfield, PA) overnight. Cells were permeabilized (1X perm buffer (Biolegend, San Diego, CA)) and stained with metal-conjugated antibodies (Table S7) intracellularly. Then cells were washed, metal barcoded and resuspended in water. EQ Four Element Calibration Beads (Fluidigm Corporation, South San Francisco, CA) were added at a concentration of 1% prior to acquisition. Cells were acquired and analyzed using a CyTOF1 Mass cytometer. The data was normalized and events with parameters having zero values were replaced with a randomized value between 0 and −1.
Mass cytometry data analysis
Normalized MC data was exported in .fcs format and pre-processed in FlowJo Version 9.9.4 (Tree Star Inc). Pre-processing included removal of cell debris (Iridium-191/193 DNA intercalating fraction) and dead cells (cisplatin+). Myeloid CD45+lin(CD3+CD7+CD15+CD19+CD20+CD56)- cells were exported and used for analysis. Downstream analysis was performed using the Cytofkit R package. For comparison of the different differentiating conditions data from 1000 myeloid cells was randomly sampled per donor and condition (3 donors; 12,000 cells in total) and autoLgcl-transformed including expression values for 36 surface markers (CD45, CD14, CD5, CD62L, CD48, CD68, CD66ace, CLA, HLA-DR, CD115, CD64, CD1c, FceR1, SEPP1, CD123, CD163, CXCR3, CD226, CD169, SIRP1a, Dectin1a, CD1a, CD141, MARCO, CD86, CX3CR1, CD206, VSIG4, CD88, CD34, MerTK, CD39, CD26, CD11c, CD11b, CD16). Detailed analysis on the Mo-GM-CSFIL-4 condition was based on 3500 cells from one individual. To define clusters of cell subpopulations, PhenoGraph was used. Points representing individual cells in the t-SNE plots were color-coded to illustrate amount of protein expression or affiliations to clusters, treatment conditions or donors, respectively. Alternatively, the gplots R package was used to generate heatmaps of marker expression of individual cells or mean values over identified cell clusters. Dendrograms represent hierarchical clustering results based on the Euclidean distance measure.
Uptake of fluorescent microbeads or yeast
Cells were incubated either with fluorescent monodispersed polystyrene microspheres (1 μm diameter, Fluoresbrite YG Beads, Polysciences) or yeast (GFP-expressing Pichia Pastoris) in a cell-to-bead ratio of 1/10 for 4h or 60 min at 37°C, respectively. Afterward, cells were harvested, washed and bead/yeast uptake was analyzed by flow cytometry using an LSR II (BD). Negative control samples were kept at 4°C. Data analysis was performed using FlowJo software (Tree Star).
Immuno blot
The Nuclear Extract Kit (Active Motif) was used to fractionate the cytosolic and nuclear proteins, fractions were separated by SDS-PAGE and transferred onto a nitrocellulose membrane (Amersham) by wet blotting. Probing was performed using hNCOR2 (Abcam), β-Tubulin (Cell Signaling) and Lamin A/C (Active Motif) antibodies. For MMP12 protein detection, the cytosolic whole protein fractions (Dignam extraction) were separated by SDS-PAGE and transferred onto a nitrocellulose membrane (Amersham) by semi-dry blotting. Probing was performed by using MMP-12 (R&D) and β-actin antibodies. Signal detection and analysis was performed on the LI-COR Odyssey system. Cell compartment separation efficiency was validated by enrichment of cytosolic proteins, such as β-Tubulin or nuclear proteins such as Lamin A/C. Signal expression values of hNCOR2 and MMP-12 were calculated in semiquantitative relation to the signal expression values of β-Tubulin, Lamin A/C and β-actin following the equation target/reference.
Migration assay
Migration was analyzed in μ-Slide 8 well chambered coverslips (Ibidi) coated with 50μg/ml humanfibronectin (Alfa Aesar). 0.7x105 cells in 300μl VLE-RPMI (Biochrom) were seeded in each well. Live cell imaging of adherent cells was performed at 37°C and 5% CO2 using a fully automated inverted Olympus Fluoview 1000 confocal microscope equipped with motorized xyz stage. Cell motility was monitored over a period of 3h by capture of differential interference contrast images every 5min with a 0.40 UPLAPO 10x Objective (Olympus). Migration parameters were calculated using the Manual Tracking and Chemotaxis Tool plugins in ImageJ.
Confocal microscopy and data analysis
Cover slides were coated with fibronectin (Sigma Aldrich) prior to the seeding of Mo-GM-CSF or Mo-GM-CSFIL-4(0-144h). After an incubation time of 3h at 37°C, adherent cells were washed with PBS and fixed in pre-cooled methanol at −20°C for 10min. Subsequently slides were washed twice in PBS and blocked with 5% BSA in PBS for 30min. Staining for NCOR2 and tubulin was done using αNCOR2 antibody (ab24551, Abcam) at 4°C overnight or using tubulin antibody (MAB1864 (Chemicon)) at room temperature for 45min. For visualization slides were stained with anti-rabbit-Alexa488 (A11034, Life Technologies) or anti-ratAlexa555 (A21434, Invitrogen) and DAPI (Sigma Aldrich) at room temperature for 45min. Slides were mounted with fluoroshield (ImmunoBioScience) plus 1% DABCO and confocal image acquisition of the medial nuclear region was performed using an inverted Olympus Fluoview 1000 confocal microscope equipped with a Plan Apochromat 60x, NA1.4 oil immersion objective (Olympus) and 405nm/488nm/543nm laser lines. Quantification of mean intensity of the green fluorescence in the nuclear region was performed using Imaris 7.6.5 (Bitplane).
NCOR2 siRNA silencing
24h prior to the experiment, the nanostraw cargo chamber was washed 3x with 10μl of 0.5% PEG 3500 (P3640, sigma) in PBS and equilibrated; chambers were equilibrated with 100μl RPMI 1640 media. CD14+ human MO were resuspended in RPMI 1640 with supplements (10% FCS, 1% Pen/Step, 1% GlutaMax and 1% NaPyruvat) and activated. siRNA solutions were prepared in 1x siRNA Buffer (Dharmacon). hNCOR2 siRNA (Dharmacon) ON Target Plus pool was used for silencing (siRNA sequences targeting hNCOR2: 1. UGGUUUACAUGUCGACUAA, 2. UGGUUUACAUGUUGUGUGA, 3. UGGUUUACAUGUUUUCUGA, 4. UGGUUUACAUGUUUUCCUA; non-targeting control siRNA sequences: 1. UGGUUUACAUGUCGACUAA, 2. UGGUUUACAUGUUGUGUGA, 3. UGGUUUACAUGUUUUCUGA, 4. UGGUUUACAUGUUUUCCUA). After removing the equilibration media the tubing system was filled with the siRNA solution or PBS. Subsequently, the cell suspension was filled into the chamber and incubated for 72h at 37°C and 5%CO2. After 72h, cells were directly lysed within the chambers by adding Trizol. qRT-PCR was performed to check silencing efficiency.
Real time PCR
Total RNA was isolated with the miRNeasy Micro kit (QIAGEN) and analyzed with the High Sensitivity RNA Tapestation system (Agilent). cDNA synthesis was prepared using the Transcriptor kit (Roche). qRT-PCR was performed using a Lightcycler 480 system (Roche). Primer sequences: NCOR2 (F, GCGAGGTCTCCCTGAGTCTT; R, CCAGTCCTCGTCATCAGCTC); CD209 (F, CCAGGTGAAGCGGTTACTTC; R, GCTCGTCGTAATCAAAAGTGC); GAPDH (F, AGCCACATCGCTCAGACAC; R, GCCCAATACGACCAAATCC). Expression values were calculated by the comparative CT method. GAPDH served as internal control gene. FC is presented relative to amount of expression in MO
Cell stimulation
Mo-M-CSF, Mo-GM-CSF, Mo-GM-CSFIL-4(0-72h) or Mo-GM-CSFIL-4(0-144h) cells were stimulated overnight under the following conditions: Media, 100ng/ml LPS ultrapure, 100ng/ml LPS ultrapure + 1000U/ml IFNγ, 1μg/ml CL097, 100ng/ml flagellin. Supernatants were harvested 19h after stimulation and stored at −80°C.
Cytokine measurement
Cytokines were measured using LEGENDplex (Biolegend, USA). Diluted cell culture supernatants were incubated for 2 hours with the beads and detection antibodies, followed by 30min incubation with SA-PE. After washing, beads were resuspended in washing buffer and acquired using a LSRII flow cytometer (BD). Data were analyzed with the LEGENDplex Data Analysis Software; concentration values were exported to Excel and visualized in R.
Oxygen consumption rate (OCR) and ECAR measurements
OCR and ECAR were determined using a XF-96 Extracellular Flux Analyzer (Seahorse Bioscience). For ECAR analysis the media was changed to bicarbonate-free RPMI supplemented 10mM glucose, 1mM pyruvate & 2 mM glutamine 1h prior to the assay and the plate was kept in a non-carbonated incubator. Measurements were performed under basal conditions and after the sequential addition of final 1μM oligomycin A, 1.5μM FCCP (fluoro-carbonyl cyanide phenylhydrazone) and 0.5μM rotenone & antimycin each. For ECAR analysis the media was changed to bicarbonate free RPMI supplemented 2 mM glutamine 1h prior to the assay and the plate was kept in a non-carbonated incubator. Measurements were performed under basal conditions and after the sequential addition of final 10 mM glucose, 1μM oligomycin and 100 mM 2-Deoxyglucose. OXPHOS was calculated as (basal OCR – OCR after rotenone & antimycin treatment), ATP production was calculated as (basal OCR – OCR after oligomycin A treatment), maximal respiration was calculated as (OCR after FCCP treatment – OCR after rotenone & antimycin treatment), glycolysis was calculated as (basal ECAR – ECAR after 2-Deoxyglucose treatment), glycolytic capacity was calculated as (ECAR after oligomycin A treatment – ECAR after 2-Deoxyglucose treatment. All reagents were purchased from Sigma, except FCCP was purchased from Tocris. ECAR and OCR raw data was normalized to DNA content using the CyQuant Assay kit (Thermo Fisher).
Microarray data generation
Up to 5 × 106 cells were harvested and lysed in TRIzol (Invitrogen) and RNA was isolated and concentration and purity was assessed using a NanoDrop 1000 UV-Vis Spectrophotometer (Thermo Scientific). Subsequently, the TargetAmp-Nano Labeling Kit for Illumina Expression BeadChip (Epicenter) was utilized to generate biotin labeled anti-sense RNA (cRNA) according to the manufacturer’s protocol. As a quality control, 100 ng cRNA were reverse transcribed to cDNA and a PCR specific for ACTB amplification was performed. For expression profiling, 750 ng cRNA were hybridized to Human HT-12v3 BeadChip arrays (Illumina), stained and imaged on an Illumina iScan system.
RNA-sequencing
Total RNA was converted into libraries of double stranded cDNA molecules as a template for high throughput sequencing following the manufacturer’s recommendations using the Illumina TruSeq RNA Sample Preparation Kit v2. Shortly, mRNA was purified from 100 ng of total RNA using poly-T oligo-attached magnetic beads. Fragmentation was carried out using divalent cations under elevated temperature in Illumina proprietary fragmentation buffer. First strand cDNA was synthesized using random oligonucleotides and SuperScript II. Second strand cDNA synthesis was subsequently performed using DNA Polymerase I and RNase H. Remaining overhangs were converted into blunt ends via exonuclease and polymerase activities and enzymes were removed. After adenylation of 3′ ends of DNA fragments, Illumina PE adaptor oligonucleotides were ligated to prepare for hybridization. DNA fragments with ligated adaptor molecules were selectively enriched using Illumina PCR primer PE1.0 and PE2.0 in a 15 cycle PCR reaction. Size-selection and purification of cDNA fragments with preferentially 200 bp in length was performed using SPRIBeads (Beckman-Coulter). Size-distribution of cDNA libraries was measured using the Agilent high sensitivity DNA assay on a Bioanalyzer 2100 system (Agilent). cDNA libraries were quantified using KAPA Library Quantification Kits (Kapa Biosystems). After cluster generation on a cBot, a 75 bp single-read run was performed on a HiSeq1500.
From a publicly available dataset (GEO: GSE35457, (Haniffa et al., 2012)) the raw data of CD14+ and CD14+CD16+ monocyte samples, as well as of CD1c+ and CD141+ DC samples from blood of up to four donors (I, II, V and VI) were used. The data was generated on Human HT12-V4 BeadChips and contained 47,323 probesets, where values for missing bead types were imputed. Since all other in-house datasets contained CD14+ MO only, the raw expression values for the CD14+ and CD14+CD16+ monocyte samples were averaged always for the same donor. In the end, the designation “CD14+ monocyte” was kept. The resulting raw expression text file was imported into PGS as well.To compare the MO-derived cells to inflammatory conditions, a dataset containing inflammatory Mac and dendritic cells, isolated from the ascites of ovarian cancerpatients, as well as primary BDCA1+ dendritic cells and CD14+CD16- and CD14dimCD16+ MO was downloaded (GEO: GSE40484, (Segura et al., 2013)). The raw Affymetrix Human Gene 1.1 ST Array data was imported into PGS using RMA background correction, quantile normalization and log2-transformation. For the inflammatory cells, only donors 1 to 4 were considered. Afterward, 12,067 present probes were determined requiring at least one group mean to be larger than 7 on log2 scale. Of those, only 9,358 probes being annotated with a gene symbol were further considered. Finally, the probes were reduced to 8,612 present genes, by considering only the probe with the highest expression for a given gene. In order to link the Affymetrix data to the Illumina datasets of MO-derived cells, three different approaches were applied.
Probe matching between platforms
Using the Illumina Probe ID (“ILMN”), 39,423 common probesets were identified on both HT12-V3 and HT12-V4 arrays. All in-house datasets (pre-datasets 1-3) and the dataset from Haniffa et al. were filtered down to those probes to enable dataset assembly and comparability. Since for some of the probes different gene name synonyms were used in the two different platform files, the annotation presented in the newer file (HT12-V4) was used for all filtered datasets later on.
Assembly of datasets and pre-processing
Based on the identified common probes, four different datasets were generated. The first one (dataset 1) was assembled by combining the complete pre-datasets 1 and 2, as well as the dataset from Haniffa et al. (GEO: GSE35457, (Haniffa et al., 2012)). The second (dataset 2) contains the pre-dataset 1 without the B, T and NK cells, plus pre-dataset 2 and the dataset from Haniffa et al.. The raw expression values of both datasets were first quantile normalized and log2-transformed independently. Afterward, using the CD14+ blood MO as bridging samples, the batch effect introduced by the combination of the three different datasets was removed, again separately for datasets 1 and 2. The third dataset (dataset 3) was composed of CD14+ MO, Mo-M-CSF, Mo-GM-CSF and Mo-GM-CSFIL-4(0-72h) from pre-dataset 1. The fourth (dataset 4) equals the filtered pre-dataset 3. Both datasets were quantile-normalized and log2-transformed independently.
Primary data analysis and visualization
Each of the four datasets was filtered down to those probes being expressed across the corresponding dataset. To do so, the maximum of the group mean expression values was determined. Probes with a maximum expression lower than 7, again 7, 6.95 and 6.75 (all on log2 scale) were excluded as not expressed from datasets 1, 2, 3 and 4, resulting in 21,250, 23,592, 18,318 and 18,857 present probes, respectively. The samples of each dataset were always visualized via principal component analysis (PCA) using the determined present probes as input. For all datasets, principal components (PCs) 1 (x axis) versus 2 (y axis) were displayed, and for datasets 1, 2 and 4 also other combinations of PCs are shown. Furthermore, the top 1000 most variable probes of datasets 1, 3 and 4 were visualized in the form of a heatmap, where rows and columns were ordered based on hierarchical clustering using Euclidean distance and average linkage. Based on the revealed row-wise clusters, the genes were grouped together for datasets 3 and 4, dependent on the cell type(s) they appeared as highly expressed. A color code was added to indicate the corresponding group-related cell types. Additionally, a few genes of each cluster were highlighted next to the heatmap. For datasets 1, 3 and 4, the top 1000 most variable probes were also used as input to calculate a pearson correlation coefficient matrix (PCCM), which were displayed as heatmaps as well. Also there, rows and columns were ordered based on hierarchical clustering using Euclidean distance and average linkage.
Identification of signature genes
Using the normalized and batch corrected dataset 2, differentially expressed probes (Fold-Change ≤ −2 or > = 2 and FDR-adjusted p value < 0.05 using a two-way ANOVA including the batch as a random effect) were determined for the three comparisons Mo-GM-CSFIL-4(0-72h) and Mo-GM-CSFIL-4(0-144h) versus CD14+ MO, Mo-GM-CSF versus CD14+ MO and Mo-M-CSF versus CD14+ MO. The three lists of probes being upregulated against CD14+ MO were intersected, and the same was done for the probes being downregulated compared to CD14+ MO. The probes were then reduced to single genes, by keeping the probe for a corresponding gene with the highest mean expression across the dataset. This resulted in 287 up- and 361 downregulated genes in MO-derived cells compared to CD14+ MO. Both lists were further filtered by keeping only those genes, for which the mean expression values of the four MO-derived cell populations did not differ by a Fold-Change > 1.7 for all pairwise comparisons. This ensured similar expression values in the MO-derived cells, and resulted in 184 up- and 279 downregulated genes compared to CD14+ MO. The expression patterns of those genes were displayed in the form of a heatmap, where rows and columns were ordered based on hierarchical clustering using Euclidean distance and average linkage. A few gene names are depicted next to the heatmap.
Linear support vector regression (SVR)
The 8,612 present genes were used as input for CIBERSORT (Newman et al., 2015) to create a signature matrix, which contained 230 genes representing the five different cell types of the inflammatory Affymetrix dataset. After excluding all genes, which did not appear within the Illumina dataset 2, 197 genes remained. In addition, the 23,592 present probes of dataset 2 were reduced to 16,618 present genes by keeping the probe with the highest mean expression for a given gene. Then, SVR in CIBERSORT was applied to predict the fractions of the five cell types (in the form of the generated and filtered signature matrix composed of 197 genes) within the gene expression profiles (present genes) of dataset 2. To reduce complexity, the Mo-GM-CSFIL-4(0-72h) and Mo-GM-CSFIL-4(0-144h) were considered as a single group (named Mo-GM-CSFIL-4(0-72h/144h)) by calculating the mean expression across all replicates of the two groups. The resulting fractions were finally visualized in the form of stacked bar plots.
GSEA analysis
Using again the 8,612 present genes, DEG were determined between inflammatory DCs and CD14+CD16- MO, inflammatory DCs and CD14dimCD16+ MO, inflammatory Mac and CD14+CD16- MO, as well as between inflammatory Mac and CD14dimCD16+ MO by requiring a Fold-Change > 4 and an FDR-adjusted p value < 0.05 for all four comparisons using one way ANOVA. For each of the two comparisons inflammatory DCs versus CD14+CD16- MO and inflammatory Mac versus CD14+CD16- MO, specific genes were determined by excluding all genes appearing as differentially expressed in at least one of the other three comparisons. This resulted in 27 specific (not overlapping) genes for each of the two comparisons, respectively. Those two gene sets were then used as input for Gene Set Enrichment Analysis within PGS, using 1000 permutations and restricting the analysis to functional groups with more than 2 and fewer than 5000 genes. The enrichments were performed on the 16,618 present genes of dataset 2, for the pairwise comparisons CD141+ DC versus CD1c+ DC, Mo-GM-CSFIL-4(0-72h/144h) versus CD14+ monocyte, Mo-GM-CSF versus CD14+ monocyte, Mo-M-CSF versus CD14+ monocyte and CD45+lin-MHCIIhi versus CD14+ monocyte. Again here, all replicates from Mo-GM-CSFIL-4(0-72h) and Mo-GM-CSFIL-4(0-144h) were considered as a single group Mo-GM-CSFIL-4(0-72h/144h). The enrichment results were visualized within a scatterplot, with the normalized enrichment score (NES) on the x axis versus the FDR-corrected p value on the y axis.
Correlation of gene expression
As a third approach, the expression patterns of genes in the dataset of MO-derived cells were directly correlated to those in the dataset of inflammatory DCs and Mac. In both datasets, 6,959 common present genes (according to the gene symbols) were identified. By considering the results observed after the application of linear SVR, both datasets were divided into four comparable groups. Within the Illumina dataset 2 containing MO-derived cells, the four generated groups were: 1) Mo-GM-CSFIL-4(0-72h/144h), 2) CD1c+ and CD141+ DCs, 3) Mo-GM-CSF and Mo-M-CSF as well as 4) CD14+ MO and CD45+lin-MHCIIhi cells. The corresponding four groups in the Affymetrix inflammatory cell dataset were: 1) inflammatory DCs, 2) BDCA1+ DCs, 3) inflammatory Mac as well as 4) CD14+CD16- and CD14dimCD16+ MO. Finally, for each of the 6,959 genes a Pearson correlation value was calculated between the profile consisting of the four group mean expression values of Illumina dataset 2 and the corresponding profile of the Affymetrix dataset. As a first filtering step, only those genes were kept, which had a correlation value > 0.4. Within a second filtering step, certain Fold-Changes were required. For the Illumina MO-derived cell dataset, Fold-Changes were calculated for Mo-GM-CSFIL-4(0-72h/144h) compared to all other cell types of this dataset, and for the Affymetrix inflammatory cell dataset, Fold-Changes were determined for inflammatory dendritic cells compared to all other cell types within the dataset. For both datasets, only those genes were kept, for which all Fold-Changes were > 1.3. This resulted in 25 genes, which were specifically expressed in both Mo-GM-CSFIL-4(0-72h/144h) and inflammatory dendritic cells, but overall lower expressed in all other investigated cell types. The same approach was applied to find commonly expressed genes in inflammatory Mac as well as in Mo-GM-CSF and Mo-M-CSF. Fold-Changes were calculated for inflammatory Mac compared to all other cell types of the Affymetrix dataset. For each gene of the Illumina dataset, Fold-Changes were determined between the minimum of the two group mean expressions of Mo-GM-CSF and Mo-M-CSF and all other cell types (Mo-GM-CSFIL-4(0-72h/144h), CD1c+ and CD141+ DCs, CD14+ MO) of this dataset. Of the genes having a correlation score > 0.4, only 44 specific ones remained after excluding those, which did not have a Fold-Change > 1.3 for all comparisons within both datasets. Additionally, genes appearing in the final annotation file of 39,423 common probes of dataset 2 with multiple transcript variants - and therefore with potentially different expression patterns depending on the variant, of which only one was chosen based on the highest expression - were removed. The expression patterns of the remaining genes of both lists were visualized in the form of heatmaps, and the corresponding Pearson correlation values were displayed within bar plots next to the heatmaps.
Identification of subset specific genes
Using dataset 3, the 18,318 present probes were further reduced to 13,452 present genes by keeping the probe with the highest expression across the dataset for a given gene. To enable the identification of important surface markers, a comprehensive collection of important surface marker genes was assembled. For this, corresponding lists were downloaded from two different sources. First, a list of genes encoding membrane proteins was downloaded from the Human Protein Atlas, available at http://www.proteinatlas.org/search/protein_class:Predicted+membrane+proteins+AND+NOT+protein_class:Predicted+secreted+proteins. Second, a list of genes encoding surface proteins was obtained from SurfaceDB (de Souza et al., 2012), available at http://www.bioinformatics-brazil.org/surfaceome/home-. The union of both lists was generated and contained 5,582 gene symbols. Of those, 2,726 genes were found among the 13,452 present genes, which were then considered to represent the present surfaceome of dataset 3. To find surface markers being specifically expressed in one of the three MO-derived cell types (Mo-M-CSF, Mo-GM-CSF and Mo-GM-CSFIL-4(0-72h)), Fold-Changes for all pairwise comparisons including the ones against CD14+ MO were calculated. A marker was defined as condition specific, if all Fold-Changes against the three other cell types were larger than 2. This resulted in 19, 19 and 54 specific genes for Mo-M-CSF, Mo-GM-CSF and Mo-GM-CSFIL-4(0-72h), respectively. For visualization, those lists were further filtered. Genes belonging to the MHCI and II complexes (HLA-) were excluded since their specificity might be different in other individuals. Additionally, again here genes appearing in the final annotation file of 39,423 common probes with multiple transcript variants were removed. The remaining genes were displayed in the form of a heatmap, where the rows were ordered decreasingly based on the expression patterns in the three-different MO-derived cell types.
Co-expression network analysis
The union of genes being differentially expressed (Fold-Change > 2 or < −2 and FDR-adjusted p value < 0.05 using one way ANOVA) between each of the three types of MO-derived cells compared to CD14+ MO of dataset 3 were imported into BioLayout Express 3D version 3.3. Requiring a correlation of at least 0.93 to at least one other gene (based on the anti-log2-expression profiles), a co-expression network was generated. Smaller networks containing less than 5 genes were excluded, which resulted in a single network containing 2,086 genes. For each of the four cell types, the Fold-Change of the respective cell type compared to the overall mean was mapped onto the network and displayed in colors ranging from blue (negative Fold-Change) over white to red (positive Fold-Change) using the App Enrichment Map v2.0.1 in Cytoscape v3.2.0. Based on the Fold-Change patterns, the networks were divided into four main clusters, where each cluster represents one of the four cell types, respectively. Based on dataset 4, another co-expression network was generated to investigate the relationship of Mo-GM-CSFIL-4(72-144h) to the other MO-derived cell types and CD14+ MO. As input, 13,691 present genes were used, which were obtained by filtering the 18,857 probes in the same ways as previously described. The network was generated in BioLayout and finally visualized using Cytoscape. Two samples (nodes) are connected, if the correlation value calculated between their expression profiles was at least 0.9.
Prediction of master regulators
To predict master regulators for each of the four cell types, the Cytoscape App iRegulon v1.3 was applied on those genes of the cell type specific clusters identified within the co-expression network, which were upregulated in the corresponding cell type compared to the overall mean with a Fold-Change > 1.5. The prediction was performed by considering ”20kb centered around TSS” as ”Putative regulatory region,” by using” No motif collection” at ”Motif collection,” and by setting the “Enrichment Score threshold” to 1.5. The resulting lists were filtered afterward to keep only those predicted regulators, which reached a normalized enrichment score (NES) > = 3 and which were identified in only one of the four cell types after the previous filtering step.
Transcriptional regulator Co-expression network
Among the 13,452 present genes of dataset 3, 587 transcriptional regulators including transcription factors, co-factors and histone modifiers were identified (Fulton et al., 2009). Those were imported into BioLayout, and by requiring a correlation of at least 0.83 to at least one other regulator, a co-expression network was generated. After excluding smaller networks containing less than 4 regulators, a single network composed of 411 transcriptional regulators remained. For each of the four cell types, the Fold-Change of the respective cell type compared to the overall mean was mapped onto the network using Enrichment Map in Cytoscape. According to the resulting patterns, cell type specific clusters of upregulated regulators were generated. Additionally, those of the previously predicted master regulators were marked in these clusters, which were part of the network. Among the 13,691 present genes of dataset 4, 585 transcriptional regulators were found (Fulton et al., 2009). The expression profiles of those were imported into BioLayout, and the final co-expression network contained 267 regulators after requiring a correlation of at least 0.85 between a gene and at least one other gene and excluding resulting smaller networks containing less than 3 regulators. Within Cytoscape using EnrichmentMap, the Fold-Changes comparing Mo-GM-CSF versus CD14+ MO as well as Mo-GM-CSFIL-4(0-72h) compared to CD14+ MO were mapped onto the network, ranging from blue (Fold-Change < -2) over white to red (Fold-Change > 2). According to the resulting patterns, a Mo-GM-CSFIL-4(0-72h) specific cluster of upregulated regulators was identified. Additionally, regulators among the 585 present ones being specifically upregulated in Mo-GM-CSFIL-4(0-72/144h) were determined by requiring all of the following conditions: Fold-Change Mo-GM-CSFIL-4(0-72h) versus CD14+ Monocyte > 1.5, Fold-Change Mo-GM-CSFIL-4(0-144h) versus CD14+ Monocyte > 1.5, Fold-Change Mo-GM-CSFIL-4(0-72h) versus Mo-GM-CSF > 1.5, Fold-Change Mo-GM-CSFIL-4(0-144h) versus Mo-GM-CSF > 1.5 and Fold-Change Mo-GM-CSF versus CD14+ Monocyte > 1. The expression profiles of the remaining 7 transcriptional regulators were visualized in the form of a heatmap, were the rows were ordered decreasingly according to expression values in Mo-GM-CSFIL-4(0-72h). Those of the seven, which were also part of the co-expression network of transcriptional regulators, were highlighted within the network.
Expression profiles of pattern recognition receptors
Lists were collected of human gene symbols encoding Toll-like, NOD-like, RIG-I-like and C-type lectin receptors, as well as proteins belonging to the inflammasome. Only those of them were displayed in the form of barplots, for which the Fold-Change for at least one of the three MO-derived cell types compared to CD14+ MO was > = 2 or ≤ −2. The bars represent the Fold-Change between the mean values of each MO-derived cell type compared to the mean of CD14+ MO. Additionally, the Fold-Changes between the single replicate values compared to the mean of CD14+ MO are added as dots into the corresponding bar.
Comparison of Mo-M-CSF and Mo-GM-CSF cells
Mo-M-CSF and Mo-GM-CSF cells were compared both directly and indirectly. First, differentially expressed probes (Fold-Change > 2 or < −2 and FDR-adjusted p value < 0.05 using one way ANOVA) were determined between Mo-M-CSF and CD14+ MO as well as between Mo-GM-CSF and CD14+ MO. Probes were reduced to single genes by keeping the probe with the highest expression. The remaining lists of genes were intersected, which identified 348 or 379 genes being commonly up- or downregulated in Mo-M-CSF and Mo-GM-CSF compared to CD14+ MO, respectively. The intersections were displayed in the form of Venn Diagrams. Additionally, the Fold-Changes against CD14+ MO for the union of the 348 and 379 genes (y-axes) were displayed versus the ranks of the corresponding Fold-Changes (x-axes), once ranked according to the values for Mo-M-CSF and once according to the values for Mo-GM-CSF. Second, DEG were determined directly between Mo-M-CSF and Mo-GM-CSF (Fold-Change > 2 or < −2 and FDR-adjusted p value < 0.05 using one way ANOVA). After reducing the results to a single probe per gene as described before, 124 genes being upregulated in Mo-GM-CSF compared to Mo-M-CSF, and 73 genes being upregulated in Mo-M-CSF compared to Mo-GM-CSF, were identified. Both lists were used as input to perform Gene Ontology Enrichment Analysis (GOEA) in PGS, where 646 as well as 444 significantly enriched GO terms (Enrichment p value < 0.05) were identified for the Mo-GM-CSF as well as Mo-M-CSF specific genes, respectively. Both lists of significant GO terms were used as input for REVIGO, which clusters semantically similar GO terms together and generates a representative GO term name to facilitate the interpretation of long GO term lists. The identified representative GO term names were then visualized as Treemaps, where each rectangle represents one of the representative GO term names. Those rectangles were further joined into “superclusters” by REVIGO, indicated by different colors and a common super-representative GO term name. For simplicity, only the super-names were displayed.
Expression of genes linked to YG bead uptake, motility or OCR/ECAR
To link the experimental data back to the transcriptome data, lists of genes were generated, which had expression patterns being analogous to the observed functional outcomes between Mo-M-CSF, Mo-GM-CSF and Mo-GM-CSFIL-4(0-72h). For YG bead uptake, the pattern was generated by requiring a Fold-Change > 1.3 between Mo-M-CSF and Mo-GM-CSF > 1.3 as well as between Mo-GM-CSF and Mo-GM-CSFIL-4(0-72h). For migration, Fold-Changes > 1.3 were required between Mo-GM-CSFIL-4(0-72h) and Mo-GM-CSF, as well as between Mo-GM-CSF and Mo-M-CSF. For OCR/ECAR, Fold-Changes > 1.3 were required between Mo-GM-CSFIL-4(0-72h) and Mo-GM-CSF, and also between Mo-GM-CSFIL-4(0-72h) and Mo-M-CSF > 1.3, as well as a Fold-Change < 1.3 between Mo-GM-CSF and Mo-M-CSF. Using literature (references see manuscript), those lists were screened for genes, which were previously linked to either YG bead uptake, motility or OCR/ECAR, respectively. The resulting gene sets were displayed in the form of heatmaps, where rows were ordered based on hierarchical clustering using Euclidean distance and average linkage.
Visualization of cytokine responses
Protein expression values of 4 representative donors were averaged and displayed in the form of a heatmap using the function heatmap.2 of the R package gplots. Not detected values were set to 0. Rows were ordered alphabetically.
Handling of the RNA-Seq time kinetics dataset
The raw fastq-files were aligned against the human genome hg19 using TopHat v2.1.0 with default options. The resulting BAM files were imported into PGS, and mRNA quantification was performed against the hg19 RefSeq Transcript database version 2015-08-04 to obtain read counts for each individual RefSeq gene. This resulted in a read count table containing 20,359 genes. The dataset was then normalized using DESeq2 and the table of normalized read counts was imported back into PGS. Then, the technical variation introduced by the different donors was removed as batch effect. Finally, only 12,794 expressed genes were kept, which had a normalized mean read count of at least 10 in at least one of the eight investigated groups.
Analysis and visualization of RNA-Seq time kinetics dataset
Based on the 12,974 present genes, a PCA (PC1 versus PC2) as well as a self-organizing map (SOM) clustering were generated. Additionally, the top 1000 most variable genes across the dataset were visualized in the form of a heatmap, where rows and columns were ordered based on hierarchical clustering using Euclidean distance and average linkage. Based on the row-wise structure, the genes were grouped together into 5 clusters, of which a few genes were highlighted next to the heatmap.
Co-expression network for the time kinetics dataset
At first, all normalized read counts < 1 were set to 1 to avoid spurious Fold-Changes. Then, DEG (Fold-Change > 1.5 or < −1.5 and FDR-corrected p value < 0.05 using a two-way ANOVA including the donor as a random effect) were determined for each of the 7 conditions compared to the mean over all 7 conditions. The union of the resulting 7 lists comprised 3,114 genes, whose expression profiles were imported into BioLayout. A co-expression network was generated by requiring for each gene a correlation of at least 0.87 to at least one other gene. The network was finally composed of 2,775 genes, since smaller networks containing less than 4 genes were excluded. Using Cytoscape, for each cell type the Fold-Change of the respective cell type compared to the overall mean was mapped onto the network and displayed in blue (Fold-Change ≤ 1.5) or red (Fold-Change > = 1.5). Based on the Fold-Change patterns, for each condition a condition-specific cluster was generated by hand. Considering only those genes of the clusters, which had a Fold-Change > 1.5 for the corresponding condition compared to the overall mean, examples of genes being solely present in only one of the groups, or of genes being shared between the clusters of two consecutive time points, were highlighted.
IL-4 signature generation
To generate an IL-4 signature, three different publicly available datasets containing myeloid cells treated with IL-4 were downloaded and processed. From the first, CD14+ MO (Mo1-3) and Mo-GM-CSFIL4(0-120h) cells (Sample: 5d veh DC1-3) were used (GEO: GSE13762). The raw Affymetrix data of those samples was imported into PGS using RMA background correction, quantile normalization and log2-transformation. From the second, CD14+ MO (MO at time point 0h of donors # D1 - # D3) and MoIL4 (0-48h) (MO at time point 48h cultured in presence of IL-4 of donor # D1 - #D3) were considered (GEO: GSE35304). The raw Agilent data of those samples was imported into PGS using quantile normalization and log2-transformation. From the third, Mo-M-CSF(0-168h) (resting fully differentiated MAC (7 days of culture), rep 1-3) and Mo-M-CSFIL4(168-192h) (alternative activated MAC (M2) rep 1-3) were used (GEO: GSE32164). The raw Affymetrix data of those samples were imported into PGS using RMA background correction, quantile normalization and log2-transformation. Within each of the three datasets, the technical variation introduced by the donor was removed as batch effect. Then, differentially expressed probes (Fold-Change > 1.5 or < −1.5 and p value < 0.05 using a two-way ANOVA including the donor as a random effect) were determined within each dataset independently. Finally, after identifying genes being commonly upregulated in the IL-4 conditions of at least two datasets, and similarly being commonly upregulated in the untreated conditions of at least two datasets, 457 genes being induced by IL-4, and 498 genes being repressed by IL-4 were determined, respectively.
Handling of the siRNA αNCoR2 dataset
The raw fastq-files were also aligned against the human genome hg19 using TopHat v2.1.0 with default options. The resulting BAM files were imported into PGS, and mRNA quantification was performed against the hg19 RefSeq Transcript database version 2015-08-04 to obtain read counts for each individual RefSeq gene. This resulted in a read count table containing 19,384 genes. The dataset was then normalized using DESeq2 and the table of normalized read counts was imported back into PGS. Finally, only 12,816 expressed genes were kept, which had a normalized mean read count of at least 10 in at least one of the two investigated groups. All normalized read counts < 1 were set to 1 to avoid spurious Fold-Changes.
Analysis and visualization of the siRNA αNCoR2 dataset
Based on the 12,816 present genes, a PCA (PC1 versus PC2) was generated. Additionally, the expression profiles of the 1,834 genes being variable (p value < 0.05) across the dataset were visualized in the form of heatmap, where both rows and columns were ordered based on hierarchical clustering using Euclidean distance and average linkage. Based on the resulting row-wise structure, the genes were grouped into clusters, of which a few genes were highlighted next to the heatmap. The variable genes were also displayed within a scatterplot, displaying the log2-mean expression of the scrambled RNA samples (x axis) versus the log2-mean expression of the siRNA αNCoR2 samples. In red and blue were specifically highlighted, which of the genes were also part of the previously established signatures of IL-4 induced and repressed genes, respectively.
GSEA
To test for humanIL-4 signature enrichment in samples treated with either anti-NCOR2 or scrambled shRNA we performed GSEA. GSEA is a computational tool that determines whether a set of genes show statistically significant, concordant differences between two conditions (http://www.broadinstitute.org/gsea/index.jsp). We have used the normalized data table as input and either genes up- or downregulated in the IL-4 signature as gene sets for the GSEA. The ‘Signal2Noise’ ranking metric was used and the gene sets were chosen showing significant change at FDRc < c0.25, nominal P valuec < c0.05 and 1,000 gene set permutations.
Quantification and Statistical Analysis
Statistical analysis was performed using SigmaPlot (Systat Software Inc.). Statistical tests used are described in the according figure legend. Mean (±SEM) was indicated as horizontal lines. If not otherwise specified n represents number of biological replicates.
Data and Software Availabiltiy
All expression data related to this manuscript can be found at Gene Expression Omnibus under the accession number GSE96719.
Authors: Ping Yu; Yang Wang; Robert K Chin; Luisa Martinez-Pomares; Siamon Gordon; Marie H Kosco-Vibois; Jason Cyster; Yang-Xin Fu Journal: J Immunol Date: 2002-05-15 Impact factor: 5.422
Authors: Elodie Segura; Maxime Touzot; Armelle Bohineust; Antonio Cappuccio; Gilles Chiocchia; Anne Hosmalin; Marc Dalod; Vassili Soumelis; Sebastian Amigorena Journal: Immunity Date: 2013-01-24 Impact factor: 31.745
Authors: James Curtis; Yang Luo; Helen L Zenner; Delphine Cuchet-Lourenço; Changxin Wu; Kitty Lo; Mailis Maes; Ali Alisaac; Emma Stebbings; Jimmy Z Liu; Liliya Kopanitsa; Olga Ignatyeva; Yanina Balabanova; Vladyslav Nikolayevskyy; Ingelore Baessmann; Thorsten Thye; Christian G Meyer; Peter Nürnberg; Rolf D Horstmann; Francis Drobniewski; Vincent Plagnol; Jeffrey C Barrett; Sergey Nejentsev Journal: Nat Genet Date: 2015-03-16 Impact factor: 38.330
Authors: Chuan Li; Antoine Menoret; Cullen Farragher; Zhengqing Ouyang; Christopher Bonin; Paul Holvoet; Anthony T Vella; Beiyan Zhou Journal: JCI Insight Date: 2019-04-16
Authors: Inmaculada Serrano; Ana Luque; Francesca Mitjavila; Anna M Blom; Santiago Rodríguez de Córdoba; M Cristina Vega; Joan Torras; Josep M Aran Journal: Front Immunol Date: 2022-04-25 Impact factor: 8.786
Authors: Amy B Hont; Allison B Powell; Danielle K Sohai; Izabella K Valdez; Maja Stanojevic; Ashley E Geiger; Kajal Chaudhary; Ehsan Dowlati; Catherine M Bollard; Conrad Russell Y Cruz Journal: Mol Ther Date: 2022-02-09 Impact factor: 12.910
Authors: Jerome C Martin; Christie Chang; Gilles Boschetti; Ryan Ungaro; Mamta Giri; John A Grout; Kyle Gettler; Ling-Shiang Chuang; Shikha Nayar; Alexander J Greenstein; Marla Dubinsky; Laura Walker; Andrew Leader; Jay S Fine; Charles E Whitehurst; M Lamine Mbow; Subra Kugathasan; Lee A Denson; Jeffrey S Hyams; Joshua R Friedman; Prerak T Desai; Huaibin M Ko; Ilaria Laface; Guray Akturk; Eric E Schadt; Helene Salmon; Sacha Gnjatic; Adeeb H Rahman; Miriam Merad; Judy H Cho; Ephraim Kenigsberg Journal: Cell Date: 2019-08-29 Impact factor: 41.582
Authors: Xiaoying Zhou; Wong Yu; Elizabeth D Mellins; Kari C Nadeau; Shu-Chen Lyu; Claudia Macaubas; Bryan Bunning; Ziyuan He Journal: J Exp Med Date: 2021-05-04 Impact factor: 14.307
Authors: Rodolfo Nazitto; Lynn M Amon; Fred D Mast; John D Aitchison; Alan Aderem; Jarrod S Johnson; Alan H Diercks Journal: J Immunol Date: 2021-05-24 Impact factor: 5.426