Single-cell DNA methylome and 3D genome atlas of human subcutaneous adipose tissue

Abstract

The cell-type-level epigenomic landscape of human subcutaneous adipose tissue (SAT) is not well characterized. Here, we elucidate the epigenomic landscape across SAT cell types using snm3C-seq. We find that SAT CG methylation (mCG) displays pronounced hypermethylation in myeloid cells and hypomethylation in adipocytes and adipose stem and progenitor cells, driving nearly half of the 705,063 differentially methylated regions (DMRs). Moreover, TET1 and DNMT3A are identified as plausible regulators of the cell-type-level mCG profiles. Both global mCG profiles and chromosomal compartmentalization reflect SAT cell-type lineage. Notably, adipocytes display more short-range chromosomal interactions, forming complex local 3D genomic structures that regulate transcriptional functions, including adipogenesis. Furthermore, adipocyte DMRs and A compartments are enriched for abdominal obesity genome-wide association study (GWAS) variants and polygenic risk, while myeloid A compartments are enriched for inflammation. Together, we characterize the SAT single-cell-level epigenomic landscape and link GWAS variants and partitioned polygenic risk of abdominal obesity and inflammation to the SAT epigenome.

Fig. 1. Schematic overview of the study design using snm3C-seq and snRNA-seq to profile cell-type-level DNA methylation, chromatin conformation and gene expression in human SAT and partition the genetic risk of abdominal obesity.

a, Illustration of snm3C-seq and snRNA-seq on nuclei isolated from SAT biopsies from Finnish females. b–g, Comprehensive analyses of DNA methylation, chromatin conformation and gene expression profiles across SAT cell types to identify cell-type-level differences in DNA methylation patterns (b) and chromatin conformation dynamics (c). Subsequently, we used the cell-type-level SAT expression data (d) to determine whether methylation pathway genes contribute to the observed differences in methylation patterns in SAT cell types and longitudinally cluster with adipogenesis pathway genes (e), identify cell-type-level TF binding motifs associated with hypomethylated regions in SAT cell types (f) and test the contribution of variants in cell-type-level DMRs and A and B compartments to the genetic risk of abdominal obesity (g).

Fig. 2. Single-nucleus-level multiomic profiles of SAT by jointly profiling methylation and chromatin conformation with snm3C-seq, followed by an integrative analysis with transcriptomic profiles, generated using SAT snRNA-seq.

a, Dimension reduction of cells using 5 kb bin mCG (top left), 100 kb bin chromatin conformation (top right) and jointly integrating mCG and chromatin conformation (bottom), profiled by snm3C-seq and visualized with uniform manifold approximation and projection (UMAP). Cells are colored by SAT cell type. b, Sankey diagram showcases the high consistency among the SAT cell-type annotations derived from the 5 kb bin mCG (left), 100 kb bin chromatin conformation (right) and joint profiling of mCG and chromatin conformation (middle), with the exception of the transition cell-type cluster that was annotated as perivascular cells by mCG and adipocytes by chromatin conformation. c–f, Integrative analysis with snRNA-seq, evaluating the concordance of cell-type cluster annotations and cell-type marker genes across modalities. c, Comparison of gene body mCG and gene expression profiles of cell-type marker genes across the matching SAT cell types, independently identified within the respective modalities, excluding the expression profiles of the transition cell-type cluster that was not identified in the SAT snRNA-seq data. Dashed lines stratify the marker genes by cell type. Dot colors represent the average gene body mCG ratio normalized per cell (left), and the average log-transformed counts per million normalized gene expression (right). d, Co-embedding of cells profiled by snm3C-seq and snRNA-seq. Cells are colored by annotations as in c (top) and modalities (bottom). e, Concordance matrix comparing the snm3C-seq and snRNA-seq derived annotations, colored by the overlapping scores between the pairs of the SAT cell types evaluated in the co-embedding space. f, UMAP visualization of the gene body mCG ratio (left) and gene expression (right) for one adipocyte marker gene, GPAM, colored per cell similarly as in c. The dashed line highlights the group of cells annotated as adipocytes.

Fig. 3. Functional pathways and gene regulatory potential of cell-type-level gene body mCG markers and DMRs.

a, Dot plots of PPAR signaling pathway genes (ACSL1, ADIPOQ, LPL, PCK1, PLIN1 and PLIN4) that are shared adipocyte marker genes between the gene body mCG and gene expression modalities, showing their gene body mCG (left) and gene expression profiles (right) across SAT cell types. The color of the dot represents the mean percentage of mCG (left) and average expression of genes (right), while the size of the dot represents the percentage of cells in which the gene is expressed (right). b, Horizontal stacked bar plot (left) showing the marginal proportions of assigned methylation states across DMRs for each SAT cell type (n.s., non-significant) and upset plot (right) displaying the top 20 most prevalent methylation state combinations across all detected DMRs, ordered by decreasing frequency, along with their respective percentages. c, Circular plot summarizing the cell-type-level TF binding motifs associated with hypomethylated regions in SAT cell types. Statistical significance was determined using HOMER (binomial test) on each cell type separately. The false-positive rate was calibrated by a stringent cutoff on the unadjusted P values (that is, P < 1 × 10⁻¹²). Track 1 shows the negative log P values (green lollipops) and track 2 shows the enrichment scores (yellow lollipops). d, Violin plots reflecting the empirical null distribution (n = 1,000) of the percentages of DMRs overlapping ChIP-seq peaks for three TFs highlighted in c (see Methods). Diamonds mark the observed overlapping proportion of the hypomethylated regions in myeloid cells in the corresponding ChIP-seq experiment. Asterisks (*) indicate statistical significance, evaluated based on a one-tailed hypergeometric test for overrepresentation of peaks in the corresponding cell-type-level hypomethylated regions (−log₁₀P = 646, 1,860 and 1,618 for IRF4, MEF2B and CEBPG, respectively). GM12878, lymphoblastoid cell line; K562, chronic myelogenous leukemia cell line.

Fig. 4. Analysis of chromatin conformation profiles in SAT reveals cell-type-level diversity in compartments, domains and loops.

a,b, Frequency of contacts per cell against genomic distance (y axis in log scale), grouped by cell types (a) and ordered by the median short-range to long-range interaction ratios (b). Dashed lines in a mark short-range and long-range contacts. The center of the box in b represents the median; the bounds of the box indicate the 25^th and 75^th percentiles and the whiskers show the minimum and maximum values within 1.5 times the interquartile range. Statistical significance was evaluated based on pairwise one-tailed Wilcoxon rank-sum tests for a higher short-range interaction ratio in adipocytes. Asterisks (***) indicate unadjusted −log₁₀P > 50. Specifically, unadjusted −log₁₀P = 197, 116, 221, 258, 128 and 91 for ASPCs, perivascular, endothelial, myeloid, lymphoid and mast cells; n = 63, 1,212, 1,814, 482, 1,206, 1,387, 316 and 172 for transition, adipocytes, ASPCs, perivascular, endothelial, myeloid, lymphoid and mast cells, respectively. c,d, UMAP visualization of low-dimensional embeddings of cells using domains (c) and loops (d) as features, colored by the snm3C-seq annotation; adjusted Rand index (ARI) evaluates the clustering concordance against snm3C-seq annotation. e, Heatmap visualization of the normalized interaction contact map on chromosome 12 and its corresponding compartment scores. f, Upset plot (left) visualizing a subset of the differential 100 kb bins and their corresponding percentages; horizontal stacked bar plot (right) showing the marginal A compartment enrichment of differential 100 kb bins. The vertical dotted line separates cell-type-level A compartment flips, B compartment flips and differential bins detected within homogeneous A or B compartments. g, Dendrogram of the five most abundant SAT cell types constructed with compartment scores on differential 100 kb bins. h, Similar to g, except on all annotated SAT cell types, constructed with mCG fractions across DMRs.

Fig. 5. Analysis of mean gene expression and DMRs across SAT cell types reveals the potential involvement of DNA methylation pathway genes in regulating cell-type-level hypermethylation and hypomethylation in SAT.

a, A schematic representation of basic mechanisms and key players in DNA methylation and demethylation. b, Dot plot of TET1 and DNMT3A showing their expression profiles across the SAT cell types. The size of the dot represents the percentage of cells in which a gene is expressed within a cell type, and the color represents the average expression of each gene across all cells within a cell type (blue indicates higher expression). c, Proportions of assigned hypomethylated (left) and hypermethylated states (right) across DMRs. d, UMAP visualization of the average global mCG ratio in a cell. The dashed line highlights those annotated as myeloid cells. e, Bar plot reflecting the distribution of normalized mCG fractions across genes that co-clustered with TET1 in f for ASPCs and adipocytes. Asterisk (*) indicates statistical significance from a paired one-tailed Wilcoxon rank-sum test, comparing the median cluster expression across n = 5 samples, showing higher expression in ASPCs than in adipocytes (unadjusted P = 0.031). The center of the bar represents the median, and the error bars represent the highest and lowest expression across n = 5 samples. f, Longitudinal expression of TET1 is plotted across the 14-day SAT preadipocyte differentiation. The shaded ribbon behind the trajectory of TET1 reflects the mean and standard deviation of the genes that clustered into similar trajectory patterns as TET1 using DPGP.

Fig. 6. Partitioned PRS and GWAS enrichment results for four key cardiometabolic traits, stratified by cell-type-level DMRs and compartments.

a,b, Lollipop plots depict the incremental variance explained of each cell-type-level PRS for abdominal obesity (using WHRadjBMI as a proxy) from the (a) A and B compartments and (b) hypomethylated and hypermethylated regions. Each lollipop represents a WHRadjBMI PRS, in which dot size corresponds to the incremental variance explained of the PRS. Significant enrichment of incremental variance explained was determined empirically by comparing to n = 10,000 permutated PRSs (see Methods). The gray vertical dashed line indicates the significance cutoff based on unadjusted one-tailed permutation P values (that is, P_perm10,000 < 0.05). Bars and lollipops corresponding to cell types with significantly enriched PRSs are colored by cell type, whereas those without are outlined in gray without a filling. c, Summary table visualizing the overall PRS and GWAS enrichment results across four key cardiometabolic traits (that is, WHRadjBMI, MASLD, BMI and CRP), stratified by cell-type-level hypomethylated regions and A compartments. Orange indicates nominal significance (P < 0.05) on the one-tailed permutation P values for PRSs and the one-tailed hypergeometric P values for GWAS enrichment. The bolded black border highlights the following consistently significant cell-type-cardiometabolic trait pair: adipocyte–WHRadjBMI.

Extended Data Fig. 1. Integrative analysis between subcutaneous adipose tissue (SAT) cells profiled by single-nucleus methyl-3C sequencing (snm3C-seq) and single-nucleus RNA sequencing (snRNA-seq).

a, Dimension reduction of cells (n = 29,423) profiled by snRNA-seq and visualized with uniform manifold approximation and projection (UMAP). b, The total number of cells profiled by snm3C-seq and snRNA-seq stratified by the SAT cell-types. c, Co-embedding of snm3C-seq gene body mCG and snRNA-seq gene expression, visualized with UMAP, highlighting the transition cell-type in red and other SAT cell-types in grey. d, Confusion matrix comparing the concordance between the de novo snm3C-seq annotations (row) and the snRNA-seq-derived annotations (column). The confusion fraction is calculated as the multi-class confusion matrix normalized by the cell counts per row. ASPC, adipose stem and progenitor cell.

Extended Data Fig. 2. Gene body mCG and RNA expression profiles across SAT cell-type marker genes and clustering analysis of the transition cell-type.

a-f, Uniform manifold approximation and projection (UMAP) visualization of the gene body mCG ratio, normalized per cell (left) and log-transformed counts per million normalized gene expression (right) for perivascular marker gene NOTCH3 (a), ASPC marker gene COL5A1 (b), endothelial cell marker gene EGFL7 (c), lymphoid cell marker gene CD2 (d), mast cell marker gene SLC18A2 (e), and myeloid cell marker gene CSF1R (f). The dashed line highlights the group of cells annotated as the corresponding cell-type. g, Gene body hypo-methylation of adipocyte marker genes (left 10 columns) and perivascular cell marker genes (right 10 columns) across adipocytes, perivascular cells, and the transition cell-type. Dot colors represent the average gene body mCG ratio normalized per cell. h, Dimension reduction of cells profiled by snm3C-seq and restricted to adipocytes, perivascular cells, and the transition cell-type, using exclusively the 5-kb bin mCG profiles and visualized with UMAP. ASPC, adipose stem and progenitor cell.

Extended Data Fig. 3. Comparisons of unique cell-type marker genes in SAT cell-types, and biological processes and functional pathways enriched among the adipocyte marker genes between the gene body mCG and gene expression modalities, as well as both known and less-known transcription factors (TFs) in SAT.

a, Venn diagrams showing the number of shared and modality-specific unique SAT cell-type marker genes between the gene body mCG and gene expression modalities. b-c, Dot plots showing significantly enriched biological processes (b) and KEGG functional pathways (c) using unique adipocyte marker genes in gene body mCG and gene expression modalities, evaluated using WebGestalt (see Methods). P values were adjusted for multiple hypothesis testing by the Benjamini-Hochberg procedure. Statistical significance was determined as FDR < 0.05. The size of the dot represents the enrichment ratio for biological processes (b) and KEGG functional pathways (c), while the color of the dot indicates FDR. d, Dot plots of fat cell differentiation biological process genes that are shared adipocyte marker genes between the mCG and gene expression modalities, showing their average gene body mCG (left) and average gene expression profiles (right) across the SAT cell-types. Dot size (right) represents the percentage of cells where the gene is expressed. e-f, Dot plots of (e) adipocyte marker genes encoding both known and less known TFs that are shared between the gene body mCG and gene expression modalities and (f) adipocyte marker genes encoding both known and less known SAT TFs that are unique to gene body mCG, showing their gene body mCG (left) and gene expression profiles (right). The color of the dot represents the mean percentage of mCG (red is high) and the average expression of a gene. ASPC, adipose stem and progenitor cell; and FDR, false discovery rate.

Extended Data Fig. 4. Cell-type level hypo-methylated regions are enriched for specific transcription factor (TF) binding motifs.

We show the top five cell-type level TF binding motifs (sorted by P) that are enriched among the hypo-methylated regions of the SAT cell-types. Statistical significance was determined using HOMER (binomial test) on each cell-type separately (see Methods). False positive rate was calibrated by a stringent cutoff on the unadjusted P values (that is, P < 1×10⁻¹²). ASPC, adipose stem and progenitor cell.

Extended Data Fig. 5. Validation of myeloid-specific transcription factors (TFs), IRF4, MEF2B, and CEBPG, using external ENCODE ChIP-seq data.

a-c, Violin plots reflecting the empirical null distribution (n = 1,000) of the percentages of DMRs overlapping ChIP-seq peaks in each experiment, stratified by the SAT cell-types. Diamonds mark the observed overlapping proportions of the corresponding cell-type level hypo-methylated regions. *Indicates the statistical significance after multiple hypothesis correction (Bonferroni) on the empirical one-tailed P values for the over-representation of TF peaks. IRF4 in GM12878 cell line (empirical unadjusted P=0.001 for myeloid, mast, and lymphoid cells) (a); MEF2B in GM12878 cell line (empirical unadjusted P=0.001 for perivascular, myeloid, mast, and lymphoid cells) (b); and CEBPG in K562 cell line (empirical unadjusted P=0.001 for adipocytes, myeloid, and mast cells) (c). ChIP-seq, Chromatin immunoprecipitation sequencing; ASPC, adipose stem and progenitor cell; GM12878, lymphoblastoid cell line; and K562, chronic myelogenous leukemia cell line.

Extended Data Fig. 6. Cell-type level differences in chromatin conformation of subcutaneous adipose tissue (SAT).

a-b, Uniform manifold approximation and projection (UMAP) visualization of low dimensional embeddings of cells using compartments (a) and insulation scores (d) as features, colored by the snm3C-seq annotation. Adjusted rand index (ARI) evaluates the clustering concordance against the snm3C-seq annotation. c, Heatmap visualization of the normalized interaction contact map on chromosome 6 and its corresponding compartment scores across the SAT cell-types. d, Box plots visualizing the distributions of the number of consecutive 100-kb bins within segments that exhibit consistent A or B compartment annotations, stratified by the SAT cell-types. Lone compartments (that is, singular bins, the annotation of which differs from the ones of both neighbors) are excluded. The center of the box represents the median; the bounds of the box indicate the 25% and 75% percentiles, while the whiskers show the minimum and maximum values within 1.5 times the interquartile range. *Indicates unadjusted –log10P=25, one-tailed Wilcoxon rank-sum test between the two groups marked by the horizontal line: consistent compartment segments in adipocytes, ASPCs, and perivascular cells (n = 2,890) versus those in myeloid and endothelial cells (n = 2,891). e, Horizontal stacked bar plot (left) showing the marginal proportions of differential 100-kb bins stratified by their annotated A and B compartments in the 5 most abundant SAT cell-types and upset plot (right) showing all compartment combinations across differential 100-kb bins in decreasing order with their corresponding percentages. The vertical dotted line separates the following three major categories: ‘Homogeneous’, ‘Cell-type Enriched’, and ‘Heterogeneous’, which correspond to unique A or B compartments in 0, 1, or more than 1 cell-type, respectively. f, Sankey diagram breaking down the numbers of differential 100-kb bins annotated as A (red) and B (blue) compartments belonging to ASPCs (left), adipocytes (middle), and myeloid cells (right). g, Similar to (f), except on perivascular cells (left), adipocytes (middle), and endothelial cells (right). ASPC, adipose stem and progenitor cell.

Extended Data Fig. 7. Cell-type specificity in interaction domains and loops.

a-c, Box plots visualizing the distribution of the number of interaction domains (a), the total number (b), and the average span (c) of interaction domains detected in each cell, stratified by cell-types. The center of the box represents the median; the bounds of the box indicate the 25% and 75% percentiles, while the whiskers show the minimum and maximum values within 1.5 times the interquartile range. Statistical significance was evaluated using pairwise one-tailed Wilcoxon rank-sum tests against adipocytes. *** Indicates Bonferroni multiple hypothesis adjusted P < 0.05 and n.s. denotes non-significant. Specifically, the unadjusted -log₁₀P=133, 45, 92, 132, 70, and 50 for ASPCs, perivascular, endothelial, myeloid, lymphoid, and mast cells, respectively in (a); the unadjusted -log₁₀P=5, 7, and 4 for ASPCs, myeloid, and lymphoid cells, respectively in (b), and the unadjusted -log₁₀P=132, 45, 90, 132, 71, 50 for ASPCs, perivascular, endothelial, myeloid, lymphoid, and mast cells, respectively in (c). Each cell is treated as an independent replicate; thus n = 1,212, 1,814, 63, 1206, 482, 1,387, 172, and 316 for adipocytes, ASPCs, transition, endothelial, perivascular, myeloid, mast, and lymphoid cells, respectively. d, Scatter plot showing the short to long-range interaction ratio per cell against the number of interaction domains detected. Cells are colored by their snm3C-seq annotation. e-f, Scatter plots showing the aggregated cell-type level median number of unique molecular identifiers (UMIs) from snRNA-seq against the median number of interaction domains (e) and the ratio of short to long-range interaction contacts (f) from snm3C-seq, colored similarly as in d. g-h, Bar plots showing the median distance (g) and the total number (h) of loop summits detected across the SAT cell-types (x-axis is ordered by the abundance of each cell-type in snm3C-seq). ASPC, adipose stem and progenitor cell.

Extended Data Fig. 8. Browser-based stacked views of multimodal data for adipocyte marker genes in five major SAT cell-types.

a-b, From top to bottom: chromatin conformation contact heatmap, imputed at a 10-kb resolution and colored at log-normalized scale; boundary probabilities calculated at a 25-kb resolution; and chromosome compartment scores evaluated at a 100-kb resolution, where the A compartments are colored red and the B compartments blue. The shaded trapezoid connects the top panels to a zoomed-in view of the CG methylation patterns surrounding the gene body, extending 5-kb up- and downstream. Vertical dashed lines mark the gene start and end positions. Black dots represent the detected chromatin loop summits. The black dashed triangles in the contact heatmap highlight differential domains that are present only in adipocytes: one encapsulates the ADIPOQ gene (a), and two reside within the TENM3 gene body (b). SAT, subcutaneous adipose tissue; and ASPC, adipose stem and progenitor cell.

Extended Data Fig. 9. Mean gene expression of DNA methylation- and demethylation-related genes across the cell-types in subcutaneous adipose tissue (SAT).

a-b, Dot plot showing expression of (a) DNA methylation genes (DNMT1, DNMT3B, and UHRF1) and (b) DNA demethylation genes (TET2, TET3, and TDG) across SAT cell-types. The size of the dot represents the percentage of cells in which a gene is expressed within a cell-type, while the color represents the average expression of each gene across all cells within a cell-type (blue indicates a higher expression). ASPC, adipose stem and progenitor cell.

Extended Data Fig. 10. Abdominal obesity -associated variants are enriched for the adipocyte A compartment.

a, The clumped and thresholded variants (R² < 0.2, unadjusted P < 0.05) used for constructing the adipocyte compartment PRSs for abdominal obesity (employing WHRadjBMI as a proxy) are plotted by genomic position against the unadjusted -log₁₀P from the UK Biobank WHRadjBMI GWAS (n = 195,863 unrelated Europeans). SNPs landing in the adipocyte A compartments are colored red, while SNPs landing in the adipocyte B compartments are colored blue. The horizontal dashed line indicates genome-wide significance (unadjusted P < 5 × 10⁻⁸). b, Bar plot showing the number of WHRadjBMI-associated variants in (a), stratified by the adipocyte compartment assignment. c, Similar to (b) except showing the number of independent (R² < 0.2) WHRadjBMI GWAS variants, passing genome-wide significance (unadjusted P < 5 × 10⁻⁸). PRS, polygenic risk score; WHRadjBMI, waist-hip-ratio adjusted for body mass index; and SNP, single nucleotide polymorphism.