所有动物工作均由劳伦斯·伯克利国家实验室动物福利和研究委员会或加利福尼亚大学戴维斯机构动物护理和使用委员会审查和批准。

　　从雌性C57BL/6N Mus musculus的不同发​​育阶段的胚胎中剖析了小鼠胎儿组织。用于在E14.5和P0上获得组织样品的小鼠购自Charles River Laboratories（C57BL/6NCRL菌株）和Taconic Biosciences（C57BL/6NTAC菌株）。用于在剩余发育阶段获得组织样品的小鼠购自Charles River Laboratories（C57BL/6NCRL菌株）。收集的胚胎或P0幼崽的数量取决于材料是否足以进行基因组测定，并且不是基于统计考虑的。在每个阶段，每个组织的每个复制物收集15至120个胚胎或幼犬之间。

　　有关详细信息，请参见补充文件1、2。

　　如前所述8，构建了Methyc-Seq库，并提供详细的协议50。所有WGB都使用100或130台单端读数用于所有WGB。

　　对于本研究的所有分析，我们使用MM10作为参考基因组，其中包括19个常染色体和两个性别染色体（对应于编码门户中的“ MMM10-Minimal”参考，https://wwww.encodeproject.org/）。MM10的FASTA文件是从UCSC基因组浏览器（2013年6月9日）下载的51。

　　如前所述52，所有WGBS数据均映射到MM10小鼠参考基因组。WGB的处理包括映射硫酸含量的噬菌体lambda基因组尖峰，以估计硫酸钠非转换率的对照。该管道（称为Methylpy）可在GitHub（https://github.com/yupenghe/methylpy）上获得。简而言之，首先将WGBS读数中的细胞选择转换为胸腺素。然后将转换后的读数通过Bowtie（1.0.0）对齐到C – T转换的参考基因组的正向链，并分别将G -A转换的参考基因组的反向链分别对齐。我们滤除了并非唯一映射或映射到两个计算转换基因组的读取。接下来，删除了PCR重复读取。最后，在相应的参考基因组序列（MM10或LAMBDA）中，甲基属于每个胞嘧啶位置的甲基化底膜（细胞选择）和未甲基化的底膜（胸腺素）。

　　计算甲基化水平以测量单细胞或较大基因组区域的DNA甲基化的强度和程度。甲基化水平定义为甲基化的基本计数的总和与给定区域中的一个胞嘧啶或跨部位的甲基化和未甲基化的基本计数的总和，从而减去了双硫酸钠非转换率。亚硫酸钠非转化率定义为硫酸氢钠处理的lambda基因组的甲基化水平。

　　我们在CG上下文和CH上下文（H = A，C或T）中计算了该指标。前者称为CG甲基化（MCG）水平或MCG水平，后者称为CH甲基化（MCH）水平或MCH水平。

　　我们计算了所有WGB数据的几个质量控制指标，结果在补充表1中列出。对于每个组织样本，我们计算了胞质覆盖率，Bisulfite钠的转化率以及生物复制之间的可重复性。胞质覆盖范围是覆盖胞嘧啶的平均读数。在计算中，我们组合了两个链的数据。亚硫酸钠转化率测量了硫酸钠转化效率，并计算为一个减去未甲基化Lambda基因组的甲基化水平。生物重复的可重复性定义为生物学重复之间的Pearson相关系数，至少十个读数涵盖的位点。

　　所有WGBS数据传递的编码标准（https://www.encodeproject.org/data-standards/wgbs/），并被Encode Consortium接受。组织样品的几乎所有生物学重复都至少具有30×胞嘧啶覆盖率。所有生物学重复都至少含有99.5％的硫酸钠转化率。所有非肝组织样品的可重复性大于0.8。肝样品的可重复性略低，但仍大于0.7。降低的可重复性是由于采样变异的增加，这是肝脏基因组中全基因组低甲基化的结果。

　　使用编码统一的处理管道处理芯片序列的芯片数据。简而言之，首先，使用BWA54（版本0.7.10）将Illumina读取映射到MM10参考文献，该参数为“ -Q 5 -L 32 -K 2”。接下来，使用PICARD工具（http://broadinstitute.github.io/picard/，版本1.92）用于使用以下参数删除PCR重复项：'remove_duplicates = true'。

　　我们将每个组蛋白修饰标记表示为整个基因组100 bp垃圾箱的连续富集值。减去芯片输入后，富集定义为RPKM。使用DeepTools255（2.3.1）的BAMCompare使用选项使用“ –binsize 100 normaligeSustrrpkm – ExtendReads 300 – ratatio减法”来计算整个基因组的富集。对于转录共激活器EP300（E1A相关蛋白p300）的芯片 - 隔离数据，我们使用MacS56（1.4.2）使用默认参数调用峰。

　　从所有阶段的所有胎儿组织处理的RNA-seq数据都从编码门户（https://www.encodeproject.org/;补充表2）下载。

　　To further validate our findings regarding transcriptomes generated across the Wold and Ecker laboratories, we generated an additional two replicates of RNA-seq data for fetal forebrain, midbrain, hindbrain and liver tissues. We first extracted total RNA using the RNeasy Lipid tissue mini kit from Qiagen (cat no. 74804). Then, we used the Truseq Stranded mRNA LT kit (Illumina, RS-122-2101 and RS-122-2102) to construct stranded RNA-seq libraries on 4 μg of the extracted total RNA. An Illumina HiSeq 2500 was used to sequence the libraries and generate 130-base single-ended reads.

　　RNA-seq data were processed using the ENCODE RNA-seq uniform processing pipeline. In brief, RNA-seq reads were mapped to the mm10 mouse reference using STAR57 aligner (version 2.4.0k) with GENCODE M4 annotation58. We quantified gene expression levels using RSEM (version 1.2.23)59, expressed as TPM. For all downstream analyses, we filtered out non-expressed genes and only retained genes that showed non-zero TPM in at least 10% of samples.

　　ATAC–seq data for all fetal tissues from all stages were downloaded from the ENCODE portal (https://www.encodeproject.org/; Supplementary Table 2). ATAC–seq reads were mapped to the mm10 genome using bowtie (1.1.2) with flag ‘-X 2000–no-mixed–no-discordant’. Then, we removed PCR duplicates using samtools60 and mitochondrial reads. Next, we converted read ends to account for Tn5 insertion position by moving the read end position by 4 bp towards the centre of the fragment. We converted paired-end read ends to single-ended read ends. Last, we used MACS2 (2.1.1.20160309) with flags ‘—nomodel —shift 37 —ext 73 —pval 1e-2 -B —SPMR —call-sumits’ to generate signal track files in bigwig format. MACS2 calculated ATAC–seq read fold enrichment over the background MACS2 moving window model. This fold enrichment is used as the intensity/signal of chromatin accessibility.

　　在这项研究中，我们使用了Gencode M458基因注释。从UCSC基因组浏览器下载了CGI注释（2016年9月5日）51。CGI海岸定义为上游2 kb，沿CGIS下游2 KB区域。启动子被定义为TSS周围-2.5 kb至+2.5 kb的区域。CGI启动子被定义为与CGI重叠的启动子，其余启动子称为非CGI启动子。

　　我们还使用以下过程获得了可映射的可载元元素（TE）的列表。从UCSC基因组浏览器下载了MM10鼠标基因组的重复邮件器注释（2016年9月12日）51。注释包括5,138,231重复。我们通过仅选择属于以下重复类之一的重复序列（replclass）：“ DNA”，“ SINE”，“ LTR”或“ LINE”来获得转座子注释。然后，我们将任何重复元素排除在其名称（repname），class（replclass）或家人（repfamily）中的任何重复元素。对于剩下的3,643,962个转座子，我们进一步过滤了含有少于两个CG位点或情况的元素，其中当将两个重复数据的数据组合在一起时，所有样品中至少读取了少于60％的CG位点。最后，我们使用了剩余的1,688,189个可映射的转座子进行本研究的分析。

　　如前所述52，我们使用甲基菌（https://github.com/yupenghe/methylpy）鉴定了CG-DMR。简而言之，我们首先称为DMS，然后将它们合并成块，如果它们都显示出相似的样品特异性甲基化模式并且在250 bp之内。最后，我们过滤了少于三个DMS的块。在此过程中，我们合并了来自所有组织的两个生物学重复的数据，不包括由于基因组的全局甲基化而导致的肝样品。

　　我们使用BedTools61（v2.27.1）的“相互构成”（v2.27.1）与先前鉴定的CG-DMR与先前鉴定的CG-DMR重叠。参考文献的CG-DMR的MM9坐标。首先使用带有默认参数的LiftOver51将11映射到MM10。在比较列表之间的基因组坐标时，CG-DMR的重叠定义为CG-DMR，至少一个基础与另一个CG-DMR重叠的CG-DMR。

　　对于每种胎儿组织类型，我们将组织特异性的CG-DMR定义为在任何胎儿阶段（E10.5至P0）中表现出在组织样本中降低甲基化的CG-DMR。仅在基线的基线中，低甲基化才有意义，因此我们使用跨组织样品的每个CG-DMR的基线MCG水平在组织样品跨组织样品中定义了基线MCG水平，该平均数被定义为包括所有样品一半的最窄MCG水平的值。具体而言，是组织样品S中的CG-DMR I（i = 1，…，m）的MCG水平（s = 1，…，n）。假设对样品进行排序，以便将基线定义为，其中a是样本索引，以至于最小化，即。定义为大于或等于N/2的最小整数。最后，我们将次甲基化样品定义为样品，其中CG-DMR I处的MCG水平比基线BI小0.3，即。然后，CG-DMR I特定于这些组织。该分析中未包括肝数据，我们排除了在任何非肝脏样品中覆盖范围为零的CG-DMR。总共只过滤了402个CG-DMR（约0.02％）。

　　我们根据基因组距离将CG-DMR与其推定的靶基因联系起来。首先，我们仅考虑表达的基因，这些基因在至少10％的所有胎儿组织样品中显示出非零TPM。接下来，我们使用BedTools61中的“盖帽”获得了表达基因TSS的坐标，并将每个CG-DMR与最接近的TSS配对。这样，我们推断了每个CG-DMR的靶基因。这些基因– TSS关联在本研究的所有后续分析中使用。

　　Reptile43算法用于识别显示出增强子样染色质特征的CG-DMR。我们称这些FedMrs。爬行动物使用随机的森林分类器来学习，然后区分增强子和基因组背景的表观基因组特征。爬行动物的一个独特特征是，通过合并其他样本的数据（作为外组/参考），它可以使用表观基因组变异信息来改善增强子预测。在这项研究中，使用来自CG甲基化（MCG），染色质可及性（ATAC -SEQ）和六个组蛋白标记的输入数据（H3K4ME1，H3K4ME2，H3K4ME3，H3K27AC，H3K27AC，H3K27ME3和H3K9AC）运行。

　　先前以类似的方式对爬行动物增强剂模型进行了训练。43。简而言之，在小鼠胚胎干细胞（MES细胞）和所有八个E11.5小鼠组织的甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基甲基疗法中进行简而言之，简而言之。CG-DMR必须至少包含两个DMS，并且它们在每个方向（5'和3'）上延伸150 bp。使用E11.5小鼠组织在MES细胞数据上训练爬行动物模型作为外组。这些样品可用于MCG和六个组蛋白修饰的数据。培训数据集由5,000个积极实例（推定的已知增强剂）和35,000个负面实例组成。阳性是在MES细胞中前5,000个EP300峰的山顶的2-KB区域。负面因素包括5,000个随机选择的启动子和30,000个随机选择的2-KB基因组垃圾箱。垃圾箱与任何积极因素或启动子没有重叠。爬行动物学会了将积极实例与负面实例区分开的染色质特征。

　　接下来，使用此增强器模型，我们将爬行动物应用于从所有非肝组织中鉴定的1,808,810 CG-DMR中划定FedMrs。根据MCG的数据和六个核心组蛋白标记，预测了每个样本的FEDMR，而其余的非肝样品则被用作外组。在爬行动物中，CG-DMR的随机森林分类器在每个样本中为每个CG-DMR分配置信度范围从0.0到1.0。该分数对应于随机森林模型中决策树的比例，该模型投票支持CG-DMR是增强剂。先前的基准测试表明，得分越高，CG-DMR显示增强剂活性43。我们命名了这个信心分数增强子得分。对于每个组织样品，FEDMR被定义为增强子得分大于0.3的CG-DMR。在每种组织类型中，FEDMR也被定义为CG-DMR，在该组织类型的至少一个组织样品中被识别为FEDMR。例如，如果仅在E14.5前脑中预测CG-DMR为FEDMR，则将其归类为前脑特异性FEDMR。

　　我们与Ref的假定成人增强剂重叠。26。我们使用一组坐标来从http://mouseencode.org/publication/mcp00/识别每种组织和细胞类型的假定增强剂的中心基础位置。接下来，我们将推定的增强剂定义为中心周围的±1-kb区域。如果它们重叠，则将来自不同组织和细胞类型的推定增​​强子组合在一起。然后使用LiftOver51将合并的推定增强剂（MM9）映射到MM10参考。最后，使用BedTools61的“相交”被用来与这些推定的增强剂重叠。

　　我们使用了Vista增强器浏览器28的增强器数据来估计在体内显示增强子活性的FEDMR的比例。具体而言，我们计算了经过实验验证为增强剂的FEDMR重叠远景元素的比例，我们称其为真正的正速率。我们评估了FEDMRS的真正正速率（前脑，中脑，后脑，心脏，肢体和神经管），其中至少30个Vista元素已在实验中被用作增强剂（阳性）。

　　但是，远景元素的选择是有偏见的。与随机选择的序列相比，它们对增强子更加丰富，这将导致对真正的正速率的高估。为了减少选择偏见的效果，我们需要首先估计如果选择最小的选择偏见，在给定组织中呈阳性（正率）的远景元素的分数。我们称这一小部分为真正的正率。详细信息可以在补充注释4中找到。然后，我们可以采样当前的Vista数据集，以构建与真实正率相匹配的正速率的数据集。由于在构建的数据集中没有膨胀正率，因此它将允许对我们的增强器预测方法进行公平评估（另请参见补充说明4有关详细信息）。

　　使用偏置控制的数据集，我们计算了每个E11.5组织的FEDMR的真实正速率。首先，我们通过其增强子得分（从最高到最低）对FEDMR进行了排名。然后，我们将带有Vista元素的给定E11.5组织的前2,500（或前2,501-5,000个顶部2,501-5,000）重叠，要求将至少一个FEDMR完全包含，以使Vista元素被计算为重叠。最后，我们计算了通过实验验证的给定组织中实验验证的增强子（即，真正的正速率）的vista元素的比例。

　　为了更好地解释FEDMR的真正正速率，我们还评估了5,000个随机选择的基因组垃圾箱，其GC含量和进化节约程度（直属评分）与前5,000名FEDMR相匹配。我们将此方法用作基线。对于每个E11.5组织，我们重复了此随机选择过程十次，并生成了十组随机区域。接下来，我们计算了偏置控制数据集中每组随机区域的真实正率。作为另一种基线方法，我们还计算了与任何FEDMRS或H3K27AC峰重叠的远景元件的正速率。

　　染色质状态的假定增强子是基因组区域，在非肝组织样品中通过Chromhmm63标记为增强子态（状态5、6和7）（参考文献23）。为了将其验证率与FedMR的验证率进行比较，我们需要选择前2500个推定的增强器。Chromhmm没有分配分数，因此我们使用H3K27AC信号对这些元素进行排名。然后，我们计算了与FEDMR重叠的前2500个推定增强器的比例。

　　为了测试FEDMR是否比染色质状态捕获更多的增强子，我们计算了非重叠FEDMR的验证率。此外，我们通过与Vista元素重叠来计算Chromhmm增强子的验证率。这用作评估FEDMRS的附加基线。

　　为了识别富含FEDMR的TF基序，我们扫描了基因组以描绘TF基序的发生，如前所述33。In brief, we used TF binding position weight matrices (PWMs) from the MEME motif database (v11, 2014 Jan 23. motif sets chen2008, hallikas2006, homeodomain, JASPAR_CORE_2014_vertebrates, jolma2010, jolma2013, macisaac_theme.v1, uniprobe_mouse,WEI2010_MOUSE_MW，WEI2010_MOUSE_PBM，ZHAO2011）。然后，使用FIMO64来扫描基因组，以使用选项“ –Output-Pthresh 1E-5-玛克斯存储得分500000”识别TF基序。

　　接下来，我们进行了超几何测试，以识别明显的基序富集。对于每种组织类型，我们计算了该组织（前景）FEDMR的基序富集，以识别出针对与前景组织列表重叠的其他组织识别的FEDMR列表。对于此分析，我们将前景和背景的平均尺寸扩展到400 bp，以避免由于尺寸差异而造成偏差。对于给定的组织t，前景和背景FEDMR的总数分别为nf，t和nb，t，nt = nf，t+nb，t是FEDMR的总数。对于给定的TF结合基序M，TF基序的发生与NF，T，M前景和NB，T，M背景FEDMR重叠，而NT，M = NF，T，T，M+NB，T，M是重叠的FEDMRS的总数。观察NF，T，M或更多重叠前景FEDMR（P）的概率定义为：

　　对于每种组织类型，我们对所有基序进行了此测试（n = 532）。然后，使用Benjamini -Hochberg方法调整每个组织的P值，如果通过1％的错误发现率（FDR）截止，则将基序称为显着。最后，我们排除了其TF表达水平小于10 tpm的任何TF结合基序。结果在补充表7中列出。

　　对于每个组织阶段，我们使用Great65找到了在该组织中鉴定出的FEDMR附近的基因的富集途径和生物过程。对于每个组织阶段，在10,000个增强子得分最高的FEDMR的“单一基因”关联策略下进行了伟大。补充表8列出了出色的分析结果。

　　我们应用了分层的LD评分回归47来测试FEDMR中不同性状的遗传力富集。LD分数回归的代码来自https://github.com/bulik/ldsc（2018年3月2日）。LD分数回归是在从https://data.broadinstitute.org/alkesgroup/ldscore/weights_hm3_no_hla.tgz下载的HAPMAP366 SNP上进行的。然后，将SNP列表进一步过滤到预处理的基线模型中使用的SNP（https://data.broadinstitute.org/alkesgroup/ldscore/1000G_PHASE3_baseLINELD_V1.V1.1.1_LDSCORES.TGZ）。使用1000个基因组Project67（https://data.broadinstitute.org/alkesgroup/ldscore/1000g_phase3_phase3_plinkfiles.tgz）中的欧洲人群的数据计算LD评分。https://data.broadinstitute.org/alkesgroup/ldscore/1000g_phase3_frq.tgz。从https://data.broadinstitute.org/alkesgroup/sumstats_formatted/下载了27个特征的摘要统计信息。“ pass_years_of_education1.sumstats”被忽略了，因为可以使用多年教育的最新研究的摘要统计数据。

　　为了获得CG-DMR的人类直系同源区域，我们使用升力将小鼠CG-DMR（MM10）映射到HG19，要求CG-DMR中至少50％的碱基可以分配给HG19（使用选项-minmatch = 0.5）。在1,880,810个小鼠DMR区（55％）中，总共1,034,801个可以与人类基因组保持一致。

　　然后，对于每个组织样品，我们用1000个基因组SNP中的FEDMR的人类直系同源区域重叠，并使用1000个基因组数据计算了LD评分。但是，仅报告了预告片基线模型中SNP的LD评分，并用于以后进行分析。LD分数是使用选项“ –LD-Wind-CM 1”计算得出的。

　　最后，我们对每个组织样本的每个特征进行了LD评分回归，并使用选项为“ - 跨lap-annot”。测试中使用的回归模型包括FEDMRS和验证基线模型中的注释，如47。后者用于控制通用调节元件中的非组织特异性富集，例如所有启动子47。我们总共进行了1,953次测试（27个特征×59个组织样本）。使用R函数PNOR和参数为'lower.tail = f'的R函数PNorm使用报告的系数z得分（系数_z得分）计算P值。系数_z得分是基于200块折刀重新采样的200个重复序列，因此该统计测试的样本量为200个。要纠正由多个比较导致的P值通胀，我们分别应用了本杰米尼 - 霍赫伯格（Benjamini-Hochberg）的方法，在每个组织样品的FEDMR上的测试中，对P值进行了Benjamini-Hochberg方法。给定5％FDR的P值截止被用来调用大量富集。

　　为了更好地了解CG-DMR的潜在功能，我们根据基因组位置和染色质特征将它们分为各个类别。首先，我们将CG-DMR与启动子，CGI和CGI岸重叠，并定义了与这些位置重叠为近端CG-DMR的CG-DMR。在153,019个近端CG-DMR中，分别与CGI启动子，非CGI启动子，CGIS和CGI海岸重叠。我们避免将近端CG-DMR分配为多种类别，通过将四个基因组特征作为CGI启动子，非CGI启动子，CGI和CGI海岸（以降低优先级订购）。每个CG-DMR被分配给具有最高优先级的类别。

　　我们进一步对剩余的1,655,791个远端CG-DMR分类如下：（1）其中397,320个被预测为远端FEDMR（CG-DMRS（CG-DMRS，它们显示增强子样染色质标志）44,68）如上所述。（2）接下来，我们将侧面FEDMRS定义为远端FEDMR 1 kb的CG，但并未被预测为增强子（FEDMR）。我们总共发现了212,620此类CG-DMR。（3）然后，在其余未分类的CG-DMR中，在至少一个组织中，将159,347个CG-DMR鉴定为组织特异性CG-DMR，因为它们显示出强组织特异性的低甲基化模式（MCG差异≥0.3）。通过检查其低甲基化组织中组蛋白标记的富集，我们发现它们富含H3K4ME1，但没有其他组蛋白标记，并且这些染色质特征类似于Primended Enhancers69。因此，我们将这些CG-DMR定义为引发远端FEDMR。（4）最后，我们将其余的CG-DMR定义为无法解释的CG-DMR（UNXDMR），因为它们的功能角色尚未分配。我们发现UNXDMR与转座有很强的重叠，我们将它们分为两个类别：TE-UNXDMR（n = 449,623）和NTE-UNXDMRS（n = 436,881）。TE-UNXDMR是与转座重叠的UNXDMR，其余为NTE-UNXDMR。

　　使用UCSC基因组浏览器51（http://hgdownload.cse.ucsc.edu/goldenpath/mmm10/phylop60way/mm10.60.60.60way.phylop60way.bw）中测量CG-DMRS的进化保护。接下来，使用DeepTools255来生成CG-DMR中心进化保护的曲线，使用选项“参考点– Referencepoint =中心-a 5000 -b 5000”。

　　为了找到在进化上保守的CG-DMR的比例，我们与不同类别的CG-DMR重叠，其中小鼠基因组中具有保守的DNA元素。从UCSC基因组浏览器51（MM10鼠标参考中的PhastConselements60way）下载了保守元素的列表。

　　我们将CG-DMR的效果大小定义为最甲基化组织样本和大容量样品平均水平之间MCG水平的绝对差异。大量样品中某些CG-DMR的平均MCG水平估计该基因组区域的基线MCG水平。将大量样品选择为所有样品中的50％，以使其MCG水平的范围最窄（有关详细信息，请参见“组织特异性CG-DMR”的识别）。在此定义中，效果大小表示CG-DMR的甲基化程度。DMS的效果大小以相同的方式定义。

　　为了鉴定相对于FEDMRS富含侧面FEDMR的TF结合基序，我们使用以前的前景和后者作为背景进行了基序分析。具体而言，对于每个组织，将组织特异性的FedMR用作背景，而侧翼远端FEDMR则在这些组织特异性FEDMR的1 kb之内用作前景。为了避免尺寸分布差异引起的潜在偏差，前景和背景区域均从两侧（5'和3'）延伸，使得两者的平均大小为400 bp。接下来，进行了超几何测试，以找到在前景中显着富集的TF结合基序。该测试与FEDMR中TF结合基序的鉴定相同。

　　我们还进行了基序分析，以识别富含底漆远端FEDMR的TF结合基序。该程序类似于FEDMR的基元富集分析。对于每个组织，在该组织中脱甲基化的底漆远端FEDMR被视为前景，而其余的底漆远端FEDMR则被视为背景。然后，进行超几何测试以识别明显的基序富集。

　　接下来，对于每种组织类型，我们比较了富含引发远端FEDMR和组织特异性FEDMR的TF结合基序。超几何测试用于测试重叠的重要性 - 如果这两个列表是基于TF表达水平大于10 tpm的TF结合基序的随机抽样（无替换），则获得观察到的重叠的机会。

　　为了估计UNXDMR和可转座元素（TES）之间重叠的重要性，我们使用默认设置的BedTools61中使用“洗牌”工具将UNXDMR的位置改组，并重新计算了重叠。重复此步骤1000次后，如果UNXDMR随机分布在基因组中，我们获得了重叠的经验估计。让观察到的TE叠加的UNXDMR的数量为XOB，并且在置换中，我是XOB的te叠加unxdmr的数量。然后，我们将p值计算为

　　其中。

　　使用与先前所述的33相同的程序来称呼大型CG-DMR。对于每种组织类型，如果组织特异性的CG-DMR相互合并，则它们彼此1 kb之内。然后，我们滤除了合并的CG-DMRS的长度小于2 kb。

　　我们具有大型cg-dMR的基因重叠，然后滤除了以“ rik”或“ gm [0-9]开头的名称开头的任何基因，其中[0-9]代表一个数字，因为这些基因的本体论是不确定的。

　　使用Rose 36,71管道鉴定出超级增强剂。首先，使用MACS256 CallPeak模块调用H3K27AC峰，并带有选项“ - extsize 300 -Q 0.05 – NOMODEL -G MM”。控制数据用于峰值呼叫步骤。接下来，Rose的运行方式是“ -s 12500 -T 2500”，H3K27AC峰，映射的H3K27AC芯片奇数读取和映射的控件读取为输入。为每个组织样品生成了超增强剂调用。然后，我们通过在胎儿发育的每个阶段融合了一种超级增强剂（E10.5至P0），获得了一种组织类型的超级增强剂。最后，我们通过合并超级增强剂的呼吁来生成一份合并的超级增强剂列表，除了肝脏以外的所有组织类型。

　　为了量化MCG动力学，我们定义并计算了MCG的损失和MCG事件。在一个阶段间隔中，MCG水平的MCG或MCG事件损失分别是MCG水平的减少或增加至少0.1。例如，如果在E11.5和E12.5处的一个CG-DMR的MCG水平分别为0.8和0.7，则在阶段E11.5 – E12.5时被认为是MCG事件的损失。阶段间隔定义为两个采样相邻阶段之间的过渡（例如，E15.5和E16.5）。

　　我们使用K-均值聚类来根据MCG和H3K27AC动力学识别前脑特异性CG-DMR的亚组。首先，对于每个前脑特异性CG-DMR，我们计算了MCG水平和H3K27AC在E10.5到成人阶段的H3K27AC富集。在这里，我们使用了已发表的甲基化数据来用于产后1-，2周和6周的额叶皮质9来近似成人前脑的DNA甲基化景观。我们还将H3K27AC数据纳入了产后1-，3和7周的前脑样品。接下来，为了使H3K27AC富集值的范围与MCH水平相当，对于每个前脑特异性的CG-DMR，负H3K27AC富集值被阈值阈值为零，然后每个值除以最大值。如果某些前脑特异性CG-DMR的最大值为零，则将所有值设置为零。进行了亚组的K-均值聚类，但未观察到新的模式。最后，我们使用了使用“最近的基因”关联策略的Great65来识别每个亚组的CG-DMR附近基因的丰富基因本体学术语。

　　To investigate the association between mCG and H3K27ac, for each tissue and each developmental stage, we first divided the tissue-specific CG-DMRs into three categories on the basis of mCG methylation levels: H (high CG methylation; mCG level > 0.6), M (moderate CG methylation; 0.2 < mCG level ≤ 0.6) and L (low CG methylation; mCG level ≤ 0.2).然后，我们通过计算四个H3K27AC中的每个级别的CG-DMR数量：[0,2]，（2，4]，（4，6]和（6，（6，∞）），我们检查了不同组CG-DMR中H3K27AC富集的分布。

　　我们确定了DMV，如前所述37。首先，将基因组分为1 kb的非重叠垃圾箱。然后，对于每个组织样品（重复），将MCG水平小于0.15的连续垃圾箱合并为块。没有数据（没有CG站点或没有读取）的垃圾箱。接下来，将至少五个与数据箱合并的任何块称为DMV。对于每个组织样品，我们首先选择一个重复中鉴定的DMV在两个重复中重复复制的DMV过滤，以将其他重复中调用的所有DMV重叠，然后合并重叠的DMV。使用此策略，我们从每个发育阶段获得了每个组织的DMV呼吁。最后，我们通过合并任何从任何发育阶段中的组织中鉴定出的所有DMV，生成了所有组织样品合并的DMV列表。

　　我们与DMV重叠，然后以“ Rik”或“ GM [0-9]开头的名称滤除了所有基因，其中[0-9]代表一个数字，因为这些基因的本体论不明显。

　　使用随机森林分类器鉴定为PMD，如前所述8。为了训练分类器，我们首先在视觉上选择了19号染色体的区域，作为E14.5肝样品中PMD或非PMD的强候选者。具体而言，我们手动注释了五个PMD，与邻近的基因组区域相比显示出明显较低的MCG水平（CHR19：46110000-46240000，CHR19：45820000-45960000，CHR19：47140000-477340000，以及CHR19：48060000-529110000-PM4713800–4928700, chr19: 7420700–7541100, chr19: 8738100–8967000, chr19: 18633300–18713800, chr19: 53315500–53390000, chr19: 55256600–55633900 and chr19:59281600–59329200）。

　　Next, these regions were divided into 10-kb non-overlapping bins and we calculated the percentiles of the methylation levels at the CG sites within each bin. CG sites that were within CGIs, DMVs37 or any of four Hox loci (see below) were excluded as these regions are typically hypomethylated which may result in incorrect PMD calling. Additionally, sites with fewer than five reads covered were also excluded. We trained the random forest classifier using data from E14.5 liver (combining the two replicates) and we then predicted whether a 10-kb bin was a PMD or non-PMD in all liver samples (considering replicates separately). We chose a large bin size (10 kb) to reduce the effect of smaller-scale variations in methylation (such as DMRs) as PMDs were first discovered as large (mean length 153 kb) regions with intermediate methylation level (<70%)7. Furthermore, the features (the distribution of methylation level of CG sites, which measured the fraction of CG sites that showed methylation levels at various methylation level ranges) used in the classifier required enough CG sites within each bin to robustly estimate the distribution, which necessitated a relatively large bin. Also, we excluded any 10-kb bins containing fewer than ten CG sites for the same reason. These percentiles were used as features for the random forest. The random forest implement was from scikit-learn (version 0.17.1)72 python module and the following arguments were supplied to the Python function RandomForestClassifier from scikit-learn: n_estimators = 10000, max_features = None, oob_score = True, compute_importances = True.

　　最后，我们将连续的10-KB垃圾箱合并为PMD，并将其滤出小于100 kb的块。我们进一步排除了在MM10基因组中与空白重叠的块（从UCSC基因组浏览器下载，2013年9月21日）。为了获得两种重复中可再现的一组PMD，我们仅考虑了大于100 kb的基因组区域，并且在两个重复中都被PMD调用覆盖。这些区域是用于以后分析的最终集合。由于成人肝脏只有一个复制品，因此我们在此阶段使用单个复制来称呼PMD。

　　最初使用上述过程称为PMD，而无需排除HOX基因簇中的CG位点。但是，由于这些HOX基因座更可能被认为是大型DMVS37，因此我们删除了与四个HOX簇重叠的所有PMD（CHR11：9625739-96358516，CHR15：102896908-10308-10303806452146273–52277140）。

　　为了检查正常小鼠肝细胞中PMD和LAD之间的关系（AML12肝细胞），我们使用了参考文献2的补充表2中的LAD数据。73。使用带有默认设置的升降机将LAD的MM9坐标转换为MM10。然后，我们使用蒙特卡洛测试来检查PMD和LAD之间的重叠的重要性。类似于检查TES和UNXDMR之间重叠的过程，我们（1,000次）PMD的基因组位置（1000次），并记录了混乱的PMD和LAD之间的重叠碱基数（用于置换i）。然后，我们将PMD和LAD之间观察到的重叠碱基（XOB）的数量进行比较，并计算出P值为：

　　其中。

　　从replicationDomain74使用了三种鼠标单元类型的复制正时数据（构建MM10）。这些分析使用的细胞类型是MES细胞（ID：1967902＆4177902_TT2ESMOCKCGHRT），神经祖细胞（ID：4180202＆4181802和4181802_TT2NSMOCKCGHRT）和小鼠胚胎成纤维纤维细胞（ID：3040067-1A）。

　　我们使用“相交”获得了有关PMD重叠蛋白编码基因的信息。使用类似的方法来鉴定与PMD侧翼区域（PMD上游和下游100 kb）重叠的蛋白质编码基因；与PMD重叠的基因已从此列表中删除。最后，我们比较了PMD重叠基因（n = 5,748）和基因（n = 2,555）的表达，这些基因与侧翼区域重叠。

　　如先前所述8，我们首先确定了CH位点的审问，该CH位点明显高于低水平噪声（在甲基化水平上约为0.005），这是由Bisuim bisuim bisulfite非转化引起的。对于每个CH站点，我们计算了支持甲基化的读数和不支持甲基化的读数。接下来，我们进行了一项二项式测试，其成功概率等于亚硫酸钠非转换率。FDR（1％）使用Benjamini – Hochberg方法控制75。该分析是针对每个三核苷酸上下文进行独立执行的（例如，计算了CAG胞质的P值截止）。最后，我们计算了甲基化MCH位点的三核苷酸上下文周围±5bp的序列基序，并使用Seqlogo76可视化序列偏好。

　　我们使用了一个迭代过程来调用MCH结构域，MCH域是与侧翼区域相比富含MCH的基因组区域。首先，我们选择了一组未显示MCH证据的样品。在以下步骤中使用了来自这些样品的数据，以滤除容易发生对齐的基因组区域，并显示出可疑的MCH丰度。对全球MCH水平和MCH基序的分析表明，E10.5和E11.5组织（不包括心脏样本）的MCH极低，并且显着甲基化的非CG位点几乎没有CA的偏好。因此，我们假设这些位点不包含MCH结构域，并且算法在对照样品中调用的任何MCH结构域都可能是人工制品。通过滤除对照样品中调用的域，我们能够排除容易绘制误差并避免处理管道中其他潜在缺陷的基因组区域。

　　为了确定在MCH水平发生急剧变化的基因组区域，我们应用了一个变化点检测算法，其中所有5-KB非重叠箱的MCH水平跨基因组作为输入。我们仅包括至少包含500个CH位点的垃圾箱，其中至少50％的CH站点被10个或更多读覆盖。确定的区域定义了将MCH结构域与显示背景MCH水平的基因组区域分开的边界。我们使用r软件包“更改点”中的函数cpt.mean实现了此步骤，并带有选项“方法=” pelt'，pen.value = 0.05，nunty =“渐近paticotic”和minseglen = 2'。为了符合所选罚款的范围，我们将MCH水平扩大了1,000倍。

　　迭代过程如下：1）创建了排除区域的空列表。2）对于每个对照样品，将变化点检测算法应用于5-kb非重叠箱的缩放MCH水平。忽略了重叠的排除区域的垃圾箱。3）基于确定的变更点将基因组分为块。4）每个块的MCH水平计算为没有重叠排除区域的重叠的5-KB箱的平均MCH水平。5）MCH域被确定为MCH水平至少比上游和下游块的MCH水平高50％的块。使用0.001的伪MCH水平避免除以零。6）将MCH域添加到排除区域的列表中。7）重复步骤2至6，直到排除区域的列表停止扩展为止。8）然后将步骤2到5应用于所有样品。9）对于每个组织或器官，仅保留了两个重复中的MCH结构域（一部分）的区域，并且将长度小于15 kb的区域滤除；MCH域必须至少跨三个垃圾箱。上述标准用于定义每个组织或器官的MCH结构域。10）将每个组织和器官的单个MCH结构域合并，以获得384 MCH结构域的单个组合列表。

　　我们根据每个MCH域的归一化MCH积累曲线和相应的侧面区域（下游100 kb和100 kb下游），将k均值聚类用于将384个确定的MCH域分组为5个簇。具体而言，1）在每个组织样品中，一个MCH结构域的MCH积累曲线表示为长度为50的向量：MCH域上游的20个5-KB垃圾箱的MCH水平，10个箱子，将MCH域的10个垃圾箱同样划分为MCH结构域和20 5-KB bin。2）然后，我们通过侧翼区域的箱的平均MCH水平（上游的20个5-kb箱和20个5-kb bin的MCH垃圾箱）对所有值进行了归一化。3）我们接下来计算了六种组织类型的样品（中脑，后脑，心脏，肠，胃和肾脏）的样品，这些组织显示出胎儿发育中最明显的MCH积累。4）使用这些组织样品的轮廓，使用K-均值（R V3.3.1）来群集具有K = 5的MCH域。我们还尝试了更高的群集数（例如，6），但没有确定任何新模式。即使使用当前的K设置（K = 5），簇1（C1）和3（C3）中的MCH结构域也具有相似的MCH累积模式。

　　我们通过使用BedTools61中的“相交”与MCH结构域重叠的基因体与MCH结构域重叠，从而获得了每个MCH结构域的重叠基因信息。仅考虑蛋白质编码基因。我们进一步滤除了以“ rik”或“ gm [0-9]开头的名称的基因，其中[0-9]代表一个数字，因为这些基因的本体论已不确定。对于每个MCH域群集的重叠基因，我们使用富集49,77来找到富集的基因本体论术语（“ go_biological_process_2015”）。

　　接下来，我们询问确定的重叠基因是否富含TF编码基因。为此，使用了来自AnimalTFDB78（2017年2月27日）的小鼠TF列表。然后，我们进行了蒙特卡洛测试，以估计发现的重要性。具体而言，XOB是所有重叠基因中TF编码基因的数量。我们随机选择（1,000次）与编码TFS的随机选择基因相同的基因数量。最后，p值计算为

　　其中。

　　为了评估MCH丰度和基因表达之间的关联，我们追踪了MCH结构域内基因的表达动力学。对于每个群集中的MCH域，我们首先计算每个重叠基因的TPM z得分。具体而言，对于每种组织类型和每种重叠的基因，我们将该组织类型样品中的TPM值归一化为Z分数。Z分数显示了动态表达的轨迹，其中消除了表达的能力信息。如果未表达基因，我们没有执行归一化。接下来，我们计算了与任何MCH域没有重叠的所有基因的Z分数。最后，我们从MCH结构域外的所有基因的z得分中减去了重叠基因的Z分数。所得值表明相对于MCH结构域而非基因，基因在MCH结构域中的表达水平。

　　我们使用WGCNA79（一种无监督的方法）来检测样品跨样品中具有相似表达曲线的基因集（R包，“ WGCNA”版本1.51）。简而言之，首先将tpm值转换为log2（伪计数1×10-5）。然后，将所有样品中每个基因的TPM值与所有其他基因的表达谱进行了比较，并获得了相关矩阵。为了获得任何两个基因之间的连接强度，我们使用功率邻接函数将此矩阵转换为邻接矩阵。为了选择功率相邻函数的参数（软阈值），我们使用了无尺度拓扑（SFT）标准，其中需要构造的网络至少近似于无标度的拓扑。SFT标准建议使用第一个阈值参数值，只要它高于0.8，就可以使用模型拟合饱和度。在这项研究中，达到了5个阈值。

　　接下来，将邻接矩阵进一步转换为拓扑重叠矩阵（TOM），该基质根据连接强度找到了每个基因迭代的“社区”。根据使用签名的混合网络类型，Biweight Mid相关性和WGCNA中tomsimilityfromexpr模块的签名tomepe参数来计算TOM。使用平均方法使用FlashClust模块进行了TOM的分层聚类。接下来，我们使用混合方法使用碎屑动物模块，深插图= 3，而minclustersize = 30个参数来识别至少30个基因的模块。使用特征代表的给定模块中基因的表达创建了汇总的模块特异性表达曲线。特征力定义为模块中所有基因的Log2转换TPM值的第一个主要成分。换句话说，这是一个虚拟基因，代表给定模块中所有基因的表达曲线。接下来，在距离阈值为0.15的所有模块的特征元素的层次聚类之后，合并了非常相似的模块。最后，合并后的所有模块重新计算了特征根。

　　为了更好地了解每个CEM中基因的生物学过程，我们使用了富集49,77（http://amp.pharm.mssm.edu/enrichr/）来识别go_biological_process_2015类别中的富集基因本学术语。

　　我们研究了CEM中调节元件的基因表达与表观基因组特征之间的关联。首先，对于每个CEM，我们都使用征元表达来总结模块中所有基因的转录模式。然后，我们计算了与CEM中基因相关的所有FEDMR的归一化平均增强子评分和标准化的平均MCG水平。具体而言，为了降低每个组织和每个阶段的潜在批次效应，我们通过所有FEDMR的平均增强子得分将每个FEDMR的增强子得分标准化。除了所有DMR的数据用于计算每个组织和每个阶段的平均MCG水平，以相似的方式对MCG FEDMR进行了标准化。接下来，对于每个CEM，其本元素的TPM，归一化的平均增强子评分和链接的FEDMR的MCG水平在每种组织类型的所有胎儿阶段转换为z评分（用于分析组织特异性表达）或跨组织类型的每个组织类型（用于每个组织类型）（用于分析时间表达）。最后，对于每个CEM，我们计算了每个模块的Z分数和归一化增强子得分（或MCG水平）的Z得分之间的Pearson相关系数（R 3.3.1）。计算了两个不同设置的相关系数：1）对于每种组织类型，使用归一化的征元表达值和增强子得分（或MCG水平）跨不同发育阶段计算的相关性；或2）对于每个发育阶段，相关性是在不同组织类型上计算的。前者分析的系数表明时间基因表达与增强子评分或调节元件的MCG水平如何相关，而后者则测量与组织特异性基因表达的关联。

　　然后，我们测试了我们观察到的相关性是否与基于洗牌数据的相关性进行比较，是否显着。在分析给定组织类型中组织特异性表达的过程中，我们将一个CEM的特征性表达映射到与随机选择的CEM中基因相关的FEDMR的增强子评分（或MCG水平）。例如，在洗牌设置中，当给定的组织类型是心脏时，我们计算了CEM14的征元表达与与CEM6基因相关的FEDMR的增强子得分之间的相关性。在对时间表的分析中，如果特定的发育阶段，我们进行了相似的排列。接下来，我们计算了此置换设置的Pearson相关系数。最后，使用两尾Mann-Whitney测试，我们比较了观察到的相关系数的中位数和基于洗牌数据的中位数。

　　有关研究设计的更多信息可在与本文有关的自然研究报告摘要中获得。