Accelerated epigenetic aging in Huntington's disease involves polycomb repressive complex 1

Abstract

Loss of epigenetic information during physiological aging compromises cellular identity, leading to de-repression of developmental genes. Here, we assessed the epigenomic landscape of vulnerable neurons in two reference mouse models of Huntington neurodegenerative disease (HD), using cell-type-specific multi-omics, including temporal analysis at three disease stages via FANS-CUT&Tag. We show accelerated de-repression of developmental genes in HD striatal neurons, involving histone re-acetylation and depletion of H2AK119 ubiquitination and H3K27 trimethylation marks, which are catalyzed by polycomb repressive complexes 1 and 2 (PRC1 and PRC2), respectively. We further identify a PRC1-dependent subcluster of bivalent developmental transcription factors that is re-activated in HD striatal neurons. This mechanism likely involves progressive paralog switching between PRC1-CBX genes, which promotes the upregulation of normally low-expressed PRC1-CBX2/4/8 isoforms in striatal neurons, alongside the down-regulation of predominant PRC1-CBX isoforms in these cells (e.g., CBX6/7). Collectively, our data provide evidence for PRC1-dependent accelerated epigenetic aging in HD vulnerable neurons.

Fig. 1. Cellular identity-associated histone marks are specifically altered in striatal neurons of HD mice.

a Scheme illustrating FANS-ChIPseq and FANS-CUT&Tag experiments conducted on 15–20 week-old R6/1 and WT mice (Created in BioRender. Boutillier, A. (2025) https://BioRender.com/e61g067). H3K27ac and H3K27me3 FANS-ChIPseq, N = 2 biological replicates; H3K9ac FANS-CUT&Tag, N = 5 biological replicates (WT) and N = 4 biological replicates (R6/1); H3K18ac FANS-CUT&Tag, N = 2 biological replicates. b IGV genome browser capture showing H3K27ac and H3K27me3 signals in R6/1 and WT NeuN+ and NeuN- samples at super-enhancer-regulated SPN identity gene locus (e.g., Drd2). The grey box highlights Drd2 super-enhancer. Rep1, biological replicate 1; rep2, biological replicate 2. c Top, metaprofiles showing H3K27ac and H3K27me3 mean signals in WT and R6/1 NeuN+ and NeuN- nuclei, along SPN identity genes (left) and all genes (right) in mouse striatum. SPN identity genes correspond to module M2 in WGCNA study by ref. . Bottom, boxplots showing H3K27ac and H3K27me3 signal variations, expressed as percentages, in R6/1 vs WT at SPN identity genes and all genes. Statistical analysis was performed using Kruskal–Wallis test and Bonferroni correction for multiple testing. Boxplots show median, first quartile (Q1), third quartile (Q3) and range (min, Q1-1.5*(Q3-Q1); max, Q3+1.5*(Q3-Q1). Median values are indicated in blue. H3K27ac, N = 2 biological replicates in each group; H3K27me3, N = 2 biological replicates in each group d Volcano plots showing H3K27ac (left) and H3K27me3 (right) differentially enriched regions in R6/1 vs WT for NeuN+ (top) and NeuN- (bottom) samples. Decreased and increased regions in R6/1 vs WT are represented in blue and red, respectively (DESeq2 method, adj. P-val < 0.1 using the Benjamini–Hochberg method for multiple testing correction). e R6/1 NeuN+ H3K27ac (blue) and H3K27me3 (red) ChIPseq signals in downregulated (HD<WT), Non-significantly changed (n.s.) or upregulated (HD>WT) genes in dSPN (left) and iSPN (right) of HD R6/2 mice. Integration with RNAseq data by ref. (dSPN, N = 7 biological replicates in each group; iPSN, N = 4 biological replicates in each group). Boxplots show median, first quartile (Q1), third quartile (Q3) and range (min, Q1-1.5*(Q3-Q1); max, Q3+1.5*(Q3-Q1). Statistical analysis was performed using Mann–Whitney test (two-sided). H3K27ac, N = 2 biological replicates in each group; H3K27me3, N = 2 biological replicates in each group. Source data are provided as a source data file.

Fig. 2. Developmental genes are epigenetically de-repressed in striatal neurons of HD mice.

a Top 10 gene ontology (GO) processes (biological processes, BP) enriched in H3K27me3-depleted regions in R6/1 vs WT NeuN+ samples. BP terms are shown as a function of gene ratio, gene count and adj. Pval. Adj. Pval were calculated using the Benjamini–Hochberg method for multiple testing correction b Top 10 HD-related signatures (HDsigDB database) enriched in H3K27me3-depleted regions in R6/1 vs WT NeuN+ samples. HDsigDB terms are shown as a function of gene ratio, gene count and adj. Pval. Adj. Pval were calculated using the Benjamini–Hochberg method for multiple testing correction c Overlap between H3K27me3-depleted (HD < WT, blue) and H3K27ac-enriched (HD>WT, red) genes in R6/1 NeuN+ samples. 125 genes, representing 42% of H3K27me3-depleted genes, are enriched in H3K27ac. Overlapping genes include Onecut1, Pax6, Runx2, Cbx4 and Cbx8. d Top 10 gene ontology (GO) processes (biological processes, BP) enriched in 125 overlapping genes. BP terms are shown as a function of gene ratio, gene count and adj. Pval. Adj. Pval were calculated using the Benjamini–Hochberg method for multiple testing correction e Boxplots showing z-scores of expression values for the 125 overlapping genes in dSPN (left) and iSPN (right) in WT and HD R6/2 mice. RNAseq data in ref. . Boxplots show median, first quartile (Q1), third quartile (Q3) and range (min, Q1-1.5*(Q3-Q1); max, Q3+1.5*(Q3-Q1). Statistical analysis was performed using Kruskal–Wallis test and Bonferroni correction for multiple testing. dSPN, N = 7 biological replicates in each group; iPSN, N = 4 biological replicates in each group. f Heatmap representing R6/1 vs WT log2FC in H3K27me3, H3K27ac, H3K9ac and H3K18ac NeuN+ samples at the 125 H3K27me3-depleted and H3K27ac-enriched overlapping genes. g Top 10 gene ontology (GO) processes (biological processes, BP) enriched in H3K27ac, H3K9ac and H3K18ac-increased regions in R6/1 vs WT NeuN+ samples. BP terms are shown as a function of gene count and adj. Pval. Adj. Pval were calculated using the Benjamini–Hochberg method for multiple testing correction h IGV genome browser capture showing H3K27me3, H3K27ac, H3K9ac and H3K18ac signals in R6/1 and WT NeuN+ and NeuN- samples at Onecut1 gene locus. The grey box highlights Onecut1 promoter. i Protein-protein interaction network using STRING on 125 overlapping genes. The major subnetwork, which includes Onecut1, Pax6, Runx2 and Cbx4 is highlighted. PPI (protein protein interaction), statistics showing the strength of the network. j mRNA levels of Onecut1 in bulk striatal tissue of R6/1 and WT mice (RNAseq data by ref. ; N = 3 biological replicates in each group), in dSPN and iSPN of R6/2 and WT mice (RNAseq data by ref. ; dSPN, N = 7 biological replicates in each group; iPSN, N = 4 biological replicates in each group) and in striatal astroglia, dSPN and iSPN of Q175 knockin (HD-KI-Q175) and control (CT) mice (RNAseq data by ref. ; Astroglia, N = 10 biological replicates in each group; dSPN, N = 10 biological replicates in each group; iPSN, N = 10 biological replicates in each group). mRNA levels, reads per kilobases (RPK). Mean values +/- sem are shown. Statistics show adj. Pval, multiple testing correction was performed using the Benjamini–Hochberg method upon analysis of RNAseq data. Source data are provided as a source data file.

Fig. 3. Bivalent promoters are de-repressed in striatal neurons of HD mice.

a Heatmap of the 36,873 annotated mm10 RefSeq gene transcripts, integrating H3K4me3 and H3K27me3 gene profiles from striatal ChIPseq data, showing five distinct epigenetic profiles generated by kmeans clustering. TSS, transcription start site; TTS, transcription termination site. The red box highlights the cluster containing bivalent promoters, including Onecut1, Pax6 and Runx2. Top, Metaprofiles showing H3K27me3 (b) and H3K27ac (c) mean signals in WT (blue) and R6/1 (orange) NeuN+ and NeuN- samples at bivalent and all promoters. Bottom, boxplots showing H3K27ac and H3K27me3 signal variations, expressed as percentage, in R6/1 vs WT at bivalent and all promoters. Statistical analysis was performed using Kruskal–Wallis test and Bonferroni correction for multiple testing. Boxplots show median, first quartile (Q1), third quartile (Q3) and range (min, Q1-1.5*(Q3-Q1); max, Q3+1.5*(Q3-Q1). Median values are indicated in blue. H3K27me3, N = 2 biological replicates in each group; H3K27ac, N = 2 biological replicates in each group. d Left, Heatmap of Log2FC expression values of bivalent genes in R6/2 vs WT dSPN and iSPN (RNAseq data by ref. ). Right, overlap between upregulated bivalent genes in R6/2 vs WT in iSPN and dSPN. 49% of increased bivalent genes in R6/2 dSPN are increased in R6/2 iSPN. e Boxplot showing z-scores of expression values for bivalent genes in dSPN and iSPN in R6/2 and WT. RNAseq data by ref. ; dSPN, N = 7 biological replicates in each group; iPSN, N = 4 biological replicates in each group. Boxplots show median, first quartile (Q1), third quartile (Q3) and range (min, Q1-1.5*(Q3-Q1); max, Q3+1.5*(Q3-Q1). Statistical analysis was performed using Kruskal–Wallis test and Bonferroni correction for multiple testing. f Left, heatmap of Log2FC expression values of bivalent genes of HD patients vs control individuals in dSPN and iSPN (RNAseq data by ref. ). Right, zoom showing ONECUT1, RUNX2 and PAX6. Source data are provided as a source data file.

Fig. 4. PRC1-CBX proteins undergo paralog switch in striatal neurons of HD mice.

a IGV genome browser capture showing H3K27ac, H3K27me3, H3K9ac and H3K18ac signals in R6/1 and WT NeuN+ and NeuN- samples at Cbx2/4/8 genomic locus. The grey boxes highlight significantly H3K27me3-depleted and H3K27ac-enriched regions in R6/1 vs WT NeuN+ samples. The enhancer region between Cbx4 and Cbx8 was top significantly H3K27me3-depleted and H3K27ac-enriched region. b mRNA levels of Cbx2/4/6/7/8 in bulk striatal tissue of R6/1 and WT mice (RNAseq data and statistics by ref. ; N = 3 biological replicates in each group). mRNA levels, reads per kilobases (RPK). Mean values +/- sem are shown. Statistics show adj. Pval, multiple testing correction was performed using the Benjamini–Hochberg method upon analysis of RNAseq data. c Heatmap of Log2FC (R6/2 / WT) expression values of PRC1 Cbx paralogs in dSPN and iSPN (RNAseq data by ref. ). Statistics show adj. Pval, multiple testing correction was performed using the Benjamini–Hochberg method upon analysis of RNAseq data. *, adj. Pval < 0.05; **, adj. Pval < 0.01; ***, adj. Pval < 0.001. d Left, Immunoblots showing CBX4, CBX6, CBX8 and total proteins levels in the striatum of R6/1 and WT mice. Male and female samples are specified. MW, molecular weight. Right, bargraphs showing CBX4, CBX6 and CBX8 protein levels. Values were normalized to total protein. Mean values +/- sem are shown. Statistical analysis was performed using Mann–Whitney test (two-sided). N = 6 biological replicates in each group. e Heatmap of Log2FC PRC1-CBX protein level values in striatal tissue of HD-R6/1 vs WT mice of 1, 2 and 3 months (proteomic data, https://www.ebi.ac.uk/pride/archive/projects/PXD013771 and https://www.hdinhd.org/). f Representative images of CBX8 co-stained with neurons (NeuN) and DNA (DAPI) in the striatum region of 15 week-old WT and R6/1 mice (left panel), arrows indicate double-positive cells (CBX8/NeuN), with quantification of CBX8 intensity in neuron (right panel, N = 5 animals per group; WT, 4920 neurons count; R6/1, 5906 neurons count). Source data are provided as a source data file.

Fig. 5. H2AK119ub is depleted at subcluster of bivalent promoters in the striatum of HD mice.

a Scheme illustrating bulk ChIPseq and FANS-CUT&Tag experiments conducted on 15–20 week-old R6/1 and WT mice (Created in BioRender. Boutillier, A. (2025) https://BioRender.com/x74t675). H2AK119ub bulk ChIPseq, N = 2 biological replicates; H2AK119ub and H3K27me3 FANS-CUT&Tag, N = 3 biological replicates. b Left, kmeans clustering analysis of bivalent promoters using striatal H2AK119ub ChIPseq data identified H2AK119ub-high subcluster (cluster 1) and H2AK119ub-low subcluster (cluster 2). Top right, metaprofile showing H2AK119ub signals in R6/1 and WT samples at cluster 1, cluster 2 and all promoters. Bottom right, boxplots showing H2AK119ub signal variations, expressed as percentages, in R6/1 vs WT at cluster 1, cluster 2 and all promoters. Boxplots show median, first quartile (Q1), third quartile (Q3) and range (min, Q1-1.5*(Q3-Q1); max, Q3+1.5*(Q3-Q1). Median values are indicated in blue. Statistical analysis was performed using Kruskal–Wallis test and Bonferroni correction for multiple testing. c Top 15 gene ontology (GO) processes (biological processes, BP) enriched in cluster 1 and cluster 2. d Metaprofiles showing H3K27me3 (Top left) and H3K27ac (Bottom left) mean signals in WT (blue) and R6/1 (orange) NeuN+ ChIPseq samples at cluster 1, cluster 2 and all promoters. Boxplots showing H3K27me3 (Top right) and H3K27ac (Bottom right) signal variations in R6/1 vs WT at cluster 1, cluster 2 and all promoters. Boxplots show median, first quartile (Q1), third quartile (Q3) and range (min, Q1-1.5*(Q3-Q1); max, Q3+1.5*(Q3-Q1). Median values are indicated in blue. Statistical analysis was performed using Kruskal–Wallis test and Bonferroni correction for multiple testing. e Boxplots showing mRNA levels of genes in cluster 1 and cluster 2 in iSPN and dSPN of R6/2 and WT mice (RNAseq data by ref. ). Fold change (FC) of HD/WT values are shown. mRNA levels, reads per kilobases (RPK). Boxplots show median, first quartile (Q1), third quartile (Q3) and range (min, Q1-1.5*(Q3-Q1); max, Q3+1.5*(Q3-Q1). Statistical analysis was performed using Kruskal–Wallis test and Bonferroni correction for multiple testing. Source data are provided as a source data file.

Fig. 6. De-repression of developmental genes is an aging signature accelerated in striatal neurons of HD mice.

a Scheme illustrating temporal FANS-CUT&Tag experiments conducted on HD Q140 knockin (HD KI) mice and WT mice of 2, 6 and 10 months of age (Created in BioRender. Boutillier, A. (2025) https://BioRender.com/t46w032). H3K27ac and H3K27me3 FANS-CUT&Tag at 2 months, N = 3 biological replicates; H3K27ac and H3K27me3 FANS-CUT&Tag at 6 months, N = 3 biological replicates; H3K27ac and H3K27me3 FANS-CUT&Tag at 10 months, N = 2 biological replicates. b IGV genome browser capture showing H3K27ac and H3K27me3 signals at Onecut1 locus in HD KI and WT NeuN+ samples produced using 2, 6 and 10 month-old mice. c Plots showing z-score values for H3K27ac, H3K27me3 and euchromatin score for Onecut1 gene, in 2, 6 and 10 month-old HD KI and WT NeuN+ samples. Mean values +/- sem are shown. d mRNA levels of Onecut1 in bulk striatal tissue of 2, 6 and 10 month-old HD Q140 KI and WT mice (RNAseq data by ref. ; N = 10 biological replicates in each group). mRNA levels, reads per kilobases (RPK). Mean values +/- sem are shown. Statistics show adj. Pval, multiple testing correction was performed using the Benjamini–Hochberg method upon analysis of RNAseq data. e Left, Boxplot showing euchromatin score for cluster 1 genes (i.e., developmental transcription factors -TF-) in HD KI and WT NeuN+ samples that were generated using 2, 6 and 10 month-old mice (left). Boxplots show median, first quartile (Q1), third quartile (Q3) and range (min, Q1-1.5*(Q3-Q1); max, Q3+1.5*(Q3-Q1). Statistical analysis was performed using Kruskal–Wallis test and Bonferroni correction for multiple testing; all pairwise comparisons are statistically different (p < 0.05). Right, linear regression computed using cluster 1 mean euchromatin score shows distortion between chronological and epigenetic ages in striatal neurons of HD KI mice (green area), which reflects accelerated epigenetic aging. f Boxplots showing H3K27ac (top) and H3K27me3 (bottom) z-scores for cluster 1 genes in HD KI and WT NeuN+ samples generated using 2, 6 and 10 month-old mice. Boxplots show median, first quartile (Q1), third quartile (Q3) and range (min, Q1-1.5*(Q3-Q1); max, Q3+1.5*(Q3-Q1). Statistical analysis was performed using Kruskal–Wallis test and Bonferroni correction for multiple testing; pairwise comparisons are statistically different (p < 0.05), except when indicated (ns). g Plots showing z-score values of H3K27ac, H3K27me3 and euchromatin score for Cbx4 and Cbx8 genes, in 2, 6 and 10 month-old HD KI and WT NeuN+ samples. Mean values +/- sem are shown. h Mean mRNA levels of Cbx4 and Cbx8 in bulk striatal tissue of 2, 6 and 10 month-old HD Q140 KI and WT mice (RNAseq data by ref. ; N = 10 biological replicates in each group). mRNA levels, reads per kilobases (RPK). Mean values +/- sem are shown. Statistics show adj. Pval, multiple testing correction was performed using the Benjamini–Hochberg method upon analysis of RNAseq data. i Metaprofiles showing H3K27ac and H3K27me3 mean signals at cluster 1 gene promoters in HD KI vs WT NeuN+ samples prepared using 2, 6 and 10 month-old mice. j Above H3K27ac and H3K27me3 metaprofiles were grouped according to age. k Gene set enrichment analysis showing showing module normalized enrichment score (NES) identified by CEMiTool analysis using euchromatin score. M2 and M4 modules include Onecut1 and Glis1, respectively. l. Overlap between M2 (green) and cluster1 (blue) genes. Statistics of overlap (P) was assessed using a binomial test. Source data are provided as a source data file.

Fig. 7. Epigenetic regulation of stress response during aging is abnormal in HD mouse striatal neurons.

a Top 10 predicted transcriptional regulators (ChEA) enriched in M4 module. ChEA terms are shown as a function of gene ratio, gene count and adj. Pval. Adj. Pval were calculated using the Benjamini–Hochberg method for multiple testing correction b Network representation of functional enrichment analysis of M4 genes. Green dots represent M4 genes, blue dots, enriched regulatory transcription factors identified by ChEA analysis, pink dots, enriched BP and grey dots, enriched KEGG pathways. The red-green arrow highlights Glis1. c Boxplot showing z-scores of expression values of M4 genes in the striatum of HD Q140 KI and WT mice of 2, 6 and 10 months of age (RNAseq data by ref. ; N = 8 biological replicates in each group). Boxplots show median, first quartile (Q1), third quartile (Q3) and range (min, Q1-1.5*(Q3-Q1); max, Q3+1.5*(Q3-Q1). Statistical analysis was performed using Kruskal–Wallis test and Bonferroni correction for multiple testing. d. Boxplot showing z-scores of expression values of M4 genes in astroglia, dSPN and iSPN of 6 month-old HD Q175 KI (HD KI) and control (CT) mice. RNA seq data by ref. ; astroglia, N = 10 biological replicates in each group; dSPN, N = 10 biological replicates in each group; iPSN, N = 10 biological replicates in each group. Boxplots show median, first quartile (Q1), third quartile (Q3) and range (min, Q1-1.5*(Q3-Q1); max, Q3+1.5*(Q3-Q1). Statistical analysis was performed using Kruskal–Wallis test and Bonferroni correction for multiple testing. e Plots showing H3K27me3 z-score values at Glis1 gene, in 2, 6 and 10 month-old HD KI and WT NeuN+ samples. Mean values +/- sem are shown. f mRNA levels of Glis1 in bulk striatal tissue of 2, 6 and 10 month-old HD Q140 KI (HD KI) and control (CT) mice (RNAseq data by ref. ). Mean values +/- sem are shown. g mRNA levels of Glis1 in striatal astroglia, dSPN and iSPN of Q175 knockin (HD KI) and control (CT) mice (right, RNAseq data by ref. ). Mean values +/- sem are shown. h mRNA levels of Glis1 in dSPN and iSPN of HD R6/2 (HD R6/2) and control (CT) mice (right, RNAseq data by ref. ). mRNA levels, reads per kilobases (RPK). Mean values +/- sem are shown. Statistics show adj. Pval, multiple testing correction was performed using the Benjamini–Hochberg method upon analysis of RNAseq data. i. Left, immunoblots showing Glis1 proteins in the striatum of R6/1 and WT mice. Male and female animals are indicated. MW, molecular weight. Right, bargraph showing quantifications of Glis1 protein levels. Values were normalized to total proteins. Mean values +/- sem are shown. Statistical analysis was performed using Mann–Whitney test (two-sided). Source data are provided as a source data file.