Article Figures & Data
Figures
- Figure 1. Application of parsimonious gene correlation network analysis to time course gene expression data of total peripheral blood B-cell differentiation to plasmablast state.
(A) Network representation of the modular pattern of gene expression during the transition of B-cell to plasmablast. Module designations derived from gene ontology enrichment indicated with color code and ovals. Module genes are shown in Table S1, high-resolution version shown in Fig S1, and interactive version at https://mcare.link/abctopb. (B) Heat map summary representation of gene ontology and signature separation between network modules (filtered false discovery rate [FDR] < 0.05 and ≥5 and ≤1,000 genes; selecting the top 15 most significant signatures per module). Significant enrichment or depletion illustrated on red/blue scale, x-axis (signatures) and y-axis (modules). Hierarchical clustering according to gene signature enrichment. For high-resolution version and extended data, see Fig S2 and Table S2. (C) Overlay of gene expression z-scores for all genes in the network shown in blue (low) to red (high) z-score color scale. Day 0 (D0) provides the starting reference point for the sequential expression patterns observed at the subsequent time points indicated following decimal point for samples between D3 and D4. (D) Heat map displaying the pattern of gene expression across the time course module numbers indicated on the right, z-score gene expression blue (−1.8 low)–red (+1.8 high) color scale as indicated in the right lower edge, showing the median expression across three donors per time point. Modules divided into three broad categories of kinetics on at D0 going off, transient expression between D0 and D6, up-regulated at late time points.
- Figure S1. Accompanies Fig 1.
High-resolution version of Fig 1A, network representation showing color-coded modules and genes. Refer also to https://mcare.link/abctopb for an interactive searchable version including Z-score expression overlays as in Fig 1C. For individual gene expression values and module gene lists, see Table S1.
- Figure S2. Accompanies Fig 1.
High-resolution version of Fig 1B. Heat map of gene ontology and signature term enrichments linked to the parsimonious gene correlation network analysis modules of the time course network analysis (filtered FDR < 0.05 and ≥5 and ≤1,000 genes; selecting the top 15 most significant signatures per module). For full signature enrichment lists, please see Table S2. Modules are shown along the x-axis, and selected signature terms along the y-axis. Signature terms and modules are hierarchically clustered to illustrate relationships. Enrichment (red) and depletion (blue) of signatures are shown on the color scale of z-score.
- Figure 2. Kinetics of exemplar genes.
Violin plots of individual gene expression for selected genes showing expression z-score on the y-axis and time point in days and hours along the x-axis for the indicated genes above each graph. Violin plots display the distribution (n = 3 donors) along with median (blue square) and the inter-quartile range. (A, B, C, D, E) genes from M2 enriched for signaling response/immediate early genes, (B) M4, (C) M7, (D) M16 reflecting different patterns of cell cycle gene expression, and (E) core transcriptional regulators of the plasma cell state.
- Figure S3. Accompanies Fig 1.
Proliferation and matching phenotypic maturation during resting B-cells to plasmablast transition. (A) Left panels example density plots for CD19 versus CD20 and CD38 versus CD138 expression. Right panels corresponding CFSE dilution time course, cells were labelled at day 0 and proliferation and phenotypical maturation was assessed at day 0, day 3, and then monitored every 24 h up to day 6 in standard treatment conditions. Data representative of three biological replicates. (A, B) Summary of percentage cells with CD20 low (light grey) and CD38 high (dark grey) phenotype based on gates shown in (A). Data represent mean and SD of three independent donors. (C) Expression of CD30 at day 0, day 3, and day 6 against isotype control. Upper panel example histogram, middle panel ΔMFI, and bottom panel percentage CD30 high (light grey fill) versus CD30 low (dark grey fill) at each time point. Data representative of five independent donors. (D) ELISpot assessment of IgM (left) and IgG (right) secretory activity in day 6 plasmablasts. 2,000 cells were seeded per well, well diameter 6 mm (∼1 mm scale bar indicated). Data representative of two biological replicates.
- Figure 3. Application of parsimonious gene correlation network analysis to time course gene expression data of memory B-cell differentiation from activated B-cell to plasmablast state.
(A) Network representation of the modular pattern of gene expression during the transition of memory B-cell–derived activated B-cells to plasmablast state. Module designations derived from gene ontology enrichment indicated with color code and ovals. Module genes shown in Table S3, high-resolution version shown in Fig S4 and interactive version at https://mcare.link/abctopb. (B) Heat map summary representation of gene ontology and signature separation between network modules (filtered FDR < 0.05 and ≥5 and ≤1,000 genes; selecting the top 15 most significant signatures per module). Significant enrichment or depletion illustrated on red/blue scale, x-axis (signatures), and y-axis (modules). Hierarchical clustering according to gene signature enrichment. For high-resolution version and extended data, see Fig S5 and Table S4. (C) Overlay of gene expression z-scores for all genes in the network shown in blue (low) to red (high) z-score color scale. Day 3 (D3) provides the starting reference point for the sequential expression patterns observed at the subsequent time points indicated following decimal point for samples between D3 and D4. (D) Heat map displaying the pattern of gene expression across the time course module numbers indicated on the right, z-score gene expression blue (−1.8 low)–red (+1.8 high) color scale as indicated at the right lower edge, showing the median expression across three donors per time point. Modules divided into three broad categories of kinetics: (left) on at D3 going off, transient expression between D3 and D6, up-regulated at D6.
- Figure S4. Accompanies Fig 3.
High-resolution version of Fig 3A, network representation showing color-coded modules and genes. Refer also to https://mcare.link/abctopb for an interactive searchable version, including z-score expression overlays as in Fig 3C. For individual gene expression values and module gene lists, see Table S3.
- Figure S5. Accompanies Fig 3.
High-resolution version of Fig 3B. Heat map of gene ontology and signature term enrichments linked to the parsimonious gene correlation network analysis modules of the time course network analysis of memory B-cell differentiation (filtered FDR < 0.05 and ≥5 and ≤1,000 genes; selecting the top 15 most significant signatures per module). For full signature enrichment lists, please see Table S4. Modules are shown along the x-axis and selected signature terms along the y-axis. Signature terms and modules are hierarchically clustered to illustrate relationships. Enrichment (red) and depletion (blue) of signatures are shown on color scale of z-score.
- Figure 4. IRF4 and BLIMP1 occupancy in human plasmablasts.
(A) Relative distribution of BLIMP1 (left) and IRF4 (right) peaks identified in human memory B-cell derived plasmablasts divided according to genomic distribution, transcription termination site (TTS), promoter, exonic, intronic and intergenic as indicated in the color code to the right of the stacked bar graph (Promoter: −1 kb–100 bp, TTS: −100 bp–1 kb, Exonic/Intronic: > 100 bp from Promoter/TTS within gene, Intergenic: >1 kb from Promoter/TTS outside gene). ChIP data derive from individual samples for day 6. (B) Venn diagram depiction of BLIMP1 and IRF4 binding site overlap genome wide. (C) Relative genomic distribution and de novo motifs discovered at sites of BLIMP1-only, IRF4-only, and BLIMP1/IRF4 overlapping occupancy. (A) Shown is the genomic distribution as stacked bar graph color-coded as in (A) and the most significantly enriched motifs with percentage of peak regions with match to represented motif variant to the right. For each motif a summary designation is provided to the left, relating to a known motif match.
- Figure S6. Accompanies Figs 4–8.
Example data for each pattern of TF occupancy. For each instance, the designation of binding type is followed by relative distribution of peaks identified divided according to the genomic distribution (transcription termination site, promoter, exonic, intronic, and intergenic as indicated in the color code to the right of the figure), example tracks including all TFs (CTCF, H3K4me3, and H3K27ac) and de novo motif data. (A, B, C, D, E, F, G) Data are shown for the following: (A) BLIMP1 Alone (all BLIMP1 peaks); (B) D6.IRF4 Alone (all IRF4 peaks at day 6 plasmablasts); (C) BLIMP1 Only, BLIMP1 & IRF4 and IRF4 Only peaks (as in Fig 4); (D) CTCF; (E) XBP1; (F) D3.IRF4 Alone (all IRF4 peaks in day 3 activated B-cells); (G) D3.IRF4 Only, D3.IRF4 & D6.IRF4 and D6.IRF4 Only (as in Fig 8).
- Figure 5. Epigenetic patterns associate with core TF occupancy at the plasmablast state.
(A) deepTools heat map representation of K-means–clustered integrated ChIP-seq data from the plasmablast state. Data are clustered across the union of peaks for IRF4, BLIMP1, and XBP1 and encompassing data for CTCF, H3K4me3, and H3K27ac from equivalent cells. Six regulatory clusters are derived designated as U.K1-K6 on the left. (A, B) Relative distribution of K-means clusters U.K1-K6 derived from (A) according to the genomic distribution, transcription termination site, promoter, exonic, intronic, and intergenic as indicated in the color code to the right of the stacked bar graph. (A, C) Percentage occupancy of individual TF binding across the K-means clusters derived (A) for each of the TFs indicated by the color code to the right of the figure (BLIMP1-red, IRF4-green, and XBP1-blue).
- Figure S7. Accompanies Fig 5.
(A) deepTools heat map representation of K-means–clustered integrated ChIP-seq data from the plasmablast state. Data are clustered across the set of XBP1 peaks and encompassing data for CTCF, H3K4me3, and H3K27ac from equivalent cells. Six regulatory clusters are derived designated X.K1-K6 on the left. (A, B) Relative distribution of XBP1 K-means clusters X.K1-K6 derived from (A) according to genomic distribution, transcription termination site, promoter, exonic, intronic, and intergenic as indicated in the color code to the right of the bar graph. (A, C) Percentage occupancy of individual TF binding across the K-means clusters derived from (A) for each of the other TFs indicated by the color code to the right of the figure (BLIMP1-red and IRF4-green).
- Figure S8. Accompanies Fig 5.
(A) deepTools heat map representation of K-means clustered integrated ChIP-seq data from the plasmablast state. Data are clustered across the set of BLIMP1 peaks and encompassing data for CTCF, H3K4me3, and H3K27ac from equivalent cells. Six regulatory clusters are derived designated B.K1-K6 on the left. (A, B) Relative distribution of BLIMP1 K-means clusters B.K1-K6 derived from (A) according to the genomic distribution, transcription termination site, promoter, exonic, intronic, and intergenic as indicated in the color code to the right of the bar graph. (A, C) Percentage occupancy of individual TF binding across the K-means clusters derived (A) for each of the other TFs indicated by the color code to the right of the figure (IRF4-green and XBP1-blue).
- Figure S9. Accompanies Fig 5.
(A) deepTools heat map representation of K-means–clustered integrated ChIP-seq data from the plasmablast state. Data are clustered across the set of IRF4 peaks and encompassing data for CTCF, H3K4me3, and H3K27ac from equivalent cells. Six regulatory clusters are derived designated d6I.K1-K6 on the left. (A, B) Relative distribution of IRF4 K-means clusters d6I.K1-K6 derived from (A) according to the genomic distribution, transcription termination site, promoter, exonic, intronic, and intergenic as indicated in the color code to the right of the bar graph. (A, C) Percentage occupancy of individual TF binding across the K-means clusters derived (A) for each of the other TFs indicated by the color code to the right of the figure (BLIMP1-red and XBP1-blue).
- Figure 6. Integration of gene regulatory modules of the activated B-cell to plasmablast transition with TF occupancy patterns and epigenetic state.
Signature enrichment heat map displaying the enrichment/depletion of the memory B-cell parsimonious gene correlation network analysis expression modules (Fig 3) for genes associated with the TF peaks in the K-means clusters of epigenetic state (Fig 5). Significance of association between TF occupancy and genes belonging to a network module is shown as a z-score color scale (−5:blue to +5:red) divided according to hierarchical clustering of K-means modules (top) (z-scores with a P-value > 0.05 were set to 0). Results ordered left to right: BLIMP1 K-means clusters (B.All & B.K1-6), IRF4 K-means clusters (I.All & I.K1-6), XBP1 K-means clusters (X.All & X.K1-6), and Union K-means clusters (U.All & U.K1-6). Module identity is indicated to the right. Median expression pattern of the module across the time course illustrated as a z-score (−1.8: dark blue to 1.8: yellow) (left).
- Figure S10. Accompanies Fig 6.
Functional analysis of G9A inhibition with UNC0638 (A) WB of total H3K9me2 in activated B-cells (ABC) differentiating under standard conditions (-), with vehicle control DMSO, or after treatment at day 3 with UNC0638 at 2 μM concentration. Cells were sampled at the indicated time points. Upper blot: H3K9me2 and lower blot: total H3 loading control. (A, B) WB of LC3A/B-I and LC3A/B-II for samples treated as in (A) as indicated above the blot. Lower blot shows total H3 loading control. For each blot, one of five biological replicates is shown. (C) Flow cytometric analysis of autophagosomes with CYTO-ID assay (Enzo) at day 3 +24, +48, and +72 h in ABCs differentiating under standard conditions or in the presence of DMSO (vehicle) or UNC0638 (2 μM). Autophagy in standard, DMSO, and UNC0638 treatment conditions are represented in different shades of grey, from darker to lighter, respectively. Data are representative of three biological replicates from a single experiment. (D) Flow cytometric analysis of surface phenotype of differentiating populations at day 6 generated under standard conditions (standard) or in the presence of DMSO (vehicle) or UNC0638 (2 μM) plots showing CD19 versus CD20 (left) and CD38 versus CD138 (right). Percentages of total cell populations are shown following single cells and viability discrimination with gates established using matching isotype controls. One of six replicates from three independent experiments is shown. Summary mean flourescent intensity (MFI) data for all six replicates are shown to the right. (E) Representative density flow plots showing gating strategy for cell division generations determined from CFSE dilution at 24, 48, and 72 h of the ABC to plasmablast transition under standard conditions (left panel). Matching bar charts showing the frequency (% of cells/generation) at 24, 48, and 72 h of the ABC to plasmablast transition for standard, DMSO-, and UNC0638-treated conditions. Bars represent mean ± SD of three replicates from a single experiment. *P < 0.05 (paired t test) (right panel). (F) ELIspot results for IgM and IgG secretion for day 6 cells generated under standard, DMSO- (vehicle control), or UNC0638-treated conditions; 2,000 cells seeded per well, well diameter 6 mm (∼1 mm scale bar indicated). Representative wells are shown on the left and matching bar chart showing the spot counts on the right. Bars are mean ± SD of three replicates from two independent experiments.
- Figure 7. G9A inhibition with UNC0638 produces a focused impact on gene expression during plasmablast differentiation.
(A) Graphical representation of differential gene expression across a range of fold-change thresholds from 1.2+ to 2.0+ across a time course after UNC0638 treatment at D3 of memory B-cell differentiation. Number of differentially expressed genes across each fold-change threshold indicated by the color-coded graphical representation, according to the color code on the right of the figure (y-axis: number of significantly differentially expressed genes, x-axis: time point in hours and days). Left graph: genes up-regulated in the absence of inhibitors, right graph: genes up-regulated in the presence of inhibitors. Data derived from three independent donors. (B) The overlap of BLIMP1 occupancy in memory-derived plasmablasts in the absence (blue) or presence (brown) of UNC0638 treatment. (C) deepTools heat maps of K-means clusters derived from the union of BLIMP1 binding sites for standard and UNC0638 conditions and considering associated epigenetic marks as indicated for H3K4me3, and H3K27ac. (D) Heat map of genes up-regulated upon G9A inhibition (fold change > 1.8, FDR < 0.05) showing patterns of gene expression as z-scores (−1.8: dark blue to 1.8: yellow) across the differentiation in the absence or presence of UNC0638 treatment. To the right are the genes identified as BLIMP1 bound highlighted with red bars. (A, E) Dumbbell graph of the relative enrichment or depletion of differentially expressed genes shown in (A) against the modules of gene expression derived from the memory B-cell parsimonious gene correlation network analysis network in Fig 3. Y-axis shows the order of modules ranked from most significantly enriched for genes up-regulated in the presence of UNC0638 through most significantly enriched in the standard conditions. For each module, the enrichments or depletion are shown for the genes up-regulated in the presence of UNC0638 (brown circles) and genes up-regulated in standard conditions (blue circles). These are plotted against the x-axis displaying Z-score of enrichment/depletion with the vertical dotted red lines indicating the point of FDR corrected significance (P-value < 0.05). (D, F) Representative tracks for BLIMP1 ChIP-seq in standard and UNC0638-treated samples, as indicated for representative genes selected from (D) whose expression is increased in the presence of UNC0638.
- Figure 8. Differential IRF4 occupancy and gene network associations between the activated B-cell (ABC) and plasmablast states.
(A) Venn diagram of overlap in IRF4 occupancy between ABCs (D3.IRF4 light green) and plasmablasts (D6.IRF4 darker green). (A, B) Comparison of known motif enrichments at sites bound by IRF4 as shown in part (A) either in ABCs only (D3.IRF4), ABCs and plasmablasts (D3.D6.IRF4), or plasmablasts only (D6.IRF4) and compared with IRF4 occupancy in cell lines representative of transformation at the ABC state (OCILY_Comb derived from OCI-LY3 and OCI-LY10 ABC-DLBCL lines) and plasma cell myeloma (MM_Comb derived from H929 and U266 malignant myeloma cell lines). Enrichment of known motifs indicated across the top is illustrated as percentage of sites with motif match (circle diameter–top right) and heat map color code (dark blue to yellow, ceiling set at −logP ≥ 2,000, bottom right). (C) deepTools heat maps of K-means clusters derived from D3.IRF4-bound regions alongside D6.IRF4 and sites bound by BLIMP1, CTCF, H3K27ac, and H3K4me3 in day 6 plasmablasts. (C, D) Integration of gene regulatory modules of the ABC to plasmablast transition with TF occupancy patterns and epigenetic state as in Fig 6 but with the added inclusion of K-means clusters for D3.IRF4-occupied regions alone, and the union of all occupied regions including D3.IRF4 (Union 2) as in (C). Heat map displays the enrichment/depletion of TF peaks in the K-means clusters relative to the memory B-cell parsimonious gene correlation network analysis expression modules. Significance of TF occupancy versus genes belonging to a network module is shown as a z-score color scale (−5: blue to +5: red) divided according to the hierarchical clustering of K-means modules (top) (z-scores with a P-value > 0.05 were set to 0). Module identity is indicated on the right, and median expression pattern of the module is shown across the time course as a z-score on the left (−1.8: dark blue to 1.8: yellow).
- Figure S11. A dichotomy in gene regulatory associations during the activated B-cell to plasmablast transition.
The upper panel briefly addresses the differentiation window from ABC to plasmablast. This is initiated by activation of B-cells with a combination of CD40L, B-cell Receptor (BCR) ligation and IL2/IL21 stimulation. Release from CD40L and BCR stimuli, and continued IL2/IL21, supports the transition to the plasmablast state. The associated shift in gene expression is shown in the context of PGCNA networks. The underlying regulatory transition is summarised in the lower panel. In ABCs IRF4 occupancy in the context of AICE and EICE regulatory elements dominates with presumptive partner factors being BATF at the former and SPIB and PU.1 at the latter. In plasmablasts IRF4 occupancy is associated with ISRE and IRF elements as well as EICEs. Loss of AICE occupancy parallels loss of BATF expression and is linked with BLIMP1 occupancy at the BATF promoter. BLIMP1 occupancy at the SPIB promoter and loss of SPIB expression is also established at this time-point. Thus, BLIMP1 may contribute to the shift in IRF4 occupancy by supporting repression of two of its key partner factors. At the same time in plasmablasts, and in myeloma cell lines, IRF4 acquires a new association with the occupancy pattern of CTCF, suggesting a potential role in chromatin looping. The third key transcription factor of the plasmablast state, XBP1, is highly enriched at regulatory elements of the secretory pathway genes induced at the plasmablast state. While XBP1 and IRF4 occupancy links to modules of genes characteristic of the plasmablast/plasma cell state and secretory pathway, BLIMP1 is associated with genes modules that are repressed in plasmablasts. Thus, a dichotomy in gene regulatory associations is established with BLIMP1 linking to repression of genes defining the activated B-cell state, while IRF4 and XBP1 link to expression of characteristic features of human plasmablasts.
Supplementary Materials
Supplemental Data 1.
A supplemental methods document providing further details for gene expression analysis, PGCNA networks, ChIP-seq analysis, data processing and data and software availability.[LSA-2020-00654_Supplemental_Data_1.pdf]
In this Issue
Subjects
Related Articles
Cited By...
- No citing articles found.