Tissue slide preparation<\/p>\n
Mouse C57 embryo sagittal frozen sections (MF-104-11-C57 and MF-104-13-C57) were purchased from Zyagen. Freshly collected E11 or E13 mouse embryos were snap frozen in optimal cutting temperature (O.C.T.) compounds and sectioned at 7\u201310\u2009\u03bcm thickness. Tissue sections were collected on poly-l-lysine-coated glass slides (Electron Microscopy Sciences, 63478-AS).<\/p>\n
Juvenile mouse brain tissue (P21) was obtained from the C57BL\/6 mice housed in the University of Pennsylvania Animal Care Facilities under pathogen-free conditions. All procedures used were approved by the Institutional Animal Care and Use Committee.<\/p>\n
Mice were euthanized at P21 using CO2 inhalation, followed by transcranial perfusion with cold Dulbecco\u2019s PBS (DPBS). After isolation, brains were embedded in Tissue-Tek O.C.T. compound and snap frozen on dry ice and a 2-methylbutane bath. Coronal cryosections of 8\u201310\u2009\u03bcm were mounted on the back of Superfrost Plus microscope slides (Fisher Scientific, 12-550-15).<\/p>\n
Preparation of transposome<\/p>\n
Unloaded Tn5 transposome (C01070010) was purchased from Diagenode and the transposome was assembled following the manufacturer\u2019s guidelines. The oligonucleotides used for transposome assembly were: Tn5ME-B, 5\u2032-\/5Phos\/CATCGGCGTACGACTAGATGTGTATAAGAGACAG-3\u2032; Tn5MErev, 5\u2032-\/5Phos\/CTGTCTCTTATACACATCT-3\u2032.<\/p>\n
DNA barcode sequences, DNA oligonucleotides and other key reagents<\/p>\n
DNA oligonucleotides used for PCR and library construction are shown in Supplementary Table 1<\/a>. All DNA barcode sequences are provided in Supplementary Tables 2<\/a> (barcode A) and 3<\/a> (barcode B) and all other chemicals and reagents are listed in Supplementary Table 4<\/a>.<\/p>\n Fabrication of the polydimethylsiloxane microfluidic device<\/p>\n Chrome photomasks were purchased from Front Range Photomasks, with a channel width of either 20 or 50\u2009\u03bcm. The moulds for polydimethylsiloxane (PDMS) microfluidic devices were fabricated using standard photolithography. The manufacturer\u2019s guidelines were followed to spin-coat SU-8-negative photoresist (Microchem, SU-2025 and SU-2010) onto a silicon wafer (WaferPro, C04004). The heights of the features were about 20 and 50\u2009\u03bcm for 20- and 50-\u03bcm-wide devices, respectively. PDMS microfluidic devices were fabricated using the SU-8 moulds. We mixed the curing and base agents in a 1:10 ratio and poured the mixture onto the moulds. After degassing for 30\u2009min, the mixture was cured at 66\u201370\u2009\u00b0C for 2\u201316\u2009h. Solidified PDMS was extracted from the moulds for further use. The detailed protocol for the fabrication and preparation of the PDMS device can be found in our previous research24<\/a>.<\/p>\n Spatial joint profiling of DNA methylation and RNA transcription<\/p>\n Frozen tissue slides were quickly thawed for 1\u2009min in a 37\u2009\u00b0C incubator. The tissue was fixed with 1% formaldehyde in PBS containing 0.05\u2009U\u2009ml\u22121 RNase inhibitor (Enzymatics) for 10\u2009min and quenched with 1.25\u2009M glycine for another 5\u2009min at room temperature. After fixation, tissue was washed twice with 1\u2009ml of DPBS\u2013RNase inhibitor and cleaned with deionized H2O.<\/p>\n The tissue was subsequently permeabilized with 100\u2009\u03bcl of 0.5% Triton X-100 plus 0.05\u2009U\u2009ml\u22121 RNase inhibitor for 30\u2009min at room temperature, then washed twice with 200\u2009\u03bcl DPBS\u2013RNase inhibitor for 5\u2009min each. After permeabilization, the tissue was treated with 100\u2009\u03bcl of 0.1\u2009N HCl for 5\u2009min at room temperature to disrupt histones from the chromatin, then washed twice with 200\u2009\u03bcl of wash buffer (10\u2009mM Tris-HCl pH\u20097.4, 10\u2009mM NaCl, 3\u2009mM MgCl2, 1% BSA and 0.1% Tween\u200920) plus 0.05\u2009U\u2009ml\u22121 RNase inhibitor for 5\u2009min at room temperature. Next, 50\u2009\u03bcl of transposition mixture (5\u2009\u03bcl of assembled transposome, 16.5\u2009\u03bcl of 1\u00d7 DPBS, 25\u2009\u03bcl of 2\u00d7 Tagmentation buffer, 0.5\u2009\u03bcl of 1% digitonin, 0.5\u2009\u03bcl of 10% Tween\u200920, 0.05\u2009U\u2009ml\u22121 RNase inhibitor (Enzymatics) and 1.87\u2009\u03bcl nuclease-free water) was added and incubated at 37\u2009\u00b0C for 60\u2009min. After 60\u2009min incubation, the first round of transposition mixture was removed and a second round of 50\u2009\u03bcl of fresh transposition mixture was added and incubated for another 60\u2009min at 37\u2009\u00b0C. To stop the transposition, 200\u2009\u03bcl of 40\u2009mM EDTA with 0.05\u2009U\u2009ml\u22121 RNase inhibitor was added with incubation for 5\u2009min at room temperature. After that, 200\u2009\u03bcl 1\u00d7 NEB3.1 buffer plus 1% RNase inhibitor was used to wash the tissue for 5\u2009min at room temperature. The tissue was then washed again with 200\u2009\u03bcl of DPBS\u2013RNase inhibitor for 5\u2009min at room temperature before proceeding with the in situ reverse transcription reaction.<\/p>\n In situ reverse transcription<\/p>\n For the in situ reverse transcription, the following mixture was added: 12.5\u2009\u03bcl 5\u00d7 reverse transcription buffer, 4.05\u2009\u03bcl RNase-free water, 0.4\u2009\u03bcl RNase inhibitor (Enzymatics), 1.25\u2009\u03bcl 50% PEG-8000, 3.1\u2009\u03bcl 10\u2009mM dNTPs, 6.2\u2009\u03bcl 200\u2009U\u2009\u03bcl\u22121 Maxima H Minus Reverse Transcriptase, 25\u2009\u03bcl 0.5\u00d7 DPBS\u2013RNase inhibitor and 10\u2009\u03bcl 100\u2009\u03bcM reverse transcription primer (biotinylated-dT oligo). The tissue was incubated for 30\u2009min at room temperature, then at 45\u2009\u00b0C for 90\u2009min in a humidified container. After the reverse transcription reaction, tissue was washed with 1\u00d7 NEB3.1 buffer plus 1% RNase inhibitor for 5\u2009min at room temperature.<\/p>\n In situ barcoding<\/p>\n For in situ ligation with the first barcode (barcode\u2009A), the first PDMS chip was covered at the tissue ROI. For alignment purposes, a 10\u00d7 objective (KEYENCE BZ-X800 fluorescence microscope, BZ-X800 Viewer Software) was used to take the bright-field image. The PDMS device and tissue slide were clamped tightly with a custom acrylic clamp. Barcode\u2009A was first annealed with ligation linker 1 by mixing 10\u2009\u03bcl of 100\u2009\u03bcM ligation linker, 10\u2009\u03bcl of 100\u2009\u03bcM individual barcode\u2009A and 20\u2009\u03bcl of 2\u00d7 annealing buffer (20\u2009mM Tris-HCl pH\u20097.5\u20138.0, 100\u2009mM NaCl2 and 2\u2009mM EDTA). For each channel, 5\u2009\u03bcl of ligation master mixture was prepared with 4\u2009\u03bcl of ligation mixture (27\u2009\u03bcl T4 DNA ligase buffer, 0.9\u2009\u03bcl RNase inhibitor (Enzymatics), 5.4\u2009\u03bcl 5% Triton X-100, 11\u2009\u03bcl T4 DNA ligase and 71.43\u2009\u03bcl RNase-free water) and 1\u2009\u03bcl of each annealed DNA barcode\u2009A (A1\u2013A50, 25\u2009\u03bcM). Vacuum was applied to flow the ligation master mixture into the 50 channels of the device and cover the ROI of the tissue, followed by incubation at 37\u2009\u00b0C for 30\u2009min in a humidified container. Then the PDMS chip and clamp were removed after washing the tissue with 1\u00d7 NEB 3.1 buffer for 5\u2009min. The slide was then washed with deionized water and dried using compressed air.<\/p>\n For in situ ligation with the second barcode (barcode\u2009B), the second PDMS chip was covered at the ROI and a bright-field image was taken with the 10\u00d7 objective. An acrylic clamp was applied to clamp the PDMS and tissue slide together. Annealing of barcode\u2009B (B1\u2013B50, 25\u2009\u03bcM) and preparation of the ligation mixture are the same as barcode\u2009A. The whole device was incubated at 37\u2009\u00b0C for 30\u2009min in a humidified container. The PDMS chip and clamp were then removed, and the slide was washed with deionized water and dried using compressed air. A bright-field image was then taken for further alignment.<\/p>\n Reverse crosslinking<\/p>\n For tissue lysis, the ROI was digested with 100\u2009\u03bcl of the reverse crosslinking mixture (0.4\u2009mg\u2009ml\u22121 proteinase\u2009K, 1\u2009mM EDTA, 50\u2009mM Tris-HCl pH\u20098.0, 200\u2009mM NaCl and 1% SDS) at 58\u201360\u2009\u00b0C for 2\u2009h in a humidified container. The lysate was then collected in a 0.2-ml PCR tube and incubated on a 60\u2009\u00b0C shaker overnight.<\/p>\n gDNA and cDNA separation<\/p>\n For gDNA and cDNA separation, the lysate was purified with Zymo DNA Clean and Concentrator kit and eluted with 100\u2009\u03bcl nuclease-free water. Next, 40\u2009\u03bcl of Dynabeads MyOne Streptavidin C1 beads were used and washed three times with 1\u00d7 B&W buffer containing 0.05% Tween\u200920 (50\u2009\u03bcl 1\u2009M Tris-HCl pH 8.0, 2,000\u2009\u03bcl 5\u2009M NaCl, 10\u2009\u03bcl 0.5\u2009M EDTA, 50\u2009\u03bcl 10% Tween\u200920 and 7,890\u2009\u03bcl nuclease-free water). After washing, beads were resuspended in 100\u2009\u03bcl of 2\u00d7 B&W buffer (50\u2009\u03bcl 1\u2009M Tris-HCl pH\u20098.0, 2,000\u2009\u03bcl 5\u2009M NaCl, 10\u2009\u03bcl 0.5\u2009M EDTA and 2,940\u2009\u03bcl nuclease-free water) containing 2\u2009\u03bcl of SUPERase In RNase inhibitor, then mixed with the gDNA\u2013cDNA lysate and allowed to bind for 1\u2009h with agitation at room temperature. A magnet was then used to separate the beads, which bind to the cDNA that contains dT, from the supernatant that contains the gDNA.<\/p>\n gDNA library generation<\/p>\n Supernatant (200\u2009\u03bcl) was collected from the above separation process for further methylated gDNA detection and library construction. Next, 1\u2009ml of DNA binding buffer was added to the 200\u2009\u03bcl supernatant and purified with the Zymo DNA Clean and Concentrator kit again, then eluted in 84\u2009\u03bcl (3\u2009\u00d7\u200928 \u03bcl) nuclease-free water. The NEBNext enzymatic methyl-seq conversion module (EM-seq) was then used to detect methylated DNA in the sample by converting unmethylated cytosines to uracil; the manufacturer\u2019s guidelines were followed. Then, 28\u2009\u03bcl of DNA sample was aliquoted into each PCR tube, TET2 reaction mixture (10\u2009\u03bcl TET2 reaction buffer containing reconstituted TET2 reaction buffer supplement, 1\u2009\u03bcl oxidation supplement, 1\u2009\u03bcl DTT, 1\u2009\u03bcl oxidation enhancer and 4\u2009\u03bcl TET2) was added to the DNA sample on ice. In brief, 5\u2009\u03bcl of diluted 1:1,300 of 500\u2009mM Fe (II) solution was added to the mixture and incubated for 1\u2009h at 37\u2009\u00b0C in a thermocycler. After the reaction, the sample was transferred to ice and 1\u2009\u03bcl of stop reagent from the kit was added. The sample was then incubated for another 30\u2009min at 37\u2009\u00b0C. TET2 converted DNA was then purified with 90\u2009\u03bcl of solid-phase reversible immobilization (SPRI) beads and eluted with 16\u2009\u03bcl nuclease-free water. The thermocycler was preheated to 85\u2009\u00b0C, 4\u2009\u03bcl formamide was added to the converted DNA and incubated for 10\u2009min at 85\u2009\u00b0C in the preheated thermocycler. After the reaction, the heated sample was immediately placed on ice to maintain the open chromatin structure, then reagents from the kit were added (68\u2009\u03bcl nuclease-free water, 10\u2009\u03bcl APOBEC reaction buffer, 1\u2009\u03bcl BSA and 1\u2009\u03bcl APOBEC) to deaminate unmethylated cytosines to uracil for 3\u2009h at 37\u2009\u00b0C in a thermocycler. Deaminated DNA was then cleaned up using 100\u2009\u03bcl (1:1 ratio) of SPRI beads and eluted in 20\u2009\u03bcl nuclease-free water.<\/p>\n Splint ligation<\/p>\n The gDNA tube was heat-shocked for 3\u2009min at 95\u2009\u00b0C and immediately put on ice for 2\u2009min. Then, 10\u2009\u03bcl of 0.75\u2009\u03bcM pre-annealed Splint Ligate P5 (SLP5) adapter was added. This adapter was diluted from a 12\u2009\u03bcM stock, which contained 6\u2009\u03bcl of 100\u2009\u03bc\u039c SLP5RC oligo, 8.4\u2009\u03bcl of 100\u2009\u03bc\u039c SLS5ME-A-H10 oligo, 5\u2009\u03bcl of 10\u00d7 T4 RNA Ligase Buffer and 30.6\u2009\u03bcl nuclease-free water in a PCR tube that was incubated at 95\u2009\u00b0C for 1\u2009min, then gradually cooled by \u22120.1\u2009\u00b0C\u2009s\u22121 to 10\u2009\u00b0C on a thermocycler. Next, 80\u2009\u03bcl of ligation master mixture was added to the gDNA tube at room temperature. The mixture contained 40\u2009\u03bcl preheated 50% PEG-8000, 12.5\u2009\u03bcl SCR buffer (666\u2009mM Tris-HCl pH\u20098.0 and 132\u2009mM MgCl2 in nuclease-free water), 10\u2009\u03bcl of 100\u2009mM DTT, 10\u2009\u03bcl of 10\u2009mM ATP, 1.25\u2009\u03bcl of 10,000\u2009U\u2009ml\u22121 T4 PNK and 6.25\u2009\u03bcl of 400,000\u2009U\u2009ml\u22121 T4 ligase. The splint ligation mixture was then splinted into five 0.2-ml PCR tubes, 20\u2009\u03bcl per tube. The tubes were shaken at 1,000\u2009rpm for 10\u2009s and spun down, then incubated for 45\u2009min at 37\u2009\u00b0C, followed by 20\u2009min at 65\u2009\u00b0C to inactivate the ligase. For splint ligation indexing PCR, 80\u2009\u03bcl of the PCR reaction mixture was mixed in each splint ligated tube. The mixture contained 20\u2009\u03bcl 5\u00d7 VeraSeq GC Buffer, 4\u2009\u03bcl 10\u2009mM dNTPs, 3\u2009\u03bcl VeraSeq Ultra Enzyme, 5\u2009\u03bcl 20\u00d7 EvaGreen dye, 2\u2009\u03bcl of 10\u2009\u03bcM N501 primer and 2\u2009\u03bcl of 10\u2009\u03bcM N70X-HT primer (Supplementary Table 1<\/a>). The mixture was then aliquoted into a clean PCR tube with 50\u2009\u03bcl volume and run on a thermocycler with the setting below, 98\u2009\u00b0C for 1\u2009min, then cycling at 98\u2009\u00b0C for 10\u2009s, 57\u2009\u00b0C for 20\u2009s and 72\u2009\u00b0C for 30\u2009s, for 13\u201319 cycles, followed by 72\u2009\u00b0C for 10\u2009s. The reaction was removed once the quantitative PCR (qPCR) signal began to plateau. The amplified PCR products were pooled and purified with a 0.8\u00d7 volume ratio of SPRI beads (bead-to-sample ratio)\u00a0and the completed DNA library was eluted in 15\u2009\u03bcl nuclease-free water.<\/p>\n cDNA library generation<\/p>\n The separated beads containing cDNA were used for cDNA library generation. In brief, 400\u2009\u03bcl of 1\u00d7 B&W buffer with 0.05% Tween\u200920 was used to wash the beads twice. Then, the beads were washed once with 400\u2009\u03bcl of 10\u2009mM Tris-HCl pH\u20098.0 containing 0.1% Tween\u200920 for 5\u2009min at room temperature. Streptavidin beads with bound cDNA molecules were placed on a magnetic rack and washed once with 250\u2009\u03bcl nuclease-free water before being resuspended in a template switching oligonucleotide solution (44\u2009\u03bcl 5\u00d7 Maxima reverse transcription buffer, 44\u2009\u03bcl of 20% Ficoll PM-400 solution, 22\u2009\u03bcl of dNTPs, 5.5\u2009\u03bcl of 100\u2009mM template switch oligo, 11\u2009\u03bcl Maxima H Minus reverse transcriptase, 5.5\u2009\u03bcl of RNase inhibitor (Enzymatics) and 88\u2009\u03bcl nuclease-free water). Resuspended beads were then incubated for 30\u2009min with agitation at room temperature and for 90\u2009min at 42\u2009\u00b0C, with gentle agitation. After the reaction, beads were washed with 400\u2009\u03bcl of 10\u2009mM Tris pH\u20098.0 containing 0.1% Tween\u200920 and washed without resuspension in 250\u2009\u03bcl nuclease-free water. Water was removed on the magnetic rack and the beads were resuspended in the PCR solution (100\u2009\u03bcl of 2\u00d7 Kappa Master mix, 8.8\u2009\u03bcl of 10\u2009\u03bcM primers 1 and 2 and 92.4\u2009\u03bcl nuclease-free water). Next, the beads were mixed well and 50 \u03bcl of the PCR mixture was split into four 0.2-ml PCR tubes. The PCR programme was run as follows: 95\u2009\u00b0C for 3\u2009min and cycling at 98\u2009\u00b0C for 20\u2009s, 65\u2009\u00b0C for 45\u2009s and 72\u2009\u00b0C for 3\u2009min, for a total of five cycles, followed by 4\u2009\u00b0C on hold. After five cycles of PCR reaction, four PCR tubes were placed on a magnetic rack and 50\u2009\u03bcl of the clear PCR solution was transferred to four optical-grade qPCR tubes, adding 2.5\u2009\u03bcl of 20\u00d7 Evagreen dye to each tube. The sample was run on a qPCR machine with the following conditions: 95\u2009\u00b0C for 3\u2009min, cycling at 98\u2009\u00b0C for 20\u2009s, 65\u2009\u00b0C for 20\u2009s and 72\u2009\u00b0C for 3\u2009min, for 13\u201317 cycles, followed by 72\u2009\u00b0C for 5\u2009min. The reaction was removed once the qPCR signal began to plateau. The amplified PCR product was purified with a 0.8\u00d7 volume ratio of SPRI beads and eluted in 20\u2009\u03bcl nuclease-free water.<\/p>\n A Nextera XT DNA Library Prep Kit was used for cDNA library preparation. In brief, 2\u2009\u03bcl (2\u2009ng) of purified cDNA (1\u2009ng\u2009\u03bcl\u22121), 10\u2009\u03bcl Tagment DNA buffer, 5\u2009\u03bcl Amplicon Tagment mix and 3\u2009\u03bcl nuclease-free water were mixed and incubated at 55\u2009\u00b0C for 5\u2009min. Then, 5\u2009\u03bcl NT buffer was added to stop the reaction with incubation at room temperature for 5\u2009min. PCR master mix (15\u2009\u03bcl 2\u00d7 N.P.M. Master mix, 1\u2009\u03bcl of 10\u2009\u03bcM P5 primer (N501) and 1\u2009\u03bcl of 10\u2009\u03bcM indexed P7 primer (N70X) and 8\u2009\u03bcl nuclease-free water) was added. The PCR reaction was run with the following programme: 95\u2009\u00b0C for 30\u2009s, cycling at 95\u2009\u00b0C for 10\u2009s, 55\u2009\u00b0C for 30\u2009s, 72\u2009\u00b0C for 30\u2009s and 72\u2009\u00b0C for 5\u2009min, for a total of 12\u2009cycles. The PCR product was then purified with a 0.7\u00d7 ratio of SPRI beads and eluted in 15\u2009\u03bcl nuclease-free water to obtain the cDNA library.<\/p>\n Library quality check and next-generation sequencing<\/p>\n An Agilent Bioanalyzer D5000 ScreenTape was used to determine the size distribution and concentration of the library before sequencing. Next-generation sequencing was conducted on an Illumina NovaSeq 6000 sequencer and NovaSeq X Plus system (150\u2009bp paired-end mode).<\/p>\n Data preprocessing<\/p>\n For RNA-sequencing data, Read\u20092 was processed to extract barcode\u2009A, barcode\u2009B and the UMIs. Using the STARsolo pipeline56<\/a> (v.2.7.10b), these processed data were mapped to the mouse genome reference (mm10). This step generated a gene matrix that captures both gene-expression and spatial-positioning information, encoded through the combination of barcodes\u2009A and B. The gene matrix was then imported into R for downstream spatial transcriptomic analysis using Seurat package (v.5.1.0)57<\/a>.<\/p>\n For DNA-methylation data, adaptor sequences were trimmed before demultiplexing the FASTQ files using the combination of barcodes\u2009A and B. We used the BISulfite-seq CUI Toolkit (BISCUIT) (v.0.3.14)58<\/a> to align the DNA sequences to the mouse reference genome (mm10). Methylation levels at individual CG and CH sites were stored as continuous values between 0 and 1, representing the fraction of methylated reads after quality filtering. These processed CG\u2013CH files were then analysed independently using the MethSCAn pipeline18<\/a> to identify VMRs59<\/a>, defined as fused genome intervals with methylation-level variance in the top 2%. We used default parameter settings when running MethSCAn, and the MethSCAn filter min-sites parameter was determined from the read coverage knee plot (Extended Data Fig. 3a<\/a>). The methylation levels and residuals of VMRs were then imported into R for downstream DNA-methylation analysis.<\/p>\n Clustering and data visualization<\/p>\n We mapped the exact location of pixels on the bright-field tissue image using a custom Python script (https:\/\/github.com\/zhou-lab\/Spatial-DMT-2024\/tree\/main\/Data_preprocess\/Image<\/a>), before removing additional empty barcodes on the basis of read-count thresholds determined by the knee plot (Extended Data Fig. 3a<\/a>). Clustering and data visualization were conducted using R in RStudio.<\/p>\n For RNA data, we used the SCTtransform function in the Seurat package (v.5.1.0), built using a regularized negative binomial model, for normalization and variance stabilization. Dimensionality reduction was performed using RunPCA function with the SCTtransformed assay. We then constructed the nearest-neighbour graph using the first 30 principal components with the FindNeighbors function and identified clusters with the default Leiden method in the FindClusters function. Finally, a UMAP embedding was computed using the same principal components with RunUMAP function.<\/p>\n Owing to the inherent sparsity of DNA-methylation data, it is impractical to analyse methylation status solely at the individual CpG level. Binary information at sparse loci cannot be used directly to construct a feature matrix suitable for downstream analysis. In our study, we adopted the VMR framework, which divides the genome into variable-sized tiles and calculates the average methylation level across CpGs in each tile for each pixel18<\/a>. This approach results in a continuous-valued matrix, in which rows correspond to pixels and columns represent genomic tiles, with values ranging from 0 to 1. VMR methylation levels and residuals were then imputed using the iterative principal component analysis approach as suggested in the MethSCAn instructions. Initially, missing residual values were replaced with zero and missing methylation levels were replaced with the average values for that VMR interval. The principal component analysis approach was iteratively applied until updated values stabilized to a threshold. The imputed residual matrix for VMRs was then imported into the existing Seurat object as another modality. Similar to the RNA-clustering pipeline, dimensionality was reduced using the RunPCA function. The first ten principal components from the residual matrix were used for clustering and UMAP embedding.<\/p>\n To visualize clusters in their spatial locations, the SpatialDimPlot function was used after clustering on the basis of gene expression or VMR residuals. UMAP embedding was visualized with the DimPlot function. The FindMarkers function was applied to select genes and VMRs that were differentially expressed or methylated for each cluster. For spatial mapping of individual VMR methylation levels or gene expression, we applied the smoothScoresNN function from the FigR package60<\/a>. The SpatialFeaturePlot function was then used to visualize VMR methylation levels and gene expression across all pixels. To illustrate the relationships between clustering results from different modalities, we generated the confusion matrix and alluvial diagram using the pheatmap and ggalluvial R package61<\/a>.<\/p>\n Integrative analysis of DNA methylation and RNA data<\/p>\n To integrate spatial DNA methylation and RNA data, WNN analysis in Seurat was applied using the FindMultiModalNeighbors function19<\/a>. On the basis of the WNN graph, clustering, UMAP embedding and spatial mapping of identified clusters were performed for integrated visualization.<\/p>\n For the integration of spatial transcriptomics data of E11 and E13 mouse embryos, the top 3,000 integration features were selected, followed by the use of PrepSCTIntegration and IntegrateData functions to generate an integrated dataset. Similarly, to integrate with public single-cell transcriptomic data25<\/a>,44<\/a>, we first identified anchors using the FindIntegrationAnchors function in Seurat, followed by data integration using the IntegrateData function. To integrate DNA-methylation data, common VMRs between both developmental stages were obtained and the integrated CCA method from the IntegrateLayers function was used to join the methylation data from the two developmental stages. A Wilcoxon signed-rank test was performed to compare the methylation levels and gene-expression differences between the two time points.<\/p>\n TF motif enrichment<\/p>\n To perform TF motif enrichment, we first used the MethSCAn diff function on distinct groups of cells to identify differentially methylated VMRs on the basis of the clustering assignment. The HOMER62<\/a> findMotifsGenome function was then applied to analyse the enrichment of known TF motifs using its default database. We followed the same parameter settings used in MethSCAn, with motif lengths of 5, 6, 7, 8, 9, 10, 11 and 12.<\/p>\n CpGs enrichment analysis<\/p>\n Enrichment analysis of individual CpGs in the differential regions (Fig. 5a,b<\/a>) was performed using knowYourCG (https:\/\/github.com\/zhou-lab\/knowYourCG<\/a>), which provides a comprehensive annotation database for each CpG, including chromatin states, TF binding sites, motif occurrences, PMD annotations and more. To avoid inflated odds ratios for high-coverage data, genomic uniformity was quantified using fold enrichment, defined as the ratio of observed overlaps to expected overlaps. The expected number of overlaps was calculated as: (number of CpGs sequenced \u00d7 number of CpGs in the chromatin state feature)\/total number of CpGs in the genome.<\/p>\n Correlation and GO enrichment analysis<\/p>\n Correlation analysis was performed for different clusters. We first used the findOverlaps function in GenomicRanges package (v.4.4)63<\/a> to map VMRs to overlapped genes. Then, the Pearson correlation test was applied to obtain the correlation between mapped genes and corresponding VMRs. The Benjamini\u2013Hochberg procedure was used to adjust all P\u2009values.<\/p>\n GO enrichment analysis was conducted using the enrichGO function from clusterProfiler package64<\/a> (v.4.2). For GO enrichment in the comparative analysis of E11 and E13 mouse embryos, the FindMarkers function in Seurat package was used to find differential genes and VMRs in the same cluster from integrated data across two developmental stages. Differentially upregulated genes (false discovery rate\u2009\u2264\u20090.05) with demethylated VMRs were used for the GO analysis.<\/p>\n Reporting summary<\/p>\n