Dynamic profiling and functional interpretation of histone lysine crotonylation and lactylation during neural development

The N-terminal tails of histones are subjected to various post-translational modifications (PTMs), which shape the local chromatin landscape to regulate transcription, replication, DNA repair and higher-order chromatin structure (Diehl and Muir, 2020; Klemm et al., 2019). Some metabolites including propionyl-CoA, butyl-CoA, crotonyl-CoA, lactyl-CoA, 2-hydroxyisobutyl-CoA, β-hydroxybutyryl-CoA, succinyl-CoA, benzoyl-CoA, malonyl-CoA and glutaryl-CoA have recently been identified as donors for histone PTMs to influence gene expression (Diehl and Muir, 2020; Li et al., 2018). Among these, crotonyl-CoA and lactyl-CoA act as substrates for histone lysine crotonylation (Kcr) and lysine lactylation (Kla), respectively, catalyzed by CBP (also known as CREBBP) or p300 (or EP300) (Sabari et al., 2015; Zhang et al., 2019). Histone Kcr and Kla have attracted great attention due to their emerging roles in linking cellular metabolism and epigenetic regulation to orchestrate gene expression, cell state transition and disease progression (Li and Wang, 2021; Li et al., 2020; Yu et al., 2021; Zhang et al., 2019).

Kcr and Kla can occur on all core histone proteins, and they share the most common modification sites with histone lysine acetylation (Kac) (Tan et al., 2011; Zhang et al., 2019). Over the last five years, regulatory mechanisms for histone Kcr have been thoroughly elucidated, including its ‘writers’ (enzymes that catalyze Kcr) (Kollenstart et al., 2019; Liu et al., 2017c; Sabari et al., 2015), ‘erasers’ (enzymes responsible for removing Kcr) (Fellows et al., 2018; Kelly et al., 2018; Wei et al., 2017), ‘readers’ (enzymes capable of recognizing Kcr) (Andrews et al., 2016; Li et al., 2016; Xiong et al., 2016; Zhao et al., 2016) and metabolic regulators (enzymes involving in crotonyl-CoA metabolism) (Fang et al., 2021; Jiang et al., 2018; Liu et al., 2017b, 2019; Sabari et al., 2015; Tang et al., 2021), which extend our knowledge about the versatile functions of histone Kcr. Recent studies have revealed that p300 and glycolysis are major pathways to regulate histone Kla levels, and that H3K18la is involved in M1 macrophage polarization, cell reprograming and tumorigenesis (Li et al., 2020; Yu et al., 2021; Zhang et al., 2019). Both histone Kcr and Kla can directly stimulate transcription in a cell-free system (Sabari et al., 2015; Zhang et al., 2019). Specifically, recognition of H3K9cr, H3K18cr and H3K27cr by the YEATS domain of AF9 (MLLT3) may regulate the transcription process by recruiting the super elongation complex (SEC) and the Dot1L complex (DOT1LC) (Li et al., 2016, 2014; Luo et al., 2012). However, it is unclear whether histone Kcr and Kla are causal or consequential for transcription under physiological conditions.

Genome-wide associations between histone Kcr and gene expression have been investigated in various organisms and consistently suggest that histone Kcr is a widely distributed mark of active chromatin (Crespo et al., 2020; Fellows et al., 2018; Gowans et al., 2019; Kelly et al., 2018; Lu et al., 2018; Sabari et al., 2015; Tan et al., 2011). In addition, histone Kla is associated with inflammation and oncogenesis-related transcriptional programs (Yu et al., 2021; Zhang et al., 2019). However, it remains unknown whether and how histone Kcr and Kla are involved in transcriptional regulation during mammalian development. Considering the important roles of histone Kac in neurodevelopment and neurological diseases (Li et al., 2019; Tang et al., 2019; Wang et al., 2010; Yao et al., 2016), the associations of histone Kcr and Kla with gene expression and their biological functions merit further investigation. Moreover, there is growing evidence that histone deacetylases (HDACs) regulate multiple histone deacylation reactions as well as deacetylation (Huang et al., 2018, 2021; Wei et al., 2017), but little is known about how different histone acylations are coordinated to regulate gene expression and biological processes when HDACs are perturbed.

In this study, we systematically investigated the dynamic changes of histone Kcr and Kla and their associations with gene expression, as well as the interplays among histone acylations during neural development. We observed the global existence and significant changes of histone Kcr and Kla in the developing telencephalon. Further multi-omics profiling uncovered that H3K9cr and H3K18la mark transcriptionally active chromatin both in vivo and in vitro. There are distinct patterns for genome-wide alterations in H3K9ac, H3K9cr and H3K18la, but integrated analysis of chromatin immunoprecipitation-sequencing (ChIP-seq) and RNA-sequencing (RNA-seq) data across different developmental stages reveals that these histone marks are extensively involved in transcriptome remodeling to favor neural differentiation. In addition, inhibition of transcription does not influence global histone Kcr and Kla levels. Finally, we discovered that inhibition of HDAC1-3 promotes neuronal transcriptional programs by multiple histone acylations.

We first determined the presence of Kcr and Kla in the H3 N-terminal tails during neural development, and found that histone Kcr and Kla were extensively present in the embryonic cortex (Fig. 1A). Kcr and Kla could be observed in TUJ1⁺ (TUBB3) neurons in the intermediate zone and cortical plate regions, and TUJ1⁻ neural progenitor cells in the ventricular/subventricular zone (Fig. 1A). By quantifying the enrichment level of these modifications across the layers of the cortex, we found that H3K9cr and H3K18la levels were significantly higher in the cortical plate regions than in the ventricular/subventricular zone (Fig. 1B; Fig. S1A). Next, we examined the levels of these histone marks across different developmental stages (Fig. 1C) and found that H3K9cr levels were increased but H3K18la and overall H3-Kla levels were reduced at the late stage of neurogenesis compared with those at the early stage of neurogenesis; however, there were no apparent changes in H3K18cr, H3K14la and overall H3-Kcr levels (Fig. 1D-F).

We then checked changes in these histone marks at the cellular level (Fig. S1B) and observed that the levels of H3K18cr and H3K18la, but not those of H3K9cr and H3K14la, declined sharply during neural differentiation of P19 embryonic carcinoma (EC) cells (Fig. S1C,D), indicating different dynamics of site-specific histone acylations during cell-fate transitions. Altogether, these results suggest that histone Kcr and Kla are widely distributed and undergo global changes during neural development.

In order to explore the chromatin function of histone Kcr and Kla during neural development, we performed specific anti-H3K9ac, anti-H3K9cr and anti-H3K18la (lactylation is not found at the H3K9 site in mouse cells) ChIP-seq assays by using highly specific antibodies (Fig. S2) to analyze the telencephalon at embryonic stage (E) 13.5 (Zhang et al., 2019), and obtained highly reproducible sequencing data (Fig. S3A). First, H3K9ac and H3K9cr were mainly located in promoter regions, whereas H3K18la was mostly enriched in introns and intergenic regions (Fig. 2A). Next, k-means clustering divided annotated genes into distinct chromatin environments by the H3K9ac signal (Fig. 2B). Genes with high levels of H3K9ac were usually accompanied by high levels of H3K9cr, H3K18la, chromatin accessibility and gene expression (Fig. 2B). In addition, H3K9ac-marked genes were also mostly regulated by H3K9cr and H3K18la (Fig. S3B), and genes with all three marks displayed higher chromatin accessibility in their promoter regions and higher gene expression levels than that of the singly marked genes (Fig. 2C). These results suggested that multiple histone acylations might co-regulate transcriptionally active chromatin.

In our analysis, we found that there were many genes annotated in the histone acylation peaks (Fig. S4A), so we sought to determine the predominant biological processes involved in histone acylations. We divided the H3K9ac, H3K9cr and H3K18la peaks into three clusters based on their modification levels in order to explore the biological functions of these histone marks (Fig. S4A). Interestingly, genes in cluster 1 with the strongest H3K9ac and H3K9cr intensity were involved in RNA and protein metabolism, and in neurodegenerative disorder-related pathways (Fig. S4B; Fig. 2D). By contrast, genes with the highest H3K18la levels were mainly involved in cell proliferation-related pathways (Fig. 2D). De novo motif analysis for peak regions of these modifications revealed several transcription factors, such as ELK4 and SOX10, which are important regulators of neural development (Fig. S4C-E). Consistent with in vivo data, we found the wide distributions (Fig. S5A), tight interplays (Fig. S5B,C) and participated functions (Fig. S5D) of these modifications in P19 EC cell-derived neural stem and progenitor cells (NSPCs) by ChIP-seq assays. Collectively, these results indicate that H3K9ac, H3K9cr and H3K18la are tightly correlated with chromatin state and gene expression, but they are involved in different signaling pathways both in vivo and in vitro.

Histone PTMs undergo dynamic changes and play important roles in embryonic neurogenesis, during which NSPCs exit the cell cycle and differentiate into neurons (Yao et al., 2016). However, it is still not clear how histone Kcr and Kla participate in this process. Thus, we systemically analyzed genome-wide changes in H3K9ac, H3K9cr and H3K18la at early (E13.5) and late (E16.5) stages of embryonic neurogenesis by ChIP-seq (Fig. S3A) and RNA-seq (Fig. S6) analyses. We found positive correlations between changes in H3K9ac, H3K9cr and H3K18la and altered gene expression (Fig. S7A-C). Compared with H3K9ac and H3K9cr, we discovered that H3K18la enrichment declined significantly in the promoter regions of genes encoding histones, which play important roles during cell proliferation (Fig. S7D) (Marzluff and Koreski, 2017). Importantly, genes with enrichments of all three marks displayed higher expression changes than those with alterations of single marks (Fig. S7E,F).

We next assessed the peaks with differential histone acylations during neural development. There were different patterns of changes in differential peaks and regulated genes among these histone marks (Fig. 3A,B; Fig. S8A,B). Specifically, compared with H3K9cr and H3K18la, the numbers of peaks and genes with increased H3K9ac enrichment were the highest (Fig. 3A; Fig. S8A), whereas the numbers of peaks and genes with reduced H3K9ac enrichment were the lowest (Fig. 3B; Fig. S8B). However, the predicted biological processes of genes with differential histone marks were highly consistent with neural development (Fig. S8C,D). To correlate the alterations in histone acylations with gene expression changes, we performed binding and expression target analysis (BETA) (Wang et al., 2013) and found that peaks with increased H3K9cr and H3K18la had significant effects as gene activators (Fig. 3C-F), whereas peaks with reduced H3K9cr and H3K18la had significant effects as gene repressors (Fig. S9). Further enrichment analysis indicated that upregulated genes with increased H3K9cr or H3K18la enrichment were mainly involved in neuronal differentiation and maturation (Fig. 3G,H; Fig. S10), whereas downregulated genes with reduced H3K9cr or H3K18la enrichment were tightly associated with cell proliferation (Fig. S9E,F). Collectively, these data support that H3K9ac, H3K9cr and H3K18la might undergo genome-wide alterations and are extensively involved in gene expression changes associated with neural differentiation in vivo.

A recent study showed that the inhibition of transcription results in the loss of histone Kac in bulk histones, and that the majority of histone Kac is dependent on transcription (Martin et al., 2021). However, the associations between non-acetylation histone acylations and transcription under physiological conditions are still unclear. To investigate whether global histone Kcr and Kla levels were influenced by transcription, we treated P19 EC cell-derived NSPCs with the well-studied transcription inhibitors triptolide (inhibiting transcription initiation) and flavopiridol (inhibiting transcription elongation) at high concentrations (Bensaude, 2011). Transcription was significantly inhibited as indicated by the loss of RNA polymerase II (RNAP2) C-terminal domain (CTD) repeat serine 5 phosphorylation (S5p) and RNAP2 CTD repeat serine 2 phosphorylation (S2p) (Fig. 4A) and reduced nascent transcripts from the Actb and Gapdh genes (Fig. 4B). Surprisingly, we did not detect any significant changes in H3K9cr, H3K18cr, H3K14la and H3K18la in P19 EC cell-derived NSPCs and neurons (Fig. 4C,D), indicating that the majority of histone Kcr and Kla was not a consequence of transcription.

Given the fact that some HDACs have been demonstrated to ‘erase’ several histone acylations except Kac (Huang et al., 2018; Wei et al., 2017), we next asked whether HDACs could regulate histone Kla. To this end, we applied two pan-HDAC inhibitors [suberoylanilide hydroxamic acid (SAHA) and valproic acid (VPA)] and found that SAHA treatment increased H3K18la levels, whereas VPA treatment stimulated both H3K14la and H3K18la levels in P19 EC cells (Fig. 5A,B). Furthermore, we found that MS-275, a selective inhibitor of HDAC1-3, robustly promoted H3K14la and H3K18la levels compared with EX 527, a selective SIRT1 inhibitor (Fig. 5C). We also observed an increase of H3K9cr and H3K18cr levels after SAHA, VPA or MS-275 treatment, which was consistent with previous studies (Fig. 5A-C) (Fellows et al., 2018; Kelly et al., 2018; Wei et al., 2017).

To confirm our findings above in vivo, we injected pregnant mice with MS-275 and detected a significant increase in H3K14la and H3K18la levels, as well as an elevation in H3K9cr and H3K18cr levels (Fig. 5D-F). We next investigated whether the dynamic changes in histone Kla levels persisted after interfering with the of expression of HDACs and found that simultaneous knockdown of Hdac1, Hdac2 and Hdac3 significantly increased H3K14la and H3K18la levels in P19 EC cells (Fig. 5G,H). Given the growing literature on regulation of histone Kcr, we explored whether knockdown of the writer (GCN5, also known as KAT2A) and erasers (HDAC1-3) for H3K9cr affect histone methylation and transcription. The results showed that H3K9ac, H3K9cr and S2p (a marker for transcriptional elongation; Chen et al., 2018) levels were significantly reduced after the knockdown of GCN5, whereas H3K9ac, H3K9cr and S2p levels were clearly increased after the knockdown of HDAC1-3 (Fig. S11). However, there were no evident changes in histone methylations, such as H3K4me3, H3K9me3 and H3K27me3, and in S5p (a marker for transcriptional initiation; Chen et al., 2018) levels (Fig. S11). Overall, these results suggest that HDAC1-3 are novel histone lysine delactylases both in cells and in vivo.

H3K18la has recently been reported to regulate M1 macrophage polarization, cell reprograming and tumorigenesis (Li et al., 2020; Yu et al., 2021; Zhang et al., 2019), however, we still do not know how H3K18la is regulated and co-operates with other histone acylations during neural stem cell-fate decisions. To this end, we treated differentiating P19 EC cells with MS-275 and performed anti-H3K9ac, anti-H3K9cr and anti-H3K18la ChIP-seq assays. Inhibition of HDAC1-3 by MS-275 induced genome-wide enhancement of H3K9ac, H3K9cr and H3K18la (Fig. S12A-C). The increased enrichment of these histone marks suggested that they acted as gene activators but not as gene repressors (Fig. S12D-F), and genes under the regulation of these modifications were involved in neuronal differentiation and maturation (Fig. S12G-I). Genes with an increase in all three marks, such as Neurod2, were tightly associated with neural cell-fate decisions (Fig. 6A,B; Fig. S13A) and displayed more prominent expression alterations compared with genes with changes in single marks (Fig. 6C).

Taking gene expression changes into account, there were 635 upregulated genes with an increase in all three histone marks after MS-275 treatment (Fig. S13B), and could be divided into three clusters based on their patterns of expression changes across different timepoints of differentiation (Fig. 6D). Genes in cluster 1 were gradually upregulated during differentiation, were tightly associated with neuronal differentiation and were pre-activated by simultaneous increases in multiple histone acylations (Fig. 6D; Fig. S13C-E). In contrast, genes in cluster 2 might be associated with the development of other germ layers, because their expression levels and histone acylations enrichment were unchanged during whole differentiation processes (Fig. 6D; Fig. S13F-H). Furthermore, genes in cluster 3 were highly expressed in P19 EC cells but were downregulated during differentiation, and they might be involved in maintaining identity of P19 EC cells (Fig. 6D).

In cultured NSPCs isolated from E13.5 mice forebrain, we also validated that the genes with increased histone Kac and Kcr enrichment were upregulated after MS-275 treatment, and that these genes were tightly associated with neuronal differentiation and maturation (Fig. S14, Fig. S15A). Functionally, MS-275-treated NSPCs displayed increased H3K9ac, H3K9cr and H3K18la levels (Fig. S15B) and differentiated into more TUJ1⁺ neurons than those in the control group (Fig. S15C,D). Collectively, these results suggest that MS-275 treatment promotes transcriptional programs associated with neuronal differentiation in differentiating P19 EC cells and NSPCs by the simultaneous activation of multiple histone acylations.

In this study, we profiled dynamic changes and performed functional analyses of histone Kcr and Kla in the developing telencephalon. We demonstrated for the first time that HDAC1-3 act as ‘erasers’ of H3K18la and furthermore linked multiple histone acylations with neural differentiation. First, histone Kcr and Kla exist extensively and undergo significant changes during neural development. Next, H3K9ac, H3K9cr and H3K18la are histone PTMs marking active chromatin regions and are tightly correlated with gene expression. Systematic dissection of alterations in these histone marks across different developmental stages revealed that genome-wide dynamic changes in H3K9cr and H3K18la are extensively involved in neuronal differentiation and cell proliferation processes, which highlights the remodeling of histone acylations to orchestrate gene expression changes as well as cell-fate transitions (Fig. 7A). Additionally, we found that global levels of Kcr and Kla are not influenced by transcriptional inhibition. Importantly, we uncovered that HDAC1-3 influence H3K18la levels and regulate neural differentiation of P19 EC cells by multiple histone acylations (Fig. 7B).

The extensive existence of histone Kcr and Kla that mark many genes in the developing telencephalon indicates that histone acylations are intrinsic pathways to regulate gene expression during mammalian development. In addition, the difference in patterns of changes in H3K9ac, H3K9cr and H3K18la is to be expected as the levels of donors for various histone acylations may alter during neural development. For instance, there is a transition from glycolytic to mitochondrial oxidative phosphorylation during neurogenesis, so it is possible that reduction in the production of lactate and lactyl-CoA may influence histone Kla levels during development (Khacho et al., 2019).

NSPCs in the embryonic forebrain with a deletion of HDAC1, HDAC2 or HDAC3 are prone to differentiate into neurons (Li et al., 2019; Tang et al., 2019). Similarly, we found pre-activation of neuronal transcriptional programs by multiple histone acylations when HDAC1-3 were inhibited during neural differentiation of P19 EC cells. These findings indicate that HDAC1-3 are repressors of neuronal differentiation during neural development. In contrast, CBP knockdown in cultured NSPCs inhibits their neuronal differentiation (Wang et al., 2010). H3K9ac and H3K14ac are enriched at the promoter of Tuj1 and are bound by CBP at E13; their enrichment peaks at E16 and then decreases postnatally in the developing cortex (Wang et al., 2010). In the adult brain, CBP and p300 maintain the excitatory neuron identity through the regulation of histone Kac in cell-type-specific promoter and enhancer regions (Lipinski et al., 2020). These results suggest important roles of CBP and p300 during the establishment and maintenance of neuron identity. At the early stages of neurogenesis, neural differentiation-related genes are kept in a primed state, and both histone acetyltransferases (HATs) and HDACs are simultaneously targeted to these development-related primed genes (Liu et al., 2017a; Wang et al., 2009). Considering the dynamic changes in histone Kac, Kcr and Kla during neural development and the capacity of CBP/p300 to catalyze multiple histone acylations (Liu et al., 2017c; Sabari et al., 2015; Zhang et al., 2019), we propose the existence of a switch for histone acylations at the regulatory elements of genes regulating neural fate decisions, which is triggered by the balance between the activity of CBP/p300 and HDACs (Fig. 7A).

Almost all studies of the causal relationship between histone Kcr/Kla and gene transcription make conclusions based mainly on the consistency between changes in histone acylations and gene expression. However, there is no evidence to support direct regulation of transcription by histone Kcr and Kla at the cellular level. In addition, whether the link between histone Kcr and Kla and transcription under physiological condition is a causal or consequential event has not been addressed previously. Although the present study cannot rule out the possibility of changes of these histone marks in the local sites of the genome, our data clearly demonstrate that transcriptional inhibition has no effect on global histone Kcr and Kla levels, indicating that the majority of histone Kcr and Kla are not dependent on transcription. Consistently, a recent study reports that there is no global effect on histone Kac within a short time frame of transcriptional inhibition (Vaid et al., 2020).

At the cellular level, some HDACs have been shown to ‘erase’ histone Kcr (HDACs 1, 2, 3 and 8, and SIRT1) (Wei et al., 2017), histone lysine benzoylation (SIRT2) (Huang et al., 2018) and histone lysine β-hydroxybutyrylation (HDACs 1 and 2) (Huang et al., 2021). Interestingly, we find that inhibition of HDAC1-3 significantly stimulates H3K14la and H3K18la levels both in vitro and in vivo, which suggests that HDAC1-3 are novel ‘erasers’ of histone Kla. A recent study by Olsen and colleagues reports the histone delactylation activity of HDAC1-3 along with unaltered H3K18la levels after knockdown of HADC1-3 in HeLa cells under physiological conditions (Moreno-Yruela et al., 2022). In the present study, however, we observed a remarkable increase in H3K14la and H3K18la levels after simultaneous knockdown of Hdac1-3 in P19 EC cells. We speculate that this discrepancy may be due to the knockdown efficiency of HDAC1-3 and specific cell types.

Convincingly, our in vivo and ChIP-seq results support the conclusion that HDAC1-3 can ‘erase’ H3K18la at the cellular level. Functionally, genes that gain H3K18la after MS-275 treatment are directly involved in neural differentiation (Fig. S12I). Considering that the associations between histone acylations and their biological functions are less understood (Fellows et al., 2018; Kelly et al., 2018), our study identifies a novel link between HDACs and multiple histone acylations to regulate cell-fate transitions associated with the transcriptional program (Fig. 7B). Future studies should focus on the crosstalk between Kcr and/or Kla with other histone modifications or regulatory elements, such as the co-occurrence of these modifications with super-enhancers. However, the present study has several limitations. We did not directly determine the regulatory relationships of Kcr/Kla with neuronal fate, due to the limited specific pathways that stimulate or block these modifications. Moreover, due to technical limitations, we still need to confirm the results from P19 EC cell-derived NSPCs in other physiological models, such as embryonic NSPCs and in vivo.

Overall, our study provides novel insights into the dynamic changes, epigenetic interactions, transcriptional regulation and biological functions of histone acylations during neural development. We also link HDACs and multiple histone acylations with NSPC-fate decisions and highlight the close relationships and crucial roles of these modifications in the epigenetic regulation of neural development and disease.

All mice were from the ICR/CD-1 background. Pregnant mice were purchased from SiPeiFu Biotechnology (Beijing, China). All mouse experiments were approved by the Animal Committee of the Institute of Zoology, Chinese Academy of Sciences, Beijing, China.

P19 EC cells were obtained from the American Type Culture Collection (CRL-1825) and maintained in the proliferation medium [minimum essential medium α (MEMα; Gibco) containing 10% fetal bovine serum (FBS; Gibco)] in an incubator at 37°C and 5% CO₂. For induced neural differentiation, 1×10⁶ cells were cultured on a 9.6 cm bacterial-grade Petri dish in the induction medium [proliferation medium containing 1 μM all-trans retinoic acid (Sigma-Aldrich)] for 4 days, with the medium being exchanged every 2 days. To induce neuronal differentiation, the aggregated neurospheres were digested to single cells by Accutase (STEMCELL Technologies) and plated at a density of 1.2×10⁷ cells per 10 cm culture dish (Corning) pre-coated with 0.1 mg/ml poly-D-lysine (PDL; Sigma-Aldrich) in the differentiation medium [Dulbecco's modified Eagle medium (DMEM)/F12 (Gibco) containing 2% B27 (Gibco)] for 4 days, with the medium being exchanged every 2 days. 293T cells were cultured in the DMEM medium (Gibco) containing 10% fetal bovine serum. NSPCs isolated from E13.5 mice forebrain were cultured and treated with MS-275 as previously described (Dai et al., 2021). The differentiation protocols for NSPCs were also followed as per our published methods (Dai et al., 2021), but 5% FBS was employed in the present study instead of forskolin and retinoic acid.

All inhibitors were purchased from Selleck Chemicals. P19 EC cells treatments were performed in the proliferation medium at the following concentrations for 12 h or 24 h: 0.5 µM SAHA (S1047, CAS: 149647-78-9), 1 mM VPA (S3944, CAS: 99-66-1), 1 µM MS-275 (S1053, CAS: 209783-80-2) and 1 µM EX 527 (S1541, CAS: 49843-98-3). For in vivo experiments, 20 mg/kg MS-275 was injected every 24 h for two consecutive days in pregnant mice at day 13.5 of pregnancy. P19 EC cell-derived NSPCs (day 1) were treated with 1 µM MS-275 for 1 day and then collected for ChIP-seq and RNA-seq assays. P19 EC cell-derived NSPC (day 2) treatments were performed in the induced medium at the following concentration for 1 h: 5 µM triptolide (S3604, CAS: 38748-32-2) and 5 µM flavopiridol HCl (S2679, CAS: 131740-09-5).

P19 EC cells were transfected with double-stranded siRNA using Lipofectamine RNAiMAX transfection reagent (Invitrogen) at a final concentration of 50 nM for 3 days and then collected for western blotting. 293T cells were transfected with siRNA using Lipofectamine 3000 transfection regent (Invitrogen) at a final concentration of 50 nM for 2 days and then collected for western blotting. The negative control siRNA (siR NC #1, siN0000001-1-5) and siRNA targeting Hdac1 (si-m-Hdac1_004, siB140729172936), Hdac2 (si-m-Hdac2_004, siB140729171225), Hdac3 (si-m-Hdac3_001, siG141217112523) and GCN5/KAT2A (si-h-KAT2A_004, siB14815190908) were purchased from Ribobio (Guangzhou, China).

DAPI (Sigma-Aldrich) was used to label cell nuclei. For immunocytochemistry analysis, cells on coverslips were fixed in 4% formaldehyde for 15 min at room temperature and then washed three times with phosphate-buffered saline (PBS). Next, cells were incubated with blocking solution (2% bovine serum albumin and 0.3% Triton X-100 in PBS) at room temperature for 1 h. Coverslips were then incubated with primary antibodies diluted in blocking solution at 4°C overnight and labeled using appropriate secondary antibodies at room temperature for 2 h. For immunostaining, embryonic mouse brains were cut into 12 µm-thick cryosections using Leica CM 1950. Staining procedures for brain slices were the same as those used for immunocytochemistry analysis. Detailed information on the primary and secondary antibodies used for immunocytochemistry and immunostaining analysis is summarized in Table S1. All images were acquired using Zeiss LSM710 and LSM880 confocal microscopes. For quantitative analysis of the enrichment levels of histone acylations in Fig. 1B, three 40 µm×40 µm bins were randomly selected across each layer of cortex, then ImageJ was used to perform quantitative analysis of mean fluorescence intensity.

Total proteins were extracted using the RIPA lysis buffer (Beyotime Biotechnology, Shanghai, China) containing 1× protease inhibitor cocktail (Bimake). Nuclear proteins were acquired as previously described (Dai et al., 2021) and protein concentrations were measured using a bicinchoninic acid (BCA) protein assay kit (Beyotime). Histones were isolated from cells and tissues following standard acid extraction protocols (Shechter et al., 2007) and their concentrations were determined using the Bradford protein assay kit (Beyotime). In brief, cells and tissues were incubated with extraction buffer (PBS with 0.5% Triton X-100, 1 mM phenylmethylsulfonyl fluoride and 5 mM sodium butyrate) with gentle rotation at 4°C for 10 min, and nuclei were isolated by centrifugation (1500 g for 10 min at 4°C). Next, 0.4 M HCl was used to resuspend nuclei, followed by incubation on ice for 30 min, and soluble histones were isolated by centrifugation (20,000 g for 10 min at 4°C).

Proteins were separated using 6-15% SDS-PAGE and transferred onto polyvinylidene fluoride (PVDF) membranes (Millipore). Blotted membranes were blocked in 5% skim milk (BD Biosciences) in Tris-buffered saline with 0.1% Tween 20 (TBS-T) and incubated with primary antibodies at 4°C overnight. Membranes were then washed three times with TBS-T and incubated with secondary antibodies at room temperature for 1-2 h. The immunoreactive products were detected using SuperSignal West Pico PLUS Chemiluminescent Substrate (Thermo Fisher Scientific) and a Tanon-5200 Chemiluminescent Imaging System (Tanon, Shanghai, China). ImageJ was used to perform quantitative analysis of western blot results as the ratio of the band intensities of target histone marks to the band intensities of reference proteins (H3, β-actin or β-tubulin). For dot blot assays, peptide samples were spotted onto PVDF membranes, which were first blocked at room temperature for 1 h and then washed in TBS-T. The membranes were incubated with primary and secondary antibodies at room temperature for 2 h and 45 min, respectively, followed by three washes with TBS-T and signal detection. Detailed information on the primary and secondary antibodies used for western blot assay is summarized in Table S1.

RNA was extracted using Trizol reagent (Invitrogen), then RNA samples were subjected to quantitative reverse transcription PCR (qRT-PCR) using a Roche LightCycler 480 II. In brief, total RNA was transcribed into cDNA using the TransScript One-Step gDNA Removal and cDNA Synthesis Kit (TransGen Biotech, Beijing, China) with random primers. For qRT-PCR analysis, cDNA was quantified by using Hieff qPCR SYBR Green Master Mix (YEASEN, Shanghai, China) in a 20 μl reaction system according to the manufacturer's instructions and PCR steps were performed as follows: 5 min initial pre-denaturation at 95°C, followed by 45 cycles of 10 s at 95°C, 30 s at 60°C and 30 s at 72°C. Samples were run in triplicate. The relative expression was calculated using the 2^−ΔΔCt method. The PCR primers were designed to target both exon and intronic regions to capture the nascent transcripts of Actb and Gapdh, and their sequences are provided in Table S2.

We used the TruePrep DNA Library Prep Kit V2 for Illumina (Vazyme, Nanjing, China) to construct assay for transposase-accessible chromatin using sequencing (ATAC-seq) libraries in two biological replicates. In brief, a total of 5×10⁴ cells (P19 EC cell-derived NSPCs and digested telencephalon cells at E13.5) were lysed in 50 µl lysis buffer [10 mM Tris-HCl (pH 7.4), 10 mM NaCl, 3 mM MgCl₂ and 0.5% NP-40] on ice for 10 min. The suspension was then centrifuged at 500 g for 5 min at 4°C, followed by the addition of 50 µl transposition reaction mix containing 5×TTBL and TTE Mix V50 (Vazyme) and incubation at 37°C for 30 min. After tagmentation, fragmented DNA was purified with 100 μl VAHTS DNA clean beads (Vazyme) and amplified for 12 cycles using the following PCR conditions: 72°C for 3 min; 98°C for 30 s; and thermocycling at 98°C for 15 s, 60°C for 30 s and 72°C for 30 s; followed by 72°C for 5 min. After PCR amplification, libraries were purified and selected using DNA clean beads. Libraries were then sequenced on Illumina platforms by Annoroad Gene Technology (Beijing, China).

mRNA-seq was performed in two biological replicates. Briefly, RNA was extracted using Trizol reagent (Invitrogen), and mRNA was purified from total RNA using poly-T oligo-attached magnetic beads (New England BioLabs). Libraries were then constructed using the NEBNext Ultra RNA Library Prep Kit for Illumina (New England Biolabs), following the manufacturer's instructions, and sequenced on Illumina platforms by Annoroad Gene Technology.

Native ChIP was carried out in two biological replicates following the published protocol with minor changes (Kelly et al., 2018). In brief, cell pellets (P19 EC cell-derived NSPCs and digested telencephalon cells at E13.5 and E16.5) were resuspended in 500 μl swollen buffer [0.3 M sucrose, 60 mM KCl, 15 mM NaCl, 5 mM MgCl₂, 0.1 mM EDTA, 15 mM Tris-HCl (pH 7.5), 5 mM sodium butyrate and 1× protease inhibitor cocktail]. Nuclei were isolated using an equal volume extraction buffer [0.3 M sucrose, 60 mM KCl, 15 mM NaCl, 5 mM MgCl₂, 0.1 mM EDTA, 15 mM Tris-HCl (pH 7.5), 0.4% NP-40, 5 mM sodium butyrate and 1× protease inhibitor cocktail] with gentle rotation (30 rpm) at 4°C for 10 min. The cell suspensions were transferred into new tubes containing 500 μl sucrose cushion buffer [0.6 M sucrose, 60 mM KCl, 15 mM NaCl, 5 mM MgCl₂, 0.1 mM EDTA, 15 mM Tris-HCl (pH 7.5), 5 mM sodium butyrate and 1× protease inhibitor cocktail] and centrifuged at 1500 g for 10 min at 4°C. Nuclei were softly resuspended by MNase digestion buffer [0.32 M sucrose, 4 mM MgCl₂, 1 mM CaCl₂, 50 mM Tris-HCl (pH 7.5), 5 mM sodium butyrate and 1× protease inhibitor cocktail]. The amount of chromatin present in the sample was assessed by spectrophotometry in the presence of 0.1% SDS and chromatin was then digested using Micrococcal nuclease (Takara) (1 U per 10 µg chromatin) at 600 rpm in a thermomixer for 5 min at 37°C. The reaction was stopped by adding EDTA to a final concentration of 5 mM and incubated on ice for 5 min, then cleared at 12,000 g for 10 min at 4°C. The supernatant (S1 fraction) was stored at 4°C, and the pellets were resuspended in dialysis buffer [0.2 mM EDTA, 1 mM Tris-HCl (pH 7.5) and 5 mM sodium butyrate] and then dialyzed against 1.5 L dialysis buffer by Slide-A-Lyzer Dialysis Cassettes, 10 K MWCO, 3 ml (Thermo Fisher Scientific) overnight at 4°C. The dialyzed pellets (S2 fraction) were cleared at 20,000 g for 10 min at 4°C and combined with the S1 fraction for immunoprecipitation.

Additionally, we measured the double-stranded DNA (dsDNA) amount of the combined S1 and S2 fractions using Qubit 4.0 fluorometer (Invitrogen) by Equalbit 1× dsDNA HS Assay Kit (Vazyme) to enable an accurate calculation of the chromatin amount used for immunoprecipitation experiments. Finally, 25 µg chromatin (dsDNA) was immunoprecipitated using 1.25 μg anti-H3K9ac (PTM Biolabs, Hangzhou, China) and anti-H3K9cr (PTM Biolabs) antibodies, and 150 µg chromatin was immunoprecipitated using 3 µg anti-H3K18la antibody (PTM Biolabs) overnight at 4°C with gentle rotation (20 rpm) in ChIP dilution buffer [1% Triton X-100, 2 mM EDTA, 150 mM NaCl, 50 mM Tri-HCl (pH 7.5), 5 mM sodium butyrate and 1× protease inhibitor cocktail]. Antibody-chromatin complexes were collected using 20 µl Magna ChIP protein A+G magnetic beads (Millipore) by incubation with rotation (20 rpm) for 2 h at 4°C, and the beads were sequentially washed once with low-salt wash buffer [0.1% SDS, 1% Triton X-100, 2 mM EDTA, 150 mM NaCl and 20 mM Tris-HCl (pH 8.1)], once with high-salt wash buffer [0.1% SDS, 1% Triton X-100, 2 mM EDTA, 500 mM NaCl and 20 mM Tris-HCl (pH 8.1)], once with LiCl wash buffer [250 mM LiCl, 1% NP-40, 10 mM Tris-HCl (pH 8.1), 1% sodium deoxycholate and 1 mM EDTA] and twice with TE buffer [10 mM Tris-HCl (pH 8.1) and 1 mM EDTA] for 10 min each at 4°C and 30 rpm. Immune complexes were eluted from the beads twice for 15 min at 65°C and 1000 rpm in a thermomixer using 250 µl elution buffer (1% SDS and 100 mM NaHCO₃). The eluted immune complexes were then treated with RNase A (100 µg/ml) for 30 min at 37°C, followed by Proteinase K (50 µg/ml) for 2 h at 58°C and 600 rpm in a thermomixer. Then, the immunoprecipitated DNA was purified using phenol-chloroform-isoamyl alcohol (25:24:1) and precipitated with two volumes of 100% ethanol and 10 µg linear acrylamide (Invitrogen) at −20°C for 2 h. Libraries were then constructed with the TruSeq ChIP Library Preparation Kit (Illumina) following the manufacturer's instructions and sequenced on Illumina platforms by Annoroad Gene Technology.

The original image data obtained after sequencing was converted into the sequence data stored in FASTQ file format by base calling. The raw reads of the datasets were filtered by Trimmomatic (v0.36) to remove contaminations and low-quality reads with the parameters ‘ILLUMINACLIP:Adapter.fa:2:30:10 SLIDINGWINDOW:4:15 LEADING:3 TRAILING:3 MINLEN:36’ and quality-controlled using FastQC (v0.11.7) (Bolger et al., 2014).

The reads of ATAC-seq data were aligned using Bowtie 2 (v2.4.1) to the mouse reference genome with parameters ‘-X 1000 --mm –local’ (Langmead and Salzberg, 2012). Samtools (v1.10) was then used to convert files into BAM format and filter-mapped reads with the parameters ‘-F 1804 -f 2 -q 30’ (Li et al., 2009). Next, we removed PCR duplicates using the Mark Duplicates function in Picard (v2.22.0) and reads mapped to the mitochondrial genome. Finally, MACS2 (v2.2.7.1) was used to call peaks with the parameters ‘--keep-dup all -q 0.01 -f BAMPE’ (Li et al., 2009). The quality control of ATAC-seq data is shown in Table S3.

High-quality ChIP-seq reads were aligned using Bowtie 2 to the mouse reference genome using default parameters. Samtools was used to convert files into BAM format and filter-mapped reads with the parameters ‘-F 1804 -f 2 -q 30’. PCR duplicates were removed using the Mark Duplicates function in Picard, and MACS2 was used to perform peak calling (-q 0.01) for the immunoprecipitated samples relative to the input samples. The quality control of ChIP-seq data is shown in Table S3.

High-quality reads of RNA-seq data were quantified using Salmon (v1.1.0) with the parameters ‘-i -g --gcBias –validateMappings’ (Patro et al., 2017). Differential gene expression analysis was conducted using DESeq2 (v1.26.0) and differentially expressed genes were defined by P-adjust<0.05 and absolute fold change >0 (Love et al., 2014). The quality control of RNA-seq data is shown in Table S3.

Peak annotation was performed using ChIPseeker (v1.22.1) at the gene level and promoter regions were defined as ±2 kb of the transcription start site (TSS) (Yu et al., 2015). Genes with at least one peak annotated in their promoter regions or gene body regions were considered to be under the regulation of histone marks. Gene enrichment analysis was performed using clusterProfiler (v3.14.3) with default parameters (Yu et al., 2012).

To analyze dynamic changes in histone acylations at promoter regions, the BEDTools (v2.29.0) ‘coverage’ function was used for read counting within ±2 kb of TSS (Quinlan and Hall, 2010), and DESeq2 was then used to perform quantitative analysis of genes with differential histone acylations, which were defined by P-adjust<0.05 and absolute fold change >1.5. DiffBind (v2.14.0) and MAnorm (v1.2.0) were adopted to perform analysis of dynamic changes in histone acylations at peak regions, and peaks with differential histone acylations were defined by P-adjust<0.05 and absolute fold change >1.5 (Shao et al., 2012). Genes with differential histone acylations were defined as at least one peak annotated in their promoter regions or gene body regions. For comparison analysis of histone acylations enrichment among different groups of NSPCs, the BEDTools ‘coverage’ function was used for read counting within both ±2 kb of TSS (promoter regions) and gene body regions, then the raw read counts were used to calculate reads per kilobase per million mapped reads (RPKM) values to represent ChIP signals as follows: [(read counts)/(region length in kb)]/(total mapped and filtered reads in Mb).

BAM files were converted into bigwig files using the deepTools (v3.4.0) ‘bamCoverage’ function with the parameters ‘--binSize 50 --normalizeUsing RPKM’, and the ‘computeMatrix’, ‘plotHeatmap’, ‘multiBigwigSummary’ and ‘plotCorrelation’ functions were used to generate heatmaps (Ramírez et al., 2016). For genome browser representation, data in bigwig files generated by deepTools were visualized using IGV (v. 2.4.10) (Thorvaldsdóttir et al., 2013). The BETA ‘basic’ function (--d.f. 0.05 --da 1 -c 0.001) was used for activating and repressive function prediction of peak regions with differential histone acylations (Wang et al., 2013). Peaks between replicates were combined by the BEDtools ‘intersect’ function, and bigwig files between replicates were merged by the deepTools ‘bigwigCompare’ functions with the parameters ‘--operation mean’. The mouse reference genome sequence (vM24) and gene annotation (vM24) were downloaded from GENCODE (https://www.gencodegenes.org/).

The results shown are the mean±s.e.m. One-way ANOVA with post hoc Tukey's test, unpaired two-tailed Student's t-test, Kruskal–Wallis with post hoc Dunn's test, unpaired two-tailed Student's t-test with Welch's correction and two-sided Fisher's exact test were used to analyze statistical significance. Differences were regarded as significant by their P-values: NS, no significance; *P<0.05; **P<0.01; and ***P<0.001. All statistical analyses were performed and diagrams generated using GraphPad Prism 8.0 and R (v3.6.3). Statistical details and methods used in each experiment can be found in the respective figure legends.

We are grateful to Shi-Wen Li and Xi-Li Zhu for their help with confocal imaging.

The peer review history is available online at https://journals.biologists.com/dev/article-lookup/doi/10.1242/dev.200049.