Main

Development and neuroinflammation in the mammalian brain are meticulously regulated; for example, mouse corticogenesis proceeds through an inside-out mechanism. Cells of the deeper cortical layers V and VI, such as corticothalamic projection neurons (PNs; layer VI CThPNs, embryonic day 12.5 (E12.5)), subcerebral projection neurons (layer V SCPNs, E13.5) and some callosal PNs (CPNs; E12.5–13.5) that populate deeper layers are defined earlier than those that populate the outer layers (IV and II/III), including granular (spiny stellate) neurons (E14.5) in layer IV and CPNs (E15.5) in layer II/III2,3. Oligodendrocyte and astrocyte genesis occurs simultaneously, populating the cerebral cortex (CTX) and the underlying corpus callosum (CC)4,5. Neurodegeneration reactivates many of the developmental programs.

In recent years, disease-associated glial states have been a primary focus, driven by large-scale sequencing strategies to delineate cell states by examining RNA, protein expression and chromatin accessibility6,7,8,9,[10](#ref-CR10 “Kirby, L. & Castelo-Branco, G. Crossing boundaries: Interplay between the immune system and oligodendrocyte lineage cells. Semin. Cell Dev. Biol. https://doi.org/10.1016/j.semcdb.2020.10.013

(2020).“),11,12,13. Spatial omics analyses offer an important view of the molecular architecture of biological systems. One spatial proteomics approach, CODEX, provides single-cell-resolution cell type and morphological information14. We developed a platform, deterministic barcoding in tissue (DBiT), that enables spatial co-profiling of transcriptome and epigenome1, as well as transcriptome and proteome15. Extending this approach to spatially map the epigenome, transcriptome and proteome within the same tissue section would enable a comprehensive investigation of the molecular mechanisms across all layers of the central dogma16. Here, we developed and applied DBiT-based spatial tri-omic technologies, including spatial assay for transposase-accessible chromatin (ATAC)–RNA–protein-seq (spatial ARP-seq) and CUT&Tag–RNA–protein-seq (spatial CTRP-seq), to simultaneous profile genome-wide chromatin accessibility or histone modifications (H3K27me3), the whole transcriptome and the proteome (around 150 proteins) within the same tissue section at the cellular level. To elucidate the cellular and molecular processes during brain development and neuroinflammation, we applied these technologies to mouse brains from newborn to juvenile stages (benchmarking to the developing human brain), and to a focal demyelination mouse model induced by lysolecithin (LPC) (Fig. 1a, Supplementary Fig. 1a,b and Supplementary Tables 7–9; online resource: https://spatial-omics.yale.edu/).

Fig. 1: Spatial tri-omics mapping of the developing mouse and human brains at the onset of myelination.

a, Schematic of the workflow. Some elements of the diagram were generated using BioRender and adapted from ref. 1, Springer Nature Limited, under a Creative Commons Attribution 4.0 licence. APC, antigen-presenting cells; astro, astrocytes; Gran/Neu, granulocytes and neutrophils; Mono/mac, monocytes and macrophages; NP, neural progenitors; prolif, proliferating cells; vasc, vasculature. b, CODEX images of the mouse brain coronal sections from P0 to P21 with two replicates (S1 and S2). Scale bar, 500 μm. c, CODEX images of the human brain V1 region at the second trimester (154 days after conception (d.p.c.)), the third trimester (224 days after conception) and infancy (507 days after conception). Scale bars, 1 mm. d, Uniform manifold approximation and projection (UMAP) and spatial distribution of spatial RNA clusters of the human brain samples at different stages. e, UMAPs of the mouse brain spatial RNA clusters at different ages from P0 to P21, with sample information. f, UMAPs of the mouse brain spatial ATAC clusters at different ages from P0 to P21, with sample information. g, The spatial distribution of RNA clusters in e. h, The spatial distribution of ATAC clusters in g. i, The spatial domain of the mouse brains generated by the integration of RNA and ATAC data in spatial ARP-seq (S1). Scale bars, 500 μm. j,k, Spatial mapping of gene expression (j) and the GAS (k) of Mbp and Foxp4 for RNA and ATAC in spatial ARP-seq. l, MBP expression from ADT protein data in spatial ARP-seq (S1). m, Label transfer of MFOLs from scRNA-seq28,29,30 to P10 and P21 mouse brains.

Multiplex CODEX imaging of the developing brain

We designed a 23-plex CODEX antibody panel for mouse (Fig. 1a) and a 19-plex panel for human (Supplementary Fig. 1c), including markers for major cell types in brain development and neuroinflammation. We imaged the mouse brain coronal sections (0.5 mm relative to bregma for P21, and closely analogous positions for P0–10), and the human brain primary visual (V1) cortical regions[17](https://www.nature.com/articles/s41586-025-09663-y#ref-CR17 “Wang, L. et al. Molecular and cellular dynamics of the developing human neocortex. Nature https://doi.org/10.1038/s41586-024-08351-7

(2025).“) at the second trimester (estimated 154 days after conception), third trimester (224 days after conception) and at infancy (year 0–1, 507 days after conception). Segmentation, Seurat clustering and cell typing of the CODEX data revealed the arealization of the mouse brain with expected marker distribution, specificity and morphology (Fig. 1b,c, Extended Data Fig. 1 and Supplementary Figs. 2–9 and 11–13).

The ramified astrocytic/neural progenitor marker GFAP was expressed in the mouse glial limitans and medial ventricular zones at P0, and began to expand at P5 to specific areas in the parenchyma and CC, a region where it was particularly abundant at P21 (Fig. 1b and Supplementary Figs. 2–7). We found Ki-67+ proliferating cells densely concentrated in the ventricular region, where neuroblasts, precursor cells and glia reside18,19. The Ki-67+ periventricular and white-matter (WM) density contracted during later development, as progenitor cells became terminally differentiated (Fig. 1b, Extended Data Fig. 1 and Supplementary Figs. 2–9).

Neuronal markers populated mostly the grey matter regions. Each marker showed cortical layer specificity with CUX1/2 (layers II/III/IV) in the upper layers and TBR1 (IV/V/VI) and CTIP2 (V/VI and striatum) in deeper cortical layers (Figs. 1b and 2a,b, Extended Data Fig. 1 and Supplementary Figs. 2–7). While these layer-defining transcription factors (TFs) (CUX1/2, CTIP2 and SATB2) were expressed at P21, their protein expression and the density of cells expressing these markers was reduced compared with at early postnatal stages (Figs. 1b and 2a,b and Supplementary Figs. 2–7 and 10), consistent with the Allen Brain Atlas profiling20 (Supplementary Fig. 10).

Fig. 2: Spatiotemporal dynamics of transcriptome and chromatin accessibility of the mouse and human brain cortical layers.

a, CODEX images of the mouse brain coronal sections from P0 to P21. Scale bar, 500 μm. b, CODEX staining of delineated cortical layers in the region of interest indicated by dashed rectangles in a. SSp, primary somatosensory area. Scale bar, 250 μm. c, Schematic of the spatiotemporal regression model on mouse brain cortical layers. d, Different cell types in cortical layers (S2) from the spatial RNA data. Epend, ependymal; Ex, excitatory; In, inhibitory; Mac, macrophages; MG, microglia; neu, neurons; NIPC, neuronal intermediate progenitor cell; MSN, medium spiny neurons; Cajal R, Cajal–Retzius cells; Mig Int Neu, migrating interneurons. e, The 15 RNA clusters generated from the regression model (S2). f, The 10 ATAC clusters generated from the regression model (S2). g, Representative genes in the cortical layers. h, RNA gene expression (top) and the ATAC GAS (bottom) calculated on the basis of the regression model for specific genes in mouse brains, segregated into groups based on variance of ATAC GAS in comparison to RNA expression. i, RNA gene expression (left) and the ATAC GAS (right) calculated on the basis of the regression model for specific genes in mouse brains. j, RNA gene expression (top) and the ATAC GAS (bottom) for myelination-related genes from the regression model. k, Representative gene expression (top) and ATAC GAS (bottom) for human brains, segregated into groups based on variance of ATAC GAS (bottom) in comparison to RNA expression (top). 2nd, second trimester; 3rd, third trimester. l, Representative RNA gene expression (top) and ATAC GAS (bottom) for oligodendrocyte-lineage-associated and myelination-associated genes through human V1 cortical development.

OLIG2, a pan oligodendrocyte marker was expressed in cells mainly in the CC but also throughout the brain (Figs. 1b and 3a,b and Supplementary Figs. 2–9). We found MBP and MOG (markers of mature oligodendrocytes, MOLs) were absent between P0 and P5, with expression starting at P7 and P10, respectively, as previously reported4 (Figs. 1b and 3a and Supplementary Figs. 2–9). Notably, the expression of MBP and MOG was initially limited to the lateral part of the CC (P7–10), beneath the primary somatosensory area and only spread throughout the entire CC at P21, suggesting a lateral-to-medial progression of myelination. Myelin in the CTX was also identified at P21 (Figs. 1b and 3a and Supplementary Fig. 8).

Fig. 3: Spatiotemporal dynamics of the mouse CC during development and myelination.

a, CODEX images of the mouse brain CC from P0 to P21. Scale bar, 500 μm. Individual channels are shown in Supplementary Figs. 2–9. b, Magnified CODEX images of the medial CC regions in a from P5, P7 and P10, showing staining for OLIG2, PDGFRα and Ki-67. Scale bar, 250 μm. c, Schematic of the spatiotemporal regression model on the mouse brain CC from P0 to P21. T&S, time and space. d, The six RNA clusters generated from the regression model. e, The 4 ATAC clusters generated from the regression model. f, Spatial RNA patterns for RNA clusters (R1–6) from the regression model. g, RNA gene expression calculated on the basis of the regression model for specific genes. h, Schematic of the injection sites and tracing pathways for CPNs and CThPNs (top) using retro-AAV-eGFP at the terminal location for both tracts separately. Bottom, the segmentation strategy using confocal imaging at high magnification (×40) for MBP staining and retro-AAV-eGFP labelling quantification throughout the CC. Thal, thalamus. i, Confocal images (×20) stitched together to capture the CC and overlying cortex of one hemisphere for CPN tracing in P10 mouse after retro-AAV-eGFP injection at P1. Representative retro-AAV-eGFP tracing from the contralateral cortical layer II/III injected terminal. j, Representative retro-AAV-eGFP tracing from the ipsilateral thalamic injected terminal. For i and j, scale bar, 250 μm. k, The percentage area of MBP staining per bin across the defined CC in each image moving medial to lateral (n = 6; 3 CPN traced + 3 CThPN traced). Data are mean ± s.d. for the region. l, The percentage area of retro-AAV-eGFP-labelled axons of the area of the CC defined per image. Data are mean ±s.d. for the binned region. n = 3 (CTX injected) and n = 3 (thalamus injected). m, The percentage area of co-labelled MBP and retro-AAV-eGFP-labelled axons of the area of the CC defined per image. Data are mean ± s.d. for the binned region. n = 3 (CTX injected) and n = 3 (thalamus injected). Significance was determined by calculating the area under the curve (AUC) formed by all the 15 datapoints per group, then running a two-tailed unpaired t-test (P = 0.0001) on the mean total area values (CTX, 18.79; thalamus, 143.2) with the standard error values (CTX, 2.770; thalamus, 7.942). ***P <0.001. n, Summary of how myelination progresses in the CC and CTX.

We included markers for immune cells: granulocytes/neutrophils (Ly6G+), monocytes/macrophages (IBA1+, CD11b+, CD169+), dendritic cells (CD11c+, MHC class II+), T cells (CD3+, CD4+C, CD8+) and B cells (CD19+, CD45R+) (Fig. 1a). As expected, most of the immune subtypes were sparse and primarily localized to the meningeal compartment. However, CD169 was abundant in the meninges and present in the choroid plexus, identifying a developmental border-associated macrophage counterpart to central nervous system (CNS)-resident CD169+ cells that increase after ischaemia21,22 (Fig. 1b and Supplementary Figs. 2–7).

Imaging with the human CODEX panel revealed distinct arealization (Fig. 1c and Supplementary Fig. 1c). Although SATB2, CTIP2 and TBR1 protein levels appeared to be reduced in infancy compared with at earlier stages, adjusting for cellular density indicates a domain expansion aligning with increased brain volume during development (Supplementary Figs. 11–14). Myelination marked by MOG expression intensified during infancy, consistent with previous reports23 (Fig. 1c and Supplementary Figs. 12 and 13). Overall, our CODEX profiling provides a comprehensive single-cell protein benchmark and spatial cell typing across mouse and human brain development.

DBiT-based spatial tri-omics profiling

To further investigate the molecular mechanisms across all levels of the central dogma, we developed an all-encompassing spatial assay for co-profiling of the epigenome, transcriptome and a large protein panel simultaneously within the same tissue section—spatial ARP-seq (Fig. 1a and Supplementary Fig. 1a). A frozen tissue section was fixed with formaldehyde, incubated with a cocktail of antibody-derived DNA tags (ADTs) for mouse (Supplementary Table 4) or human (Supplementary Table 5), and treated with Tn5 transposase loaded with a universal ligation linker to insert the adapter at transposase accessible genomic loci. The tissue was then incubated with a biotinylated poly(T) adapter to bind to the poly(A) tail of mRNAs and ADTs for in-tissue reverse transcription. Next, two microfluidic channel array chips with perpendicular channels were applied to introduce spatial barcodes Ai (i = 1–100 or 220) and Bj (j = 1–100 or 220) into the tissue, forming a 2D grid of spatially barcoded tissue pixels. Each pixel, defined by a unique combination of Ai and Bj barcodes, measured 20 μm (100 barcodes; data S1) or 15 μm (220 barcodes; data S2), resulting in a total of 10,000 or 48,400 barcoded pixels. After releasing the barcoded cDNAs (from both mRNAs and ADTs) and genomic DNA (gDNA) fragments, separate libraries for gDNA and cDNA were constructed for next-generation sequencing (Supplementary Fig. 1a). We also developed spatial tri-omic profiling method that simultaneously measures genome-wide histone modifications, the transcriptome and proteins—termed spatial CTRP-seq. The procedure followed a workflow similar to that described above, except that an antibody against H3K27me3 was applied to the tissue section, followed by a protein-A-tethered Tn5–DNA complex to perform cleavage under targets and tagmentation (CUT&Tag)24,25 (Fig. 1a and Supplementary Fig. 1b; details of quality control and metrics are provided in the Methods).

Spatial deconvolution of the brain around birth

To investigate the spatial molecular dynamics of mammalian cortical and WM development, we used spatial ARP-seq on tissue sections that were either adjacent to or from regions anatomically matched to those processed with CODEX. We clustered RNA and ATAC data from postnatal mouse brains separately, integrating across all timepoints with the same pixel size (Supplementary Table 7). For the coronal mouse brains from P0–21 processed with spatial ARP-seq (replicate S1), we identified 22 major RNA (R0–21) and 15 major ATAC clusters (A0–14) (Fig. 1e–h). We further integrated our RNA and ATAC data using SpatialGlue, leading to 18 refined spatial domains (D1–18) (Fig. 1i and Supplementary Fig. 18a,b). The spatial distribution of these domains aligned with tissue histology from the Allen Brain Atlas20 and provided better arealization when compared to CODEX (Fig. 1b, Extended Data Fig. 1 and Supplementary Figs. 2–7), with the division of cortical layers II/III and V, subdivision between primary and secondary motor cortical areas and eventually the anterior cingulate at each layer. Medial D4 and D14 clusters segregated from clusters D9 (layer II/III) and D2 (layer V), respectively (Fig. 1i). Similar spatial domains were observed in the replicate (S2) (Extended Data Fig. 2a–e).

Consistent with CODEX, chromatin accessibility and RNA expression of Mbp and Mog in the CC were not present at P0–5 but were rather induced at P7–10 (Fig. 1j,k, Extended Data Fig. 2f,g and Supplementary Fig. 18c–e). Moreover, spatial ARP-seq data revealed that MBP and MOG proteins arose also at later stages (Fig. 1l and Extended Data Fig. 3c). While the striatum was characterized by markers of medium spiny neurons such as Bcl11b (clusters D5, D12 and D17) (Fig. 1i and Supplementary Fig. 19), the expression and ATAC gene activity score (GAS) of Foxp4, a TF that is known to regulate morphogenesis, disappeared gradually, while Foxp1 and Foxp2 retained dominance in the postnatal striatum, consistent with previous studies24,25,26,27 (Fig. 1j,k, Extended Data Fig. 2f,g and Supplementary Fig. 18c–e). Spiked-in spatial ARP-seq TF ADTs also showed general concordance in the positional signal for CTIP2, CUX1/2, NEUN, SATB2 and TBR1 as compared to RNA expression, ATAC GAS and the CODEX imaging at P0 (Extended Data Fig. 5). Most of the antibodies in the ADT panel target immune proteins that were not detectable during development. Nevertheless, we were able to detect markers for monocytes, macrophages and microglia (CD11b, CD68 and CD86), dendritic cells (CD11c) and natural killer cells (CD49b) at P0 (Extended Data Fig. 5a).

After mapping cell states defined by previous scRNA-seq studies28,29,30 to our RNA data (Extended Data Figs. 3a and 4), we identified that the cluster localized to the CC corresponds to myelin-forming oligodendrocytes (MFOLs) (Fig. 1m and Extended Data Fig. 3b) and is represented only in the P10 and P21 samples, consistent with our CODEX data. Integrating single-cell ATAC-seq mouse brain atlas data31 with our P10 ATAC-seq data also identified all major cell types (Extended Data Fig. 3e,g). Furthermore, label transfer based on epigenetic states showed consistent identification of dominant cell types between ATAC-seq and RNA-seq (Extended Data Fig. 3d–g).

Spatial ARP-seq analysis of three human brain V1 region tissue sections, adjacent to those analysed using CODEX, identified 16 RNA and 9 ATAC clusters (Fig. 1d, Supplementary Fig. 20a and Supplementary Table 9). These data suggest that transcriptional profiles strongly define cell types during prenatal and neonatal brain development, with spatial patterns closely matching MERFISH data from adjacent tissue sections[17](https://www.nature.com/articles/s41586-025-09663-y#ref-CR17 “Wang, L. et al. Molecular and cellular dynamics of the developing human neocortex. Nature https://doi.org/10.1038/s41586-024-08351-7

(2025).“) (Supplementary Fig. 20d). RNA clusters R1, R2 and R9 (second trimester), R6 and R9 (third trimester) and R5, R3, R4 and R10 (infancy) represent cortical layers. Moreover, the ADT protein expression of CD56 (encoded by NCAM1) correlated with NCAM1 RNA expression and chromatin accessibility mainly in the second trimester, but to a lesser extent later in development (Supplementary Fig. 20b,c). RNA expression and ATAC GAS for RBFOX3, primarily in the cortical layers, and MBP, which began to dominate WM starting from the third trimester, were consistent with observations from CODEX (Fig. 1c and Supplementary Figs. 11–14 and 20b,c). Thus, our spatial ARP-seq and CODEX analyses in mouse and human brains provide a comprehensive dataset, capturing spatial conservation and temporal evolution across all three layers of the central dogma of molecular biology.

Cortical chromatin accessibility persistence

We developed a computational framework to systematically examine RNA and ATAC patterns crossing three dimensions, that is, space, time and modality (Fig. 2c and Methods). The framework starts with applying generalized additive regression to model the spatial (cortical layers II/III, IV, V and VI) and temporal (P0–21) effects, and their interaction on GAS for ATAC and gene expression for RNA. We undertook a two-step procedure to categorize patterns of genes that exhibit significant spatial/temporal changes in either RNA or ATAC. First, we performed joint clustering on concatenated RNA and ATAC data to simultaneously capture the spatiotemporal patterns across both modalities. Subsequently, we combined patterns that exhibited similar RNA and ATAC profiles. Using this framework, applied to two postnatal datasets separately (S1 and S2), we identified 13–15 distinct RNA and 7–10 ATAC patterns, resulting in 27–30 unique combinatorial gene sets that display varied patterns of changes across time and space (Fig. 2e,f, Extended Data Fig. 6, Supplementary Fig. 21 and Supplementary Table 7). Although each cortical dataset (S1 and S2) has robust similarity (Supplementary Fig. 17), subtle variation does occur, probably due to the pixel size difference and the addition of P7 timepoint in the second dataset.

Previous single-cell RNA-seq and ATAC–seq mapping of PNs through development from embryonic (E10.5) to early postnatal (P4) has identified gene expression programs that define progenitors (Eomes, Hes1, Hmgb2, Id4, Jun and Neurog2) and lineage-spanning signatures (Mn1, Mef2c, Neurod6, Nfib, Sox4, Sox11 and Zeb2)3. We verified the expression of expected genes and confirmed that early progenitor-restricted genes were not detected in any of our postnatal gene sets. We then categorized the 15 RNA clusters identified in the S2 datasets into three major groups based on shared spatial and/or temporal expression patterns. Layer-specific markers were generally found within the spatial and/or temporal categories, but also in the binary cluster exhibiting both temporal and spatial similarities (Fig. 2e and Supplementary Fig. 21).

We next investigated the distribution of the different PN subtypes. Expression profiles were also previously defined for PN subclasses, including CThPNs (Chgb,* Elavl2*,* Ndn*, Tbr1, Tcf4, Tle4, Zfp428), SCPNs (Bcl11b, Crym,* Etv1*,* Fezf2*, Ldb2, Pex5l, Thy1) and CPNs (Cux1, Cux2, Dok5,* Lhx2*,* Plxna4*, Ptprk, Satb2)3. We mapped these markers to determine whether their spatial segregation was accurate. Overall, we found that the upper cortical layer clusters (R1, R7, R8 and R13) were enriched for CPN markers. The deeper layers contained SCPN markers populating R2, R4/R5 and R15, and CThPNs that had a slightly broader cluster representation in R4/R5, R8, R12 and R13 (Supplementary Table 10). Each subclass of PNs had marker enrichment in RNA clusters related to space, in particular Ptprk (R1) for CPNs and Rorb (R2) for layer IV stellate neurons*, Etv1* (R5), Fezf2 (R5) and Ldb2 (R4) for SCPNs in layer V, and Elavl2 (R5), Ndn (R4) and Tbr1 (R5) for CThPNs in layer VI. By contrast, for R12–15 clusters, time was the main driver of RNA expression changes, with consistent temporal expression variation observed across a subset of layers. For example, Satb2 (R13), Plxna4 (R13) and Bcl11b (R13) presented such a pattern. The third group (R7–R11) exhibited a binary pattern influenced by both time and space; for example, Cux1, Cux2, Dok5 and Lhx2 were present in the binary RNA cluster R7 (Fig. 2e,g,h, Supplementary Figs. 22 and 23 (for S1) and Supplementary Table 10).

The ATAC GAS patterns exhibited a more diffuse spatial/temporal distribution, and the identified patterns did not always reflect the variations observed in the RNA clusters (Fig. 2e,f,h, Extended Data Fig. 6 and Supplementary Fig. 21). While RNA expression of many cortical layer-specific TFs (Bcl11b,* Cux1*, Cux2, Lhx2, Satb2 and Tle4) aligned with chromatin accessibility, we observed a misalignment for a subset of key TFs (Fig. 2g–i and Supplementary Figs. 22–25). Specific subsets of genes presented chromatin accessibility trailing in time and/or more diffuse in space, when compared to their RNA expression. ATAC signal for lineage spanning genes like Sox4, Sox11 and Mef2c and for layer-specific genes like Sox5 and Plxna4 trailed temporally when compared to RNA (Fig. 2g–i and Supplementary Figs. 22–25). We observed few instances where chromatin accessibility preceded gene expression, such as Lmo4 (Fig. 2h and Supplementary Fig. 25 (for S1)). Notably, Lmo4 has reported roles in astrocytes and may therefore also be relevant postnatally3,32,33. Thus, although epigenetic priming occurs during cortical layer development3, in the postnatal development, we observe a decrease in the expression of a subset of cortical genes—at both the RNA and protein levels—that is not concurrent with the ATAC signals. Thus, these genes may retain an epigenetic memory of their previous transcriptional states of an earlier developmental stage.

Cortical open chromatin spatial spreading

Cortical-layer specific TFs/cofactors such as Bhlhe22, Fezf2, Ldb2, Tshz2, Etv1, Foxp2 and Tbr1 exhibited chromatin accessibility spreading across layers (Fig. 2g–i and Supplementary Figs. 22–25). For example, Fezf2 and Tbr1 were quite restricted at the RNA expression level to layer V and VI, respectively. However, their gene accessibility score extends to layer VI (for Fezf2) and IV (for Tbr1) (Fig. 2h). The dynamics between these two factors has been recently documented in which the CThPNs (Tbr1+) in Tbr1-knockout mice are not distinct from their corticofugal counterpart, SCPN (Fezf2*+*)34. The reduced ATAC coverage of Tbr1 in layer V might be relevant, considering that Tbr1 is a direct transcriptional repressor of Fezf234. Conversely, Fezf2 accessibility in layer VI may be permitting, as Fezf2 is not a direct repressor of Tbr134. We have previously shown that the repressive histone mark H3K27me3 is deposited at the TSS of genes that are not expressed in specific cortical layers, and absent from genes that are expressed in those layers1,35. Thus, despite the spatially diffuse chromatin accessibility observed in postnatal cortical layers, Polycomb-mediated repression mechanisms may prevent ectopic expression of layer-specific genes.

Spatiotemporal RNA–ATAC cluster GO analysis

We also performed Gene Ontology (GO) analysis on each of the 27 combined RNA/ATAC patterns. General GO terms affiliated with neuronal biology were described across the patterns associated with CPNs, CThPNs and SCPNs. However, other terms were more specific to PN subtype associated patterns. Glutamatergic synaptic transmission, including NMDA and AMPA receptor activity terms, were found in four clusters of R2-A3/A5/A7 (IV stellate neuron and SCPN-associated) and R5-A2 (CThPN and SCPN-associated) (Supplementary Fig. 26). Moreover, other neurotransmission terms (calcium-ion-regulated exocytosis of neurotransmitter, neurotransmitter secretion and positive regulation of synaptic transmission) were also identified in the CThPN and SCPN-associated R5-A2/A5 and R15-A5/A8 clusters. CPN-associated clusters (R7-A3, R8-A5 and R13-A6/A7) had the terms integrin signalling, axonogenesis, dendrite morphogenesis and neuron migration. Thus, all clusters included biological processes involved in general developmental changes, with CPNs being less mature due to the ongoing axonal guidance programs, while CThPNs and SCPNs were more mature with features of exocytosis, and glutamatergic synaptic transmission already underway. The maturity regarding postnatal development is consistent with embryonic development order3.

Notably, we also identified clusters exhibiting GOs related to oligodendrogenesis13,36,37 in R4/R5 (OPCs/committed oligodendrocyte precursors (COPs)/newly formed oligodendrocytes (NFOLs): Fxyd6, Grin3a, Spon1, Pcdh15,* 2610035D17Rik*), R8 (OPCs/COPs/NFOLs: Ptprz1, Sox6, Frmd4a), R9 (OPCs/COPs/NFOLs: Lhfpl3,* Sdc3*; MOLs: Marcksl1), R13 (OPCs/COPs/NFOLs: Nnat,* Midn*,* Vcan*) and R15 (OPCs/COPs: Bcas1, Tnr,* Kcnip3*; MOLs: Cnp, Nrgn, S100b). NMDA receptor subunits indicating primed OPCs38 were also segregated to R4 (Grin3a) and R9 (Grin2d) (Fig. 2j, Supplementary Figs. 26–28 and Supplementary Table 10). As clusters R4 and R9 are CThPN and SCPN-related, our data suggest an overlap between the time and space of OPCs and corticofugal PN maturation. Myelin-associated genes identified across the gene sets were exclusively clustered to R15 (Mbp, Mag, Mobp, Plp1 and Mal), which is the final temporal program to arise beginning at P10, increasing at P21, and is associated with CThPNs and SCPNs (Fig. 2e,j and Supplementary Figs. 22 and 23).

There was a greater intensity of the spatial segregation of myelin-associated gene expression, where layer VI, but also V, had the highest levels for all the myelin genes such as Cnp, Mag and Mbp at P10 and P21. Chromatin accessibility for these myelination-related genes was already increased in several cortical layers as early as P5, suggesting chromatin priming for myelination during postnatal cortical development. Layer II/III was delayed in OPC differentiation and myelination programs (Fig. 2j and Supplementary Figs. 22 and 23), suggesting that myelination of upper-layer CPN neurons occurs later than for CThPN and SCPN neurons. Thus, our data indicate that OPCs and MOLs follow a similar spatial and temporal trajectory to that of CThPNs and SCPNs in early postnatal development while, at P10 and P21, the spatial myelination takes over, which is seen more profoundly in the deep layers where the SCPNs and CThPNs heavily populate.

Human cortical development benchmarking

We used our human V1 cortical maps