The single-cell epigenomic and transcriptional landscape of immunity to influenza vaccination

Abstract

Emerging evidence indicates a fundamental role for the epigenome in immunity. Here, we mapped the epigenomic and transcriptional landscape of immunity to influenza vaccination in humans at the single-cell level. Vaccination against seasonal influenza induced persistently diminished H3K27ac in monocytes and myeloid dendritic cells (mDCs), which was associated with impaired cytokine responses to Toll-like receptor stimulation. Single-cell ATAC-seq analysis revealed an epigenomically distinct subcluster of monocytes with reduced chromatin accessibility at AP-1-targeted loci after vaccination. Similar effects were observed in response to vaccination with the AS03-adjuvanted H5N1 pandemic influenza vaccine. However, this vaccine also stimulated persistently increased chromatin accessibility at interferon response factor (IRF) loci in monocytes and mDCs. This was associated with elevated expression of antiviral genes and heightened resistance to the unrelated Zika and Dengue viruses. These results demonstrate that vaccination stimulates persistent epigenomic remodeling of the innate immune system and reveal AS03's potential as an epigenetic adjuvant.

Trial registration: ClinicalTrials.gov NCT01910519 NCT02154061.

Keywords: adjuvant; antiviral immunity; epigenome; influenza; innate memory; monocyte; single cell; systems biology; trained immunity; vaccines.

Figure 1 -. TIV alters the global histone modification profile of immune cells.

(A) Study overview. (B) UMAP was used to create a dimensionality-reduced representation of the global histone mark profiles of all immune cell subset. (C) UMAP was used to visualize epigenomic profiles at the sample level. (D, E) The effect size of vaccine-induced changes to the global histone modification profiles at day 30 vs day 0 were calculated. D) Top-10 most significantly increased and reduced histone modifications. E) Heatmap showing histone modification changes in innate immune cells. Only changes with an FDR <= 20% are shown. (F) Change in histone modification levels relative to day 0 before vaccination for a set of highly reduced histone modifications in C monos and mDCs. Dots and lines indicate average modification levels, error bars indicate the standard error of mean. (G) UMAP representation of single monocytes and mDCs using H2BK5ac, H3K37ac, H3K9ac, H4K5ac, and PADI4. Left panel: cell density at each time point, right panel: H3K27ac levels in each single cell. Red ellipses indicate high-density areas corresponding to bright areas in left panel. See also Figure S1, S2, S3

Figure 2 -. TIV-induced histone modification changes correlate with cytokine production.

(A) Schematic overview of experiment. (B) Heatmap showing the relative change in cytokine levels at indicated time points compared to day 0. (C) Cytokine levels before (d0) and after (d30) vaccination for each investigated subject. Wilcoxon signed rank test was used for hypothesis testing. * p <= 0.05, ** p <= 0.01, *** p <= 0.001, **** p <= 0.0001, n = 18–19 (D) Change in cytokine levels relative to day 0 for cytokines in C. Dots and lines indicate average. (E, F) Pearson correlation of the cytokine levels of the 10 cytokines in C) with histone modification levels in C monos as well as C mono frequency in PBMCs as determined by EpiTOF and sample viability. (n = 87 samples from all time points) E) Boxplots of correlation coefficients for each cytokine after stimulation with either viral or bacterial cocktail. F) Scatter plots for the indicated histone modifications and cytokines. G) Gating scheme showing the production of IL-1b and TNFa in C monos after indicated treatment. H) Boxplot summary of the fraction of IL-1b+ or TNFa+ cells in multiple donors. Wilcoxon rank sum test, * p <= 0.05, ** p <= 0.01, *** p <= 0.001, **** p <= 0.0001, n = 4–11. See also Figure S4

Figure 3 -. TIV induces reduced chromatin accessibility in immune response genes and AP-1 controlled regions.

(A) Schematic overview of the experiment. (B) Differentially accessible chromatin regions (DARs) at day 30 vs day 0 were identified using DESeq2. Pval <= 0.05. (C) Heatmap representation of the normalized accessibility at the top 200 as well as cytokine-associated DARs in C monos for each analyzed sample. p: promoter −2000 bp to +500 bp; d: distal −10kbp to +10kbp – promoter; t: trans < −10kbp or > +10kbp. (D) Network representation of gene set enrichment analysis of DARs in C monos using the Reactome database. Only significantly enriched terms (p <= 0.05) are shown. Color indicates whether majority of enriched regions showed enhanced (red) or reduced (blue) accessibility. Heatmaps show signed –log10(pval) for significantly enriched terms in highlighted clusters. (E) Motif-based overrepresentation analysis of transcription factor binding sites in DARs at day 30 vs day 0. (F) Scatter plot showing the change in TF gene expression (x-axis) plotted against the enrichment in DARs for selected transcription factors in the Encode database. Blue color indicates AP-1 members with significantly reduced expression. (G) Change in gene expression of AP-1 family members using bulk transcriptomics data from 3–9 independent flu vaccine trials previously conducted. Heatmap indicates average log2 fold change in gene expression over all trials. N indicates subject and study number at each time point. Wilcoxon signed rank test, * p <= 0.05, ** p <= 0.01, *** p <= 0.001. (H) DARs in indicated cell type were correlated with H3K27ac levels as measured by EpiTOF and DARs with significant correlation (p <= .05) were analyzed for transcription factor target gene enrichment using the Encode database. Blue color indicates significantly changed AP-1 members. (I) Histogram showing the level of phospho-c-Jun in C monos in the indicated conditions. (J) Box plot summary of the fraction of phospho-c-Jun positive cells in classical monocytes. Wilcoxon rank sum test, * p <= 0.05, ** p <= 0.01, *** p <= 0.001, **** p <= 0.0001, n = 4–11 See also Figure S5

Figure 4 -. Heterogeneity within monocyte population drives TIV induced epigenomic changes.

(A) Schematic overview of the experiment. (B) UMAP representation of scATAC-seq landscape after pre-processing and QC filtering. (C) Heatmap showing differences in chromatin accessibility at indicated time points for the top5 transcription factors per subset. (D) UMAP representation of epigenomic subclusters within the classical monocyte population. (E) Density plot showing the relative contribution of different epigenomic subclusters to the total monocyte population at a given vaccine time point. (F) Variability in TF accessibility within the monocyte population. Value indicates range of accessibility values in all single monocytes. (G) Heatmap showing differences in chromatin accessibility between monocyte subclusters subset. (H) UMAP representation of monocyte subclusters showing differences in AP-1 accessibility. (I) UMAP representation of monocyte subclusters showing difference in accessibility at Hotspot module 2,3 gene loci. (J) Enrichment analysis of genes associated with loci in Hotspot module 2,3. (K) UMAP representation of the transcriptional landscape of single monocytes. Color indicates expression of genes associated with Hotspot modules 2,3.

Figure 5 -. H5N1+AS03 induces repressive epigenomic state akin to TIV.

(A) Schematic overview of experiment. (B) UMAP representation of EpiTOF landscape. (C) Histone modification levels in classical monocytes at day 0 vs day 42 as measured by EpiTOF. (D) Cytokine levels in supernatant of TLR-stimulated PBMCs at day 0 and day 42 after vaccination with H5N1+AS03. (C, D) Wilcox signed rank test; *p <= 0.05, ** p <= 0.01; EpiTOF: n = 9/9, Luminex: n = 13. (E) UMAP representation of scATAC-seq (left) and scRNA-seq (right) landscape after pre-processing and QC filtering. (F) Change in accessibility of detected AP-1 family TFs in classical monocytes. Color indicates whether cells are derived from subjects vaccinated with H5N1 (green) or H5N1+AS03 (orange). (n = 2/2) (G) Overrepresentation analysis of significantly different DARs in classical monocytes using the Reactome database. Color indicates whether enriched genes were predominantly up- or down-regulated. (H) Volcano plot showing changes in expression of AP-1 TF genes in classical monocytes at D42 compared to D0. See also Figure S6

Figure 6 -. H5N1+AS03 induces epigenomic state of enhanced antiviral immunity.

(A) Heatmap showing the change in chromatin accessibility at day 42 vs day 0 for the top5 transcription factors per subset. Color indicates the difference in accessibility, grey fields indicate non-significant changes (fdr > 0.05). (B) Line graph showing the difference in transcription factor (TF) accessibility during vaccination. Each line represents a separate TF within the indicated family. (C) Volcano plot showing the change in gene expression for IRF/STAT TF genes. (D) Scatter plot showing chromatin accessibility values for IRF1 (x-axis) and FOS (y-axis) in single cells. Indicated statistics are based on Pearson correlation. (E) MA plot showing the average accessibility and log2(FC) accessibility for genomic regions containing an IRF1 binding motif. Red color indicates regions with significantly changed accessibility (P <= 0.05). (F) Gene set enrichment analysis of significantly changed regions in E) occurring in at least 5% of C monos using the Reactome database. (G) Interferon gamma and IP10 levels in plasma of vaccinated subjects. Dots and lines indicate average, ribbons indicate standard error of mean. (H5N1/H5N1+AS03: IFNG, n = 7/14; IP10, n = 16/34) (H) Scatter plot showing changes in chromatin accessibility (x-axis) and changes in gene expression (y-axis) at day 21 vs day 0 for C monos (scATAC P <= 0.05 and occurring in at least 5% of cells). Indicated statistics are based on Pearson correlation analysis and Chi-square test. (I) Change in gene expression for selected antiviral and interferon-related BTMs in bulk RNA-seq analysis for subjects vaccinated with H1N1 (green) and H1N1+AS03 (orange) at indicated time points. (H1N1: n = 16, H1N1+AS03: n = 34). (J) Scatter plot showing the change in chromatin accessibility at day 21 vs day 0 in C monos (x-axis) and the significant change (p <= 0.05, log2(FC) > +/−0.03) in vaccine-induced gene expression at the booster vaccination compared to the prime vaccination (y-axis, Day22day21 vs Day1day0). Chi-square test was used to determine whether both variables were related. (K) Bubble plot showing enrichment results using the Encode TF target gene database. Color indicates the origin of the analyzed genes in J). See also Figure S7

Figure 7 -. H1N1+AS03 induces enhanced resistance to in-vitro infection with heterologous viruses.

(A) Schematic overview of the experiment. (B) Boxplot showing viral titers in Dengue-, Zika-, and mock-infected samples. (C) Line graph showing the viral growth curve for Dengue virus (red) and Zika virus (blue). Dots and lines indicate average, error bars indicate standard error of mean. n > 21 samples (D) Log2 fold change in viral titers relative to day 0 before vaccination. Wilcoxon signed rank test was used to compare changes within group; ** p <= 0.01, *p <= 0.05, n = 8–9. (E) Boxplot showing the concentration of IFNa, IFNg, and IP10 in Dengue-, Zika-, and mock-infected cultures at 24h after incubation. Wilcoxon rank-sum test was used to compare groups. (F, G) Pearson correlation analysis of the change in viral titers (d0 vs d21) with change in vaccine-induced, in-vivo expression of enhanced antiviral genes at prime (d0 vs d1) and boost (d21 vs d22) (red genes Figure 6G). F) Boxplot showing correlation coefficient per viral condition. G) Scatter plot showing change in vaccine-induced expression of IRF1 (x-axis) and viral titers (y-axis). (H) Model of bi-directional epigenomic reprogramming. B, E) Wilcoxon rank sum test was used to compare groups.