所有使用活动物的实验程序均由SALK Institute动物护理和使用委员会批准 ，该协议编号为18-00006。成人（p56）C57BL/6J雄性小鼠购自杰克逊实验室，并在12小时的深色灯光周期中维持在Salk Animal屏障设施中
，最多10天 。在冰冷的解剖缓冲液中，从额头上从额头上萃取大脑并切成冠状 ，并在含有2.5 mm kcl
，0.5 mm CaCl2，7 mm mgcl2 ，7 mm mgcl2
，1.25 mm naH2pO4，1.25 mm naH2pO4 ，110 mm soscrose
，110 mm铜，10 mm glucose和10 mm glucose和10 mm glucose和10 mmmmmmmmn的Co中提取冠状。在开始切片程序之前 ，将溶液保持冰冷
，并用95％O2燃烧至少15分钟
。将切片保存在包含冰冷解剖缓冲液的12孔板中（最多20分钟），直到解剖由配备SDF PLAPO 1XPF物镜的SZX16 Olympus显微镜帮助。Olympus Cellens维度1.8用于图像采集 。根据Allen Brain参考Atlas CCFV316 ，将每个大脑区域从切片中沿前后轴进行解剖（有关每个冠状切片的后视图，请参见扩展数据图1）
。在解剖过程中将切片保存在冰冷的解剖培养基中
，并立即将其冷冻在干冰中，以进行后层池和核产生。对于细胞核分离 ，将每个解剖区域从6-30只动物中汇总
，并为每个切片处理两个生物复制品 。 　　如先前所述1,6分离核
。通过与Alexa Fluor 488偶联的抗Neun抗体（MAB377X，Millipore）和1：1,000在4°C下在4°C下在4°C稀释1 h的抗体488浓度的抗体抗体（MAB3777X ，Millipore）通过1：1,000稀释的Alexa Fluor 488偶联抗体稀释1,000个稀释型
，在4°C下孵育1小时，将分离出的核标记。使用85μm喷嘴的BD涌入分把在22.5 psi鞘压力下进行单核的粉丝 。将单核分类的384孔板的每个孔都用2μL蛋白酶K消化缓冲液（1μLM-Digestion缓冲液（Zymo ，d5021-9）
，0.1μL20μgμgμl-1蛋白酶K和0.9μLH2O）进行了预加载
。通过将护套流动到空板的孔中并根据液滴位置进行调整来进行接收384孔板的比对。选择单细胞（单滴单）模式以确保排序的严格度 。对于每个384孔板，用neun+（488+）栅极和Neun-（488-）栅极对1-22列进行排序 ，达到Neun+与Neun-nuclei的比率为11：1
。BD涌入软件v1.2.0.142用于选择细胞群体。 　　先前针对SNMC-Seq21,2进行了详细的亚硫酸盐转化和库制备方法 。使用150 bp配对的末端模式
，使用Illumina Novaseq 6000仪器测序了由小鼠脑组织产生的SNMC-SEQ2和SN-M3C-SEQ（见下文）。Freedom Evoware v2.7用于库准备，Illumina Miseq控制软件v3.1.0.13和Novaseq 6000控制软件v1.6.0/实时分析（RTA）v3.4.4用于测序
。 　　如前所述4进行了单核甲基-3C测序（SN-M3C-Seq）。简而言之 ，从背侧齿状回（DG-1和DG-2，补充表2）
，腹侧齿状回（DG-3和DG-4）（DG-3和DG-4），背hip（Ca-1和ca-2） ，以及腹侧hip（Ca-3和Ca-4）中的相同批次的解剖组织样品
，是液体nitro的 。然后使用砂浆和杵冷冻时将样品粉碎，然后在DPBS中立即用2％甲醛固定10分钟
。将样品用0.2 M甘氨酸淬灭 ，并在-80°C下储存
，直到准备进一步加工。如前所述4分解核后，将核与Nlaiii消化过夜 ，并连接4 h 。然后用HOECHST 33342（但不用Neun抗体染色）将核染色
，并通过0.2-μm滤波器过滤，并与SNMC-SEQ2样品相似
。使用SNMC-Seq2方法生成库。 　　这项研究中的小鼠脑解剖和解剖结构的命名遵循艾伦小鼠脑公共坐标框架（CCF）16 。根据Allen CCF的层次结构 ，我们使用了一个三级空间区域组织来促进描述：（1）例如
，主要区域，例如等皮层 ，髋关节；（2）等级内的子区域，例如MOP
，SSP；（3）在MOP内进行解剖区域，例如MOP-1和MOP-2
。补充表1包含本研究中使用的所有缩写的全名。所有大脑图像均基于Wang等人的16和©2017 Allen Brain Science Institute创建 。艾伦大脑参考地图集
。可从：http：//atlas.brain-map.org/获得。 　　以下方法部分分为三个阶段 。第一阶段“映射和特征生成”描述了单细胞甲基化特异性数据格式的映射和生成文件。第二阶段“与聚类相关 ” ，描述了聚类
，识别DMG或集成其他数据集，这些数据集发生在单细胞级别
。第三阶段是“细胞类型特异性调节元件” ，描述了使用簇合并甲基甲基甲基甲基组的推定细胞类型特异性调节元素的鉴定。其他特定于图的分析主题可能会结合一个以上阶段的结果 。 　　我们为我们组开发的所有基于1,2,37的基于单一的基于 - 甲基 - 甲基群的技术实施了一条多功能映射管道YAP（https://hq-1.gitbook.io/mc/）
。该管道的主要步骤包括：（1）将FASTQ文件删除到单个单元格中；（2）读取水平质量控制（QC）；（3）映射；（4）BAM文件处理和QC；（5）最终分子剖面生成。先前描述了SNMC-SEQ2的五个步骤的详细信息2 。我们将所有读数映射到小鼠MM10基因组
。我们计算了映射后每个细胞中两组基因组区域的甲基胞嘧啶计数和总胞嘧啶计数。MM10基因组的非重叠染色体100-KB仓（由“ BedTools MakeWindows -W 100000
”产生）进行聚类分析和ANN模型训练
，并且使用小鼠Gencode VM22定义的基因身体区域±2 Kb用于聚类和集成其他模态 。 　　如先前所述4，在金牛座-MH管道后映射甲基甲基测序读数4
。具体而言 ，对Illumina适配器进行了修剪
，然后在两侧修剪了另外10 bp。然后，使用Bowtie的Bismark将R1和R2读数分别映射到MM10基因组 。收集了未封闭的读数 ，并将其分为较短的读数
，代表前40个bp，最后40 bp和原始读数的中间部分（如果修剪后读取长度> 80 bp）
。使用Bowtie的Bismark再次映射了分裂读数。用mapq读取 <10 were removed. The filtered bam files from split and unsplit R1 and R2 reads were deduplicated with Picard and merged into a single bam file to generate the methylation data. Methylpy (v1.4.2)38 was used to generate an ALLC file (base-level methylation counts) from the bam file for every single cell. We paired the R1 and R2 bam files where each read-pair represents a potential contact to generate the Hi-C contact map. For generating contact files, read pairs where the two ends mapped within 1 kbp of each other, were removed. 　　Cell filtering. We filtered the cells on the basis of these main mapping metrics: (1) mCCC level <0.03; (2) overall mCG level >0.5; (3) overall mCH level <0.2; (4) total final reads >500,000; and (5) Bismark mapping rate >0.5 。在过滤过程中还产生和评估了其他指标 ，例如基因组覆盖率，PCR重复速率和指数比。但是
，在以1-5的主要指标删除异常值后，发现其他异常值很少
。注意 ，MCCC水平用作硫酸硫酸硫酸非转化率上限的估计。 　　功能过滤 。通过去除平均总胞嘧啶基呼叫的垃圾箱过滤100 kb的基因组箱特征 <250 (low coverage) or >3,000（异常高覆盖区域）
。与编码BlackList39重叠的区域也被排除在进一步分析之外。 　　甲基化水平的计算和归一化 。对于CG和CH甲基化
，来自甲基胞嘧啶和总胞嘧啶基质的甲基化水平计算包含两个步骤：（1）（1）β-二比分布的先前估计，以及（2）每个细胞的后验计算和归一化
。 　　步骤1：对于每个单元 ，我们计算了原始甲基胞嘧啶水平（MC/COV）的样品平均值和方差
，其中COV是总胞嘧啶碱基覆盖率，而MC是每个序列上下文（CG或CH）的甲基胞嘧啶碱基覆盖率。然后 ，使用矩的方法估算β分布的形状参数（α
，β）： 　　该方法为每个细胞使用不同的先验，用于不同的甲基化类型 ，并使用更多信息（较高的原始方差）将较弱的先验使用 。 　　步骤2：然后
，我们计算后验：对于每个单元格中的所有垃圾箱
。就像单细胞RNA-seq分析中的百万读数（CPM）归一化一样，我们通过细胞的全局平均甲基化 ，M =α/（α+β）将这种后甲基化归一化。因此，归一化后
，所有具有0 COV的后值将其恒定值为1 。所得的归一化MC级矩阵不包含NA（不可用）值，并且COV较低的特征往往具有接近1的平均值
。 　　选择高度可变功能。使用scanpy.pp.highly_variable_genes函数从scanpy 1.4.3 package40选择了高度可变的甲基化特征 。简而言之 ，scanpy.pp.highly_variable_genes功能通过与落入给定垃圾箱的基因平均基因表达的基因分散体的平均分散体缩放来归一化基因的分散。在我们的修改方法中
，我们认为平均甲基化水平和特征的平均COV（100 kb箱或基因）都可能影响MC水平分散
。我们将落入平均值和COV的组合箱中的特征分组。然后，我们将每个均值-COV组内的色散标准化 。分散归一化后 ，我们选择了基于标准化分散的前3,000个特征进行聚类分析
。 　　不同MC类型的尺寸降低和组合。对于每个选定的功能
，将MC级别缩放到单位方差和零均值 。然后，我们对缩放的MC级矩阵进行了主成分分析（PCA）
。通过使用肘方法检查每个PC的方差比 ，选择主组件（PC）的数量（PC）。然后将CH和CG PC连接在一起
，以在聚类和多种学习中进行进一步分析（补充表6，用于PCA和聚类分析的参数） 。 　　在串联PC上共识聚类
。我们使用基于k-nearest邻居（KNN）图的多个leiden clustering41的共识聚类方法来说明莱顿聚类算法的随机性。在MCH和MCG矩阵中选择了PCA的主要PC之后 ，我们将PC串联在一起
，使用euclidean距离的Scanpy.pp.neighbours构造KNN图 。给定的固定分辨率参数，我们在KNN图上重复了Leiden聚类的300次 ，并将这些群集分配组合为新功能矩阵，其中每个单个Leiden结果都是一个功能
。然后
，我们使用Scikit-Learn软件包中使用的异常值DBSCAN算法来使用锤子距离在Leiden特征矩阵上进行共识聚类。DBSCAN的不同EPSILON参数遍历以生成共识群集版本，其簇数量从最小值到最大群集数量 ，在多个Leiden运行中观察到的群集数量 。每个版本都包含一些离群值；这些通常分为三类：（1）位于两个簇之间的细胞具有梯度差而不是透明边界
，例如，IT层的边界；（2）读取数量少的细胞潜在地缺乏基本特征的信息来确定特定群集；（3）具有大量读数的细胞。类型1和2异常值的数量取决于分辨率参数 ，并在“分辨率参数”部分的选择中进行了讨论
。细胞过滤后
，3型离群值非常罕见。然后，下面的监督模型评估确定了最终共识群集版本 。 　　对聚类分配的监督模型评估
。我们从Scikit-Learn软件包中进行了交叉验证（RFECV）42个过程 ，以评估每个共识聚类版本的群集可重复性。我们首先从此过程中删除了异常值
，然后将10％的单元格作为最终测试数据集 。对于其余90％的单元格，我们使用十倍的交叉验证使用输入PC作为特征和Sklearn.Metrics.balanced_accuracy_score43作为评估得分来训练多类预测模型
。多类预测模型基于Imblearn软件包中的BalancedRandomForestClassifier ，该软件包说明了分类不平衡问题44。训练后
，我们使用10％测试数据集使用Balanced_Accuracy_score的分数来测试模型性能 。我们将最佳模型和相应的聚类分配保留为最终聚类版本
。最后，我们使用此预测模型来预测离群值的集群分配。我们以预测概率> 0.3挽救了异常值 ，否则将其标记为异常值 。 　　可视化的流形学习
。在每一轮聚类分析中，使用Scanpy40软件包的实现
，在PC矩阵上运行T-SNE45,46和UMAP17嵌入。来自两种算法的坐标在补充表5中 。 　　分辨率参数的选择。选择莱顿算法的分辨率参数对于确定簇的最终数至关重要
。我们通过三个标准选择了分辨率参数：（1）。异常值的一部分 <0.05 in the final consensus clustering version; (2) the ultimate prediction model accuracy >0.9;（3）每个群集≥30的平均单元，它控制着簇的大小以达到进一步的表观基因组分析所需的最小覆盖范围 ，例如DMR调用 。这三个标准都阻止了集群的过度分裂
。因此
，我们选择了使用网格搜索符合标准的最大分辨率参数。 　　单元格（1级聚类）注释 。我们根据Neun-栅极起源和低全局MCH馏分对非神经元细胞注释
。鉴于CH甲基化和基因表达之间的强抗相关性，我们在基因体上使用了pan兴奋性标记的基因±2 kb（例如SLC17A7和SV2B）的低ch-甲基化 ，以及诸如GAD1和GAD2之类的泛抑制标记
，以分别注释兴奋性和抑制性细胞类别。 　　主要类型（2级）和亚型（第3级）注释 。我们使用了众所周知的标记基因的基因体±2 kb降低甲基化（或非神经元的低-CG-甲基化）和解剖信息来注释神经元和非神经元簇
。所有群集标记基因均在补充表7中列出，以及群集名称的描述 ，对标记基因信息的引用以及数据浏览器的URL。基于先前研究中报道的众所周知的标记基因
，对主要细胞类型进行了注释 。19,31,47,48,49
。只要有可能，我们就将这些群集命名为规范名称（例如IT-L23 ，L6B）或使用反映集群特定空间位置的描述性名称（例如
，EP，CLA ，IG-CA2）。对于子类型，我们通过其母体主要类型名称命名群集
，然后是子类型标记基因名称 。 　　我们使用成对策略在同一轮分析中计算每对簇的DMG。我们使用了所有蛋白质编码和长的非编码RNA基因的基因体±2 kb区域，这些基因具有来自小鼠Gencode VM22的证据1或2的基因
。我们使用了按全局MCH水平归一化的单细胞水平MCH级分（如以上聚类步骤中的“甲基化水平的计算和归一化 ”）来计算所有神经元簇之间的标记。我们使用按全局MCG水平归一化的MCG分数分别比较了非神经元簇 。对于每个成对比较 ，我们使用Wilcoxon秩和测试来选择具有显着降低的基因（低乙基化）
。根据调整后的P选择标记基因 < 10−3 with multitest correction using the Benjamini–Hochberg procedure, delta-normalized methylation level change <−0.5 (hypo-methylation) and area under the receiver-operating curve (AUROC) >0.8. We required each cluster to have ≥5 DMGs compared to any other cluster. Otherwise, the smallest cluster that did not meet this criterion was merged to the closest cluster based on Euclidean distance between cluster centroids in the PC matrix used for clustering. Then the marker identification process was repeated until all clusters found enough marker genes. 　　On the basis of the consensus clustering steps described above, we used an iterative approach to cluster the data into three levels of categories. In the first level, termed CellClass, clustering analysis is done using all cells and then manually merged into three canonical classes: excitatory neurons, inhibitory neurons, and non-neurons based on marker genes. Within each CellClass, we performed all the preprocessing and clustering steps again to obtain clusters for the MajorType level using the same stop criteria. Furthermore, within each MajorType, we obtained clusters for the SubType level. All clusters’ annotations and relationships are presented in Supplementary Table 7. 　　To build the taxonomy tree of subtypes, we selected the top 50 genes that showed the most significant changes for each subtypes’ pairwise comparisons. We then used the union of these genes from all subtypes and obtained 2,503 unique genes. We calculated the median mCH level of these genes in each subtype and applied bootstrap resampling-based hierarchical clustering with average linkage and the correlation metric using the R package pvclust (v.2.2)50. 　　We defined the impact score (IS) to summarize pairwise comparisons for two subtype groups, where one group, A, contains M clusters and the other group, B, contains N clusters. For each gene or motif, the number of total related pairwise comparisons is M × N, the number of significant comparisons with desired change (hypo-methylation for gene or enrichment for motif) is a in group A and b in group B. The IS is then calculated as  and for the two directions. For either group, IS ranges from −1 to 1, and 0 means no impact, 1 means full impact and −1 means full impact in the other group (Extended Data Fig. 7e). 　　We explored two scenarios using the IS to describe cluster characteristics (Extended Data Fig. 7e). The first scenario is considering each pair of branches in the subtype taxonomy tree as comprising group A and group B. Thus, the IS can quantify and rank genes or motifs to the upper nodes based on the leaves’ pairwise comparisons (Fig. 3d–f). The second scenario summarizes the total impact for specific genes or motifs regarding the taxonomy tree based on the calculation in the first scenario (Extended Data Fig. 7f–k). In a subtype taxonomy tree with n subtypes, the total non-singleton node was n − 1, and each node i had a height hi and associated ISA for one of the branches (ISB = −ISA). The node-height-weighted total IS (IStotal) was then calculated as: 　　The larger total IS indicated that a gene or motif shows more cell-type-taxonomy-related significant changes. The total IS can also be calculated in a sub-tree or any combination of interests to rank genes and motifs most related to that combination (See ‘Figure-specific methods’ for Fig. 5 regarding calculating layer and region total IS from the same tree). 　　A portion of the same brain tissue sample used in this study for methylome profiling was also processed using snATAC-seq in a parallel study of chromatin accessibility3. The final high-quality snATAC-seq cells were assigned to 160 chromatin accessibility clusters (a-types). The snATAC-seq-specific data analysis steps are described in Li et al.3. Here, we performed cross-modality data integration and label-transferring to assign the 160 a-types to the 161 methylome subtypes in the following steps: 　　We used the overlap score to match a-type and methylome subtypes. The overlap score, range from 0 to 1, was defined as the sum of the minimum proportion of samples in each cluster overlapped within each co-cluster52. A higher score between one methylome subtype and one a-cluster indicates they consistently co-clustered within one or more co-clusters. Besides matching clusters in integration analysis, the overlap score was also used in two other cases: (1) to quantify replicates and region overlaps over methylome subtypes (Extended Data Fig. 2e–g); and (2) to quantify the overlap of each L5-ET subtype overlapping with ‘soma location’ and ‘projection target’ labels from epi-retro-seq cells (Extended Data Fig. 5j) through integration with the epi-retro-seq dataset. 　　After clustering analysis, we used the subtype cluster assignments to merge single-cell ALLC files into the pseudo-bulk level and then used methylpy (v1.4.2)38 DMRfind function to calculate mCG DMRs across all subtypes. The base calls of each pair of CpG sites were added before analysis. In brief, the methylpy function used a permutation-based root mean square test of goodness of fit to identify differentially methylated sites (DMS) simultaneously across all samples (subtypes in this case), and then merge the DMS within 250 bp into the DMR. We further excluded DMS calls that have low absolute mCG level differences by using a robust-mean-based approach. For each DMR merged from the DMS, we ordered all the samples by their mCG fraction and calculated the robust mean m using the samples between 25th and 75th percentiles. We then reassigned hypo-DMR and hyper-DMR to each sample when a region met two criteria: (1) the sample mCG fraction of this DMR is lower than (m − 0.3) for hypo-DMR or (m + 0.3) for hyper-DMR, and (2) the DMR is originally a significant hypo- or hyper-DMR in that sample judged by methylpy. DMRs without any hypo- or hyper-DMR assignment were excluded from further analyses. On the basis of these filtering criteria, we estimate the false discovery rate of calling DMRs is 2.7% (Supplementary Note 2, Extended Data Fig. 6). 　　We performed enhancer prediction using the REPTILE53 algorithm. REPTILE is a random-forest-based supervised method that incorporates different epigenomic profiles with base-level DNA methylation data to learn and then distinguish the epigenomic signatures of enhancers and genomic background. We trained the model in a similar way as in the previous studies8,53, using CG methylation, chromatin accessibility of each subtype and mouse embryonic stem cells (mouse ES cells). The model was first trained on mouse ES cell data and then predicted a quantitative score that we termed enhancer score for each subtype’s DMRs. The positives were 2 kb regions centred at the summits of the top 5,000 EP300 peaks in mouse ES cells. Negatives include randomly chosen 5,000 promoters and 30,000 2-kb genomic bins. The bins have no overlap with any positive region or gene promoter8. 　　Methylation and chromatin accessibility profiles in bigwig format for mouse ES cells were from the GEO database (GSM723018). The mCG fraction bigwig file was generated from subtype-merged ALLC files using the ALLCools package (https://github.com/lhqing/ALLCools). For chromatin accessibility of each subtype, we merged all fragments from snATAC-seq cells that were assigned to this subtype in the integration analysis and used deeptools bamcoverage to generate CPM normalized bigwig files. All bigwig file bin sizes were 50 bp. 　　We used 719 motif PWMs from the JASPAR 2020 CORE vertebrates database54, where each motif was able to assign corresponding mouse transcription factor genes. The specific DMR sets used in each motif-enrichment analysis are described in figure specific methods below. For each set of DMRs, we standardized the region length to the centre ±250bp and used the FIMO tool from the MEME suite55 to scan the motifs in each enhancer with the log-odds score P <10−6 as the threshold. To calculate motif enrichment, we use the adult non-neuronal mouse tissue DMRs10 as background regions unless expressly noted. We subtracted enhancers in the region set from the background and then scanned the motifs in background regions using the same approach. We then used Fisher’s exact test to find motifs enriched in the region set and the Benjamini–Hochberg procedure to correct multiple tests. We used the TFClass56 classification to group transcription factors with similar motifs. 　　To calculate DMR–DMG partial correlation, we used the mCG fraction of DMRs and the mCH fraction of DMGs in each neuronal subtype. We first used linear regression to regress out variance due to global methylation difference (using scanpy.pp.regress_out function), then use the residual matrix to calculate the Pearson correlation between DMR and DMG pairs where the DMR centre is within 1 Mb of the TSSs of the DMG. We shuffled the subtype orders in both matrices and recalculated all pairs 100 times to generate the null distribution. 　　After merging the chromatin contacts from cells belonging to the same type, we generated a .hic file of the cell-type with Juicer tools pre. HICCUPS57 was used to identify loops in each cell type. The loops from eight major cell types were concatenated and deduplicated and used as the total samples for differential loop calling. A loop-by-cell matrix was generated, in which each element represents the number of contacts supporting each loop in each cell. The matrix was used as input of EdgeR to identify differential interactions with ANOVA tests. Loops with FDR <10−5 and minimum–maximum fold change >2 were used as differential loops. Note that the abundance of cell types is highly variable, leading to different coverages of contact maps after merging all the cells from each cell type. Since HICCUPS loop calling is sensitive to the coverage, more loops were identified in the abundant cell types (for example, 12,614 loops were called in DG, containing 1,933 cells) compared to the less abundant ones (for example, 1,173 loops were called in MGE, containing 145 cells). Therefore, we do not compare the feature counts related to the loops across cell types directly in our analyses. 　　3D model of dissection regions (Fig. 1b–e). We created in silico dissection regions based on the Allen CCF16 3D model using Blender 2.8 that precisely follow our dissection plan. To ease visualization of all different regions, we modified the layout and removed some of the symmetric structures, but all the actual dissections were applied symmetrically to both hemispheres. 　　Calculating the genome feature detected ratio (Extended Data Fig. 2a). The detected ratio of chromosome 100-kb bins and gene bodies is calculated as the percentage of bins with >20总胞嘧啶覆盖率。非重叠的染色体100-kb垃圾箱是由Bedtools Makewindows-W 100000生成的；基因体由Gencode VM22定义 。 　　与Epi-Retro-Seq L5-ET细胞的整合（图2G – J
，扩展数据图5G – J）
。Epi-Retro-Seq是一种基于SNMC-Seq2的方法，结合了逆行AAV Labelling22。Epi-Retro-seq数据集收集的L5-ET细胞的非重叠染色体100-kb bin基质与本研究中的所有L5-ET细胞串联 ，用于共聚类和嵌入
，如“聚类相关方法
”中所述 。然后，我们计算了本研究中亚型与Epi-Retro-Seq细胞的“ SOMA位置”或“投影目标 ”标签之间的OS
。第一个OS有助于量化两项研究之间的空间位置的一致性。第二个操作系统使我们能够在这项研究中将亚型的投影靶标算 。 　　成对DMR和基序富集分析（图3C ，F）
。通过比较所有亚型
，如“细胞类型的调节元件
”中所述鉴定了总亚型DMR。然后，如果DMR为：（1）仅在其中一种亚型中 ，我们将DMR分配给每个亚型对 。（2）两个亚型之间的MCG分数差异> 0.4。每个亚型对与两组成对DMR相关联
。我们在每个DMR集中使用另一组作为背景的“细胞类型调节元件”中描述的基序 - 富集分析。然后使用富含任一方向的图案计算影响评分
，并与分类学上的上点相关 。 　　与基因组区域重叠的EDMR（图4B）
。在Li等人3中鉴定出簇特异性的SNATAC-SEQ峰。我们使用床托合并来汇总总非重叠峰区域，而床托相交以计算峰和EDMR之间的重叠 。在HE等人中 ，使用甲基塞克58在HE等人中鉴定了发育中的前脑和其他组织FEDMR
。图4B中使用的所有基因组特征均定义为He等人3，除了使用更新的MM10 CGI区域和ReponMaster转换元素列表（UCSC Table浏览器于2019年10月9日下载）。 　　基因 - 增强剂景观的热图（扩展数据图8E） 。通过> 0.3的EDMR – Gene相关性选择每个基因的EDMR
。扩展数据中的热图的部分通过（1）本研究中的161个亚型中的每个EDMR的MCG分数收集。（2）每个EDMR的SNATAC-SEQ子类型级片段每千座映射片段（fpkm）在相同的亚型订单中 。如“与聚类相关方法 ”中所述的集成结果合并了亚型SNATAC轮廓
。（3）从胚胎第10.5天（e10.5）到p0的十个发育时间点
，前脑组织中每个EDMR的MCG分数（来自He等人的数据）；（4）每个EDMR的H3K27AC FPKM在7个从E11.5到P0的时间点中的7个开发时间点（Gorkin等人的数据）；（5）p56额脑组织中每个EDMR的H3K27AC FPKM（Lister等人的数据）；（6）EDMR使用BedTools相交与前脑FEDMR重叠。 　　与染色体相互作用的细胞嵌入（图4E） 。Schicluster60用于生成SN-M3C-Seq细胞的T-SNE嵌入。具体而言，为每个细胞的每个染色体生成以1-MB分辨率的触点矩阵
。然后 ，通过线性卷积将矩阵平滑
，板= 1，然后随机行走 ，重新启动概率= 0.5。将平滑地图上最强相互作用的最强相互作用的前20个百分点提取
，将二进化并用于PCA 。前20个PC用于T-SNE
。 　　IT层解剖区域DMG和DMR分析（图5A – C）。为了收集足够的细胞进行解剖区域分析，我们仅使用了IT神经元的主要类型（对应于L2/3 ，L4
，L5和L6） 。我们根据层解剖区域将细胞分组分组，并保持> 50个细胞的组以进行进一步分析（扩展数据图9b）
。我们进行了成对的DMG ，DMR和基序 - 富集分析
，与图3中的亚型分析相同，但使用层解剖区域组标签。然后 ，我们为这些群体建立了空间分类法，并用它来计算影响分数 。为了分别对与层相关或解剖相关的基因和图案进行排列
，我们使用了两组分类（扩展数据图9C，图9C ，层的顶部集合
，底部为区域），并使用上述方程式计算了两个总影响得分
。 　　DG细胞组和梯度DMR分析（图5E）。根据细胞的整体MCH水平将DG细胞分为四个均匀大小的组 ，截止阈值为0.45％
，0.55％和0.69％ 。然后，我们使用“与聚类相关方法
”中描述的方法随机选择了每个组中的400个单元
。为了确保在DG组之间识别的DMR不是由于随机性引起的 ，我们还随机采样了来自所有DG细胞的400个细胞的15组
，无论其全局MCH如何，并将其称为Control-DMRs（使用相同的滤波条件）为Control-DMR（2,003）。只有0.04％的梯度DMR与对照DMR重叠；这些被从进一步分析中删除 。根据线性序列（1 、2、3
、4）计算每个梯度DMR的MCG级分的Pearson相关性（ρ） ，以量化梯度趋势。DMR具有ρ <−0.75 or ρ >0.75被认为是显着相关的
。在进一步分析中未包括弱相关的DMR（占DMR的10％）。 　　DMR-和DMS富集的基因（图5F
，G） 。为了研究特定基因体中相关的DMR或DMS富集，我们使用Fisher的精确测试将基因体内的DMS和胞嘧啶的数量与±1 MB区域中的DMS和胞嘧啶的数量进行了比较
。我们选择了通过两个标准的基因：（1）调整后的P <0.01 with multitest correction using the Benjamini–Hochberg procedure, and (2) overlap with >20个DMS。DMR和DMS富基因的基因本体分析使用Goatools61进行 。所有具有基因体长度> 5 kb的蛋白质编码基因均用作防止基因长度偏置的背景
。 　　隔室强度分析（扩展数据图10G）。我们通过Z分数在DG接触矩阵的每个1-MB箱中通过Z评分进行了总染色体接触 ，并保留具有归一化覆盖率在-1和2之间的垃圾箱进行分析 。过滤后，将全基因组骑士 - 瑞氏级别范围62接触矩阵的PC1用作隔室评分
。分数分为50个类别
，尺寸从低到高，垃圾箱被分配给类别。每组的染色体内观察/预期（OVE）矩阵用于量化车厢强度 。我们计算了每对类别中的平均OVE值 ，以生成50×50的鞍矩阵
。计算隔室强度的平均左上和右下10×10矩阵除以右上和左下10×10矩阵的平均值63。 　　域分析（扩展数据图10i） 。我们使用Arrowhead57在DG中以10-kb分辨率以10-kb分辨率确定了4,580个接触域。对于bin
，绝缘得分I由 　　其中a是骑士 - 瑞式标准级矩阵的ove，平均值是亚字符中范围内的平均值
。对于每个组 ，在所有边界上都计算并平均域边界和100 kb侧翼区域的绝缘得分。 　　与图6有关
，为了降低计算复杂性，我们在100-kb bin MCH功能的数据集上应用了PCA ，以获取前3,000个PC
，该PC保留了原始数据方差的61％ 。然后，将这3,000个PC用于训练和测试预测模型
。我们使用一个具有两个隐藏层的ANN同时预测细胞亚型及其解剖区域。输入层包含3,000个节点 ，然后包含一个带有1,000个节点的共享层 。共享层同时连接到该小节区域亚型的两个分支隐藏层
，每个层都包含200个节点
。相应的单热编码输出层遵循分支隐藏的图层。我们使用五倍的交叉验证来访问模型性能 。我们在每个隐藏层上使用辍学率P = 0.5的辍学技术64，以防止在训练过程中过度拟合
。Adam Optimization65用于训练网络具有交叉渗透损失函数。训练时期的数字和批量大小分别为10和100 。训练和测试过程是通过Tensorflow 2.066进行的
。 　　模型性能。这两个输出层分别为每个单个单元的输入生成两个概率向量 ，分别是细胞亚型和解剖区域的预测 。使用最高概率的亚型和解剖区域标签用作每个细胞计算准确性的预测结果。当计算细胞解剖区域的精度（图6C）时，我们以不同的严格度定义了两种精度：（1）使用预测标签的精确精度
，以及（2）使用预测的标签或其潜在的重叠邻居的模糊精度
。潜在的重叠邻居基于Allen CCF（扩展数据图11A，补充表2）策划了特定解剖区域的相邻区域。ANN模型的确切精度为69％ ，模糊精度为89％ 。为了通过ANN评估提高解剖区域的精度的多少
，我们仅根据剖析区域组成的每个亚型中的天真猜测来计算模糊精度（扩展数据中的灰点图11c）
。我们还使用逻辑回归和基准的随机森林训练了其他模型。ANN对亚型预测的性能与逻辑回归和随机森林相当 。相比之下，针对其他两个模型（扩展数据图11b） ，位置预测中的性能显着提高
，这表明将细胞与不同解剖区域区分开可能需要基因组区域之间的非线性关系
。我们将Scikit-Learn（V0.23）用于逻辑回归和随机森林实现以及多级分类的多项式目标函数。随机森林的n_estimators设置为1,000 。 　　生物学特征对解剖区域预测的重要性（图6E）
。为了评估哪些DNA区域存储可使用我们的模型区分的细胞空间起源的信息，我们通过检查跨单元的每个PC置换量如何影响预测准确性 ，从而评估了PC特征的重要性。我们测试了每个功能的五个排列
，并使用降低的平均精度表明PC特征的重要性 。我们检查了100 kB箱中包含的基因，其中最重要的PC特征是给定细胞类型的最重要的PC功能
。 　　有关研究设计的更多信息可在与本文有关的自然研究报告摘要中获得。