Genome-Wide DNA Methylation Profiling of Triple-Negative Breast Cancer Uncovers Epigenetic Biomarkers of Tumor Identity and the Immune Microenvironment

Genome-Wide DNA Methylation Profiling of Triple-Negative Breast Cancer Uncovers Epigenetic Biomarkers of Tumor Identity and the Immune Microenvironment | Research Square window.SnipcartSettings = { analytics: { enabled: false } }; (function() { var accessVector = localStorage.getItem('access_vector') || ''; window.dataLayer = window.dataLayer || []; if (accessVector) { window.dataLayer.push({ user: { profile: { profileInfo: { snid: accessVector } } } }); } })(); (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0],j=d.createElement(s),dl=l!='dataLayer'?'&l='+l:'';j.async=true;j.src='https://www.googletagmanager.com/gtm.js?id='+i+dl;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-K279D39R'); Browse Preprints In Review Journals COVID-19 Preprints AJE Video Bytes Research Tools Research Promotion AJE Professional Editing AJE Rubriq About Preprint Platform In Review Editorial Policies Our Team Advisory Board Help Center Sign In Submit a Preprint Cite Share Download PDF Article Genome-Wide DNA Methylation Profiling of Triple-Negative Breast Cancer Uncovers Epigenetic Biomarkers of Tumor Identity and the Immune Microenvironment Inaki Sasiain, Deborah F Nacer, Mats Jönsson, Johan Vallon-Christersson, and 3 more This is a preprint; it has not been peer reviewed by a journal. https://doi.org/ 10.21203/rs.3.rs-8250160/v1 This work is licensed under a CC BY 4.0 License Status: Under Review Version 1 posted 11 You are reading this latest preprint version Abstract DNA methylation deregulation is an essential feature of tumor biology, influencing cancer formation and pathogenesis. Analyses of DNA methylation in bulk cancer specimens are challenging due to the high data dimensionality, noise, and mixture of different cell types. Here we characterized methylation dynamics in cancer by investigating the variance structure of bulk tumor DNA methylation data through identification of groups of highly correlated CpGs, termed CpG cassettes. Using triple-negative breast cancer (TNBC) as a model system our approach identified co-occurring methylation patterns linked to different intrinsic tumor processes and pathways, gene inactivation, and composition of the tumor immune microenvironment. Our framework also demonstrated the presence of high variance DNA methylation seemingly not linked to tumor biology that could be excluded using chromatin accessibility filtering. Together, this work outlines a comprehensive approach to analyze bulk tumor DNA methylation data combining tumor purity adjustment, functional CpG filtering, stratification by CpG contexts, and identification of highly correlated CpG modules to enhance tumor-intrinsic DNA methylation patterns and our understanding of processes shaping epigenetic tumor evolution. Health sciences/Biomarkers Biological sciences/Cancer Biological sciences/Computational biology and bioinformatics Biological sciences/Genetics Health sciences/Oncology Figures Figure 1 Figure 2 Figure 3 Figure 4 Figure 5 Figure 6 Figure 7 Figure 8 INTRODUCTION Epigenetic mechanisms, such as DNA methylation and histone acetylation, regulate gene expression without modifying the genome. Deregulation of these processes is acknowledged as a major driver of cancer formation and pathogenesis 1 . DNA methylation involves the covalent addition of methyl groups in cytosine residues located at cytosine-guanine dinucleotides (CpG) across the genome and plays a vital role in transcriptional regulation, gene imprinting, X-chromosome inactivation, and suppression of repeat elements 2 . In the context of cancer, DNA methylation has been involved in several relevant molecular mechanisms, like genome repair 3,4 , immune evasion 5,6 , and treatment response 7,8 . A deeper understanding of DNA methylation dynamics in cancer and how it shapes the tumor ecosystem is therefore a highly relevant aspect of cancer research. DNA methylation profiling of bulk tumor tissue is frequently performed by high-density microarrays due to their cost-effectiveness, ease of use, and high CpG coverage 9 . However, these high dimensional data constitute a complex data modality influenced by both intrinsic confounders such as analytical noise and tumor heterogeneity, and extrinsic factors like presence of non-malignant cells with cell-type specific methylation signatures. Together, these factors cause variation in the methylation estimates of single CpGs (referred to as beta values), where tumor heterogeneity and tumor purity have been highlighted as important sources of bias in analyses if not accounted for 10 . Given the complexity of DNA methylation data, different analytic frameworks have been implemented in the context of cancer. A common approach is variance-based filtering, aiming to reduce dimensionality and noise by retaining only the N most variable CpG sites prior to unsupervised or supervised analyses 11-15 . This aims to exclude non-informative CpGs, enriching for CpGs whose variability is more likely to reflect underlying biological processes. Variance filtering, however, also has potential drawbacks as it could keep irrelevant high variance CpGs and is influenced by tumor extrinsic confounders, like the normal cell fraction. To circumvent the latter, tumor purity adjustment methods have been developed that can be applied before variance-based filtering to enrich for tumor-intrinsic signals 16-18 . An additional layer of biological complexity in the analysis of DNA methylation data arises from CpG contexts, i.e., the genomic context in which the CpGs are located. CpGs can be divided into three broad contexts based on their relative location to genes: i) promoter CpGs, affecting the promoter regions of genes, ii) proximal CpGs, located close to the promoter regions but not in the core promoter sequence, and iii) distal CpGs, located at a greater distance from gene promoters. CpG context has functional implications; while promoter and proximal CpGs are often involved in regulation of transcription factor binding to gene promoters, directly affecting gene expression, distal CpGs have a less straightforward link to genes and can be involved in regulating binding of transcription factors to enhancer regions 2,19,20 . In addition, CpGs can reside in regions with different chromatin accessibility, normally referred to as open (accessible) or closed (inaccessible) chromatin as measured by e.g. assays for transposase-accessible chromatin using sequencing (ATAC-seq). Given the complexity of bulk tumor DNA methylation data, this study aimed to deconstruct and characterize its variance structure, considering both the genomic context of CpGs and adjustment of methylation data for tumor purity. Our hypothesis is that this approach can identify groups of CpGs (referred to as CpG cassettes) that share common methylation patterns across samples and are likely involved in related or co-occurring biological processes. As a model system for these analyses, we focused on triple-negative breast cancer (TNBC) based on recent work that demonstrated the importance of both tumor purity correction and CpG context for unsupervised analysis of DNA methylation data 21 . TNBC constitutes ~10-20% of breast cancer cases and is characterized by the negativity of three molecular markers: HER2/ERBB2 and the progesterone (PR) and estrogen receptors (ER) 22 . This subtype is highly heterogeneous with an often aggressive clinical course, shows high genomic instability and has historically lacked available targeted therapies 23 . TNBC has been molecularly stratified by, e.g., mRNA profiling including both the PAM50 subtypes 24 (~70% belong to the Basal PAM50 subtype) or TNBC specific subtypes like the Lehmann et al. subtyping scheme dividing TNBC into Luminal-Androgen Receptor (LAR), Mesenchymal (M), Basal-Like 1 (BL1) and Basal-Like 2 (BL2) tumors 25,26 . In agreement with PAM50, previous studies have demonstrated that the Basal/non-Basal division is highly relevant in TNBC DNA methylation dynamics, with clearly differentiated methylation patterns between tumors, representing the two main attractor states within the TNBC methylation landscape 21 . Besides global DNA methylation patterns, gene specific methylation alterations also play an important role in TNBC. BRCA1 promoter hypermethylation has been reported as the most frequent cause of a genetic homologous recombination deficiency (HRD) phenotype (observed in nearly 60% of cases) alongside mutations in the BRCA1 or BRCA2 genes 27 , and CDH1 hypermethylation has been associated with tumor invasiveness 28 . Finally, immune infiltration introduces an additional dimension of heterogeneity in TNBC, lowering tumor cell purity in bulk specimens, while being a biomarker for better outcome and response to treatment 29 . Immune infiltration in TNBC has so far been associated with features like low clonal heterogeneity, mutation load, copy number alterations 30 and even epigenetic subtypes 21 . Despite intense research, we are still far from a mechanistic understanding of the factors driving heterogeneity in the tumor immune microenvironment (TIME), both across cancers in general and in TNBC specifically, and what role somatic DNA methylation alterations plays in the TIME. With the increasing usage of immunotherapy also for early-stage patients, understanding the TIME and how it is shaped becomes increasingly important. In this setting, TNBC provides an intriguing model system for deconstructing bulk tumor DNA methylation data in order to identify novel epigenetic features associated with both tumor-intrinsic and -extrinsic factors. In this work we demonstrate how statistical deconstruction of bulk tumor DNA methylation data can identify both biologically relevant CpG sets that mirror both gene specific alterations as well as patterns specific to molecular subtypes, as well as variance (noise) without apparent biological relevance for the disease in question that may obscure data analyses. Moreover, by integrating additional -omics layers combined with in situ methods like spatial transcriptomics we identify putative epigenetically driven tumor alterations that correlate with an altered TIME, providing support for a concept of tumor-intrinsic epigenetic immune evasion. MATERIALS AND METHODS Cohorts and Datasets Inclusion and ethics statement. Patients were enrolled in the Sweden Cancerome Analysis Network – Breast (SCAN-B) study (ClinicalTrials.gov ID NCT02306096) 31,32 approved by the Regional Ethical Review Board in Lund, Sweden (registration numbers 2009/658, 2010/383, 2012/58, 2013/459, 2014/521, 2015/277, 2016/541, 2016/742, 2016/944, 2018/267) and the Swedish Ethical Review Authority (registration numbers 2019-01252, 2024-02040-02). All patients provided written informed consent prior to enrolment, including to publish information about sex and age. All analyses were performed in accordance with patient consent and ethical regulations and decisions. SCAN-B cohorts. This work is based on two TNBC cohorts, used for discovery and validation purposes. Both cohorts belong to the Swedish SCAN-B study, and all analyzed tissue specimens were taken prior to any treatment. In Sweden, the definition of TNBC is a tumor with <10% of cells with immunohistochemistry (IHC) staining for ER and PR (thus including tumors with 1–9% stained cells) and an IHC HER2-staining score of <2, or for patients with IHC 2+, a non-amplified in situ hybridization status. All SCAN-B data for ER, PR, and HER2 status were obtained from clinical routine analyses performed by regional pathology departments and provided through the national Swedish breast cancer quality registry. The discovery cohort consists of clinical, pathological and whole genome sequencing (WGS) data from a population representative set of 235 TNBC patients originally reported by Staaf et al. 27 . Molecular data layers included DNA methylation data based on Illumina EPIC bead chips obtained from Aine et al. 21 and RNA-sequencing data summarized as fragments per kilobase million (FPKM) counts obtained from Staaf et al 33 . Different TNBC sample classifications such as PAM50 subtypes, the four refined Lehmann TNBC subtypes (TNBCtype: BL1, BL2, M, and LAR), DNA methylation epitypes (Basal1, Basal2, Basal3, nonBasal1, and nonBasal2), HRD status (HRD or HR proficient: HRP), and rank-based gene expression scores from eight biological metagene genes (including a stroma and immune response metagene 34 ) were obtained from Aine et al. 21 . Additionally, pathology estimates of tumor infiltrating lymphocytes (TILs) and PD-L1 combined positive score (PDL1-CPS) were collected from Aine et al. 35 and Sigurjonsdottir et al. 36 , respectively. CPS scores were categorized into three groups as CPS 0 (CPS score <1), CPS 1 (CPS score 1-9), and CPS 2 (CPS score ≥10). Digital cell counts of relevant immune cell type markers such as CD20, CD3, CD8, FOXP3, CD4, CD68 and PD-L1 were obtained from Roostee et al. 29 calculated with the TMArQ software. A detailed description of the inclusion and exclusion criteria of patients in the discovery cohort is available in the original publication 27 . The validation cohort consists of molecular and clinical data from 136 TNBC samples reported by Aine et al. 21 . Molecular data layers included are DNA methylation data based on Illumina EPIC bead chips obtained from Aine et al. 21 and RNA-sequencing data summarized as FPKM counts obtained from Staaf et al. 33 . Different TNBC sample classifications such as PAM50 subtypes, Lehmann TNBC subtypes and DNA methylation epitypes were obtained from Aine et al. 21 . The composition and characteristics of the discovery and validation cohorts are detailed in Table 1. Table 1 . Composition and characteristics of the discovery and validation cohorts . Age and TILs are reported as median value and interquartile range (IQR) per cohort. ND stands for unavailable data. DISCOVERY VALIDATION N 235 136 Age 62 (IQR: 51-72) 62.5 (IQR: 50-75) TILs 20 (IQR: 10-40) ND Grade 2 26 24 3 206 81 Not available 3 31 Lymph node status N0 152 87 N+ 83 39 Not available 0 10 PAM50 Basal-like 187 95 Normal-like 7 13 HER2-enriched 34 27 Luminal A 5 10 Luminal B 1 0 Unclassified 1 0 Lehmann TNBCtype LAR 40 32 M 54 29 BL1 91 29 BL2 47 20 Unspecified 0 20 Not available 3 6 HRD status HRD 139 ND HRP 96 ND Epitype Basal1 53 23 Basal2 68 37 Basal3 54 29 NonBasal1 29 22 NonBasal2 31 25 Processed and filtered tumor purity-adjusted DNA methylation data 16,17 in the form of beta values (range 0-1, representing hypomethylated to hypermethylated) of both cohorts, comprising in total 741145 CpGs in the discovery cohort and 717458 CpGs in the validation cohort, were obtained from Aine et al. 21 . Data regarding CpG characteristics was compiled from different sources. CpG context, chromosome accessibility status based on ATAC-seq, and overlap with repetitive sequences and transcription factor binding sites were obtained from Aine et al. 21 . Overlap of CpGs with enhancer regions was assessed from enhancers annotated in the GeneHancer genomic regulatory element database (v5.24) 37 . CpG density was determined based on the genomic location of each CpG using windows of 5000 base pairs centered on the CpG and calculated as the number of CpGs in the window divided by total base pairs. TCGA-BRCA external validation dataset. An external validation cohort of 645 breast cancers, representing all molecular and clinical breast cancer subtypes, from The Cancer Genome Atlas (TCGA) with DNA methylation data (Illumina Infinium 450K) and RNA-sequencing data summarized as FPKM was obtained from Aine et al. 21 . Pathologist estimated TILs values for the samples in the cohort were obtained from Cha et al. 38 . Normal breast tissue datasets. Molecular data from normal breast tissue were obtained from three different sources. Illumina Infinium 450K DNA methylation data from 96 normal breast tissue samples were obtained from Hair et al. 39 deposited in Gene Expression Omnibus (GEO) under accession number GSE67919. Illumina Infinium 450K DNA methylation profiles for 10 different flow cytometry sorted blood cell types/fractions were obtained from Reinius et al. 40 deposited in GEO under accession number GSE35069. Single cell RNA-sequencing (scRNA-seq) data from normal breast tissue together with cell type annotations were obtained from Reed et al. 41 deposited in the CellXGene platform. TNBC cell line datasets. Bulk RNA-sequencing data from 34 TNBC cell lines were obtained from Jovanovic et al. 42 summarized as FPKM data deposited under the GEO accession number GSE202770. Eight TNBC cell lines with available RNA-sequencing data, summarized as FPKM values and Illumina EPIC DNA methylation data were obtained from Aine et al. 21 . Additionally, 19 cell lines with available scRNA-seq data were obtained from Jovanovic et al. 42 deposited under the GEO accession number GSE202771 and processed as described by Aine et al 21 . TNBC scRNA-seq dataset. Processed scRNA-seq data together with cell type annotations from eight TNBC tumor samples were collected from Chen et al. 43 deposited under the GEO accession number GSE161529. Spatial transcriptomics TNBC dataset. Single cell spatial transcriptomics data from 65 tumors from the discovery cohort was used to assess OAS2 expression in situ . Spatial transcriptomics analysis was performed using the nanoString/Bruker CosMx instruments on two tissue microarray (TMA) sections, generally containing 2 cores per tumor (114 cores in total), using the Universal Gene Characterisation panel (1000 transcripts). Briefly, TMA sections 4μm thick were cut and mounted on VWR Superfrost Plus slides (company), without coverslips. After sectioning, slides were dried at 37°C overnight at an angle no greater than 45 degrees. Slides were stored at 4°C with a desiccant bag until analysis. Summarized clinicopathological and molecular characteristics for the 65 tumors are available in Supplementary Table 1. CosMx analysis was performed following the manufacturer’s instructions at the Spatial Biology Facility, Faculty of Life Sciences and Medicine, King’s College, London, UK. The image files were processed using nanoString/Bruker’s AtoMx platform and the in-house pipeline and quality assurance set up at the Spatial Biology Facility. Identification and analysis of CpG cassettes CpG cassette identification. CpG cassettes were identified using Weighted Gene-Correlation Network Analysis (WGCNA) through the WGCNA R package (v1.73) 44 . The blockwiseModules() function was applied using bicor as the correlation measure and the following parameters: networkType="unsigned", minModuleSize=3, maxBlockSize=6000, reassignThreshold=0 and mergeCutHeight=0.25. Different soft-thresholding beta values (β), the parameter that determined the stringency of the clustering approach, were used to run the analyses (5, 8, 10, 15, 20 and 25) and evaluated using the pickSoftThreshold() function. The optimal soft-threshold power was set to 10 using tumor purity-adjusted and not context stratified methylation data and was kept across runs to maximize consistency. The remaining parameters were set to the default values. The cassette detection was run on both tumor purity-adjusted and unadjusted CpG beta values, without CpG context stratification, using a variance threshold of 0.1 due to computational limitations. This approach was also applied on tumor purity-adjusted and unadjusted CpG beta values after CpG context stratification, using 0.05 as the variance threshold for promoter and proximal cassettes and 0.1 for distal cassettes. A higher value was used for distal CpGs due to computational limitations given the higher number of such CpGs. Additionally, the same approach and the same parameters were used to identify CpG cassettes from distal CpGs filtered based on chromatin accessibility, defined as CpG overlap with breast cancer ATAC-seq peaks from the study by Corces et al. 45 (see Aine et al. 21 for details). For each cassette, the first principal component (PC1) of the DNA methylation data was computed as a summary measure. The sign of PC1 was adjusted to ensure it aligned with the original data, so that higher values consistently reflected increased methylation beta values, facilitating interpretability. Cell type deconvolution from DNA methylation data. To analyze the impact of tumor purity adjustment on the methylation data and to obtain estimates of the cell composition of each sample the MethylCIBERSORT R package (v0.2.1) 46 was used to generate mixture files from adjusted and unadjusted methylation data using the Prep.CancerType() function and the reference breast cancer signatures (breast_v2). The obtained data were used to run deconvolution through the CIBERSORTx online tool 47 , using 1000 permutations and disabling quantile normalization. CpG cassette stability analysis. The stability of the detected CpG cassettes was analyzed on tumor purity-adjusted and CpG context-stratified data. After filtering CpGs based on variance as detailed previously, 20 groups of 100 samples were randomly selected from the discovery cohort, using the R sample() function with probabilities reflecting the original distribution of the DNA methylation epitypes to avoid class imbalance. Cassette detection was then performed as previously described. Association of CpG cassettes with the Basal/non-Basal division. The predictive capacity of each major CpG cassette in calling the PAM50 Basal/non-Basal division was determined through the Adjusted Rand Index (ARI) using the adjustedRandIndex() function from the mclust (v6.1.1) R package. Specifically, for each cassette the PAM50 Basal/non-Basal tumor subtype classification was compared to the resulting tow sample clusters obtained from hierarchical clustering (using Euclidean distance and Ward.D2 linkage) of CpG methylation data setting the number of clusters to two. Additionally, cassettes with ARI>0.65 were used to stratify samples in Basal/non-Basal groups using hierarchical clustering (Euclidean distance, Ward.D2 linkage) to evaluate their predictive capacity. Gene overrepresentation analyses linked to CpG cassettes. Gene overrepresentation analysis of the CpGs included in the analyzed cassettes was performed using the missMethyl R package (v1.38.0) 48 , which adjusts for the varying numbers of CpGs per gene. The gometh() function was used to perform Gene Ontology enrichment with the following parameters: collection="GO", array.type="EPIC", prior.prob=TRUE, equiv.cpg=TRUE, fract.counts=TRUE, and genomic.features="ALL". The same function and parameters were used to perform enrichment of KEGG (Kyoto Encyclopedia of Genes and Genomes) pathways but setting collection="KEGG". Cancer Hallmark enrichment was determined using the gsameth() function with the same parameters as above but using human hallmarks signatures (v7.1) obtained from the Molecular Signatures Database (MSigDB) 49 as the collection parameter. Analyses of gene expression changes linked to CpG cassettes. The association of CpG cassettes with gene expression differences was assessed by linking CpGs in individual cassettes to genes using provided UCSC annotations. For each gene, tumors were divided into two groups based on k-means clustering of CpGs connected to the gene through the kmeans() R function, using centers=2 and the default parameters, and mRNA expression differences between the two methylation groups were analyzed using Wilcoxon’s test. Assessment of proportion of CpG cassettes in N most variable CpGs. The top 1000, 5000, 10000, 20000 and 30000 most variable CpGs were determined in tumor purity adjusted and unadjusted DNA methylation data, respectively. The proportion of each detected cassette within each CpG context and CpG set was calculated by identifying the CpGs shared between the cassette and the top N groups. Differential methylation analysis. Adjusted DNA methylation data from the discovery cohort was used to identify differentially methylated CpGs per molecular subtype between PAM50 Basal/non-Basal and the Lehmann TNBC subtype division. Differential methylation was tested pairwise between each molecular subtype using Wilcoxon’s test followed by Bonferroni multiple testing correction. Finally the overlap of significant CpGs with identified distal CpG cassettes was analyzed. Identification and analyses of CpG cassettes associated with immune infiltration CpG cassette identification. To identify CpG cassettes associated with TILs we performed an association analysis between cassettes and TIL proportions. First, the association was tested using Kendall’s correlation value between TILs and PC1 of the analyzed cassettes. Additionally, aiming to capture the two potential methylation states of CpGs, we defined two methylation groups per cassette by applying k-means clustering on their PC1, setting the number of clusters to 2. Next, we compared TIL proportions between these groups using Wilcoxon’s test. To prioritize cassettes relevant for TILs we used the Boruta algorithm, using PC1 of each CpG cassette as summary value. Feature selection was performed using the Boruta() function from the Boruta R package (v8.0.0) 50 , setting doTrace=2, maxRuns=500, and remaining parameters to default. scRNA-seq analyses. Analyses of scRNA-seq data from tumor and normal breast samples were performed using the Seurat R package (v5.3.0) 51 together with cell type annotations from the original publications. The DotPlot() function with default parameters was used for gene expression analyses. Spatial transcriptomics analyses. The CosMx spatial transcriptomics data were analyzed using the python Squidpy (v1.6.5) 52 and Scanpy (v1.10.4) 53 libraries. Cells with a proportion of negative probabilities >4%, <30 transcripts detected, and <30 different genes detected were removed. Transcript counts were then normalized using the normalize_total() function and setting target_sum=1e4 followed by log1p transformation through the log1p() function. Tumor cells were identified per TMA slide through Gaussian Mixture Modelling applied on the mean panCK intensity detected per cell using the scikit-learn python library (v1.6.1) 54 . Statistical analyses, data handling and plotting Data handling, statistical analyses and plotting were run in R (v4.4.1). Data handling was performed through base R functions, and the dplyr (v1.1.4), tidyr (v1.3.1) and reshape2 (v1.4.4) packages. Pairwise comparisons were performed using two-sided Wilcoxon’s test. Comparisons between three or more groups were performed using the Kruskal-Wallis test. The Chi-squared statistical test was used for comparisons between categorical variables. The binomial test was used to compare observed and expected events. False Discovery Rate (FDR) or Bonferroni multiple testing correction was implemented when necessary to control for the expected proportions of false positives using the p.adjust() function setting the method to “fdr” or “bonferroni”, respectively. Asterisks indicate p-value significance as follows: * when p≤0.05, ** when p≤0.01, *** when p≤0.001, **** when p≤0.0001, and ns when p>0.05. Correlation values were determined using the Spearman and Kendall methods. Plotting was performed using base R functions, ggplot2 (v3.5.2), ComplexHeatmap (v2.20.2), ggsankey (v0.0.99999), corrplot (v0.95), VennDiagram (v1.6.0) and other R packages’ specific plotting functions. Column clustering of every heatmap was performed using Euclidean distance and the Ward.D2 method while the rows were kept unclustered if not specified otherwise. Data availability statement For the discovery cohort, used WGS data originally reported by Staaf et al. 27 is available from [https://data.mendeley.com/datasets/2mn4ctdpxp/3], while used RNA-sequencing data originally reported by Staaf et al. 33 is available from [https://data.mendeley.com/datasets/yzxtxn4nmd/3]. DNA methylation data used for the discovery cohort was originally reported by Aine et al. 21 and is deposited in the Gene Expression Omnibus under the GSE148748 and GSE148906 accession numbers. DNA methylation data used for the SCAN-B validation cohort was previously reported in Aine et al. 21 and it is deposited in the Gene Expression Omnibus (GEO) under accession number GSE290981. RNA-sequencing data from this cohort is available from [https://data.mendeley.com/datasets/yzxtxn4nmd/3]. TILs and PDL1-CPS estimates for the discovery cohort from Aine et al. 35 and Sigurjonsdottir et al. 36 , respectively, are available in the original publications. TMArQ determined immune cell counts per sample in the discovery cohort by Roostee et al. 29 are available from the original publication. Sample annotations such as PAM50 subtypes, the four refined Lehmann TNBC subtypes, DNA methylation epitypes, HRD status, and rank-based gene expression scores for both the discovery and validation cohorts were obtained from the original publication Aine et al. 21 Raw RNA-sequencing data used in this study and generated by Aine et al. 21 for eight TNBC cell lines are available through the Short Read Archive (SRA) archive at NCBI under BioProject accession PRJNA1189708 and study SRP547133. Corresponding DNA methylation data from Aine et al. 21 for the eight TNBC cell lines is deposited in the GEO database (accession number GSE282347). The previously reported breast cancer TCGA data used in this study is available from the GDC data portal [https://portal.gdc.cancer.gov]. TILs estimates for the TCGA external validation cohort from Cha et al. 38 are available from their original publication. The previously reported normal breast tissue DNA methylation data from Hair et al. 39 used in this study is deposited in the Gene Expression Omnibus under the accession number GSE67919. The previously reported scRNA-seq data of normal breast tissue from Reed et al. 41 is deposited in the CellXGene platform and available from [https://cellxgene.cziscience.com/collections/48259aa8-f168-4bf5-b797-af8e88da6637]. The previously reported RNA-sequencing data for TNBC cell lines from Jovanovic et al. 42 used in this study are deposited in GEO (accession number GSE202770). The previously reported TNBC scRNA-seq data from Chen et al. 43 used in this study are deposited in GEO (accession number GSE161529). The previously reported sorted immune cell DNA methylation data from Reinius et al. used in this study are deposited in GEO (accession number GSE35069). The gene enhancer dataset was retrieved from the GeneHancer genomic regulatory element database (v5.24) 37 and it is available from [https://www.genecards.org/Guide/GeneHancer]. The CosMx spatial transcriptomics data required for the analyses performed in this work is available as supplementary_file_5.csv . Code availability The code used for the analyses implemented in this paper is available in the following GitHub repository: https://github.com/StaafLab/methylation_dynamics_tnbc RESULTS CpG cassettes in tumor purity-adjusted TNBC DNA methylation data To detect CpG cassettes (groups of covarying CpGs) in TNBC we applied Weighted Gene Co-expression Network Analysis (WGCNA) on methylation data from 741145 CpGs in the 235 discovery cohort tumors both before and after adjusting beta values for tumor purity, and without initially stratifying into CpG context. WGCNA clusters CpGs based on correlation values weighted with a soft-thresholding parameter (β) controlling the stringency of the clustering. Higher β values led to more homogeneous cassettes, i.e., smaller and more tightly co-methylated CpG sets that comprised a smaller fraction of all CpGs available (Figure 1A). We selected cassettes defined with β=10 for further analyses to balance cassette homogeneity and the fraction of clustered CpGs, as well as because this value produced a network with approximate scale-free topology (Supp. Figure 1). In a scale-free network, most nodes (CpGs) have few connections, while a few nodes serve as highly connected hubs, a property commonly observed in biological networks and considered desirable for WGCNA. To assess the impact of tumor purity adjustment on the CpG cassette detection we compared WGCNA results for adjusted versus unadjusted methylation data. WGCNA on adjusted CpG data identified more distinct inter-cassette patterns and greater consistency across samples compared to unadjusted data (Figure 1B). This is presumably due to the removal of the effect of the cellular admixture from the tumor sample, as purity adjusted data showed, as expected, higher tumor fractions based on MethylCIBERSORT analysis (Figure 1C). Given our focus on tumor-intrinsic patterns, tumor purity adjusted methylation data were used for all subsequent analyses. We next focused on the impact of CpG context (promoter, proximal, and distal) on cassette detection in the discovery cohort. CpG context is not inherently captured by variation in methylation levels across samples, nor can it be addressed by analyzing only the N most variable CpGs or by correcting for tumor purity (Figure 1D). However, it plays a substantial role in methylation dynamics given its functional role in gene regulation 19,55 . Ignoring CpG context leads to identification of CpG cassettes with substantial differences in CpG context proportions (Figure 1E). For instance, the detected cassette 1 is almost exclusively comprised of distal CpGs, while cassette 2 shows high enrichment of promoter and proximal CpGs, which indicates that some of the detected methylation patterns are CpG context specific. Consistently, performing WGCNA-based cassette detection within each CpG context, respectively, identified cassettes with clearly different methylation patterns (Figure 1F). An important aspect of our CpG cassette detection approach is reproducibility and robustness. To evaluate the stability of the identified CpG context-specific cassettes in a smaller data set, we randomly subsampled the discovery cohort data down to 100 samples (while maintaining the proposed DNA methylation epitype distribution 21 ) and performed a new CpG context cassette detection. This was performed 20 times. Comparison of results to the original CpG context cassettes showed that CpGs from the original distal cassettes remained largely grouped together in subsampling cassettes with similar methylation patterns as the original ones (Figure 1G), and a similar observation was made for proximal and promoter CpGs (Supp. Figure 2). Finally, the DNA methylation patterns of the CpG context-specific cassettes in the discovery cohort were consistent in the SCAN-B validation cohort (Supp. Figure 3). Together, these results support the robustness and reproducibility of the cassette approach. The identified CpG cassettes for all, promoter, proximal and distal CpGs obtained from adjusted beta values using β=10 are available as supplementary_file_1.csv , supplementary_file_2.csv , supplementary_file_3.csv and supplementary_file_4.csv, respectively. CpG cassettes reflect major biological classes in TNBC A distinct feature of the clustered cassette heatmaps in Figure 1 is the aggregation of tumors into mainly a PAM50 Basal and a non-Basal branch (e.g., Figure 1F). To examine this observation in more detail we analyzed the association of detected CpG cassettes in the three CpG contexts with major biological divisions in TNBC, in particular the Basal/non-Basal division proposed to represent the two main global methylation states in TNBC 21 . Among the seven main cassettes identified within each CpG context, several, but not all, captured the two primary methylation states in TNBC. This was quantified using the Adjusted Rand Index (ARI) between sample clusters derived from each cassette and the PAM50 Basal/non-Basal classification (Figure 2A). Specifically, five cassettes (promoter cassette 1, n=3015 CpGs; proximal cassette 1, n=3446 CpGs; proximal cassette 3, n=562 CpGs; distal cassette 2, n=2021 CpGs; distal cassette 3, n=812 CpGs) exhibited ARI values greater than 0.65, indicating strong concordance with the Basal/non-Basal split. Hierarchical clustering of the CpGs included in high-ARI cassettes successfully differentiated Basal from non-Basal samples (Figure 2B). Moreover, CpGs included in the detected cassettes were linked to biologically relevant processes based on pathway enrichment analysis, which was restricted to promoter and proximal cassettes due to available gene mappings. The largest promoter CpG cassette linked to the Basal/non-Basal division, promoter cassette 1, showed enrichment for the Epithelial-Mesenchymal Transition (EMT) hallmark and several KEGG pathways related to cancer processes (Figure 2C). When dividing samples into two groups based on DNA methylation values of CpGs included in this cassette, changes in gene expression correlated with methylation values were observed between the groups, including for genes such as LDHB , MUCL1 , ANP32E and TFAP2B (Figure 2D). Similar findings were observed when analyzing the main proximal cassette linked to the Basal/non-Basal division, proximal cassette 1. This cassette showed enrichment for hallmarks such as EMT, KRAS signaling and Pancreatic beta cells, and for KEGG pathways such as Breast cancer, Basal cell carcinoma, various signaling pathways and glycosphingolipid biosynthesis (Figure 2E). This cassette was also associated with expression changes in genes such as SPINK8 , ANP32E and PAX2 (Figure 2F). As more detailed examples, LDHB and SPINK8 were selected to show the expression changes linked to CpG cassettes (Figure 2G). LDHB shows low mRNA expression in non-Basal TNBC consistent with the hypermethylation of involved promoter CpGs, and the opposite pattern is seen in Basal tumors. In contrast, SPINK8 expression is notably correlated with the DNA methylation status of a proximal CpG, with mRNA expression only in non-Basal tumors. Chromatin accessibility distinguishes biologically relevant variation in distal CpGs As shown in Figure 2A, two major distal cassettes, cassettes 2 and 3, showed a clear correlation with the Basal/non-Basal PAM50 division of the discovery cohort. In contrast, distal cassette 1 (n=6743 CpGs) and other smaller distal and non-distal cassettes showed no association to these PAM50 groups. A notable feature of distal cassette 1 is a high-variance, gradient-like DNA methylation pattern across tumors (Figure 3A). Moreover, distal cassette 1 showed low overlap with enhancer regions (~20%), similar to other cassettes not linked to the Basal/non-Basal division like distal cassettes 4, 5 and 7, and a relatively high fraction of CpGs overlapping with repetitive sequences (~27%) (Figure 3B). The overlap of cassette CpGs with transcription factor binding sites (TFBSs) also showed distinct patterns: cassettes 1, 4, 5 and 7 showed minimal overlap with TFBSs, mirroring the trend seen in enhancer regions (Figure 3C). Moreover, CpGs in distal cassette 1 were often located in generally sparse genomic regions (Figure 3D). Further analyses on potential biological implications of distal cassette 1, with DNA methylation values summarized through the first principal component (PC1), showed no clear difference regarding important tumor characteristics in TNBC such as the PAM50 Basal/non-Basal division, HRD status, TNBC mRNA subtypes, tumor purity, or TIL levels (Figure 3E). Concordantly, when selecting CpGs through supervised methods, such as differential methylation analysis, CpGs included in the noise cassettes were rarely selected (Supp. Figure 4). Taken together, both genomic and biological features suggest that distal cassette 1 represents a potential source of high variance, largely non-informative noise rather than biologically relevant CpGs in the context of TNBC. Importantly, distal cassette 1 could not be removed from the data using variance-based filtering (Figure 3F), even without applying tumor purity correction (Supp. Figure 5). In further support of the non-informativeness of a set of high variance distal CpGs, we found similar results to what we described above also after applying a lower soft-thresholding value that minimizes the fraction of unclustered CpGs (Supp. Figure 6). Again, the largest detected distal cassette showed low overlap with enhancer regions, lower CpG density, reduced overlap with TFBSs and could not be excluded by variance-based filtering on adjusted or unadjusted beta values (Supp. Figure 7). Consequently, if not addressed less informative, high variance cassettes like distal cassette 1 may influence unsupervised DNA methylation analyses. However, applying an alternative filtering strategy based on chromatin accessibility could successfully remove noise-associated cassettes (Figure 3G). This approach was based on overlap of CpGs with ATAC-seq peaks (indicating regions of open, accessible chromatin) from Corces et al. 45 that have a distinctive distribution per CpG context (Table 2). Specifically, the original distal cassette 1, as well as other cassettes with low overlap with TFBSs and enhancer regions, were removed, while approximately 50% of CpGs in cassettes 2 and 3 were retained (Figure 3H). Interestingly, the new cassettes derived from ATAC-seq filtered CpGs had clear counterparts among the original cassettes enriched for TFBSs and enhancer regions (Figure 3I). A similar result was observed when using a lower soft-thresholding value (Supp. Figure 8). To explore whether the identified distal cassette patterns were only a feature of TNBC, representing a limited proportion of all breast cancers, we analyzed the methylation patterns of the seven largest distal cassettes in 645 TCGA breast cancers of all molecular and clinical subtypes. Importantly, a similar apparent noise pattern of distal cassette 1,4, 5 and 7 was observed in all PAM50 subtypes, whereas distal cassette 2, 3 and 6 strongly aligned with a Basal/non-Basal division of the cohort, irrespective of ER, PR, and HER2 status (Supp. Figure 9). Together, this highlights the biological relevance and consistency of the identified patterns also in a broader breast cancer setting. Table 2. Distribution of CpGs from the discovery cohort overlapping with ATAC-seq peaks per CpG context. ATAC + overlap with Corces et al ATAC-seq peaks, and ATAC – have no overlap All Promoter Proximal Distal ATAC - 205799 44722 96285 394339 ATAC + 535346 80567 44531 80701 Specific promoter CpG cassettes are linked to TIME status in TNBC After analyzing the major cassettes reflecting broad biological patterns in TNBC, we focused on smaller cassettes that could be involved in the molecular heterogeneity seen within TNBC. We restricted the analysis to promoter CpG cassettes to facilitate a more direct biological interpretation. As a validation step, we first investigated whether known promoter methylation patterns such as BRCA1 promoter hypermethylation were identified by our cassette approach. Indeed, we identified a cassette (promoter cassette 11, n=19 CpGs) including only CpGs associated with BRCA1 that on a sample level showed excellent agreement with both HRD phenotype and pyrosequencing-based BRCA1 methylation status (representing an orthogonal methylation analysis method), consistent with previous studies 3,27 (Figure 4A). We next focused on epigenetic mechanisms potentially associated with TIME composition, considering the established significance of the TIME for prognosis and treatment prediction in TNBC 56 . A correlation analysis between PC1 of the 1224 promoter cassettes and TIL estimates (acting as a proxy for TIME status) identified a series of cassettes (n=119) that were linked to significant changes in TILs when tumors were split into two groups based on cassette methylation patterns (Figure 4B). To prioritize cassettes for further study, we used the random forest-based Boruta feature selection algorithm identifying 28 promoter cassettes associated with TIL variation (Figure 4C, Supplementary Table 2). Among them, we focused on promoter cassette 10 (Figure 4D), representing both the largest of the 28 cassettes (n=30 CpGs) and the cassette with the highest absolute Kendall correlation to TILs. Eight genes were associated with promoter cassette 10 ( GBP4 , SAMD9L , CARD16 , OAS2 , ZBP1 , APOL1 , LIPA , BATF2 ) and showed enrichment for immune response terms, particularly innate and antiviral processes (Figure 4E). Importantly, a similar methylation pattern for cassette 10 CpGs was also observed in the independent SCAN-B validation cohort (Supp. Figure 10), where tumors with a hypermethylated pattern showed significantly lower rank scores of a gene expression-based immune response metagene 34 (Figure 4F). Finally, a similar link between the methylation pattern for cassette 10 CpGs and TIME status was also in the larger 645-sample general TCGA breast cancer cohort, with hypomethylation being associated with increased immune infiltration providing evidence for generalizability (Supp. Figure 11). Features of a promoter CpG cassette associated with TIME status in TNBC To understand the biological implications of promoter cassette 10 in more depth we focused on its associated genes with more than one mapped CpG to avoid false positives, leaving five genes: GBP4 , OAS2 , CARD16 , ZBP1 and SAMD9L . Firstly, all five genes showed two distinctive methylation states (hypermethylated and hypomethylated) among the analyzed samples linked to statistically significant differences in gene expression (Figure 5A-B). The specific DNA methylation and expression patterns per gene were also confirmed in the SCAN-B validation cohort (Supp. Figure 12) and in the TNBCs from the general TCGA breast cancer validation cohort (Supp. Figure 13A-B). In line with cassette 10’s association with TIL levels, grouping samples by the methylation status of the five genes also showed significant associations with both TIL levels and in situ PDL1-CPS scores, a measure of PD-L1 expression in malignant and non-malignant cells (Figure 5C-D). Specifically, hypermethylated tumors showed a more immune cold TIME exemplified by lower TIL levels and PDL1-CPS scores. A similar trend of association between the methylation status and TILs was obtained in the TNBCs from the general TCGA breast cancer validation cohort (Supp Figure 13C). Additionally, the methylation states of the identified genes were generally linked to higher infiltration of specific immune cell types such as B cells and cytotoxic T cells, measured through both MethylCIBERSORT immune cell fractions and in situ immune cell counts, with higher infiltration in the hypomethylated samples (Supp. Figure 14). Next, we analyzed the co-occurrence of the methylation patterns for the five genes in different TNBC subtypes, finding a high co-occurrence of gene hypermethylation in tumors with a non-Basal PAM50 phenotype and/or Luminal-Androgen Receptor and Mesenchymal Lehmann TNBC subtypes (Figure 5E). In support of this finding, the expression patterns of the five genes were also highly correlated across tumors in the discovery cohort (Supp. Figure 15). Similar methylation co-occurrence patterns were also observed in the SCAN-B validation cohort (Supp. Figure 16). Finally, in the discovery cohort a trend of higher cumulative hypermethylation status of the five genes with lower immune response metagene rank scores and TILs was observed (Supp. Figure 17A). In contrast, an opposite positive trend was observed for rank scores of a gene expression stroma metagene 34 and fibroblast fraction estimates by MethylCIBERSORT (Supp. Figure 17B), together highlighting the potential association of the methylation status with the composition of the tumor microenvironment. Characterization of the promoter CpG cassette associated with TIME status across breast cancer subtypes Aiming to characterize the role of the identified CpG cassette and genes in a more general breast cancer population we analyzed their DNA methylation patterns, gene expression and association with immune infiltration in the 645-sample TCGA breast cancer validation cohort. First, we assessed the differences in immune infiltration per PAM50 breast cancer subtype, observing the highest infiltration in the Basal intrinsic group, followed by HER2-enriched tumors. In contrast, Luminal B, and particularly Luminal A and Normal-like samples showed the lowest TIL levels (Figure 6A). Next, we evaluated whether the detected promoter cassette 10 in TNBC was still present in the complete breast cancer population. Although the structure of the cassette was less binary than in TNBC (Figure 6B), clustering on the CpGs included in the cassette still revealed a strong association between general hypomethylation and increased immune infiltration (Figure 6C). Next, we analyzed GBP4 , OAS2 , ZBP1 and CARD16 associated with promoter cassette 10 in the TCGA cohort. SAMD9L was excluded due to lack of CpG coverage in the available Infinium 450K data. After grouping samples based on the methylation status of the CpGs overlapping with respective gene, similar hyper- and hypomethylation patterns to the ones observed in TNBC were identified, that, consistent with TNBC observations, were associated with gene expression levels across all PAM50 subtypes (Supp. Figure 18). Moreover, a generally higher fraction of hypermethylated samples per gene was observed in the non-basal subtypes (Figure 6D), along with a general co-occurrence of hypermethylation or hypomethylation phenotypes among the genes consistent with TNBC findings, with the co-hypomethylated samples showing the highest TILs (Figure 6E). Finally, we evaluated the association of the DNA methylation status of each gene with immune infiltration within the breast cancer PAM50 subtypes. In Basal, HER2-enriched and Luminal B tumors, hypomethylation of the genes was generally associated with higher TILs levels, whereas no difference was identified in Luminal A and Normal-like tumors, likely due to their generally low immune infiltration (Figure 6F). Methylation and expression patterns of promoter cassette 10-associated genes in normal breast cells To investigate the potential tumor extrinsic expression patterns of the five genes, we analyzed the DNA methylation levels of CpGs in promoter cassette 10 and the expression of the five associated genes in normal breast cells using a combination of reported DNA methylation and scRNA-seq data. We observed different methylation states of cassette 10 CpGs per gene in 96 normal breast tissue specimens (Figure 7A). Two genes showed clear CpG hypomethylation phenotypes in normal breast tissue, OAS2 and GBP4 , although the latter showed less extreme values. CpGs for CARD16 and ZBP1 showed beta values around 0.5, indicating a potential mixture of cells with different methylation states in the analyzed bulk tissue. SAMD9L could not be analyzed as no CpGs close to this gene’s promoter were present in the Infinium 450K methylation data used. The methylation status of the genes of interest was also assessed in specific immune cell populations using data from flow cytometry-sorted blood samples, revealing distinct hypomethylation for all CpGs in all genes ( SAMD9L excluded due to lack of CpG coverage, Supp. Figure 19). Next, we analyzed gene expression of the five genes in normal tissue using scRNA-seq data, observing heterogeneous expression in epithelial, immune and stromal cells (Figure 7B). In epithelial cells, we observed expression of CARD16 and GBP4 in specific subtypes (luminal adaptive secretory precursor cell of mammary gland and mammary gland epithelial cells), while the expression measured for the remaining genes was generally low. In immune cells, the expression was equally heterogeneous: CARD16 was mainly expressed in dendritic cells and mature natural killer (NK) cells, GBP4 mostly in different NK cell populations, SAMD9L in dendritic cells and macrophages, and ZBP1 in different NK cell populations and plasma cells. OAS2 , on the other hand, showed low expression levels in all immune cell types. Finally, in stromal cells, GBP4 and CARD16 were expressed in most cell types except fibroblasts of mammary gland and perivascular cells. OAS2 and SAMD9L showed lower expression, particularly focused on different endothelial cell subtypes, while ZBP1 showed little expression in any stromal cell type. Taken together, these analyses demonstrate a heterogenous expression of the five genes in different normal breast and immune cell types, combined with typically a hypomethylated promoter pattern, especially in immune cell types. Tumor-intrinsic expression of promoter cassette 10-associated genes To investigate the tumor-intrinsic expression of the five genes associated with promoter cassette 10 we first analyzed bulk gene expression data from 34 TNBC cell lines, representing a tumor microenvironment-free context. We identified tumor-intrinsic expression of SAMD9L , CARD16 , OAS2 , and GBP4 in some cell lines, while only very low to no expression of ZBP1 overall (Figure 8A). Expression of the five genes was correlated across the 34 cell lines (Supp. Figure 20), in line with findings in the discovery cohort. Additionally, in eight previously reported TNBC cell lines with matched RNA-sequencing data and Illumina EPIC DNA methylation, we observed concordance between the methylation state of the promoter CpGs and gene expression levels, i.e., higher expression with CpG hypomethylation (Figure 8B). Together, these findings support a hypothesis of tumor-intrinsic methylation patterns and expression of the genes by tumor cells. Furthermore, scRNA-seq from 19 TNBC cell lines confirmed tumor-intrinsic expression of the genes in a subset of cell lines, while at the same time revealing substantial cell heterogeneity in expressing cell lines, i.e., cell lines with gene expression comprised both expressing and non-expressing cancer cells (Supp. Figure 21). To substantiate the cell line observations in actual tumor tissue we evaluated the tumor-intrinsic expression of the five genes in tumor tissue using scRNA-seq data from eight TNBC patients. Despite the small sample set and substantially different expression levels, we observed co-expression of the identified genes in different samples and a trend that tumors with expression of the genes tended to have an activated TIME, determined as higher counts of immune cells in the scRNA-seq data (Figure 8C). To further assess the tumor-intrinsic expression of the five genes in situ we analyzed TMA cores from 65 tumors in the discovery cohort by spatial transcriptomics using Nanostring CosMX. Focusing on OAS2 (due to measured expression levels and CosMX gene panel content), we detected statistically significant higher mean OAS2 expression and higher proportion of tumor cells expressing OAS2 in tumors classified as hypomethylated by the DNA methylation analysis shown in Figure 5A (Figure 8D-E). Importantly, these results connect the in vitro observations with similar findings in primary early-stage treatment naïve tumor tissue, i.e., apparent tumor-intrinsic expression in a subfraction of tumor cells correlating with a tumor-intrinsic hypomethylation status of the promoter region. DISCUSSION In this work, we aimed to characterize methylation dynamics in cancer by investigating the variance structure of bulk tumor DNA methylation data based on the identification of groups of highly correlated CpGs, termed CpG cassettes. In our model system, TNBC, this approach identified CpGs with co-occurring methylation patterns linked to different intrinsic tumor processes and pathways, gene inactivation, and correlation to TIME status, but it also demonstrated the presence of high variance DNA methylation, apparently not biologically relevant for TNBC, that may influence downstream analyses. Moreover, our study further supports that tumor purity adjustment methods can reduce the confounding effects caused by non-malignant cells, thereby increasing the correlation between CpGs with tumor-intrinsic methylation patterns. Together, the current study exemplifies and highlights the utility of network-based analyses of DNA methylation data. TNBC is a clinical subgroup of breast cancer with two previously reported main methylation patterns linked to a Basal/non-Basal division 21 . Our approach successfully identified those patterns regardless of the CpG context (promoter, proximal or distal), demonstrating that even small groups of correlated CpGs (cassettes) can accurately predict this division. Notably, our work also revealed that not all high-variance methylation patterns appear biologically relevant in the context of the disease in question, an important consideration for unsupervised analysis of DNA methylation data. For instance, the largest detected distal CpG cassette (distal cassette 1) did not show any association with tumor features such as tumor subtype, tumor purity, or TIME status, and the included CpGs showed a low overlap with enhancer regions and transcription factor binding sites, in contrast to the largest promoter and proximal cassettes. Together, these observations suggest that this distal cassette may reflect stochastic noise rather than biologically relevant patterns. Similar observations were made for other distal cassettes as well, e.g., distal CpG cassettes 4,5 and 7. It is possible that these biologically non-relevant CpG cassettes may be related to assay probe design or to their location in typically non-functional genomic regions. Still, it should be noted that these CpGs are present in both the EPIC and Infinium 450K platforms despite comprehensive data pre-processing and filtering. Moreover, our analysis of the general TCGA breast cancer cohort also demonstrates that these patterns are not specific to TNBC, but instead apply to all breast cancer subtypes, highlighting the generalizability of our findings. It is likely that similar DNA methylation patterns are also present in data from other malignancies. An important observation regarding these non-relevant cassettes was that they could not be removed using variance-based CpG filtering, a common feature selection step before unsupervised analysis. In contrast, supervised analyses appeared less affected by these specific CpGs. Taken together, these non-biological patterns may represent a significant source of bias in unsupervised analyses if not properly addressed, particularly in cohorts without strong apparent DNA methylation states (like the Basal/non-Basal division in TNBC). Nevertheless, a more functionally oriented strategy for filtering CpGs based on chromatin accessibility was able to remove these high-variance non-informative cassettes, while retaining the patterns detected in the biologically relevant cassettes that were enriched for functionally active methylation regions. Based on the availability of large-scale ATAC-seq tumor data for different malignancies (e.g. from Corces et al. 45 ), our observations support the use of chromatin accessibility-based filtering approaches as an important filtering step when analyzing DNA methylation data, particularly when focusing on distal CpGs. Beyond the Basal/non-Basal DNA methylation patterns in TNBC, illustrated by several CpG cassettes irrespective of CpG context, our approach also detected cassettes that reflected other tumor specific patterns like hypermethylation of the BRCA1 tumor suppressor gene, representing alterations which may at least partially help to explain the molecular heterogeneity within the disease. In a focused exploration of promoter CpG cassettes we also identified tumor methylation patterns correlated with TIME status in TNBC. The limitation of this analysis to promoter cassettes was due to a more straightforward association of a CpG with a gene when located in a promoter region, and it should be acknowledged that careful analysis of both proximal and distal CpG methylation patterns may reveal additional findings. Supporting the latter, one of the previously reported TNBC DNA methylation epitypes based on tumor purity-adjusted and ATAC-seq filtered distal CpGs showed a notably immune warm phenotype (Basal3) 21 . Among the identified promoter CpG cassettes correlating with TIL levels we focused on the largest one, promoter cassette 10, and five of the associated genes involved in innate immune response: GBP4 , OAS2 , CARD16 , ZBP1 and SAMD9L . Promoter hypermethylation in these genes was mainly identified in non-Basal tumors and tumors with a Luminal Androgen Receptor or Mesenchymal-like TNBC mRNA subtype, suggesting an association with tumor-specific phenotypes to be further explored. Both Luminal Androgen Receptor and Mesenchymal-like TNBC subtypes have been described as the more immune-cold subtypes among TNBCs, consistent with our findings 57 . Importantly, further analyses involving cancer cell lines, scRNA-seq, and spatial transcriptomics data of both malignant and non-malignant breast cells supported that the methylation and associated gene expression patterns are likely tumor-intrinsic and not explicitly driven by contaminations from normal or immune cells in the tumor microenvironment. The five genes connected to promoter CpG cassette 10 are all involved in innate immune processes and have, in some cases, been associated with antitumor immune response. GBP4 is a GTPase belonging to the Guanylate Binding Protein family, which are central regulators of the innate immune response 58 and are involved in relevant innate immune processes such as the formation of the inflammasome 59 . Specifically, GBP4 has been widely identified as relevant for tumor biology in different cancer types; it has been recognized as a pan-cancer marker of immune warm tumors 60 , its promoter’s hypomethylation and expression have been linked to T cell infiltration and recruitment in vitro in pancreatic cancer 61 as well as a good prognosis biomarker in intrahepatic cholangiocarcinoma, where it was associated with higher activity of various innate immune pathways 62 . OAS2 is an innate immune protein involved in antiviral response 63 . ZBP1 , on the other hand, is an immune sensor of Z conformation nucleic acids and participates in the immune response against different pathogens 64 . This gene has been involved in the regulation of immune response in cancer through the activation of necroptosis or PANoptosis 65,66 . SAMD9L has a clear role promoting antiviral immunity and has been identified as a tumor suppressor 67 . Finally, CARD16 plays a less clear role in tumors. This protein belongs to the caspase-1 protein subfamily 68 , which is generally involved in proinflammatory mechanisms and the regulation of the inflammasome activation 69 . Given their tumor-intrinsic expression, potential epigenetic regulation, and functional implications for the innate immune response, the five genes associated with promoter cassette 10 may play an important role in shaping the composition of the TIME in TNBC, but potentially also in other breast cancer subtypes as suggested by our exploratory analysis of the general TCGA breast cancer cohort. A notable feature of the expression pattern of the five genes is the apparent intra-tumor heterogeneity with both expressing and non-expressing tumor cells within the same tissue specimen. Interestingly, this heterogeneity mirrors recent findings by us for PD-L1 in both TNBC cell lines and in primary tumors of a specific proposed epitype of TNBC (Basal3) 21 . Because we cannot at this point reliably deconstruct tumor methylome clonality, particularly in tumor purity adjusted bulk data due to the nature of the processing, it remains unclear whether the observed heterogeneity in expression is mirrored by a similar heterogeneity in promoter methylation status, or if it is instead determined by microenvironmental signals or stochastic expression. Furthermore, an important limitation of the current study is that causal conclusions cannot be drawn due to its observational nature. It is possible that the observed promoter methylation patterns could reflect downstream effects of immune-related signaling, particularly hypomethylation induced by interferon-γ signaling, as the specific genes are interferon-inducible. Nevertheless, several findings support the hypothesis of a causal relationship between hypermethylation of the detected genes and an “immune cold” TIME. First, the detected hypomethylation pattern is not exclusive to samples with high TIL levels. Some samples with low immune responses measured through TILs show the methylation phenotype, which could indicate that it is not a consequence of immune response-related signaling pathways and that alternative alterations may be more important in some cases. Second, CpGs associated with OAS2 and GBP4 are generally hypomethylated in normal breast tissue, implying that the TIME-associated CpG sites may become hypermethylated relative to a hypomethylated baseline, rather than undergoing demethylation. Third, prior analyses have identified a relationship between T-cell recruitment and hypomethylation of GBP4 promoter CpGs in vitro in pancreatic cancer, suggesting a causal role of hypermethylation in acquiring an immune cold phenotype 61 . Altogether, these observations support a hypothesis that hypermethylation and lack of expression of the identified genes are linked to reduced immune infiltration in TNBC, particularly in specific molecular phenotypes. However, whether these methylation patterns depend on the cell of origin of the tumor or are acquired during tumorigenesis remains unclear and cannot be assessed by the data presented in this work. Additional experimental validation in vitro and in vivo is needed to clarify the role of the detected methylation patterns in tumors and the relevance of gene expression in different normal breast cell types to establish stronger causal relationships and help identify the origin of the observed epigenetic features. In conclusion, this work outlines a new comprehensive framework to analyze bulk tumor DNA methylation data, combining tumor purity adjustment, functional CpG filtering, stratification by CpG contexts, and identification of highly correlated CpG modules to enhance the identification of tumor-intrinsic DNA methylation patterns and their biological associations. Importantly, this approach extends beyond TNBC and may help to deconstruct DNA methylation patterns and thereby molecular heterogeneity also in other malignancies. In TNBC, this methodology, in combination with additional -omics layers and in situ data, identified CpG cassettes not only associated with the major molecular subtypes of the disease, but also smaller cassettes linked to tumor suppressor gene inactivation. Moreover, it also identified cassettes correlating with TIME status that appeared enriched in distinct molecular phenotypes of TNBC as well as other tumor subgroups, lending support to the broader concept of epigenetic immunoediting, as previously proposed in glioblastoma 70 . Intriguing aspects that remain to be answered with respect to potential epigenetic immunoediting relate to the plasticity of these changes during treatment with for instance immune checkpoint inhibitors and whether treatment responses differ depending on the different mechanisms. Declarations ACKNOWLEDGEMENTS First, the authors would like to acknowledge Jelmar Quist and Rosamund Numah for their assistance with the processing the CosMx spatial transcriptomics data. Additionally, we would like to acknowledge all patients and clinicians participating in the SCAN-B study, personnel at the central SCAN-B laboratory at the Division of Oncology, Department of Clinical Sciences Lund, Lund University, the Swedish national breast cancer quality registry (NKBC), Regional Cancer Centre South, RBC Syd, and the South Sweden Breast Cancer Group (SSBCG). FUNDING Financial support for this study was provided by the Swedish Cancer Society (CAN 2021/1407 JS, 2024/3591 JS), the Mrs Berta Kamprad Foundation (FBKS-2020-5 JS and FBKS-2024-14 JS), the Swedish Research Council (2021-01800, 2025-02643 JS), the BCF-VÖS Foundation (JS), the National Society of Breast Cancer Associations in Sweden (JS), and Swedish governmental funding (ALF, grant 2022/0021 JS). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the paper. Anita Grigoriadis acknowledges support by Breast Cancer Now (147; KCL-BCN-Q3), We acknowledge use of the Spatial Biology Facility at King’s College London, supported in part by the MRC (MR/X012476/1) and the CRUK City of London Centre (CANCTA-2022/100001). AUTHORS CONTRIBUTIONS Conception and design : JS and IS Collection and assembly of data : IS, AG, MJ, DFN, MA Provision of study material or patients : JVC Data analysis and interpretation : IS with support of JS and DFN Financial support : JS Administrative support : JS Manuscript writing : IS with the support of all authors Final approval of manuscript : All authors Agree to be accountable for all aspects of the work : All authors COMPETING INTERESTS The authors declare that they have no financial or non-financial competing interests. References Hanahan, D. Hallmarks of Cancer: New Dimensions. Cancer Discov 12 , 31–46 (2022). https://doi.org/10.1158/2159-8290.CD-21-1059 Lakshminarasimhan, R. & Liang, G. The Role of DNA Methylation in Cancer. Adv Exp Med Biol 945 , 151–172 (2016). https://doi.org/10.1007/978-3-319-43624-1_7 Glodzik, D. et al. Comprehensive molecular comparison of BRCA1 hypermethylated and BRCA1 mutated triple negative breast cancers. Nat Commun 11 , 3747 (2020). https://doi.org/10.1038/s41467-020-17537-2 Kalachand, R. D. et al. BRCA1 Promoter Methylation and Clinical Outcomes in Ovarian Cancer: An Individual Patient Data Meta-Analysis. J Natl Cancer Inst 112 , 1190–1203 (2020). https://doi.org/10.1093/jnci/djaa070 Li, C. et al. Mi-2beta promotes immune evasion in melanoma by activating EZH2 methylation. Nat Commun 15 , 2163 (2024). https://doi.org/10.1038/s41467-024-46422-5 Jung, H. et al. DNA methylation loss promotes immune evasion of tumours with high mutation and copy number load. Nat Commun 10 , 4278 (2019). https://doi.org/10.1038/s41467-019-12159-9 Shigeyasu, K. et al. Clinical Significance of MLH1 Methylation and CpG Island Methylator Phenotype as Prognostic Markers in Patients with Gastric Cancer. PLoS One 10 , e0130409 (2015). https://doi.org/10.1371/journal.pone.0130409 Fackler, M. J. et al. DNA methylation markers predict recurrence-free interval in triple-negative breast cancer. NPJ Breast Cancer 6 , 3 (2020). https://doi.org/10.1038/s41523-020-0145-3 Rauluseviciute, I., Drablos, F. & Rye, M. B. DNA methylation data by sequencing: experimental approaches and recommendations for tools and pipelines for data analysis. Clin Epigenetics 11 , 193 (2019). https://doi.org/10.1186/s13148-019-0795-x Jaffe, A. E. & Irizarry, R. A. Accounting for cellular heterogeneity is critical in epigenome-wide association studies. Genome Biol 15 , R31 (2014). https://doi.org/10.1186/gb-2014-15-2-r31 Wang, Y. et al. Prognostic and Predictive Value of Three DNA Methylation Signatures in Lung Adenocarcinoma. Front Genet 10 , 349 (2019). https://doi.org/10.3389/fgene.2019.00349 DiNome, M. L. et al. Clinicopathological Features of Triple-Negative Breast Cancer Epigenetic Subtypes. Ann Surg Oncol 26 , 3344–3353 (2019). https://doi.org/10.1245/s10434-019-07565-8 Lin, L. H. et al. DNA Methylation Identifies Epigenetic Subtypes of Triple-Negative Breast Cancers With Distinct Clinicopathologic and Molecular Features. Mod Pathol 36 , 100306 (2023). https://doi.org/10.1016/j.modpat.2023.100306 Bourgon, R., Gentleman, R. & Huber, W. Independent filtering increases detection power for high-throughput experiments. Proc Natl Acad Sci U S A 107 , 9546–9551 (2010). https://doi.org/10.1073/pnas.0914005107 Sandoval, J. et al. A prognostic DNA methylation signature for stage I non-small-cell lung cancer. J Clin Oncol 31 , 4140–4147 (2013). https://doi.org/10.1200/JCO.2012.48.5516 Staaf, J. & Aine, M. Tumor purity adjusted beta values improve biological interpretability of high-dimensional DNA methylation data. PLoS One 17 , e0265557 (2022). https://doi.org/10.1371/journal.pone.0265557 Sasiain, I., Nacer, D. F., Aine, M., Veerla, S. & Staaf, J. Tumor purity estimated from bulk DNA methylation can be used for adjusting beta values of individual samples to better reflect tumor biology. NAR Genom Bioinform 6 , lqae146 (2024). https://doi.org/10.1093/nargab/lqae146 Qin, Y., Feng, H., Chen, M., Wu, H. & Zheng, X. InfiniumPurify: An R package for estimating and accounting for tumor purity in cancer methylation research. Genes Dis 5 , 43–45 (2018). https://doi.org/10.1016/j.gendis.2018.02.003 Aran, D., Sabato, S. & Hellman, A. DNA methylation of distal regulatory sites characterizes dysregulation of cancer genes. Genome Biol 14 , R21 (2013). https://doi.org/10.1186/gb-2013-14-3-r21 Aran, D. & Hellman, A. DNA methylation of transcriptional enhancers and cancer predisposition. Cell 154 , 11–13 (2013). https://doi.org/10.1016/j.cell.2013.06.018 Aine, M. et al. The DNA methylation landscape of primary triple-negative breast cancer. Nat Commun 16 , 3041 (2025). https://doi.org/10.1038/s41467-025-58158-x Harbeck, N. et al. Breast cancer. Nat Rev Dis Primers 5 , 66 (2019). https://doi.org/10.1038/s41572-019-0111-2 Zagami, P. & Carey, L. A. Triple negative breast cancer: Pitfalls and progress. NPJ Breast Cancer 8 , 95 (2022). https://doi.org/10.1038/s41523-022-00468-0 Parker, J. S. et al. Supervised risk predictor of breast cancer based on intrinsic subtypes. J Clin Oncol 27 , 1160–1167 (2009). https://doi.org/10.1200/JCO.2008.18.1370 Lehmann, B. D. et al. Identification of human triple-negative breast cancer subtypes and preclinical models for selection of targeted therapies. J Clin Invest 121 , 2750–2767 (2011). https://doi.org/10.1172/JCI45014 Lehmann, B. D. et al. Refinement of Triple-Negative Breast Cancer Molecular Subtypes: Implications for Neoadjuvant Chemotherapy Selection. PLoS One 11 , e0157368 (2016). https://doi.org/10.1371/journal.pone.0157368 Staaf, J. et al. Whole-genome sequencing of triple-negative breast cancers in a population-based clinical study. Nat Med 25 , 1526–1533 (2019). https://doi.org/10.1038/s41591-019-0582-4 Huang, R., Ding, P. & Yang, F. Clinicopathological significance and potential drug target of CDH1 in breast cancer: a meta-analysis and literature review. Drug Des Devel Ther 9 , 5277–5285 (2015). https://doi.org/10.2147/DDDT.S86929 Roostee, S. et al. Tumour immune characterisation of primary triple-negative breast cancer using automated image quantification of immunohistochemistry-stained immune cells. Sci Rep 14 , 21417 (2024). https://doi.org/10.1038/s41598-024-72306-1 Karn, T. et al. Association Between Genomic Metrics and Immune Infiltration in Triple-Negative Breast Cancer. JAMA Oncol 3 , 1707–1711 (2017). https://doi.org/10.1001/jamaoncol.2017.2140 Saal, L. H. et al. The Sweden Cancerome Analysis Network - Breast (SCAN-B) Initiative: a large-scale multicenter infrastructure towards implementation of breast cancer genomic analyses in the clinical routine. Genome Med 7 , 20 (2015). https://doi.org/10.1186/s13073-015-0131-9 Ryden, L. et al. Minimizing inequality in access to precision medicine in breast cancer by real-time population-based molecular analysis in the SCAN-B initiative. Br J Surg 105 , e158–e168 (2018). https://doi.org/10.1002/bjs.10741 Staaf, J. et al. RNA sequencing-based single sample predictors of molecular subtype and risk of recurrence for clinical assessment of early-stage breast cancer. NPJ Breast Cancer 8 , 94 (2022). https://doi.org/10.1038/s41523-022-00465-3 Fredlund, E. et al. The gene expression landscape of breast cancer is shaped by tumor protein p53 status and epithelial-mesenchymal transition. Breast Cancer Res 14 , R113 (2012). https://doi.org/10.1186/bcr3236 Aine, M. et al. Molecular analyses of triple-negative breast cancer in the young and elderly. Breast Cancer Res 23 , 20 (2021). https://doi.org/10.1186/s13058-021-01392-0 Sigurjonsdottir, G. et al. Evaluation of alternative prognostic thresholds for SP142 and 22C3 immunohistochemical PD-L1 expression in triple-negative breast cancer: results from a population-based cohort. Breast Cancer Res Treat 210 , 271–284 (2025). https://doi.org/10.1007/s10549-024-07561-x Fishilevich, S. et al. GeneHancer: genome-wide integration of enhancers and target genes in GeneCards. Database (Oxford) 2017 (2017). https://doi.org/10.1093/database/bax028 Cha, Y. J. et al. Genomic Characteristics Related to Histology-Based Immune Features in Breast Cancer. Mod Pathol 38 , 100736 (2025). https://doi.org/10.1016/j.modpat.2025.100736 Hair, B. Y. et al. Body mass index associated with genome-wide methylation in breast tissue. Breast Cancer Res Treat 151 , 453–463 (2015). https://doi.org/10.1007/s10549-015-3401-8 Reinius, L. E. et al. Differential DNA methylation in purified human blood cells: implications for cell lineage and studies on disease susceptibility. PLoS One 7 , e41361 (2012). https://doi.org/10.1371/journal.pone.0041361 Reed, A. D. et al. A single-cell atlas enables mapping of homeostatic cellular shifts in the adult human breast. Nat Genet 56 , 652–662 (2024). https://doi.org/10.1038/s41588-024-01688-9 Jovanovic, B. et al. Heterogeneity and transcriptional drivers of triple-negative breast cancer. Cell Rep 42 , 113564 (2023). https://doi.org/10.1016/j.celrep.2023.113564 Chen, Y., Pal, B., Lindeman, G. J., Visvader, J. E. & Smyth, G. K. R code and downstream analysis objects for the scRNA-seq atlas of normal and tumorigenic human breast tissue. Sci Data 9 , 96 (2022). https://doi.org/10.1038/s41597-022-01236-2 Langfelder, P. & Horvath, S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics 9 , 559 (2008). https://doi.org/10.1186/1471-2105-9-559 Corces, M. R. et al. The chromatin accessibility landscape of primary human cancers. Science 362 (2018). https://doi.org/10.1126/science.aav1898 Chakravarthy, A. et al. Author Correction: Pan-cancer deconvolution of tumour composition using DNA methylation. Nat Commun 9 , 4642 (2018). https://doi.org/10.1038/s41467-018-07155-4 Newman, A. M. et al. Determining cell type abundance and expression from bulk tissues with digital cytometry. Nat Biotechnol 37 , 773–782 (2019). https://doi.org/10.1038/s41587-019-0114-2 Phipson, B., Maksimovic, J. & Oshlack, A. missMethyl: an R package for analyzing data from Illumina's HumanMethylation450 platform. Bioinformatics 32 , 286–288 (2016). https://doi.org/10.1093/bioinformatics/btv560 Liberzon, A. et al. The Molecular Signatures Database (MSigDB) hallmark gene set collection. Cell Syst 1 , 417–425 (2015). https://doi.org/10.1016/j.cels.2015.12.004 Kursa, M. B. & Rudnicki, W. R. Feature Selection with the Boruta Package. J Stat Softw 36 , 1–13 (2010). https://doi.org/DOI 10.18637/jss.v036.i11 Butler, A., Hoffman, P., Smibert, P., Papalexi, E. & Satija, R. Integrating single-cell transcriptomic data across different conditions, technologies, and species. Nat Biotechnol 36 , 411–420 (2018). https://doi.org/10.1038/nbt.4096 Palla, G. et al. Squidpy: a scalable framework for spatial omics analysis. Nat Methods 19 , 171–178 (2022). https://doi.org/10.1038/s41592-021-01358-2 Wolf, F. A., Angerer, P. & Theis, F. J. SCANPY: large-scale single-cell gene expression data analysis. Genome Biol 19 , 15 (2018). https://doi.org/10.1186/s13059-017-1382-0 Pedregosa, F. et al. Scikit-learn: Machine Learning in Python. J Mach Learn Res 12 , 2825–2830 (2011). Jjingo, D., Conley, A. B., Yi, S. V., Lunyak, V. V. & Jordan, I. K. On the presence and role of human gene-body DNA methylation. Oncotarget 3 , 462–474 (2012). https://doi.org/10.18632/oncotarget.497 Loi, S. et al. The journey of tumor-infiltrating lymphocytes as a biomarker in breast cancer: clinical utility in an era of checkpoint inhibition. Ann Oncol 32 , 1236–1244 (2021). https://doi.org/10.1016/j.annonc.2021.07.007 Harano, K. et al. Rates of immune cell infiltration in patients with triple-negative breast cancer by molecular subtype. PLoS One 13 , e0204513 (2018). https://doi.org/10.1371/journal.pone.0204513 Tretina, K., Park, E. S., Maminska, A. & MacMicking, J. D. Interferon-induced guanylate-binding proteins: Guardians of host defense in health and disease. J Exp Med 216 , 482–500 (2019). https://doi.org/10.1084/jem.20182031 Jang, J. H., Kim, D. H. & Surh, Y. J. Dynamic roles of inflammasomes in inflammatory tumor microenvironment. NPJ Precis Oncol 5 , 18 (2021). https://doi.org/10.1038/s41698-021-00154-7 Zhou, W., Yeerkenbieke, G., Zhang, Y., Zhou, M. & Li, J. Guanylate binding protein 4 shapes an inflamed tumor microenvironment and identifies immuno-hot tumors. J Cancer Res Clin Oncol 150 , 90 (2024). https://doi.org/10.1007/s00432-024-05605-9 Tasiheng, Y. et al. DNA hypo-methylation and expression of GBP4 induces T cell exhaustion in pancreatic cancer. Cancer Immunol Immunother 73 , 208 (2024). https://doi.org/10.1007/s00262-024-03786-3 Liao, H. et al. Whole-Genome DNA Methylation Profiling of Intrahepatic Cholangiocarcinoma Reveals Prognostic Subtypes with Distinct Biological Drivers. Cancer Res 84 , 1747–1763 (2024). https://doi.org/10.1158/0008-5472.CAN-23-3298 Merold, V. et al. Structural basis for OAS2 regulation and its antiviral function. Mol Cell 85 , 2176–2193 e2113 (2025). https://doi.org/10.1016/j.molcel.2025.05.001 Xie, F. et al. ZBP1 condensate formation synergizes Z-NAs recognition and signal transduction. Cell Death Dis 15 , 487 (2024). https://doi.org/10.1038/s41419-024-06889-y Zhang, T. et al. ADAR1 masks the cancer immunotherapeutic promise of ZBP1-driven necroptosis. Nature 606 , 594–602 (2022). https://doi.org/10.1038/s41586-022-04753-7 Baik, J. Y. et al. ZBP1 not RIPK1 mediates tumor necroptosis in breast cancer. Nat Commun 12 , 2666 (2021). https://doi.org/10.1038/s41467-021-23004-3 Peng, S. et al. Structure and function of an effector domain in antiviral factors and tumor suppressors SAMD9 and SAMD9L. Proc Natl Acad Sci U S A 119 (2022). https://doi.org/10.1073/pnas.2116550119 Holland, M., Rutkowski, R. & T, C. L. Evolutionary Dynamics of Proinflammatory Caspases in Primates and Rodents. Mol Biol Evol 41 (2024). https://doi.org/10.1093/molbev/msae220 Devi, S. et al. CARD-only proteins regulate in vivo inflammasome responses and ameliorate gout. Cell Rep 42 , 112265 (2023). https://doi.org/10.1016/j.celrep.2023.112265 Gangoso, E. et al. Glioblastomas acquire myeloid-affiliated transcriptional programs via epigenetic immunoediting to elicit immune evasion. Cell 184 , 2454–2470 e2426 (2021). https://doi.org/10.1016/j.cell.2021.03.023 Additional Declarations No competing interests reported. Supplementary Files supplementarytables.pdf supplementaryfigures.pdf SupplementaryDataFiles.zip Cite Share Download PDF Status: Under Review Version 1 posted Editorial decision: Revision requested 12 Jan, 2026 Reviews received at journal 12 Jan, 2026 Reviews received at journal 05 Jan, 2026 Reviewers agreed at journal 22 Dec, 2025 Reviewers agreed at journal 21 Dec, 2025 Reviewers agreed at journal 20 Dec, 2025 Reviewers agreed at journal 09 Dec, 2025 Reviewers invited by journal 09 Dec, 2025 Editor assigned by journal 02 Dec, 2025 Submission checks completed at journal 02 Dec, 2025 First submitted to journal 01 Dec, 2025 You are reading this latest preprint version Research Square lets you share your work early, gain feedback from the community, and start making changes to your manuscript prior to peer review in a journal. As a division of Research Square Company, we’re committed to making research communication faster, fairer, and more useful. We do this by developing innovative software and high quality services for the global research community. Our growing team is made up of researchers and industry professionals working together to solve the most critical problems facing scientific publishing. Also discoverable on Platform About Our Team In Review Editorial Policies Advisory Board Help Center Resources Author Services Accessibility API Access RSS feed Manage Cookie Preferences © Research Square 2026 | ISSN 2693-5015 (online) Privacy Policy Terms of Service Do Not Sell My Personal Information {"props":{"pageProps":{"initialData":{"identity":"rs-8250160","acceptedTermsAndConditions":true,"allowDirectSubmit":false,"archivedVersions":[],"articleType":"Article","associatedPublications":[],"authors":[{"id":559008833,"identity":"331bdc88-740b-42cb-b950-1ae9f9788b7e","order_by":0,"name":"Inaki Sasiain","email":"","orcid":"","institution":"Lund University","correspondingAuthor":false,"prefix":"","firstName":"Inaki","middleName":"","lastName":"Sasiain","suffix":""},{"id":559008834,"identity":"f3caeb26-2a72-4229-a9f1-e4ed0c6cca7f","order_by":1,"name":"Deborah F Nacer","email":"","orcid":"","institution":"Lund University","correspondingAuthor":false,"prefix":"","firstName":"Deborah","middleName":"F","lastName":"Nacer","suffix":""},{"id":559008835,"identity":"438b416f-7344-4322-a87e-c2d4aebaa7d3","order_by":2,"name":"Mats Jönsson","email":"","orcid":"","institution":"Lund University","correspondingAuthor":false,"prefix":"","firstName":"Mats","middleName":"","lastName":"Jönsson","suffix":""},{"id":559008839,"identity":"a798dfb2-f762-4b0f-805a-7c1d5e35e301","order_by":3,"name":"Johan Vallon-Christersson","email":"","orcid":"","institution":"Lund University","correspondingAuthor":false,"prefix":"","firstName":"Johan","middleName":"","lastName":"Vallon-Christersson","suffix":""},{"id":559008840,"identity":"e7aba75d-e6a6-450a-a28c-ce46476dac3f","order_by":4,"name":"Anita Giorgiadis","email":"","orcid":"","institution":"King's College London","correspondingAuthor":false,"prefix":"","firstName":"Anita","middleName":"","lastName":"Giorgiadis","suffix":""},{"id":559008842,"identity":"440d55a7-eb8b-47ff-bb56-7ef57d391abe","order_by":5,"name":"Mattias Aine","email":"","orcid":"","institution":"Lund University","correspondingAuthor":false,"prefix":"","firstName":"Mattias","middleName":"","lastName":"Aine","suffix":""},{"id":559008843,"identity":"c5fe8073-5d15-43e5-948f-9aafe8a057fe","order_by":6,"name":"Johan Staaf","email":"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZAAAAAyAQMAAABI0h/eAAAABlBMVEX///8AAABVwtN+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA8ElEQVRIiWNgGAWjYBACxgYeOPPhAwYGZgY+YrVIAFUbG4C0sBG2B6HFTIIoLcztvQc/FzDU1fG3H2arulFhzcAmkcD84QM+h/WcS5aewXBYQuJMMtvtnDPpIC1skjPwaZmRYyDNw3BAguFA/rHbuW2HwVqYefBrMf7Nw1AnIX/+MVtx7j+wFubPf/BrMQPawixhcCOZjTm3AayFQRqf9xl7zphZ8xgcltx44zGzdM6xdB42nodtkj14tBi29xjf5qmo45c7n8z4OafGWo6fPfnwhx/4tDSASAOEAA8oevG5i0Eer+woGAWjYBSMAhAAANF8QUPUArT/AAAAAElFTkSuQmCC","orcid":"","institution":"Lund University","correspondingAuthor":true,"prefix":"","firstName":"Johan","middleName":"","lastName":"Staaf","suffix":""}],"badges":[],"createdAt":"2025-12-01 12:08:20","currentVersionCode":1,"declarations":"","doi":"10.21203/rs.3.rs-8250160/v1","doiUrl":"https://doi.org/10.21203/rs.3.rs-8250160/v1","draftVersion":[],"editorialEvents":[],"editorialNote":"","failedWorkflow":false,"files":[{"id":98056033,"identity":"0186fe13-c1cb-45d4-b7f5-92fdcabc2d15","added_by":"auto","created_at":"2025-12-12 10:01:25","extension":"docx","order_by":0,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":5127770,"visible":true,"origin":"","legend":"","description":"","filename":"manuscriptdraft.docx","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/bbd8542709979df1438e1faf.docx"},{"id":98056028,"identity":"dc20f408-86a9-4ee5-9b8a-cf24412de2b8","added_by":"auto","created_at":"2025-12-12 10:01:25","extension":"json","order_by":1,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":12616,"visible":true,"origin":"","legend":"","description":"","filename":"5c8ef91ac08a4af18ea05fa947da60cd.json","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/4e11b193545935bc2bd802c4.json"},{"id":98056054,"identity":"6dc624c6-b662-41c3-aa80-74b65600a71f","added_by":"auto","created_at":"2025-12-12 10:01:26","extension":"zip","order_by":2,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":55753750,"visible":true,"origin":"","legend":"","description":"","filename":"SupplementaryDataFiles.zip","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/1d8ac2966821ac3e48933295.zip"},{"id":98428354,"identity":"4ab6efe4-def4-474f-a164-ece55acdb8dc","added_by":"auto","created_at":"2025-12-17 16:41:56","extension":"pdf","order_by":3,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":10200972,"visible":true,"origin":"","legend":"","description":"","filename":"supplementaryfigures.pdf","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/4ba3b767a11d484dc38cabe9.pdf"},{"id":98427632,"identity":"37a8ebab-c5f3-46a8-9053-003becd73563","added_by":"auto","created_at":"2025-12-17 16:40:54","extension":"pdf","order_by":4,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":57685,"visible":true,"origin":"","legend":"","description":"","filename":"supplementarytables.pdf","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/1ac58e217d8ce94ac9dd0520.pdf"},{"id":98428646,"identity":"b95bc989-81a9-4f95-a3d3-ef0fcd61285b","added_by":"auto","created_at":"2025-12-17 16:42:13","extension":"jpeg","order_by":5,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":866452,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage1.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/8b422f4b95ea65529f0eca19.jpeg"},{"id":98056037,"identity":"12ae1757-7f91-4420-bcb2-d12268a05141","added_by":"auto","created_at":"2025-12-12 10:01:25","extension":"jpeg","order_by":6,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":604921,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage2.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/95b9417e11ddb0ea99f9e8da.jpeg"},{"id":98428396,"identity":"be518a57-b0e7-424a-8e26-66b256789043","added_by":"auto","created_at":"2025-12-17 16:41:59","extension":"jpeg","order_by":7,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":601325,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage3.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/34e9fe87533397ed41252cbe.jpeg"},{"id":98056038,"identity":"d325ada5-1ac6-414a-a3da-d7e1473e8909","added_by":"auto","created_at":"2025-12-12 10:01:26","extension":"jpeg","order_by":8,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":518829,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage4.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/0bf8f17f33d2ee293a6d2443.jpeg"},{"id":98428585,"identity":"729a79d2-a4e6-46c0-b7a8-bd54c37a4fd4","added_by":"auto","created_at":"2025-12-17 16:42:10","extension":"jpeg","order_by":9,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":736595,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage5.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/07d2871a52c56f5bdf177072.jpeg"},{"id":98056041,"identity":"7cb32d64-8532-4c85-bf81-e5b774457c06","added_by":"auto","created_at":"2025-12-12 10:01:26","extension":"jpeg","order_by":10,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":623934,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage6.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/81c9fa6e8b62cfe3d0efa3c6.jpeg"},{"id":98427974,"identity":"72eecce9-0fb7-4ef2-8d10-e56d46507dea","added_by":"auto","created_at":"2025-12-17 16:41:27","extension":"jpeg","order_by":11,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":476938,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage7.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/c0a02e9598eafd1f6687d758.jpeg"},{"id":98056048,"identity":"2a66d105-b817-430b-8e82-2ead45ee8071","added_by":"auto","created_at":"2025-12-12 10:01:26","extension":"jpeg","order_by":12,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":719144,"visible":true,"origin":"","legend":"","description":"","filename":"floatimage8.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/7298a7768583d6d4bf1f4caf.jpeg"},{"id":98428512,"identity":"98779331-905e-4077-802f-7b26d2034286","added_by":"auto","created_at":"2025-12-17 16:42:06","extension":"png","order_by":13,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":344277,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage1.png","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/2477b1ff86010c66a409eb3b.png"},{"id":98056043,"identity":"5ad1f7c9-2c34-407e-99cf-53f9a4d07eeb","added_by":"auto","created_at":"2025-12-12 10:01:26","extension":"png","order_by":14,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":171279,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage2.png","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/d2e69b4ca7ca097d1b2d0b3c.png"},{"id":98056045,"identity":"153132c0-106f-4448-b2b6-5b338bd4910f","added_by":"auto","created_at":"2025-12-12 10:01:26","extension":"png","order_by":15,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":226550,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage3.png","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/8a3339e8e6c743d5e5039c5e.png"},{"id":98427831,"identity":"7f2abf0c-743e-4c1e-a995-c142b38e31b5","added_by":"auto","created_at":"2025-12-17 16:41:16","extension":"png","order_by":16,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":159571,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage4.png","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/b0adbc5848997facccd9f8da.png"},{"id":98428398,"identity":"564d622b-f873-4594-b9f2-e56efc2191f0","added_by":"auto","created_at":"2025-12-17 16:41:59","extension":"png","order_by":17,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":197260,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage5.png","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/05747ce368d4d4a6af17cdcf.png"},{"id":98056053,"identity":"2a53e870-fb20-4876-9662-ff7208dddccf","added_by":"auto","created_at":"2025-12-12 10:01:26","extension":"png","order_by":18,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":167993,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage6.png","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/763530972c1ae0166f364011.png"},{"id":98429004,"identity":"c94505f9-bfc1-42c1-81d2-37d69e63c2c3","added_by":"auto","created_at":"2025-12-17 16:42:41","extension":"png","order_by":19,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":135267,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage7.png","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/0e39422eace6eb8b3e1b5e3b.png"},{"id":98056049,"identity":"286cf524-bdd9-4bba-9831-a2c744d46a76","added_by":"auto","created_at":"2025-12-12 10:01:26","extension":"png","order_by":20,"title":"","display":"","copyAsset":false,"role":"acdc-reference","size":200778,"visible":true,"origin":"","legend":"","description":"","filename":"Onlinefloatimage8.png","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/ec5d60f5a9a4ab133652e794.png"},{"id":98056026,"identity":"13c72695-dd15-466d-a09e-e6a701cf1752","added_by":"auto","created_at":"2025-12-12 10:01:25","extension":"jpeg","order_by":1,"title":"Figure 1","display":"","copyAsset":false,"role":"figure","size":866452,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eDetection of CpG cassettes.\u003c/strong\u003e\u003cem\u003e \u003c/em\u003e\u003cstrong\u003eA\u003c/strong\u003e) Impact of the soft-thresholding parameter (β) on number of detected cassettes (bars), cassette size (red line) and proportion of clustered and unclassified (Uncl.) CpGs from the 717458 CpGs of the discovery cohort. A higher β results in smaller cassettes comprising a lower portion of all CpGs, and higher proportion of unclassified CpGs. \u003cstrong\u003eB\u003c/strong\u003e) First seven cassettes detected using tumor purity-adjusted and unadjusted beta values. Columns correspond to samples and rows to CpGs. \u003cstrong\u003eC\u003c/strong\u003e) Cell type admixture estimated from unadjusted beta values using MethylCIBERSORT (top) and change in tumor fraction before and after tumor purity adjustment (bottom). \u003cem\u003eCD14\u003c/em\u003estands for macrophages and monocytes, \u003cem\u003eCD19\u003c/em\u003e for B cells, \u003cem\u003eCD4_Eff\u003c/em\u003eStands for Helper T lymphocytes, \u003cem\u003eCD56\u003c/em\u003e for NK cells, \u003cem\u003eCD8\u003c/em\u003e for cytotoxic T lymphocytes, \u003cem\u003eEnd\u003c/em\u003e for endothelial cells, \u003cem\u003eEos\u003c/em\u003e for Eosinophile, \u003cem\u003eFib\u003c/em\u003e for Fibroblasts, \u003cem\u003eNeu\u003c/em\u003e for neutrophile and \u003cem\u003eTreg\u003c/em\u003efor Regulatory T-cell. \u003cstrong\u003eD\u003c/strong\u003e) Proportion of CpG contexts after selecting the N most varying CpGs before (left) or after (right) tumor purity adjustment. \u003cstrong\u003eE\u003c/strong\u003e) Distribution of CpG contexts in the input data used for cassette detection (left), and deviations in CpG context composition of the 4 largest detected cassettes and unclassified CpGs compared to the original distribution (right). \u003cstrong\u003eF\u003c/strong\u003e) DNA Methylation heatmaps for the seven largest CpG cassettes determined in distal, promoter, and proximal CpG contexts respectively from purity-adjusted methylation data. \u003cstrong\u003eG\u003c/strong\u003e) Stability analysis of five distal cassettes detected based on 20 random samplings of 100 tumors each. Top row: for each CpG site originally assigned to a given cassette, the proportion of resamplings in which it was assigned to the most common cassette (red), to any cassette other than the most common (black) or was not assigned to any cassette (gray) is shown. A high red line indicates consistent reassignment to the same cassette across resamplings. Bottom row: proportional cassette assignments for each CpG site across resamplings, with colors representing the cassettes to which CpGs were assigned. Uniform colors within a cassette indicate strong preservation of the original group structure.\u003c/p\u003e","description":"","filename":"floatimage1.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/58e50f2d0698ed259f9d5476.jpeg"},{"id":98427732,"identity":"b0a09b08-091f-4b96-92c7-19bad5e040d1","added_by":"auto","created_at":"2025-12-17 16:41:02","extension":"jpeg","order_by":2,"title":"Figure 2","display":"","copyAsset":false,"role":"figure","size":604921,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eThe identified cassettes reflect the Basal/non-Basal division. A\u003c/strong\u003e) Adjusted Rand Index between hierarchical clustering from each of the largest CpG cassettes (Euclidean distance, Ward.D2 linkage) and PAM50 Basal/non-Basal division. The cassettes whose ARI\u0026gt;0.65 are marked with asterisks. \u003cstrong\u003eB\u003c/strong\u003e) Distribution of PAM50 subtypes in sample clusters obtained from beta values of all CpGs included in cassettes with ARI\u0026gt;0.65 (n=9856 CpGs) using a hierarchical approach (Euclidean distance, Ward.D2 linkage). \u003cstrong\u003eC\u003c/strong\u003e) Gene set overrepresentation analysis of genes linked to CpGs in promoter cassette 1. \u003cstrong\u003eD\u003c/strong\u003e) Volcano plot of genes with expression patterns driven by CpGs in promoter cassette 1. P-values were obtained from Wilcoxon’s tests between two clusters of samples obtained through hierarchical clustering (Euclidean distance, Ward.D2 linkage) of DNA methylation beta values of CpGs in promoter cassette 1, and the fold change in expression was calculated using the median between both groups. \u003cstrong\u003eE\u003c/strong\u003e) Gene set overrepresentation analysis of genes linked to CpGs in proximal cassette 1. \u003cstrong\u003eF\u003c/strong\u003e) Volcano plot of genes with expression patterns driven by CpGs in proximal cassette 1. The used values were calculated as done in D. \u003cstrong\u003eG\u003c/strong\u003e) CpGs (rows) affecting LDHB and SPINK8 genes and their FPKM gene expression per sample (columns). Samples were clustered based on promoter and proximal methylation patterns from CpGs linked to each gene, respectively, using k-means clustering and k=2. The row annotations of each heatmap indicate the CpGs included in promoter cassette 1 and proximal cassette 1, respectively.\u003c/p\u003e","description":"","filename":"floatimage2.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/96240f943da7b14f6ede27f4.jpeg"},{"id":98056056,"identity":"84cf6ee0-9880-42de-897a-12ce2820c364","added_by":"auto","created_at":"2025-12-12 10:01:26","extension":"jpeg","order_by":3,"title":"Figure 3","display":"","copyAsset":false,"role":"figure","size":601325,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eCharacteristics of distal cassettes. A\u003c/strong\u003e) DNA methylation values of CpGs belonging to distal cassette 1 (rows) and its summary PC1 value per sample (columns). \u003cstrong\u003eB\u003c/strong\u003e) Fraction of CpGs overlapping with enhancers and repetitive regions in the seven largest distal cassettes. \u003cstrong\u003eC\u003c/strong\u003e) Fraction of CpGs overlapping with transcription factor (TF) binding sites in the seven largest distal cassettes. TFs shown are the 25 with highest pairwise differences between any two of the analyzed cassettes. \u003cstrong\u003eD\u003c/strong\u003e) CpG density in the seven largest distal cassettes. \u003cstrong\u003eE\u003c/strong\u003e) Comparison between PC1 of cassette 1 and biological features. From left to right: PAM50 Basal/non-Basal, HRD, Lehmann TNBC subtypes (TNBCtype), ASCAT tumor purity estimated from WGS, and TILs. Cor. = correlation. \u003cstrong\u003eF\u003c/strong\u003e) Proportion of CpG cassettes per N top variable CpGs after variance-based filtering on tumor purity-adjusted methylation data. The top violin plot shows the measured variance from adjusted betas per group, and the bar plot the proportion of cassettes. CpGs in cassettes comprising below 1% of the total CpGs were grouped and identified as ”Others”. Uncl. = unclassified CpGs. \u003cstrong\u003eG\u003c/strong\u003e) Distal cassettes detected before and after chromatin accessibility-based filtering (considering overlap with ATAC-seq peaks). \u003cstrong\u003eH\u003c/strong\u003e) Sankey plots between cassettes determined from all CpGs (All) and only CpGs with ATAC-seq peak overlap (ATAC+). The NA category represents CpGs that were filtered out from the data, ”Other” includes cassettes with proportion below 2% and “Uncl.” unclassified CpGs. \u003cstrong\u003eI\u003c/strong\u003e) Sankey plots between cassettes determined from ATAC-seq filtered CpGs (ATAC+, left) and all CpGs (right). ”Other” class includes cassettes whose proportion is below 2%.\u003c/p\u003e","description":"","filename":"floatimage3.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/9f849b516764c28fd1199b64.jpeg"},{"id":98056031,"identity":"c1708c86-36cc-4f51-b38e-0c8827041353","added_by":"auto","created_at":"2025-12-12 10:01:25","extension":"jpeg","order_by":4,"title":"Figure 4","display":"","copyAsset":false,"role":"figure","size":518829,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eDetection of small cassettes linked to immune infiltration. A\u003c/strong\u003e) Association of \u003cem\u003eBRCA1\u003c/em\u003e promoter hypermethylation with promoter cassette 11. The row annotations represent from left to right CpGs in promoter cassette 11 and CpG context. \u003cstrong\u003eB\u003c/strong\u003e) Association between PC1 of cassettes and tumor infiltrating lymphocytes (TILs). The x axis corresponds to Kendall correlation coefficient, and the y axis to transformed Wilcoxon p-values when splitting samples based on methylation patterns using k-means clustering and k=2. \u003cstrong\u003eC\u003c/strong\u003e) Promoter cassettes selected by the Boruta algorithm, their Kendall correlation value with TILs and median importance from the random forest runs. Numbers on the y axis represent cassette identifiers. \u003cstrong\u003eD\u003c/strong\u003e) DNA methylation beta values of CpGs in promoter cassette 10 (rows) per tumor (columns). Tumors are grouped using k-means clustering and k=2 on the detected cassette’s beta values. \u003cstrong\u003eE\u003c/strong\u003e) Gene Ontology (GO) enrichment of genes associated with the 30 CpGs in promoter cassette 10. \u003cstrong\u003eF\u003c/strong\u003e) Immune metagene rank scores versus hypermethylated and hypomethylated sample groups defined by beta values of CpGs in promoter cassette 10 in the SCAN-B validation cohort. Significance was assessed using the Wilcoxon’s test. Sample groups are the same as in D.\u003c/p\u003e","description":"","filename":"floatimage4.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/e27fffaf1d855c0810301a29.jpeg"},{"id":98056032,"identity":"276e8047-1dad-4975-aa3e-1dddc7de4aba","added_by":"auto","created_at":"2025-12-12 10:01:25","extension":"jpeg","order_by":5,"title":"Figure 5","display":"","copyAsset":false,"role":"figure","size":736595,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eThe genes involved in promoter cassette 10 are linked to immune infiltration and show subtype specific methylation patterns. A\u003c/strong\u003e) Methylation state of CpGs overlapping with genes with CpGs included in promoter cassette 10. The samples are grouped using k-means clustering and k=2 on overlapping CpG’s adjusted beta values for each gene. FPKM gene expression per sample for each of the genes is detailed in the bottom annotation. The row annotations represent from left to right the CpGs belonging to promoter cassette 10 and CpG context. \u003cstrong\u003eB\u003c/strong\u003e) Association between observed methylation patterns per gene and FPKM gene expression. The p-values displayed were obtained using Wilcoxon’s test. \u003cstrong\u003eC\u003c/strong\u003e) Association between observed methylation patterns per gene and tumor infiltrating lymphocytes (TILs). The p-values displayed were obtained using Wilcoxon’s test. \u003cstrong\u003eD\u003c/strong\u003e) Methylation state per gene versus PDL1-CPS classes. \u003cstrong\u003eE\u003c/strong\u003e) Co-occurrence of the methylation patterns detected per gene in A (rows) for tumors (columns) of PAM50 Basal/non-Basal (left) and Lehmann TNBCtype (right) TNBC subtypes. Panels below heatmaps show how many of the five genes have the same methylation status for a sample. The p-value significance displayed on the left was obtained using Chi-square tests.\u003c/p\u003e","description":"","filename":"floatimage5.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/c2322c9ef6cbf14a03fc7171.jpeg"},{"id":98056029,"identity":"2813cf97-8a51-4ad5-92dd-a2106decd5dd","added_by":"auto","created_at":"2025-12-12 10:01:25","extension":"jpeg","order_by":6,"title":"Figure 6","display":"","copyAsset":false,"role":"figure","size":623934,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eCharacterization of TIME-associated CpG cassette and associated genes across breast cancer subtypes. A) \u003c/strong\u003eDistribution of TILs across PAM50 breast cancer subtypes in the TCGA breast cancer cohort. Statistical significance was assessed using the Kruskal-Wallis test. \u003cstrong\u003eB)\u003c/strong\u003e DNA methylation beta values of available CpGs in promoter cassette 10 (rows) across TCGA-BRCA tumor (columns). Tumors are grouped using k-means clustering and k=2 on the detected cassette’s beta values. \u003cstrong\u003eC)\u003c/strong\u003e Association between promoter cassette 10 and immune infiltration in breast cancer. Samples were clustered as performed in B. Statistical significance was assessed using Wilcoxon’s test. \u003cstrong\u003eD) \u003c/strong\u003eDistribution of methylation pattern per gene across breast cancer subtypes. Tumors are grouped using k-means clustering and k=2 on the CpGs overlapping with the genes of interest. Statistical significance was evaluated using the Chi-squared test. \u003cstrong\u003eE)\u003c/strong\u003eUpSet diagram of shared methylation status per gene across samples. The bottom annotation indicates measured TILs levels for each combination of methylation patterns. Statistical significance was evaluated two-sided binomial tests with FDR correction. \u003cstrong\u003eF) \u003c/strong\u003eAssociation of TIL levels with methylation status of the detected genes (rows) per breast cancer PAM50 subtype (columns). Clusters were defined as performed in D. Statistical significance was assessed using Wilcoxon’s test.\u003c/p\u003e","description":"","filename":"floatimage6.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/256f4dc21c2d2db36c66983f.jpeg"},{"id":98056055,"identity":"de84c3e9-5936-41f6-93d0-c2de141925e1","added_by":"auto","created_at":"2025-12-12 10:01:26","extension":"jpeg","order_by":7,"title":"Figure 7","display":"","copyAsset":false,"role":"figure","size":476938,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eDNA methylation and gene expression of five genes of interest in normal breast tissue.\u003c/strong\u003e \u003cstrong\u003eA\u003c/strong\u003e) Methylation state of CpGs in promoter cassette 10 affecting specific genes in 96 normal breast tissue specimens from Hair et al.\u003csup\u003e39\u003c/sup\u003e. The CpGs not included here were not available in the used dataset due to different methylation array versions. \u003cstrong\u003eB\u003c/strong\u003e) Expression of genes of interest in normal breast tissue using scRNA-seq data. Cells were grouped based on the assigned cell types in Reed et al.\u003csup\u003e41\u003c/sup\u003e. The average expression, displayed as the color of the dots, represents average scaled expression of cells belonging to each group. Values were scaled independently per epithelial, immune and stromal cells. The proportion of cells expressing a gene is represented by dot size.\u003c/p\u003e","description":"","filename":"floatimage7.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/88af8acb69c4eef43285cdff.jpeg"},{"id":98056034,"identity":"be0d31af-8c04-4cb3-a925-acba413577a2","added_by":"auto","created_at":"2025-12-12 10:01:25","extension":"jpeg","order_by":8,"title":"Figure 8","display":"","copyAsset":false,"role":"figure","size":719144,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003ePromoter cassette 10 genes are intrinsically expressed by tumor cells and correlated with tumor-intrinsic methylation patterns.\u003c/strong\u003e \u003cstrong\u003eA\u003c/strong\u003e) Bulk RNA-sequencing expression of \u003cem\u003eSAMD9L\u003c/em\u003e, \u003cem\u003eCARD16\u003c/em\u003e, \u003cem\u003eZBP1\u003c/em\u003e, \u003cem\u003eOAS2\u003c/em\u003e, and \u003cem\u003eGBP4\u003c/em\u003e in 34 TNBC cell lines. The color and size of the dots represent log1p-transformed FPKM values. \u003cstrong\u003eB\u003c/strong\u003e) DNA methylation of CpGs (rows) connected to the same five genes and their log1p FPKM gene expression (bottom annotation) in eight TNBC cell lines (columns) with matched data. The row annotation represents the CpGs belonging to promoter cassette 10. \u003cstrong\u003eC\u003c/strong\u003e) scRNA-seq gene expression of the five genes in tumor cells from eight TNBC tumor specimens as well as cell counts of malignant and non-malignant cells. For gene expression in tumor cells, the average expression represents average scaled expression of tumor cells belonging to each sample. For cell type counts, the color of the tile plot in the right represents the proportion of each cell type per sample, and the number displayed raw cell counts. Cell types were obtained from Chen et al. \u003csup\u003e43\u003c/sup\u003e. \u003cstrong\u003eD\u003c/strong\u003e) Example of TMA cores from one hypomethylated (PD35993a) and one hypermethylated (PD31097a) tumor used in the \u003cem\u003ein situ\u003c/em\u003e assessment of \u003cem\u003eOAS2\u003c/em\u003e expression in tumor cells by the Nanostring CosMx spatial transcriptomics platform. The images from left to right represent: standard H\u0026amp;E staining, immunofluorescence image showing panCK (used as a tumor marker) and DAPI stainings from CosMX analysis, Tumor/Rest masks created per cell based on panCK intensity and normalized and log1p-transformed \u003cem\u003eOAS2\u003c/em\u003e expression per cell. \u003cstrong\u003eE\u003c/strong\u003e) Comparison of mean \u003cem\u003eOAS2\u003c/em\u003e expression and the proportion of tumor cells expressing \u003cem\u003eOAS2\u003c/em\u003e in hypomethylated and hypermethylated samples, using sample labels defined in Figure 5A. The gene expression values used are normalized and log1p-transformed. N stands for number of unique samples. The p-values displayed were obtained using Wilcoxon’s test.\u003c/p\u003e","description":"","filename":"floatimage8.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/3e68792ea951f414ef289455.jpeg"},{"id":98622154,"identity":"b32937b7-8e52-4475-b729-95c1888be527","added_by":"auto","created_at":"2025-12-19 16:47:13","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"manuscript-pdf","size":7098859,"visible":true,"origin":"","legend":"","description":"","filename":"manuscript.pdf","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/50f397e5-7f8a-43cb-9ca0-567c395a3731.pdf"},{"id":98426205,"identity":"3b3fd053-b8d3-4f67-8cfc-41dd8f4e069f","added_by":"auto","created_at":"2025-12-17 16:35:51","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"supplement","size":57685,"visible":true,"origin":"","legend":"","description":"","filename":"supplementarytables.pdf","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/d725b393f1741fafca6d5717.pdf"},{"id":98426310,"identity":"cec33820-8e59-4478-8ee1-0d83f070145d","added_by":"auto","created_at":"2025-12-17 16:36:06","extension":"pdf","order_by":1,"title":"","display":"","copyAsset":false,"role":"supplement","size":10200972,"visible":true,"origin":"","legend":"","description":"","filename":"supplementaryfigures.pdf","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/95821fe7e98c58f0a5163fb7.pdf"},{"id":98056057,"identity":"08e37c62-cc9e-472e-9c81-078ae484afc5","added_by":"auto","created_at":"2025-12-12 10:01:28","extension":"zip","order_by":2,"title":"","display":"","copyAsset":false,"role":"supplement","size":55753750,"visible":true,"origin":"","legend":"","description":"","filename":"SupplementaryDataFiles.zip","url":"https://assets-eu.researchsquare.com/files/rs-8250160/v1/6e95405d112394ac6036429e.zip"}],"financialInterests":"No competing interests reported.","formattedTitle":"Genome-Wide DNA Methylation Profiling of Triple-Negative Breast Cancer Uncovers Epigenetic Biomarkers of Tumor Identity and the Immune Microenvironment","fulltext":[{"header":"INTRODUCTION","content":"\u003cp\u003eEpigenetic mechanisms, such as DNA methylation and histone acetylation, regulate gene expression without modifying the genome. Deregulation of these processes is acknowledged as a major driver of cancer formation and pathogenesis\u003csup\u003e1\u003c/sup\u003e. DNA methylation involves the covalent addition of methyl groups in cytosine residues located at cytosine-guanine dinucleotides (CpG) across the genome and plays a vital role in transcriptional regulation, gene imprinting, X-chromosome inactivation, and suppression of repeat elements\u003csup\u003e2\u003c/sup\u003e. In the context of cancer, DNA methylation has been involved in several relevant molecular mechanisms, like genome repair\u003csup\u003e3,4\u003c/sup\u003e, immune evasion\u003csup\u003e5,6\u003c/sup\u003e, and treatment response\u003csup\u003e7,8\u003c/sup\u003e. A deeper understanding of DNA methylation dynamics in cancer and how it shapes the tumor ecosystem is therefore a highly relevant aspect of cancer research.\u003c/p\u003e\n\u003cp\u003eDNA methylation profiling of bulk tumor tissue is frequently performed by high-density microarrays due to their cost-effectiveness, ease of use, and high CpG coverage\u003csup\u003e9\u003c/sup\u003e. However, these high dimensional data constitute a complex data modality influenced by both intrinsic confounders such as analytical noise and tumor heterogeneity, and extrinsic factors like presence of non-malignant cells with cell-type specific methylation signatures. Together, these factors cause variation in the methylation estimates of single CpGs (referred to as beta values), where tumor heterogeneity and tumor purity have been highlighted as important sources of bias in analyses if not accounted for\u003csup\u003e10\u003c/sup\u003e.\u0026nbsp;Given the complexity of DNA methylation data, different analytic frameworks have been implemented in the context of cancer. A common approach is variance-based filtering, aiming to reduce dimensionality and noise by retaining only\u0026nbsp;the \u003cem\u003eN\u003c/em\u003e most variable CpG sites prior\u0026nbsp;to unsupervised or supervised analyses\u003csup\u003e11-15\u003c/sup\u003e.\u0026nbsp;This aims to exclude non-informative CpGs, enriching for CpGs whose variability is more likely to reflect underlying biological processes.\u0026nbsp;Variance filtering, however, also has potential drawbacks as it could keep irrelevant high variance CpGs and is influenced by tumor extrinsic confounders, like the normal cell fraction. To circumvent the latter, tumor purity adjustment methods have been developed that can be applied before variance-based filtering to enrich for tumor-intrinsic signals\u003csup\u003e16-18\u003c/sup\u003e.\u0026nbsp;\u003c/p\u003e\n\u003cp\u003eAn additional layer of biological complexity in the analysis of DNA methylation data arises from CpG contexts, i.e., the genomic context in which the CpGs are located. CpGs can be divided into three broad contexts based on their relative location to genes: i) promoter CpGs, affecting the promoter regions of genes, ii) proximal CpGs, located close to the promoter regions but not in the core promoter sequence, and iii) distal CpGs, located at a greater distance from gene promoters. CpG context has functional implications; while promoter and proximal CpGs are often involved in regulation of transcription factor binding to gene promoters, directly affecting gene expression, distal CpGs have a less straightforward link to genes and can be involved in regulating binding of transcription factors to enhancer regions\u003csup\u003e2,19,20\u003c/sup\u003e. In addition, CpGs can reside in regions with different chromatin accessibility, normally referred to as open (accessible) or closed (inaccessible) chromatin as measured by e.g. assays for transposase-accessible chromatin using sequencing (ATAC-seq).\u003c/p\u003e\n\u003cp\u003eGiven the complexity of bulk tumor DNA methylation data, this study aimed to deconstruct and characterize its variance structure, considering both the genomic context of CpGs and adjustment of methylation data for tumor purity. Our hypothesis is that this approach can identify groups of CpGs (referred to as CpG cassettes) that share common methylation patterns across samples and are likely involved in related or co-occurring biological processes. As a model system for these analyses, we focused on triple-negative breast cancer (TNBC) based on recent work that demonstrated the importance of both tumor purity correction and CpG context for unsupervised analysis of DNA methylation data\u003csup\u003e21\u003c/sup\u003e. TNBC constitutes ~10-20% of breast cancer cases and is characterized by the negativity of three molecular markers: HER2/ERBB2 and the progesterone (PR) and estrogen receptors (ER)\u003csup\u003e22\u003c/sup\u003e. This subtype is highly heterogeneous with an often aggressive clinical course, shows high genomic instability and has historically lacked available targeted therapies\u003csup\u003e23\u003c/sup\u003e. TNBC has been molecularly stratified by, e.g., mRNA profiling including both the PAM50 subtypes\u003csup\u003e24\u003c/sup\u003e (~70% belong to the Basal PAM50 subtype) or TNBC specific subtypes like the Lehmann et al. subtyping scheme dividing TNBC into Luminal-Androgen Receptor (LAR), Mesenchymal (M), Basal-Like 1 (BL1) and Basal-Like 2 (BL2) tumors\u003csup\u003e25,26\u003c/sup\u003e. In agreement with PAM50, previous studies have demonstrated that the Basal/non-Basal division is highly relevant in TNBC DNA methylation dynamics, with clearly differentiated methylation patterns between tumors, representing the two main attractor states within the TNBC methylation landscape\u003csup\u003e21\u003c/sup\u003e. Besides global DNA methylation patterns, gene specific methylation alterations also play an important role in TNBC. \u003cem\u003eBRCA1\u003c/em\u003e promoter hypermethylation has been reported as the most frequent cause of a genetic homologous recombination deficiency (HRD) phenotype (observed in nearly 60% of cases) alongside mutations in the \u003cem\u003eBRCA1\u0026nbsp;\u003c/em\u003eor\u003cem\u003e\u0026nbsp;BRCA2\u003c/em\u003e genes\u003csup\u003e27\u003c/sup\u003e, and \u003cem\u003eCDH1\u0026nbsp;\u003c/em\u003ehypermethylation has been associated with tumor invasiveness\u003csup\u003e28\u003c/sup\u003e. Finally, immune infiltration introduces an additional dimension of heterogeneity in TNBC, lowering tumor cell purity in bulk specimens, while being a biomarker for better outcome and response to treatment\u003csup\u003e29\u003c/sup\u003e. Immune infiltration in TNBC has so far been associated with features like low clonal heterogeneity, mutation load, copy number alterations\u003csup\u003e30\u003c/sup\u003e and even epigenetic subtypes\u003csup\u003e21\u003c/sup\u003e. Despite intense research, we are still far from a mechanistic understanding of the factors driving heterogeneity in the tumor immune microenvironment (TIME), both across cancers in general and in TNBC specifically, and what role somatic DNA methylation alterations plays in the TIME. With the increasing usage of immunotherapy also for early-stage patients, understanding the TIME and how it is shaped becomes increasingly important. In this setting, TNBC provides an intriguing model system for deconstructing bulk tumor DNA methylation data in order to identify novel epigenetic features associated with both tumor-intrinsic and -extrinsic factors.\u003c/p\u003e\n\u003cp\u003eIn this work we demonstrate how statistical deconstruction of bulk tumor DNA methylation data can identify both biologically relevant CpG sets that mirror both gene specific alterations as well as patterns specific to molecular subtypes, as well as variance (noise) without apparent biological relevance for the disease in question that may obscure data analyses. Moreover, by integrating additional -omics layers combined with \u003cem\u003ein situ\u003c/em\u003e methods like spatial transcriptomics we identify putative epigenetically driven tumor alterations that correlate with an altered TIME, providing support for a concept of tumor-intrinsic epigenetic immune evasion.\u0026nbsp;\u003c/p\u003e"},{"header":"MATERIALS AND METHODS","content":"\u003cp\u003e\u003cstrong\u003eCohorts and Datasets\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eInclusion and ethics statement.\u0026nbsp;\u003c/strong\u003ePatients were enrolled in the Sweden Cancerome Analysis Network – Breast (SCAN-B) study (ClinicalTrials.gov ID NCT02306096)\u003csup\u003e31,32\u003c/sup\u003e approved by the Regional Ethical Review Board in Lund, Sweden (registration numbers 2009/658, 2010/383, 2012/58, 2013/459, 2014/521, 2015/277, 2016/541, 2016/742, 2016/944, 2018/267) and the Swedish Ethical Review Authority (registration numbers 2019-01252, 2024-02040-02). All patients provided written informed consent prior to enrolment, including to publish information about sex and age. All analyses were performed in accordance with patient consent and ethical regulations and decisions.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eSCAN-B cohorts.\u003c/strong\u003e This work is based on two TNBC cohorts, used for discovery and validation purposes. Both cohorts belong to the Swedish SCAN-B study, and all analyzed tissue specimens were taken prior to any treatment. In Sweden, the definition of TNBC is a tumor with \u0026lt;10% of cells with immunohistochemistry (IHC) staining for ER and PR (thus including tumors with 1–9% stained cells) and an IHC HER2-staining score of \u0026lt;2, or for patients with IHC 2+, a non-amplified \u003cem\u003ein situ\u003c/em\u003e hybridization status. All SCAN-B data for ER, PR, and HER2 status were obtained from clinical routine analyses performed by regional pathology departments and provided through the national Swedish breast cancer quality registry. The discovery cohort consists of clinical, pathological and whole genome sequencing (WGS) data from a population representative set of 235 TNBC patients originally reported by Staaf et al.\u003csup\u003e27\u003c/sup\u003e. Molecular data layers included DNA methylation data based on Illumina EPIC bead chips obtained from Aine et al.\u003csup\u003e21\u003c/sup\u003e and RNA-sequencing data summarized as fragments per kilobase million (FPKM) counts obtained from Staaf et al\u003csup\u003e33\u003c/sup\u003e. Different TNBC sample classifications such as PAM50 subtypes, the four refined Lehmann TNBC subtypes (TNBCtype: BL1, BL2, M, and LAR), DNA methylation epitypes (Basal1, Basal2, Basal3, nonBasal1, and nonBasal2), HRD status (HRD or HR proficient: HRP), and rank-based gene expression scores from eight biological metagene genes (including a stroma and immune response metagene\u003csup\u003e34\u003c/sup\u003e) were obtained from Aine et al.\u003csup\u003e21\u003c/sup\u003e. Additionally, pathology estimates of tumor infiltrating lymphocytes (TILs) and PD-L1 combined positive score (PDL1-CPS) were collected from Aine et al.\u003csup\u003e35\u003c/sup\u003e and Sigurjonsdottir et al.\u003csup\u003e36\u003c/sup\u003e, respectively.\u0026nbsp;CPS scores were categorized into three groups as CPS 0 (CPS score \u0026lt;1), CPS 1 (CPS score 1-9), and CPS 2 (CPS score ≥10). Digital cell counts of relevant immune cell type markers such as CD20, CD3, CD8, FOXP3, CD4, CD68 and PD-L1 were obtained from Roostee et al.\u003csup\u003e29\u003c/sup\u003e calculated with the TMArQ software. A detailed description of the inclusion and exclusion criteria of patients in the discovery cohort is available in the original publication\u003csup\u003e27\u003c/sup\u003e. The validation cohort consists of molecular and clinical data from 136 TNBC samples reported by Aine et al.\u003csup\u003e21\u003c/sup\u003e. Molecular data layers included are DNA methylation data based on Illumina EPIC bead chips obtained from Aine et al.\u003csup\u003e21\u003c/sup\u003e and RNA-sequencing data summarized as FPKM counts obtained from Staaf et al.\u003csup\u003e33\u003c/sup\u003e. Different TNBC sample classifications such as PAM50 subtypes, Lehmann TNBC subtypes and DNA methylation epitypes were obtained from Aine et al.\u003csup\u003e21\u003c/sup\u003e. The composition and characteristics of the discovery and validation cohorts are detailed in Table 1.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eTable\u0026nbsp;\u003c/strong\u003e\u003cstrong\u003e1\u003c/strong\u003e\u003cstrong\u003e. Composition and characteristics of the discovery and validation cohorts\u003c/strong\u003e\u003cstrong\u003e.\u003c/strong\u003e Age and TILs are reported as median value and interquartile range (IQR) per cohort. ND stands for unavailable data.\u003c/p\u003e\n\u003ctable border=\"1\" cellspacing=\"0\" cellpadding=\"0\"\u003e\n \u003ctbody\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e\u003cstrong\u003eDISCOVERY\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e\u003cstrong\u003eVALIDATION\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e\u003cstrong\u003eN\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e235\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e136\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e\u003cstrong\u003eAge\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e62 (IQR: 51-72)\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e62.5 (IQR: 50-75)\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e\u003cstrong\u003eTILs\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e20 (IQR: 10-40)\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eND\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd colspan=\"3\" valign=\"top\"\u003e\n \u003cp\u003e\u003cstrong\u003eGrade\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e2\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e26\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e24\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e3\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e206\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e81\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eNot available\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e3\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e31\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd colspan=\"3\" valign=\"top\"\u003e\n \u003cp\u003e\u003cstrong\u003eLymph node status\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eN0\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e152\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e87\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eN+\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e83\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e39\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eNot available\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e0\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e10\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd colspan=\"3\" valign=\"top\"\u003e\n \u003cp\u003e\u003cstrong\u003ePAM50\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eBasal-like\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e187\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e95\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eNormal-like\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e7\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e13\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eHER2-enriched\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e34\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e27\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eLuminal A\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e5\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e10\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eLuminal B\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e1\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e0\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eUnclassified\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e1\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e0\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd colspan=\"3\" valign=\"top\"\u003e\n \u003cp\u003e\u003cstrong\u003eLehmann TNBCtype\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eLAR\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e40\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e32\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eM\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e54\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e29\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eBL1\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e91\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e29\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eBL2\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e47\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e20\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eUnspecified\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e0\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e20\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eNot available\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e3\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e6\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd colspan=\"3\" valign=\"top\"\u003e\n \u003cp\u003e\u003cstrong\u003eHRD status\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eHRD\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e139\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eND\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eHRP\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e96\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eND\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd colspan=\"3\" valign=\"top\"\u003e\n \u003cp\u003e\u003cstrong\u003eEpitype\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eBasal1\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e53\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e23\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eBasal2\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e68\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e37\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eBasal3\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e54\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e29\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eNonBasal1\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e29\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e22\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003eNonBasal2\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e31\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\"\u003e\n \u003cp\u003e25\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003c/tbody\u003e\n\u003c/table\u003e\n\u003cp\u003eProcessed and filtered tumor purity-adjusted DNA methylation data\u003csup\u003e16,17\u003c/sup\u003e\u0026nbsp; in the form of beta values (range 0-1, representing hypomethylated to hypermethylated) of both cohorts, comprising in total 741145 CpGs in the discovery cohort and 717458 CpGs in the validation cohort, were obtained from Aine et al.\u003csup\u003e21\u003c/sup\u003e. Data regarding CpG characteristics was compiled from different sources. CpG context, chromosome accessibility status based on ATAC-seq, and overlap with repetitive sequences and transcription factor binding sites were obtained from Aine et al.\u003csup\u003e21\u003c/sup\u003e. Overlap of CpGs with enhancer regions was assessed from enhancers annotated in the GeneHancer genomic regulatory element database (v5.24)\u003csup\u003e37\u003c/sup\u003e. CpG density was determined based on the genomic location of each CpG using windows of 5000 base pairs centered on the CpG and calculated as the number of CpGs in the window divided by total base pairs.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eTCGA-BRCA external validation dataset.\u0026nbsp;\u003c/strong\u003eAn external validation cohort of 645 breast cancers, representing all molecular and clinical breast cancer subtypes, from The Cancer Genome Atlas (TCGA) with DNA methylation data (Illumina Infinium 450K) and RNA-sequencing data summarized as FPKM was obtained from Aine et al.\u003csup\u003e21\u003c/sup\u003e. Pathologist estimated TILs values for the samples in the cohort were obtained from Cha et al.\u003csup\u003e38\u003c/sup\u003e.\u0026nbsp;\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eNormal breast tissue datasets.\u003c/strong\u003e Molecular data from normal breast tissue were obtained from three different sources. Illumina Infinium 450K DNA methylation data from 96 normal breast tissue samples were obtained from Hair et al.\u003csup\u003e39\u003c/sup\u003e deposited in Gene Expression Omnibus (GEO) under accession number GSE67919. Illumina Infinium 450K DNA methylation profiles for 10 different flow cytometry sorted blood cell types/fractions were obtained from Reinius et al.\u003csup\u003e40\u003c/sup\u003e deposited in GEO under accession number GSE35069. Single cell RNA-sequencing (scRNA-seq) data from normal breast tissue together with cell type annotations were obtained from Reed et al.\u003csup\u003e41\u003c/sup\u003e deposited in the CellXGene platform.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eTNBC cell line datasets.\u003c/strong\u003e Bulk RNA-sequencing data from 34 TNBC cell lines were obtained from Jovanovic et al.\u003csup\u003e42\u003c/sup\u003e summarized as FPKM data deposited under the GEO accession number GSE202770. Eight TNBC cell lines with available RNA-sequencing data, summarized as FPKM values and Illumina EPIC DNA methylation data were obtained from Aine et al.\u003csup\u003e21\u003c/sup\u003e. Additionally, 19 cell lines with available scRNA-seq data were obtained from Jovanovic et al.\u003csup\u003e42\u003c/sup\u003e deposited under the GEO accession number GSE202771 and processed as described by Aine et al\u003csup\u003e21\u003c/sup\u003e.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eTNBC scRNA-seq dataset.\u003c/strong\u003e Processed scRNA-seq data together with cell type annotations from eight TNBC tumor samples were collected from Chen et al.\u003csup\u003e43\u003c/sup\u003e deposited under the GEO accession number GSE161529.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eSpatial transcriptomics TNBC dataset.\u003c/strong\u003e Single cell spatial transcriptomics data from 65 tumors from the discovery cohort was used to assess \u003cem\u003eOAS2\u0026nbsp;\u003c/em\u003eexpression \u003cem\u003ein situ\u003c/em\u003e. Spatial transcriptomics analysis was performed using the nanoString/Bruker CosMx instruments on two tissue microarray (TMA) sections, generally containing 2 cores per tumor (114 cores in total), using the Universal Gene Characterisation panel (1000 transcripts). Briefly, TMA sections 4μm thick were cut and mounted on VWR Superfrost Plus slides (company), without coverslips. After sectioning, slides were dried at 37°C overnight at an angle no greater than 45 degrees. Slides were stored at 4°C with a desiccant bag until analysis. Summarized clinicopathological and molecular characteristics for the 65 tumors are available in Supplementary Table 1. CosMx analysis was performed following the manufacturer’s instructions at the Spatial Biology Facility, Faculty of Life Sciences and Medicine, King’s College, London, UK. The image files were processed using nanoString/Bruker’s AtoMx platform and the in-house pipeline and quality assurance set up at the Spatial Biology Facility.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eIdentification and analysis of CpG cassettes\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCpG cassette identification.\u003c/strong\u003e CpG cassettes were identified using Weighted Gene-Correlation Network Analysis (WGCNA) through the WGCNA R package (v1.73)\u003csup\u003e44\u003c/sup\u003e. The blockwiseModules() function was applied using \u003cem\u003ebicor\u003c/em\u003e as the correlation measure and the following parameters: networkType=\"unsigned\", minModuleSize=3, maxBlockSize=6000, reassignThreshold=0 and mergeCutHeight=0.25. Different soft-thresholding beta values (β), the parameter that determined the stringency of the clustering approach, were used to run the analyses (5, 8, 10, 15, 20 and 25) and evaluated using the pickSoftThreshold() function. The optimal soft-threshold power was set to 10 using tumor purity-adjusted and not context stratified methylation data and was kept across runs to maximize consistency. The remaining parameters were set to the default values. The cassette detection was run on both tumor purity-adjusted and unadjusted CpG beta values, without CpG context stratification, using a variance threshold of 0.1 due to computational limitations. This approach was also applied on tumor purity-adjusted and unadjusted CpG beta values after CpG context stratification, using 0.05 as the variance threshold for promoter and proximal cassettes and 0.1 for distal cassettes. A higher value was used for distal CpGs due to computational limitations given the higher number of such CpGs. Additionally, the same approach and the same parameters were used to identify CpG cassettes from distal CpGs filtered based on chromatin accessibility, defined as CpG overlap with breast cancer ATAC-seq peaks from the study by Corces et al.\u003csup\u003e45\u003c/sup\u003e (see Aine et al.\u003csup\u003e21\u003c/sup\u003e for details). For each cassette, the first principal component (PC1) of the DNA methylation data was computed as a summary measure. The sign of PC1 was adjusted to ensure it aligned with the original data, so that higher values consistently reflected increased methylation beta values, facilitating interpretability.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCell type deconvolution from DNA methylation data.\u003c/strong\u003e To analyze the impact of tumor purity adjustment on the methylation data and to obtain estimates of the cell composition of each sample the MethylCIBERSORT R package (v0.2.1)\u003csup\u003e46\u003c/sup\u003e was used to generate mixture files from adjusted and unadjusted methylation data using the Prep.CancerType() function and the reference breast cancer signatures (breast_v2). The obtained data were used to run deconvolution through the CIBERSORTx online tool\u003csup\u003e47\u003c/sup\u003e, using 1000 permutations and disabling quantile normalization.\u0026nbsp;\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCpG cassette stability analysis.\u0026nbsp;\u003c/strong\u003eThe stability of the detected CpG cassettes was analyzed on tumor purity-adjusted and CpG context-stratified data. After filtering CpGs based on variance as detailed previously, 20 groups of 100 samples were randomly selected from the discovery cohort, using the R sample() function with probabilities reflecting the original distribution of the DNA methylation epitypes to avoid class imbalance. Cassette detection was then performed as previously described.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAssociation of CpG cassettes with the Basal/non-Basal division.\u0026nbsp;\u003c/strong\u003eThe predictive capacity of each major CpG cassette in calling the PAM50 Basal/non-Basal division was determined through the Adjusted Rand Index (ARI) using the adjustedRandIndex() function from the mclust (v6.1.1) R package. Specifically, for each cassette the PAM50 Basal/non-Basal tumor subtype classification was compared to the resulting tow sample clusters obtained from hierarchical clustering (using Euclidean distance and Ward.D2 linkage) of CpG methylation data setting the number of clusters to two. Additionally, cassettes with ARI\u0026gt;0.65 were used to stratify samples in Basal/non-Basal groups using hierarchical clustering (Euclidean distance, Ward.D2 linkage) to evaluate their predictive capacity.\u0026nbsp;\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eGene overrepresentation analyses linked to CpG cassettes.\u0026nbsp;\u003c/strong\u003eGene overrepresentation analysis of the CpGs included in the analyzed cassettes was performed using the missMethyl R package (v1.38.0)\u003csup\u003e48\u003c/sup\u003e, which adjusts for the varying numbers of CpGs per gene. The gometh() function was used to perform Gene Ontology enrichment with the following parameters: collection=\"GO\", array.type=\"EPIC\", prior.prob=TRUE, equiv.cpg=TRUE, fract.counts=TRUE, and genomic.features=\"ALL\". The same function and parameters were used to perform enrichment of KEGG (Kyoto Encyclopedia of Genes and Genomes) pathways but setting collection=\"KEGG\". Cancer Hallmark enrichment was determined using the gsameth() function with the same parameters as above but using human hallmarks signatures (v7.1) obtained from the Molecular Signatures Database (MSigDB) \u003csup\u003e49\u003c/sup\u003e as the collection parameter.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAnalyses of gene expression changes linked to CpG cassettes.\u0026nbsp;\u003c/strong\u003eThe association of CpG cassettes with gene expression differences was assessed by linking CpGs in individual cassettes to genes using provided UCSC annotations. For each gene, tumors were divided into two groups based on\u0026nbsp;k-means clustering\u0026nbsp;of CpGs connected to the gene through the kmeans() R function, using centers=2 and the default parameters, and\u0026nbsp;mRNA expression differences between the two methylation groups were analyzed using Wilcoxon’s test.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAssessment of proportion of CpG cassettes in N most variable CpGs.\u0026nbsp;\u003c/strong\u003eThe top 1000, 5000, 10000, 20000 and 30000 most variable CpGs were determined in tumor purity adjusted and unadjusted DNA methylation data, respectively. The proportion of each detected cassette within each CpG context and CpG set was calculated by identifying the CpGs shared between the cassette and the top N groups.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eDifferential methylation analysis.\u0026nbsp;\u003c/strong\u003eAdjusted DNA methylation data from the discovery cohort was used to identify differentially methylated CpGs per molecular subtype between PAM50 Basal/non-Basal and the Lehmann TNBC subtype division. Differential methylation was tested pairwise between each molecular subtype using Wilcoxon’s test followed by Bonferroni multiple testing correction. Finally the overlap of significant CpGs with identified distal CpG cassettes was analyzed.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eIdentification and analyses of CpG cassettes associated with immune infiltration\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCpG cassette identification.\u0026nbsp;\u003c/strong\u003eTo identify CpG cassettes associated with TILs we performed an association analysis between cassettes and TIL proportions. First, the association was tested using Kendall’s correlation value between TILs and PC1 of the analyzed cassettes. Additionally, aiming to capture the two potential methylation states of CpGs, we defined two methylation groups per cassette by applying k-means clustering on their PC1, setting the number of clusters to 2. Next, we compared TIL proportions between these groups using Wilcoxon’s test. To prioritize cassettes relevant for TILs we used the Boruta algorithm, using PC1 of each CpG cassette as summary value. Feature selection was performed using the Boruta() function from the Boruta R package (v8.0.0)\u003csup\u003e50\u003c/sup\u003e, setting doTrace=2, maxRuns=500, and remaining parameters to default.\u0026nbsp;\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003escRNA-seq analyses.\u0026nbsp;\u003c/strong\u003eAnalyses of scRNA-seq data from tumor and normal breast samples were performed using the Seurat R package (v5.3.0)\u003csup\u003e51\u003c/sup\u003e together with cell type annotations from the original publications. The DotPlot() function with default parameters was used for gene expression analyses.\u0026nbsp;\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eSpatial transcriptomics analyses.\u0026nbsp;\u003c/strong\u003eThe CosMx spatial transcriptomics data were analyzed using the python Squidpy (v1.6.5)\u003csup\u003e52\u003c/sup\u003e and Scanpy (v1.10.4)\u003csup\u003e53\u003c/sup\u003e libraries. Cells with a proportion of negative probabilities \u0026gt;4%, \u0026lt;30 transcripts detected, and \u0026lt;30 different genes detected were removed. Transcript counts were then normalized using the normalize_total() function and setting target_sum=1e4 followed by log1p transformation through the log1p() function. Tumor cells were identified per TMA slide through Gaussian Mixture Modelling applied on the mean panCK intensity detected per cell using the scikit-learn python library (v1.6.1)\u003csup\u003e54\u003c/sup\u003e.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eStatistical analyses, data handling and plotting\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eData handling, statistical analyses and plotting were run in R (v4.4.1). Data handling was performed through base R functions, and the dplyr (v1.1.4), tidyr (v1.3.1) and reshape2 (v1.4.4) packages. Pairwise comparisons were performed using two-sided Wilcoxon’s test. Comparisons between three or more groups were performed using the Kruskal-Wallis test. The Chi-squared statistical test was used for comparisons between categorical variables. The binomial test was used to compare observed and expected events. False Discovery Rate (FDR) or Bonferroni multiple testing correction was implemented when necessary to control for the expected proportions of false positives using the p.adjust() function setting the method to “fdr” or “bonferroni”, respectively. Asterisks indicate p-value significance as follows: * when p≤0.05, ** when p≤0.01, *** when p≤0.001, **** when p≤0.0001, and \u003cem\u003ens\u003c/em\u003e when p\u0026gt;0.05. \u0026nbsp;Correlation values were determined using the Spearman and Kendall methods. Plotting was performed using base R functions, ggplot2 (v3.5.2), ComplexHeatmap (v2.20.2), ggsankey (v0.0.99999), corrplot (v0.95), VennDiagram (v1.6.0) and other R packages’ specific plotting functions. Column clustering of every heatmap was performed using Euclidean distance and the Ward.D2 method while the rows were kept unclustered if not specified otherwise.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eData availability statement\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eFor the discovery cohort, used WGS data originally reported by Staaf et al.\u003csup\u003e27\u003c/sup\u003e\u0026nbsp; is available from [https://data.mendeley.com/datasets/2mn4ctdpxp/3], while used RNA-sequencing data originally reported by Staaf et al.\u003csup\u003e33\u003c/sup\u003e is available from [https://data.mendeley.com/datasets/yzxtxn4nmd/3]. DNA methylation data used for the discovery cohort was originally reported by Aine et al.\u003csup\u003e21\u003c/sup\u003e and is deposited in the Gene Expression Omnibus under the GSE148748 and GSE148906 accession numbers. DNA methylation data used for the SCAN-B validation cohort was previously reported in Aine et al.\u003csup\u003e21\u003c/sup\u003e and it is deposited in the Gene Expression Omnibus (GEO) under accession number GSE290981. RNA-sequencing data from this cohort is available from [https://data.mendeley.com/datasets/yzxtxn4nmd/3]. TILs and PDL1-CPS estimates for the discovery cohort from Aine et al.\u003csup\u003e35\u003c/sup\u003e and Sigurjonsdottir et al.\u003csup\u003e36\u003c/sup\u003e,\u0026nbsp;respectively, are available in the original publications. TMArQ determined immune cell counts per sample in the discovery cohort by Roostee et al.\u003csup\u003e29\u003c/sup\u003e are available from the original publication. Sample annotations such as PAM50 subtypes, the four refined Lehmann TNBC subtypes, DNA methylation epitypes, HRD status, and rank-based gene expression scores for both the discovery and validation cohorts were obtained from the original publication Aine et al.\u003csup\u003e21\u003c/sup\u003e Raw RNA-sequencing data used in this study and generated by Aine et al.\u003csup\u003e21\u003c/sup\u003e for eight TNBC cell lines are available through the Short Read Archive (SRA) archive at NCBI under BioProject accession PRJNA1189708 and study SRP547133. Corresponding DNA methylation data from Aine et al.\u003csup\u003e21\u003c/sup\u003e for the eight TNBC cell lines is deposited in the GEO database (accession number GSE282347). The previously reported breast cancer TCGA data used in this study is available from the GDC data portal [https://portal.gdc.cancer.gov]. TILs estimates for the TCGA external validation cohort from Cha et al.\u003csup\u003e38\u003c/sup\u003e are available from their original publication. The previously reported normal breast tissue DNA methylation data from Hair et al.\u003csup\u003e39\u003c/sup\u003e used in this study is deposited in the Gene Expression Omnibus under the accession number GSE67919. The previously reported scRNA-seq data of normal breast tissue from Reed et al.\u003csup\u003e41\u003c/sup\u003e is deposited in the CellXGene platform and available from [https://cellxgene.cziscience.com/collections/48259aa8-f168-4bf5-b797-af8e88da6637]. The previously reported RNA-sequencing data for TNBC cell lines from Jovanovic et al.\u003csup\u003e42\u003c/sup\u003e used in this study are deposited in GEO (accession number GSE202770). The previously reported TNBC scRNA-seq data from Chen et al.\u003csup\u003e43\u003c/sup\u003e used in this study are deposited in GEO (accession number GSE161529). The previously reported sorted immune cell DNA methylation data from Reinius et al. used in this study are deposited in GEO (accession number GSE35069). The gene enhancer dataset was retrieved from the GeneHancer genomic regulatory element database (v5.24)\u003csup\u003e37\u003c/sup\u003e and it is available from [https://www.genecards.org/Guide/GeneHancer]. The CosMx spatial transcriptomics data required for the analyses performed in this work is available as \u003cem\u003esupplementary_file_5.csv\u003c/em\u003e.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCode availability\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThe code used for the analyses implemented in this paper is available in the following GitHub repository: https://github.com/StaafLab/methylation_dynamics_tnbc\u003c/p\u003e"},{"header":"RESULTS","content":"\u003cp\u003e\u003cstrong\u003eCpG cassettes in tumor purity-adjusted TNBC DNA methylation data\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eTo detect CpG cassettes (groups of covarying CpGs) in TNBC we applied Weighted Gene Co-expression Network Analysis (WGCNA) on methylation data from 741145 CpGs in the 235 discovery cohort tumors both before and after adjusting beta values for tumor purity, and without initially stratifying into CpG context. WGCNA clusters CpGs based on correlation values weighted with a soft-thresholding parameter (\u0026beta;) controlling the stringency of the clustering. Higher \u0026beta; values led to more homogeneous cassettes, i.e., smaller and more tightly co-methylated CpG sets that comprised a smaller fraction of all CpGs available (Figure 1A).\u0026nbsp;We selected cassettes defined with \u0026beta;=10 for further analyses to balance cassette homogeneity and the fraction of clustered CpGs, as well as because this value produced a network with approximate scale-free topology (Supp. Figure 1). In a scale-free network, most nodes (CpGs) have few connections, while a few nodes serve as highly connected hubs, a property commonly observed in biological networks and considered desirable for WGCNA. To assess the impact of tumor purity adjustment on the CpG cassette detection we compared WGCNA results for adjusted versus unadjusted methylation data. WGCNA on adjusted CpG data identified more distinct inter-cassette patterns and greater consistency across samples compared to unadjusted data (Figure 1B). This is presumably due to the removal of the effect of the cellular admixture from the tumor sample, as purity adjusted data showed, as expected, higher tumor fractions based on MethylCIBERSORT analysis (Figure 1C). Given our focus on tumor-intrinsic patterns, tumor purity adjusted methylation data were used for all subsequent analyses. \u0026nbsp;\u003c/p\u003e\n\u003cp\u003eWe next focused on the impact of CpG context (promoter, proximal, and distal) on cassette detection in the discovery cohort. CpG context is not inherently captured by variation in methylation levels across samples, nor can it be addressed by analyzing only the N most variable CpGs or by correcting for tumor purity (Figure 1D). However, it plays a substantial role in methylation dynamics given its functional role in gene regulation\u003csup\u003e19,55\u003c/sup\u003e. Ignoring CpG context leads to identification of CpG cassettes with substantial differences in CpG context proportions (Figure 1E). For instance, the detected cassette 1 is almost exclusively comprised of distal CpGs, while cassette 2 shows high enrichment of promoter and proximal CpGs, which indicates that some of the detected methylation patterns are CpG context specific. Consistently, performing WGCNA-based cassette detection within each CpG context, respectively, identified cassettes with clearly different methylation patterns (Figure 1F). An important aspect of our CpG cassette detection approach is reproducibility and robustness. To evaluate the stability of the identified CpG context-specific cassettes in a smaller data set, we randomly subsampled the discovery cohort data down to 100 samples (while maintaining the proposed DNA methylation epitype distribution\u003csup\u003e21\u003c/sup\u003e) and performed a new CpG context cassette detection. This was performed 20 times. Comparison of results to the original CpG context cassettes showed that CpGs from the original distal cassettes remained largely grouped together in subsampling cassettes with similar methylation patterns as the original ones (Figure 1G), and a similar observation was made for proximal and promoter CpGs (Supp. Figure 2). Finally, the DNA methylation patterns of the CpG context-specific cassettes in the discovery cohort were consistent in the SCAN-B validation cohort (Supp. Figure 3). Together, these results support the robustness and reproducibility of the cassette approach. The identified CpG cassettes for all, promoter, proximal and distal CpGs obtained from adjusted beta values using \u0026beta;=10 are available as \u003cem\u003esupplementary_file_1.csv\u003c/em\u003e, \u003cem\u003esupplementary_file_2.csv\u003c/em\u003e, \u003cem\u003esupplementary_file_3.csv\u003c/em\u003e and \u003cem\u003esupplementary_file_4.csv,\u0026nbsp;\u003c/em\u003erespectively.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCpG cassettes reflect major biological classes in TNBC\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eA distinct feature of the clustered cassette heatmaps in Figure 1 is the aggregation of tumors into mainly a PAM50 Basal and a non-Basal branch (e.g., Figure 1F). To examine this observation in more detail we analyzed the association of detected CpG cassettes in the three CpG contexts with major biological divisions in TNBC, in particular the Basal/non-Basal division proposed to represent the two main global methylation states in TNBC\u003csup\u003e21\u003c/sup\u003e. Among the seven main cassettes identified within each CpG context, several, but not all, captured the two primary methylation states in TNBC. This was quantified using the Adjusted Rand Index (ARI) between sample clusters derived from each cassette and the PAM50 Basal/non-Basal classification (Figure 2A). Specifically, five cassettes (promoter cassette 1, n=3015 CpGs; proximal cassette 1, n=3446 CpGs; proximal cassette 3, n=562 CpGs; distal cassette 2, n=2021 CpGs; distal cassette 3, n=812 CpGs) exhibited ARI values greater than 0.65, indicating strong concordance with the Basal/non-Basal split. Hierarchical clustering of the CpGs included in high-ARI cassettes successfully differentiated Basal from non-Basal samples (Figure 2B).\u003c/p\u003e\n\u003cp\u003eMoreover, CpGs included in the detected cassettes were linked to biologically relevant processes based on pathway enrichment analysis, which was restricted to promoter and proximal cassettes due to available gene mappings. The largest promoter CpG cassette linked to the Basal/non-Basal division, promoter cassette 1, showed enrichment for the Epithelial-Mesenchymal Transition (EMT) hallmark and several KEGG pathways related to cancer processes (Figure 2C). When dividing samples into two groups based on DNA methylation values of CpGs included in this cassette, changes in gene expression correlated with methylation values were observed between the groups, including for genes such as \u003cem\u003eLDHB\u003c/em\u003e, \u003cem\u003eMUCL1\u003c/em\u003e, \u003cem\u003eANP32E\u003c/em\u003e and \u003cem\u003eTFAP2B\u003c/em\u003e (Figure 2D). Similar findings were observed when analyzing the main proximal cassette linked to the Basal/non-Basal division, proximal cassette 1. This cassette showed enrichment for hallmarks such as EMT, KRAS signaling and Pancreatic beta cells, and for KEGG pathways such as Breast cancer, Basal cell carcinoma, various signaling pathways and glycosphingolipid biosynthesis (Figure 2E). This cassette was also associated with expression changes in genes such as \u003cem\u003eSPINK8\u003c/em\u003e, \u003cem\u003eANP32E\u003c/em\u003e and \u003cem\u003ePAX2\u0026nbsp;\u003c/em\u003e(Figure 2F). As more detailed examples, \u003cem\u003eLDHB\u003c/em\u003e and \u003cem\u003eSPINK8\u003c/em\u003e were selected to show the expression changes linked to CpG cassettes (Figure 2G). \u003cem\u003eLDHB\u003c/em\u003e shows low mRNA expression in non-Basal TNBC consistent with the hypermethylation of involved promoter CpGs, and the opposite pattern is seen in Basal tumors. In contrast, \u003cem\u003eSPINK8\u003c/em\u003e expression is notably correlated with the DNA methylation status of a proximal CpG, with mRNA expression only in non-Basal tumors.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eChromatin accessibility distinguishes biologically relevant variation in distal CpGs\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eAs shown in Figure 2A, two major distal cassettes, cassettes 2 and 3, showed a clear correlation with the Basal/non-Basal PAM50 division of the discovery cohort. In contrast, distal cassette 1 (n=6743 CpGs) and other smaller distal and non-distal cassettes showed no association to these PAM50 groups. A notable feature of distal cassette 1 is a high-variance, gradient-like DNA methylation pattern across tumors (Figure 3A). Moreover, distal cassette 1 showed low overlap with enhancer regions (~20%), similar to other cassettes not linked to the Basal/non-Basal division like distal cassettes 4, 5 and 7, and a relatively high fraction of CpGs overlapping with repetitive sequences (~27%) (Figure 3B). The overlap of cassette CpGs with transcription factor binding sites (TFBSs) also showed distinct patterns: cassettes 1, 4, 5 and 7 showed minimal overlap with TFBSs, mirroring the trend seen in enhancer regions (Figure 3C). Moreover, CpGs in distal cassette 1 were often located in generally sparse genomic regions (Figure 3D). Further analyses on potential biological implications of distal cassette 1, with DNA methylation values summarized through the first principal component (PC1), showed no clear difference regarding important tumor characteristics in TNBC such as the PAM50 Basal/non-Basal division, HRD status, TNBC mRNA subtypes, tumor purity, or TIL levels (Figure 3E). Concordantly, when selecting CpGs through supervised methods, such as differential methylation analysis, CpGs included in the noise cassettes were rarely selected (Supp. Figure 4). Taken together, both genomic and biological features suggest that distal cassette 1 represents a potential source of high variance, largely non-informative noise rather than biologically relevant CpGs in the context of TNBC.\u003c/p\u003e\n\u003cp\u003eImportantly, distal cassette 1 could not be removed from the data using variance-based filtering (Figure 3F), even without applying tumor purity correction (Supp. Figure 5). In further support of the non-informativeness of a set of high variance distal CpGs, we found similar results to what we described above also after applying a lower soft-thresholding value that minimizes the fraction of unclustered CpGs (Supp. Figure 6). Again, the largest detected distal cassette showed low overlap with enhancer regions, lower CpG density, reduced overlap with TFBSs and could not be excluded by variance-based filtering on adjusted or unadjusted beta values (Supp. Figure 7). Consequently, if not addressed less informative, high variance cassettes like distal cassette 1 may influence unsupervised DNA methylation analyses. However, applying an alternative filtering strategy based on chromatin accessibility could successfully remove noise-associated cassettes (Figure 3G). This approach was based on overlap of CpGs with ATAC-seq peaks (indicating regions of open, accessible chromatin) from Corces et al.\u003csup\u003e45\u003c/sup\u003e that have a distinctive distribution per CpG context (Table 2). Specifically, the original distal cassette 1, as well as other cassettes with low overlap with TFBSs and enhancer regions, were removed, while approximately 50% of CpGs in cassettes 2 and 3 were retained (Figure 3H). Interestingly, the new cassettes derived from ATAC-seq filtered CpGs had clear counterparts among the original cassettes enriched for TFBSs and enhancer regions (Figure 3I). A similar result was observed when using a lower soft-thresholding value (Supp. Figure 8). To explore whether the identified distal cassette patterns were only a feature of TNBC, representing a limited proportion of all breast cancers, we analyzed the methylation patterns of the seven largest distal cassettes in 645 TCGA breast cancers of all molecular and clinical subtypes. Importantly, a similar apparent noise pattern of distal cassette 1,4, 5 and 7 was observed in all PAM50 subtypes, whereas distal cassette 2, 3 and 6 strongly aligned with a Basal/non-Basal division of the cohort, irrespective of ER, PR, and HER2 status (Supp. Figure 9). Together, this highlights the biological relevance and consistency of the identified patterns also in a broader breast cancer setting.\u0026nbsp;\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eTable 2. Distribution of CpGs from the discovery cohort overlapping with ATAC-seq peaks per CpG context.\u0026nbsp;\u003c/strong\u003eATAC + overlap with Corces et al ATAC-seq peaks, and ATAC \u0026ndash; have no overlap\u003c/p\u003e\n\u003ctable border=\"1\" cellspacing=\"0\" cellpadding=\"0\"\u003e\n \u003ctbody\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\" style=\"width: 66px;\"\u003e\n \u003cp\u003e\u0026nbsp;\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\" style=\"width: 76px;\"\u003e\n \u003cp\u003e\u003cstrong\u003eAll\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\" style=\"width: 94px;\"\u003e\n \u003cp\u003e\u003cstrong\u003ePromoter\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\" style=\"width: 95px;\"\u003e\n \u003cp\u003e\u003cstrong\u003eProximal\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\" style=\"width: 85px;\"\u003e\n \u003cp\u003e\u003cstrong\u003eDistal\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\" style=\"width: 66px;\"\u003e\n \u003cp\u003e\u003cstrong\u003eATAC -\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\" style=\"width: 76px;\"\u003e\n \u003cp\u003e205799\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\" style=\"width: 94px;\"\u003e\n \u003cp\u003e44722\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\" style=\"width: 95px;\"\u003e\n \u003cp\u003e96285\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\" style=\"width: 85px;\"\u003e\n \u003cp\u003e394339\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003ctr\u003e\n \u003ctd valign=\"top\" style=\"width: 66px;\"\u003e\n \u003cp\u003e\u003cstrong\u003eATAC +\u003c/strong\u003e\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\" style=\"width: 76px;\"\u003e\n \u003cp\u003e535346\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\" style=\"width: 94px;\"\u003e\n \u003cp\u003e80567\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\" style=\"width: 95px;\"\u003e\n \u003cp\u003e44531\u003c/p\u003e\n \u003c/td\u003e\n \u003ctd valign=\"top\" style=\"width: 85px;\"\u003e\n \u003cp\u003e80701\u003c/p\u003e\n \u003c/td\u003e\n \u003c/tr\u003e\n \u003c/tbody\u003e\n\u003c/table\u003e\n\u003cp\u003e\u003cstrong\u003eSpecific promoter CpG cassettes are linked to TIME status in TNBC\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eAfter analyzing the major cassettes reflecting broad biological patterns in TNBC, we focused on smaller cassettes that could be involved in the molecular heterogeneity seen within TNBC. We restricted the analysis to promoter CpG cassettes to facilitate a more direct biological interpretation. As a validation step, we first investigated whether known promoter methylation patterns such as \u003cem\u003eBRCA1\u003c/em\u003e promoter hypermethylation were identified by our cassette approach. Indeed, we identified a cassette (promoter cassette 11, n=19 CpGs) including only CpGs associated with \u003cem\u003eBRCA1\u003c/em\u003e that on a sample level showed excellent agreement with both HRD phenotype and pyrosequencing-based \u003cem\u003eBRCA1\u003c/em\u003e methylation status (representing an orthogonal methylation analysis method), consistent with previous studies\u003csup\u003e3,27\u003c/sup\u003e (Figure 4A).\u0026nbsp;\u003c/p\u003e\n\u003cp\u003eWe next focused on epigenetic mechanisms potentially associated with TIME composition, considering the established significance of the TIME for prognosis and treatment prediction in TNBC\u003csup\u003e56\u003c/sup\u003e. A correlation analysis between PC1 of the 1224 promoter cassettes and TIL estimates (acting as a proxy for TIME status) identified a series of cassettes (n=119) that were linked to significant changes in TILs when tumors were split into two groups based on cassette methylation patterns (Figure 4B). To prioritize cassettes for further study, we used the random forest-based Boruta feature selection algorithm identifying 28 promoter cassettes associated with TIL variation (Figure 4C, Supplementary Table 2). Among them, we focused on promoter cassette 10 (Figure 4D), representing both the largest of the 28 cassettes (n=30 CpGs) and the cassette with the highest absolute Kendall correlation to TILs. Eight genes were associated with promoter cassette 10 (\u003cem\u003eGBP4\u003c/em\u003e, \u003cem\u003eSAMD9L\u003c/em\u003e, \u003cem\u003eCARD16\u003c/em\u003e, \u003cem\u003eOAS2\u003c/em\u003e, \u003cem\u003eZBP1\u003c/em\u003e, \u003cem\u003eAPOL1\u003c/em\u003e, \u003cem\u003eLIPA\u003c/em\u003e, \u003cem\u003eBATF2\u003c/em\u003e) and showed enrichment for immune response terms, particularly innate and antiviral processes (Figure 4E). Importantly, a similar methylation pattern for cassette 10 CpGs was also observed in the\u0026nbsp;independent SCAN-B validation cohort\u0026nbsp;(Supp. Figure 10), where tumors with a hypermethylated pattern showed significantly lower rank scores of a gene expression-based immune response metagene\u003csup\u003e34\u003c/sup\u003e (Figure 4F). Finally, a similar link between the methylation pattern for cassette 10 CpGs and TIME status was also in the larger 645-sample general TCGA breast cancer cohort, with hypomethylation being associated with increased immune infiltration providing evidence for generalizability (Supp. Figure 11).\u0026nbsp;\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eFeatures of a promoter CpG cassette associated with TIME status in TNBC\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eTo understand the biological implications of promoter cassette 10 in more depth we focused on its associated genes with more than one mapped CpG to avoid false positives, leaving five genes: \u003cem\u003eGBP4\u003c/em\u003e, \u003cem\u003eOAS2\u003c/em\u003e, \u003cem\u003eCARD16\u003c/em\u003e, \u003cem\u003eZBP1\u003c/em\u003e and \u003cem\u003eSAMD9L\u003c/em\u003e. Firstly, all five genes showed two distinctive methylation states (hypermethylated and hypomethylated) among the analyzed samples linked to statistically significant differences in gene expression (Figure 5A-B). The specific DNA methylation and expression patterns per gene were also confirmed in the SCAN-B validation cohort (Supp. Figure 12) and in the TNBCs from the general TCGA breast cancer validation cohort (Supp. Figure 13A-B). In line with cassette 10\u0026rsquo;s association with TIL levels, grouping samples by the methylation status of the five genes also showed significant associations with both TIL levels and \u003cem\u003ein situ\u003c/em\u003e PDL1-CPS scores, a measure of PD-L1 expression in malignant and non-malignant cells (Figure 5C-D). Specifically, hypermethylated tumors showed a more immune cold TIME exemplified by lower TIL levels and PDL1-CPS scores. A similar trend of association between the methylation status and TILs was obtained in the TNBCs from the general TCGA breast cancer validation cohort (Supp Figure 13C). Additionally, the methylation states of the identified genes were generally linked to higher infiltration of specific immune cell types such as B cells and cytotoxic T cells, measured through both MethylCIBERSORT immune cell fractions and \u003cem\u003ein situ\u003c/em\u003e immune cell counts, with higher infiltration in the hypomethylated samples (Supp. Figure 14). Next, we analyzed the co-occurrence of the methylation patterns for the five genes in different TNBC subtypes, finding a high co-occurrence of gene hypermethylation in tumors with a non-Basal PAM50 phenotype and/or Luminal-Androgen Receptor and Mesenchymal Lehmann TNBC subtypes (Figure 5E). In support of this finding, the expression patterns of the five genes were also highly correlated across tumors in the discovery cohort (Supp. Figure 15). Similar methylation co-occurrence patterns were also observed in the SCAN-B validation cohort (Supp. Figure 16). Finally, in the discovery cohort a trend of higher cumulative hypermethylation status of the five genes with lower immune response metagene rank scores and TILs was observed (Supp. Figure 17A). In contrast, an opposite positive trend was observed for rank scores of a gene expression stroma metagene\u003csup\u003e34\u003c/sup\u003e and fibroblast fraction estimates by MethylCIBERSORT (Supp. Figure 17B), together highlighting the potential association of the methylation status with the composition of the tumor microenvironment.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCharacterization of the promoter CpG cassette associated with TIME status across breast cancer subtypes\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eAiming to characterize the role of the identified CpG cassette and genes in a more general breast cancer population we analyzed their DNA methylation patterns, gene expression and association with immune infiltration in the 645-sample TCGA breast cancer validation cohort. First, we assessed the differences in immune infiltration per PAM50 breast cancer subtype, observing the highest infiltration in the Basal intrinsic group, followed by HER2-enriched tumors. In contrast, Luminal B, and particularly Luminal A and Normal-like samples showed the lowest TIL levels (Figure 6A). Next, we evaluated whether the detected promoter cassette 10 in TNBC was still present in the complete breast cancer population. Although the structure of the cassette was less binary than in TNBC (Figure 6B), clustering on the CpGs included in the cassette still revealed a strong association between general hypomethylation and increased immune infiltration (Figure 6C).\u0026nbsp;\u003c/p\u003e\n\u003cp\u003eNext, we analyzed \u003cem\u003eGBP4\u003c/em\u003e, \u003cem\u003eOAS2\u003c/em\u003e, \u003cem\u003eZBP1\u003c/em\u003e and \u003cem\u003eCARD16\u003c/em\u003e associated with promoter cassette 10 in the TCGA cohort. \u003cem\u003eSAMD9L\u003c/em\u003e was excluded due to lack of CpG coverage in the available Infinium 450K data. After grouping samples based on the methylation status of the CpGs overlapping with respective gene, similar hyper- and hypomethylation patterns to the ones observed in TNBC were identified, that, consistent with TNBC observations, were associated with gene expression levels across all PAM50 subtypes (Supp. Figure 18). Moreover, a generally higher fraction of hypermethylated samples per gene was observed in the non-basal subtypes (Figure 6D), along with a general co-occurrence of hypermethylation or hypomethylation phenotypes among the genes consistent with TNBC findings, with the co-hypomethylated samples showing the highest TILs (Figure 6E). \u0026nbsp;Finally, we evaluated the association of the DNA methylation status of each gene with immune infiltration within the breast cancer PAM50 subtypes. In Basal, HER2-enriched and Luminal B tumors, hypomethylation of the genes was generally associated with higher TILs levels, whereas no difference was identified in Luminal A and Normal-like tumors, likely due to their generally low immune infiltration (Figure 6F).\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eMethylation and expression patterns of promoter cassette 10-associated genes in normal breast cells\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eTo investigate the potential tumor extrinsic expression patterns of the five genes, we analyzed the DNA methylation levels of CpGs in promoter cassette 10 and the expression of the five associated genes in normal breast cells using a combination of reported DNA methylation and scRNA-seq data. We observed different methylation states of cassette 10 CpGs per gene in 96 normal breast tissue specimens (Figure 7A). Two genes showed clear CpG hypomethylation phenotypes in normal breast tissue, \u003cem\u003eOAS2\u0026nbsp;\u003c/em\u003eand \u003cem\u003eGBP4\u003c/em\u003e, although the latter showed less extreme values. CpGs for \u003cem\u003eCARD16\u003c/em\u003e and \u003cem\u003eZBP1\u0026nbsp;\u003c/em\u003eshowed beta values around 0.5, indicating a potential mixture of cells with different methylation states in the analyzed bulk tissue. \u003cem\u003eSAMD9L\u003c/em\u003e could not be analyzed as no CpGs close to this gene\u0026rsquo;s promoter were present in the\u0026nbsp;Infinium\u0026nbsp;450K methylation data used. The methylation status of the genes of interest was also assessed in specific immune cell populations using data from flow cytometry-sorted blood samples, revealing distinct hypomethylation for all CpGs in all genes (\u003cem\u003eSAMD9L\u0026nbsp;\u003c/em\u003eexcluded due to lack of CpG coverage, Supp. Figure 19).\u003c/p\u003e\n\u003cp\u003eNext, we analyzed gene expression of the five genes in normal tissue using scRNA-seq data, observing heterogeneous expression in epithelial, immune and stromal cells (Figure 7B). \u0026nbsp;In epithelial cells, we observed expression of \u003cem\u003eCARD16\u0026nbsp;\u003c/em\u003eand \u003cem\u003eGBP4\u0026nbsp;\u003c/em\u003ein specific subtypes (luminal adaptive secretory precursor cell of mammary gland and mammary gland epithelial cells), while the expression measured for the remaining genes was generally low. In immune cells, the expression was equally heterogeneous: \u003cem\u003eCARD16\u003c/em\u003e was mainly expressed in dendritic cells and mature natural killer (NK) cells, \u003cem\u003eGBP4\u003c/em\u003e mostly in different NK cell populations, \u003cem\u003eSAMD9L\u003c/em\u003e in dendritic cells and macrophages, and \u003cem\u003eZBP1\u003c/em\u003e in different NK cell populations and plasma cells. \u003cem\u003eOAS2\u003c/em\u003e, on the other hand, showed low expression levels in all immune cell types. Finally, in stromal cells, \u003cem\u003eGBP4\u003c/em\u003e and \u003cem\u003eCARD16\u003c/em\u003e were expressed in most cell types except fibroblasts of mammary gland and perivascular cells. \u003cem\u003eOAS2\u0026nbsp;\u003c/em\u003eand\u003cem\u003e\u0026nbsp;SAMD9L\u0026nbsp;\u003c/em\u003eshowed lower expression, particularly focused on different endothelial cell subtypes, while \u003cem\u003eZBP1\u003c/em\u003e showed little expression in any stromal cell type. Taken together, these analyses demonstrate a heterogenous expression of the five genes in different normal breast and immune cell types, combined with typically a hypomethylated promoter pattern, especially in immune cell types.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eTumor-intrinsic expression of\u0026nbsp;\u003c/strong\u003e\u003cstrong\u003epromoter cassette 10-associated genes\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eTo investigate the tumor-intrinsic expression of the five genes associated with promoter cassette 10 we first analyzed bulk gene expression data from 34 TNBC cell lines, representing a tumor microenvironment-free context. We identified tumor-intrinsic expression of \u003cem\u003eSAMD9L\u003c/em\u003e, \u003cem\u003eCARD16\u003c/em\u003e, \u003cem\u003eOAS2\u003c/em\u003e, and \u003cem\u003eGBP4\u003c/em\u003e in some cell lines, while only very low to no expression of \u003cem\u003eZBP1\u003c/em\u003e overall (Figure 8A). Expression of the five genes was correlated across the 34 cell lines (Supp. Figure 20), in line with findings in the discovery cohort. Additionally, in eight previously reported TNBC cell lines with matched RNA-sequencing data and Illumina EPIC DNA methylation, we observed concordance between the methylation state of the promoter CpGs and gene expression levels, i.e., higher expression with CpG hypomethylation (Figure 8B). Together, these findings support a hypothesis of tumor-intrinsic methylation patterns and expression of the genes by tumor cells. Furthermore, scRNA-seq from 19 TNBC cell lines confirmed tumor-intrinsic expression of the genes in a subset of cell lines, while at the same time revealing substantial cell heterogeneity in expressing cell lines, i.e., cell lines with gene expression comprised both expressing and non-expressing cancer cells (Supp. Figure 21).\u003c/p\u003e\n\u003cp\u003eTo substantiate the cell line observations in actual tumor tissue we evaluated the tumor-intrinsic expression of the five genes in tumor tissue using scRNA-seq data from eight TNBC patients. Despite the small sample set and substantially different expression levels, we observed co-expression of the identified genes in different samples and a trend that tumors with expression of the genes tended to have an activated TIME, determined as higher counts of immune cells in the scRNA-seq data (Figure 8C). To further assess the tumor-intrinsic expression of the five genes \u003cem\u003ein situ\u003c/em\u003e we analyzed TMA cores from 65 tumors in the discovery cohort by spatial transcriptomics using Nanostring CosMX. Focusing on \u003cem\u003eOAS2\u003c/em\u003e (due to measured expression levels and CosMX gene panel content), we detected statistically significant higher mean \u003cem\u003eOAS2\u003c/em\u003e expression and higher proportion of tumor cells expressing \u003cem\u003eOAS2\u003c/em\u003e in tumors classified as hypomethylated by the DNA methylation analysis shown in Figure 5A (Figure 8D-E). Importantly, these results connect the \u003cem\u003ein vitro\u003c/em\u003e observations with similar findings in primary early-stage treatment na\u0026iuml;ve tumor tissue, i.e., apparent tumor-intrinsic expression in a subfraction of tumor cells correlating with a tumor-intrinsic hypomethylation status of the promoter region.\u003c/p\u003e"},{"header":"DISCUSSION","content":"\u003cp\u003eIn this work, we aimed to characterize methylation dynamics in cancer by investigating the variance structure of bulk tumor DNA methylation data based on the identification of groups of highly correlated CpGs, termed CpG cassettes. In our model system, TNBC, this approach identified CpGs with co-occurring methylation patterns linked to different intrinsic tumor processes and pathways, gene inactivation, and correlation to TIME status, but it also demonstrated the presence of high variance DNA methylation, apparently not biologically relevant for TNBC, that may influence downstream analyses. Moreover, our study further supports that tumor purity adjustment methods can reduce the confounding effects caused by non-malignant cells, thereby increasing the correlation between CpGs with tumor-intrinsic methylation patterns. Together, the current study exemplifies and highlights the utility of network-based analyses of DNA methylation data.\u003c/p\u003e\n\u003cp\u003eTNBC is a clinical subgroup of breast cancer with two previously reported main methylation patterns linked to a Basal/non-Basal division\u003csup\u003e21\u003c/sup\u003e. Our approach successfully identified those patterns regardless of the CpG context (promoter, proximal or distal), demonstrating that even small groups of correlated CpGs (cassettes) can accurately predict this division. Notably, our work also revealed that not all high-variance methylation patterns appear biologically relevant in the context of the disease in question, an important consideration for unsupervised analysis of DNA methylation data. For instance, the largest detected distal CpG cassette (distal cassette 1) did not show any association with tumor features such as tumor subtype, tumor purity, or TIME status, and the included CpGs showed a low overlap with enhancer regions and transcription factor binding sites, in contrast to the largest promoter and proximal cassettes. Together, these observations suggest that this distal cassette may reflect stochastic noise rather than biologically relevant patterns. Similar observations were made for other distal cassettes as well, e.g., distal CpG cassettes 4,5 and 7. It is possible that these biologically non-relevant CpG cassettes may be related to assay probe design or to their location in typically non-functional genomic regions. Still, it should be noted that these CpGs are present in both the EPIC and Infinium 450K platforms despite comprehensive data pre-processing and filtering. Moreover, our analysis of the general TCGA breast cancer cohort also demonstrates that these patterns are not specific to TNBC, but instead apply to all breast cancer subtypes, highlighting the generalizability of our findings. It is likely that similar DNA methylation patterns are also present in data from other malignancies. An important observation regarding these non-relevant cassettes was that they could not be removed using variance-based CpG filtering, a common feature selection step before unsupervised analysis. In contrast, supervised analyses appeared less affected by these specific CpGs. Taken together, these non-biological patterns may represent a significant source of bias in unsupervised analyses if not properly addressed, particularly in cohorts without strong apparent DNA methylation states (like the Basal/non-Basal division in TNBC). Nevertheless, a more functionally oriented strategy for filtering CpGs based on chromatin accessibility was able to remove these high-variance non-informative cassettes, while retaining the patterns detected in the biologically relevant cassettes that were enriched for functionally active methylation regions. Based on the availability of large-scale ATAC-seq tumor data for different malignancies (e.g. from Corces et al.\u003csup\u003e45\u003c/sup\u003e), our observations support the use of chromatin accessibility-based filtering approaches as an important filtering step when analyzing DNA methylation data, particularly when focusing on distal CpGs.\u0026nbsp;\u003c/p\u003e\n\u003cp\u003eBeyond the Basal/non-Basal DNA methylation patterns in TNBC, illustrated by several CpG cassettes irrespective of CpG context, our approach also detected cassettes that reflected other tumor specific patterns like hypermethylation of the \u003cem\u003eBRCA1\u003c/em\u003e tumor suppressor gene, representing alterations which may at least partially help to explain the molecular heterogeneity within the disease. In a focused exploration of promoter CpG cassettes we also identified tumor methylation patterns correlated with TIME status in TNBC. The limitation of this analysis to promoter cassettes was due to a more straightforward association of a CpG with a gene when located in a promoter region, and it should be acknowledged that careful analysis of both proximal and distal CpG methylation patterns may reveal additional findings. Supporting the latter, one of the previously reported TNBC DNA methylation epitypes based on tumor purity-adjusted and ATAC-seq filtered distal CpGs showed a notably immune warm phenotype (Basal3)\u003csup\u003e21\u003c/sup\u003e. Among the identified promoter CpG cassettes correlating with TIL levels we focused on the largest one, promoter cassette 10, and five of the associated genes involved in innate immune response: \u003cem\u003eGBP4\u003c/em\u003e, \u003cem\u003eOAS2\u003c/em\u003e, \u003cem\u003eCARD16\u003c/em\u003e, \u003cem\u003eZBP1\u003c/em\u003e and \u003cem\u003eSAMD9L\u003c/em\u003e. Promoter hypermethylation in these genes was mainly identified in non-Basal tumors and tumors with a Luminal Androgen Receptor or Mesenchymal-like TNBC mRNA subtype, suggesting an association with tumor-specific phenotypes to be further explored. Both Luminal Androgen Receptor and Mesenchymal-like TNBC subtypes have been described as the more immune-cold subtypes among TNBCs, consistent with our findings\u003csup\u003e57\u003c/sup\u003e. Importantly, further analyses involving cancer cell lines, scRNA-seq, and spatial transcriptomics data of both malignant and non-malignant breast cells supported that the methylation and associated gene expression patterns are likely tumor-intrinsic and not explicitly driven by contaminations from normal or immune cells in the tumor microenvironment.\u003c/p\u003e\n\u003cp\u003eThe five genes connected to promoter CpG cassette 10 are all involved in innate immune processes and have, in some cases, been associated with antitumor immune response. \u003cem\u003eGBP4\u003c/em\u003e is a GTPase belonging to the Guanylate Binding Protein family, which are central regulators of the innate immune response\u003csup\u003e58\u003c/sup\u003e and are involved in relevant innate immune processes such as the formation of the inflammasome\u003csup\u003e59\u003c/sup\u003e. Specifically, \u003cem\u003eGBP4\u003c/em\u003e has been widely identified as relevant for tumor biology in different cancer types; it has been recognized as a pan-cancer marker of immune warm tumors\u003csup\u003e60\u003c/sup\u003e, its promoter\u0026rsquo;s hypomethylation and expression have been linked to T cell infiltration and recruitment \u003cem\u003ein vitro\u0026nbsp;\u003c/em\u003ein pancreatic cancer\u003csup\u003e61\u003c/sup\u003e as well as a good prognosis biomarker in intrahepatic cholangiocarcinoma, where it was associated with higher activity of various innate immune pathways\u003csup\u003e62\u003c/sup\u003e. \u003cem\u003eOAS2\u003c/em\u003e is an innate immune protein involved in antiviral response\u003csup\u003e63\u003c/sup\u003e. \u003cem\u003eZBP1\u003c/em\u003e, on the other hand, is an immune sensor of Z conformation nucleic acids and participates in the immune response against different pathogens\u003csup\u003e64\u003c/sup\u003e. This gene has been involved in the regulation of immune response in cancer through the activation of necroptosis or PANoptosis\u003csup\u003e65,66\u003c/sup\u003e. \u003cem\u003eSAMD9L\u0026nbsp;\u003c/em\u003ehas a clear role promoting antiviral immunity and has been identified as a tumor suppressor\u003csup\u003e67\u003c/sup\u003e. Finally, \u003cem\u003eCARD16\u003c/em\u003e plays a less clear role in tumors. This protein belongs to the caspase-1 protein subfamily\u003csup\u003e68\u003c/sup\u003e,\u0026nbsp;which is generally involved in proinflammatory mechanisms and the regulation of the inflammasome activation\u003csup\u003e69\u003c/sup\u003e.\u0026nbsp;\u003c/p\u003e\n\u003cp\u003eGiven their tumor-intrinsic expression, potential epigenetic regulation, and functional implications for the innate immune response, the five genes associated with promoter cassette 10 may play an important role in shaping the composition of the TIME in TNBC, but potentially also in other breast cancer subtypes as suggested by our exploratory analysis of the general TCGA breast cancer cohort. A notable feature of the expression pattern of the five genes is the apparent intra-tumor heterogeneity with both expressing and non-expressing tumor cells within the same tissue specimen. Interestingly, this heterogeneity mirrors recent findings by us for PD-L1 in both TNBC cell lines and in primary tumors of a specific proposed epitype of TNBC (Basal3)\u003csup\u003e21\u003c/sup\u003e. Because we cannot at this point reliably deconstruct tumor methylome clonality, particularly in tumor purity adjusted bulk data due to the nature of the processing, it remains unclear whether the observed heterogeneity in expression is mirrored by a similar heterogeneity in promoter methylation status, or if it is instead determined by microenvironmental signals or stochastic expression. Furthermore, an important limitation of the current study is that causal conclusions cannot be drawn due to its observational nature. It is possible that\u0026nbsp;the observed promoter methylation patterns could reflect downstream effects of immune-related signaling, particularly hypomethylation induced by interferon-\u0026gamma; signaling, as the specific genes are interferon-inducible. Nevertheless, several findings support the hypothesis of a causal relationship between hypermethylation of the detected genes and an \u0026ldquo;immune cold\u0026rdquo; TIME. First, the detected hypomethylation pattern is not exclusive to samples with high TIL levels. Some samples with low immune responses measured through TILs show the methylation phenotype, which could indicate that it is not a consequence of immune response-related signaling pathways and that alternative alterations may be more important in some cases. Second, CpGs associated with \u003cem\u003eOAS2\u003c/em\u003e and \u003cem\u003eGBP4\u003c/em\u003e are generally hypomethylated in normal breast tissue, implying that the TIME-associated CpG sites may become hypermethylated relative to a hypomethylated baseline, rather than undergoing demethylation. Third, prior analyses have identified a relationship between T-cell recruitment and hypomethylation of \u003cem\u003eGBP4\u003c/em\u003e promoter CpGs \u003cem\u003ein vitro\u003c/em\u003e in pancreatic cancer, suggesting a causal role of hypermethylation in acquiring an immune cold phenotype\u003csup\u003e61\u003c/sup\u003e. Altogether, these observations support a hypothesis that hypermethylation and lack of expression of the identified genes are linked to reduced immune infiltration in TNBC, particularly in specific molecular phenotypes. However, whether these methylation patterns depend on the cell of origin of the tumor or are acquired during tumorigenesis remains unclear and cannot be assessed by the data presented in this work. Additional experimental validation \u003cem\u003ein vitro\u003c/em\u003e and \u003cem\u003ein vivo\u003c/em\u003e is needed to clarify the role of the detected methylation patterns in tumors and the relevance of gene expression in different normal breast cell types to establish stronger causal relationships and help identify the origin of the observed epigenetic features.\u0026nbsp;\u003c/p\u003e\n\u003cp\u003eIn conclusion, this work outlines a new comprehensive framework to analyze bulk tumor DNA methylation data, combining tumor purity adjustment, functional CpG filtering, stratification by CpG contexts, and identification of highly correlated CpG modules to enhance the identification of tumor-intrinsic DNA methylation patterns and their biological associations. Importantly, this approach extends beyond TNBC and may help to deconstruct DNA methylation patterns and thereby molecular heterogeneity also in other malignancies. In TNBC, this methodology, in combination with additional -omics layers and \u003cem\u003ein situ\u003c/em\u003e data, identified CpG cassettes not only associated with the major molecular subtypes of the disease, but also smaller cassettes linked to tumor suppressor gene inactivation. Moreover, it also identified cassettes correlating with TIME status that appeared enriched in distinct molecular phenotypes of TNBC as well as other tumor subgroups, lending support to the broader concept of epigenetic immunoediting, as previously proposed in glioblastoma\u003csup\u003e70\u003c/sup\u003e. Intriguing aspects that remain to be answered with respect to potential epigenetic immunoediting relate to the plasticity of these changes during treatment with for instance immune checkpoint inhibitors and whether treatment responses differ depending on the different mechanisms.\u003c/p\u003e"},{"header":"Declarations","content":"\u003cp\u003e\u003cstrong\u003eACKNOWLEDGEMENTS\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eFirst, the authors would like to acknowledge Jelmar Quist and Rosamund Numah for their assistance with the processing the CosMx spatial transcriptomics data. Additionally, we would like to acknowledge all patients and clinicians participating in the SCAN-B study, personnel at the central SCAN-B laboratory at the Division of Oncology, Department of Clinical Sciences Lund, Lund University, the Swedish national breast cancer quality registry (NKBC), Regional Cancer Centre South, RBC Syd, and the South Sweden Breast Cancer Group (SSBCG).\u0026nbsp;\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eFUNDING\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eFinancial support for this study was provided by the Swedish Cancer Society (CAN 2021/1407 JS, 2024/3591 JS), the Mrs Berta Kamprad Foundation (FBKS-2020-5 JS and FBKS-2024-14 JS),\u0026nbsp;the Swedish Research Council (2021-01800,\u0026nbsp;2025-02643\u0026nbsp;JS), the BCF-VÖS Foundation (JS), the National Society of Breast Cancer Associations in Sweden (JS), and Swedish governmental funding (ALF, grant 2022/0021 JS).\u0026nbsp;The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the paper. Anita Grigoriadis acknowledges support by Breast Cancer Now (147; KCL-BCN-Q3), We acknowledge use of the Spatial Biology Facility at King’s College London, supported in part by the MRC (MR/X012476/1) and the CRUK City of London Centre (CANCTA-2022/100001).\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAUTHORS CONTRIBUTIONS\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e\u003cem\u003eConception and design\u003c/em\u003e: JS and IS\u003c/p\u003e\n\u003cp\u003e\u003cem\u003eCollection and assembly of data\u003c/em\u003e: IS, AG, MJ, DFN, MA\u003c/p\u003e\n\u003cp\u003e\u003cem\u003eProvision of study material or patients\u003c/em\u003e: JVC\u003c/p\u003e\n\u003cp\u003e\u003cem\u003eData analysis and interpretation\u003c/em\u003e: IS with support of JS and DFN\u003c/p\u003e\n\u003cp\u003e\u003cem\u003eFinancial support\u003c/em\u003e: JS\u003c/p\u003e\n\u003cp\u003e\u003cem\u003eAdministrative support\u003c/em\u003e: JS\u0026nbsp;\u003c/p\u003e\n\u003cp\u003e\u003cem\u003eManuscript writing\u003c/em\u003e: IS with the support of all authors\u003c/p\u003e\n\u003cp\u003e\u003cem\u003eFinal approval of manuscript\u003c/em\u003e: All authors\u003c/p\u003e\n\u003cp\u003e\u003cem\u003eAgree to be accountable for all aspects of the work\u003c/em\u003e: All authors\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCOMPETING INTERESTS\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThe authors declare that they have no financial or non-financial competing interests.\u0026nbsp;\u003c/p\u003e"},{"header":"References","content":"\u003col\u003e\n\u003cli\u003eHanahan, D. Hallmarks of Cancer: New Dimensions. \u003cem\u003eCancer Discov\u003c/em\u003e \u003cstrong\u003e12\u003c/strong\u003e, 31\u0026ndash;46 (2022). https://doi.org/10.1158/2159-8290.CD-21-1059\u003c/li\u003e\n\u003cli\u003eLakshminarasimhan, R. \u0026amp; Liang, G. The Role of DNA Methylation in Cancer. \u003cem\u003eAdv Exp Med Biol\u003c/em\u003e \u003cstrong\u003e945\u003c/strong\u003e, 151\u0026ndash;172 (2016). https://doi.org/10.1007/978-3-319-43624-1_7\u003c/li\u003e\n\u003cli\u003eGlodzik, D.\u003cem\u003e et al.\u003c/em\u003e Comprehensive molecular comparison of BRCA1 hypermethylated and BRCA1 mutated triple negative breast cancers. \u003cem\u003eNat Commun\u003c/em\u003e \u003cstrong\u003e11\u003c/strong\u003e, 3747 (2020). https://doi.org/10.1038/s41467-020-17537-2\u003c/li\u003e\n\u003cli\u003eKalachand, R. D.\u003cem\u003e et al.\u003c/em\u003e BRCA1 Promoter Methylation and Clinical Outcomes in Ovarian Cancer: An Individual Patient Data Meta-Analysis. \u003cem\u003eJ Natl Cancer Inst\u003c/em\u003e \u003cstrong\u003e112\u003c/strong\u003e, 1190\u0026ndash;1203 (2020). https://doi.org/10.1093/jnci/djaa070\u003c/li\u003e\n\u003cli\u003eLi, C.\u003cem\u003e et al.\u003c/em\u003e Mi-2beta promotes immune evasion in melanoma by activating EZH2 methylation. \u003cem\u003eNat Commun\u003c/em\u003e \u003cstrong\u003e15\u003c/strong\u003e, 2163 (2024). https://doi.org/10.1038/s41467-024-46422-5\u003c/li\u003e\n\u003cli\u003eJung, H.\u003cem\u003e et al.\u003c/em\u003e DNA methylation loss promotes immune evasion of tumours with high mutation and copy number load. \u003cem\u003eNat Commun\u003c/em\u003e \u003cstrong\u003e10\u003c/strong\u003e, 4278 (2019). https://doi.org/10.1038/s41467-019-12159-9\u003c/li\u003e\n\u003cli\u003eShigeyasu, K.\u003cem\u003e et al.\u003c/em\u003e Clinical Significance of MLH1 Methylation and CpG Island Methylator Phenotype as Prognostic Markers in Patients with Gastric Cancer. \u003cem\u003ePLoS One\u003c/em\u003e \u003cstrong\u003e10\u003c/strong\u003e, e0130409 (2015). https://doi.org/10.1371/journal.pone.0130409\u003c/li\u003e\n\u003cli\u003eFackler, M. J.\u003cem\u003e et al.\u003c/em\u003e DNA methylation markers predict recurrence-free interval in triple-negative breast cancer. \u003cem\u003eNPJ Breast Cancer\u003c/em\u003e \u003cstrong\u003e6\u003c/strong\u003e, 3 (2020). https://doi.org/10.1038/s41523-020-0145-3\u003c/li\u003e\n\u003cli\u003eRauluseviciute, I., Drablos, F. \u0026amp; Rye, M. B. DNA methylation data by sequencing: experimental approaches and recommendations for tools and pipelines for data analysis. \u003cem\u003eClin Epigenetics\u003c/em\u003e \u003cstrong\u003e11\u003c/strong\u003e, 193 (2019). https://doi.org/10.1186/s13148-019-0795-x\u003c/li\u003e\n\u003cli\u003eJaffe, A. E. \u0026amp; Irizarry, R. A. Accounting for cellular heterogeneity is critical in epigenome-wide association studies. \u003cem\u003eGenome Biol\u003c/em\u003e \u003cstrong\u003e15\u003c/strong\u003e, R31 (2014). https://doi.org/10.1186/gb-2014-15-2-r31\u003c/li\u003e\n\u003cli\u003eWang, Y.\u003cem\u003e et al.\u003c/em\u003e Prognostic and Predictive Value of Three DNA Methylation Signatures in Lung Adenocarcinoma. \u003cem\u003eFront Genet\u003c/em\u003e \u003cstrong\u003e10\u003c/strong\u003e, 349 (2019). https://doi.org/10.3389/fgene.2019.00349\u003c/li\u003e\n\u003cli\u003eDiNome, M. L.\u003cem\u003e et al.\u003c/em\u003e Clinicopathological Features of Triple-Negative Breast Cancer Epigenetic Subtypes. \u003cem\u003eAnn Surg Oncol\u003c/em\u003e \u003cstrong\u003e26\u003c/strong\u003e, 3344\u0026ndash;3353 (2019). https://doi.org/10.1245/s10434-019-07565-8\u003c/li\u003e\n\u003cli\u003eLin, L. H.\u003cem\u003e et al.\u003c/em\u003e DNA Methylation Identifies Epigenetic Subtypes of Triple-Negative Breast Cancers With Distinct Clinicopathologic and Molecular Features. \u003cem\u003eMod Pathol\u003c/em\u003e \u003cstrong\u003e36\u003c/strong\u003e, 100306 (2023). https://doi.org/10.1016/j.modpat.2023.100306\u003c/li\u003e\n\u003cli\u003eBourgon, R., Gentleman, R. \u0026amp; Huber, W. Independent filtering increases detection power for high-throughput experiments. \u003cem\u003eProc Natl Acad Sci U S A\u003c/em\u003e \u003cstrong\u003e107\u003c/strong\u003e, 9546\u0026ndash;9551 (2010). https://doi.org/10.1073/pnas.0914005107\u003c/li\u003e\n\u003cli\u003eSandoval, J.\u003cem\u003e et al.\u003c/em\u003e A prognostic DNA methylation signature for stage I non-small-cell lung cancer. \u003cem\u003eJ Clin Oncol\u003c/em\u003e \u003cstrong\u003e31\u003c/strong\u003e, 4140\u0026ndash;4147 (2013). https://doi.org/10.1200/JCO.2012.48.5516\u003c/li\u003e\n\u003cli\u003eStaaf, J. \u0026amp; Aine, M. Tumor purity adjusted beta values improve biological interpretability of high-dimensional DNA methylation data. \u003cem\u003ePLoS One\u003c/em\u003e \u003cstrong\u003e17\u003c/strong\u003e, e0265557 (2022). https://doi.org/10.1371/journal.pone.0265557\u003c/li\u003e\n\u003cli\u003eSasiain, I., Nacer, D. F., Aine, M., Veerla, S. \u0026amp; Staaf, J. Tumor purity estimated from bulk DNA methylation can be used for adjusting beta values of individual samples to better reflect tumor biology. \u003cem\u003eNAR Genom Bioinform\u003c/em\u003e \u003cstrong\u003e6\u003c/strong\u003e, lqae146 (2024). https://doi.org/10.1093/nargab/lqae146\u003c/li\u003e\n\u003cli\u003eQin, Y., Feng, H., Chen, M., Wu, H. \u0026amp; Zheng, X. InfiniumPurify: An R package for estimating and accounting for tumor purity in cancer methylation research. \u003cem\u003eGenes Dis\u003c/em\u003e \u003cstrong\u003e5\u003c/strong\u003e, 43\u0026ndash;45 (2018). https://doi.org/10.1016/j.gendis.2018.02.003\u003c/li\u003e\n\u003cli\u003eAran, D., Sabato, S. \u0026amp; Hellman, A. DNA methylation of distal regulatory sites characterizes dysregulation of cancer genes. \u003cem\u003eGenome Biol\u003c/em\u003e \u003cstrong\u003e14\u003c/strong\u003e, R21 (2013). https://doi.org/10.1186/gb-2013-14-3-r21\u003c/li\u003e\n\u003cli\u003eAran, D. \u0026amp; Hellman, A. DNA methylation of transcriptional enhancers and cancer predisposition. \u003cem\u003eCell\u003c/em\u003e \u003cstrong\u003e154\u003c/strong\u003e, 11\u0026ndash;13 (2013). https://doi.org/10.1016/j.cell.2013.06.018\u003c/li\u003e\n\u003cli\u003eAine, M.\u003cem\u003e et al.\u003c/em\u003e The DNA methylation landscape of primary triple-negative breast cancer. \u003cem\u003eNat Commun\u003c/em\u003e \u003cstrong\u003e16\u003c/strong\u003e, 3041 (2025). https://doi.org/10.1038/s41467-025-58158-x\u003c/li\u003e\n\u003cli\u003eHarbeck, N.\u003cem\u003e et al.\u003c/em\u003e Breast cancer. \u003cem\u003eNat Rev Dis Primers\u003c/em\u003e \u003cstrong\u003e5\u003c/strong\u003e, 66 (2019). https://doi.org/10.1038/s41572-019-0111-2\u003c/li\u003e\n\u003cli\u003eZagami, P. \u0026amp; Carey, L. A. Triple negative breast cancer: Pitfalls and progress. \u003cem\u003eNPJ Breast Cancer\u003c/em\u003e \u003cstrong\u003e8\u003c/strong\u003e, 95 (2022). https://doi.org/10.1038/s41523-022-00468-0\u003c/li\u003e\n\u003cli\u003eParker, J. S.\u003cem\u003e et al.\u003c/em\u003e Supervised risk predictor of breast cancer based on intrinsic subtypes. \u003cem\u003eJ Clin Oncol\u003c/em\u003e \u003cstrong\u003e27\u003c/strong\u003e, 1160\u0026ndash;1167 (2009). https://doi.org/10.1200/JCO.2008.18.1370\u003c/li\u003e\n\u003cli\u003eLehmann, B. D.\u003cem\u003e et al.\u003c/em\u003e Identification of human triple-negative breast cancer subtypes and preclinical models for selection of targeted therapies. \u003cem\u003eJ Clin Invest\u003c/em\u003e \u003cstrong\u003e121\u003c/strong\u003e, 2750\u0026ndash;2767 (2011). https://doi.org/10.1172/JCI45014\u003c/li\u003e\n\u003cli\u003eLehmann, B. D.\u003cem\u003e et al.\u003c/em\u003e Refinement of Triple-Negative Breast Cancer Molecular Subtypes: Implications for Neoadjuvant Chemotherapy Selection. \u003cem\u003ePLoS One\u003c/em\u003e \u003cstrong\u003e11\u003c/strong\u003e, e0157368 (2016). https://doi.org/10.1371/journal.pone.0157368\u003c/li\u003e\n\u003cli\u003eStaaf, J.\u003cem\u003e et al.\u003c/em\u003e Whole-genome sequencing of triple-negative breast cancers in a population-based clinical study. \u003cem\u003eNat Med\u003c/em\u003e \u003cstrong\u003e25\u003c/strong\u003e, 1526\u0026ndash;1533 (2019). https://doi.org/10.1038/s41591-019-0582-4\u003c/li\u003e\n\u003cli\u003eHuang, R., Ding, P. \u0026amp; Yang, F. Clinicopathological significance and potential drug target of CDH1 in breast cancer: a meta-analysis and literature review. \u003cem\u003eDrug Des Devel Ther\u003c/em\u003e \u003cstrong\u003e9\u003c/strong\u003e, 5277\u0026ndash;5285 (2015). https://doi.org/10.2147/DDDT.S86929\u003c/li\u003e\n\u003cli\u003eRoostee, S.\u003cem\u003e et al.\u003c/em\u003e Tumour immune characterisation of primary triple-negative breast cancer using automated image quantification of immunohistochemistry-stained immune cells. \u003cem\u003eSci Rep\u003c/em\u003e \u003cstrong\u003e14\u003c/strong\u003e, 21417 (2024). https://doi.org/10.1038/s41598-024-72306-1\u003c/li\u003e\n\u003cli\u003eKarn, T.\u003cem\u003e et al.\u003c/em\u003e Association Between Genomic Metrics and Immune Infiltration in Triple-Negative Breast Cancer. \u003cem\u003eJAMA Oncol\u003c/em\u003e \u003cstrong\u003e3\u003c/strong\u003e, 1707\u0026ndash;1711 (2017). https://doi.org/10.1001/jamaoncol.2017.2140\u003c/li\u003e\n\u003cli\u003eSaal, L. H.\u003cem\u003e et al.\u003c/em\u003e The Sweden Cancerome Analysis Network - Breast (SCAN-B) Initiative: a large-scale multicenter infrastructure towards implementation of breast cancer genomic analyses in the clinical routine. \u003cem\u003eGenome Med\u003c/em\u003e \u003cstrong\u003e7\u003c/strong\u003e, 20 (2015). https://doi.org/10.1186/s13073-015-0131-9\u003c/li\u003e\n\u003cli\u003eRyden, L.\u003cem\u003e et al.\u003c/em\u003e Minimizing inequality in access to precision medicine in breast cancer by real-time population-based molecular analysis in the SCAN-B initiative. \u003cem\u003eBr J Surg\u003c/em\u003e \u003cstrong\u003e105\u003c/strong\u003e, e158\u0026ndash;e168 (2018). https://doi.org/10.1002/bjs.10741\u003c/li\u003e\n\u003cli\u003eStaaf, J.\u003cem\u003e et al.\u003c/em\u003e RNA sequencing-based single sample predictors of molecular subtype and risk of recurrence for clinical assessment of early-stage breast cancer. \u003cem\u003eNPJ Breast Cancer\u003c/em\u003e \u003cstrong\u003e8\u003c/strong\u003e, 94 (2022). https://doi.org/10.1038/s41523-022-00465-3\u003c/li\u003e\n\u003cli\u003eFredlund, E.\u003cem\u003e et al.\u003c/em\u003e The gene expression landscape of breast cancer is shaped by tumor protein p53 status and epithelial-mesenchymal transition. \u003cem\u003eBreast Cancer Res\u003c/em\u003e \u003cstrong\u003e14\u003c/strong\u003e, R113 (2012). https://doi.org/10.1186/bcr3236\u003c/li\u003e\n\u003cli\u003eAine, M.\u003cem\u003e et al.\u003c/em\u003e Molecular analyses of triple-negative breast cancer in the young and elderly. \u003cem\u003eBreast Cancer Res\u003c/em\u003e \u003cstrong\u003e23\u003c/strong\u003e, 20 (2021). https://doi.org/10.1186/s13058-021-01392-0\u003c/li\u003e\n\u003cli\u003eSigurjonsdottir, G.\u003cem\u003e et al.\u003c/em\u003e Evaluation of alternative prognostic thresholds for SP142 and 22C3 immunohistochemical PD-L1 expression in triple-negative breast cancer: results from a population-based cohort. \u003cem\u003eBreast Cancer Res Treat\u003c/em\u003e \u003cstrong\u003e210\u003c/strong\u003e, 271\u0026ndash;284 (2025). https://doi.org/10.1007/s10549-024-07561-x\u003c/li\u003e\n\u003cli\u003eFishilevich, S.\u003cem\u003e et al.\u003c/em\u003e GeneHancer: genome-wide integration of enhancers and target genes in GeneCards. \u003cem\u003eDatabase (Oxford)\u003c/em\u003e \u003cstrong\u003e2017\u003c/strong\u003e (2017). https://doi.org/10.1093/database/bax028\u003c/li\u003e\n\u003cli\u003eCha, Y. J.\u003cem\u003e et al.\u003c/em\u003e Genomic Characteristics Related to Histology-Based Immune Features in Breast Cancer. \u003cem\u003eMod Pathol\u003c/em\u003e \u003cstrong\u003e38\u003c/strong\u003e, 100736 (2025). https://doi.org/10.1016/j.modpat.2025.100736\u003c/li\u003e\n\u003cli\u003eHair, B. Y.\u003cem\u003e et al.\u003c/em\u003e Body mass index associated with genome-wide methylation in breast tissue. \u003cem\u003eBreast Cancer Res Treat\u003c/em\u003e \u003cstrong\u003e151\u003c/strong\u003e, 453\u0026ndash;463 (2015). https://doi.org/10.1007/s10549-015-3401-8\u003c/li\u003e\n\u003cli\u003eReinius, L. E.\u003cem\u003e et al.\u003c/em\u003e Differential DNA methylation in purified human blood cells: implications for cell lineage and studies on disease susceptibility. \u003cem\u003ePLoS One\u003c/em\u003e \u003cstrong\u003e7\u003c/strong\u003e, e41361 (2012). https://doi.org/10.1371/journal.pone.0041361\u003c/li\u003e\n\u003cli\u003eReed, A. D.\u003cem\u003e et al.\u003c/em\u003e A single-cell atlas enables mapping of homeostatic cellular shifts in the adult human breast. \u003cem\u003eNat Genet\u003c/em\u003e \u003cstrong\u003e56\u003c/strong\u003e, 652\u0026ndash;662 (2024). https://doi.org/10.1038/s41588-024-01688-9\u003c/li\u003e\n\u003cli\u003eJovanovic, B.\u003cem\u003e et al.\u003c/em\u003e Heterogeneity and transcriptional drivers of triple-negative breast cancer. \u003cem\u003eCell Rep\u003c/em\u003e \u003cstrong\u003e42\u003c/strong\u003e, 113564 (2023). https://doi.org/10.1016/j.celrep.2023.113564\u003c/li\u003e\n\u003cli\u003eChen, Y., Pal, B., Lindeman, G. J., Visvader, J. E. \u0026amp; Smyth, G. K. R code and downstream analysis objects for the scRNA-seq atlas of normal and tumorigenic human breast tissue. \u003cem\u003eSci Data\u003c/em\u003e \u003cstrong\u003e9\u003c/strong\u003e, 96 (2022). https://doi.org/10.1038/s41597-022-01236-2\u003c/li\u003e\n\u003cli\u003eLangfelder, P. \u0026amp; Horvath, S. WGCNA: an R package for weighted correlation network analysis. \u003cem\u003eBMC Bioinformatics\u003c/em\u003e \u003cstrong\u003e9\u003c/strong\u003e, 559 (2008). https://doi.org/10.1186/1471-2105-9-559\u003c/li\u003e\n\u003cli\u003eCorces, M. R.\u003cem\u003e et al.\u003c/em\u003e The chromatin accessibility landscape of primary human cancers. \u003cem\u003eScience\u003c/em\u003e \u003cstrong\u003e362\u003c/strong\u003e (2018). https://doi.org/10.1126/science.aav1898\u003c/li\u003e\n\u003cli\u003eChakravarthy, A.\u003cem\u003e et al.\u003c/em\u003e Author Correction: Pan-cancer deconvolution of tumour composition using DNA methylation. \u003cem\u003eNat Commun\u003c/em\u003e \u003cstrong\u003e9\u003c/strong\u003e, 4642 (2018). https://doi.org/10.1038/s41467-018-07155-4\u003c/li\u003e\n\u003cli\u003eNewman, A. M.\u003cem\u003e et al.\u003c/em\u003e Determining cell type abundance and expression from bulk tissues with digital cytometry. \u003cem\u003eNat Biotechnol\u003c/em\u003e \u003cstrong\u003e37\u003c/strong\u003e, 773\u0026ndash;782 (2019). https://doi.org/10.1038/s41587-019-0114-2\u003c/li\u003e\n\u003cli\u003ePhipson, B., Maksimovic, J. \u0026amp; Oshlack, A. missMethyl: an R package for analyzing data from Illumina\u0026apos;s HumanMethylation450 platform. \u003cem\u003eBioinformatics\u003c/em\u003e \u003cstrong\u003e32\u003c/strong\u003e, 286\u0026ndash;288 (2016). https://doi.org/10.1093/bioinformatics/btv560\u003c/li\u003e\n\u003cli\u003eLiberzon, A.\u003cem\u003e et al.\u003c/em\u003e The Molecular Signatures Database (MSigDB) hallmark gene set collection. \u003cem\u003eCell Syst\u003c/em\u003e \u003cstrong\u003e1\u003c/strong\u003e, 417\u0026ndash;425 (2015). https://doi.org/10.1016/j.cels.2015.12.004\u003c/li\u003e\n\u003cli\u003eKursa, M. B. \u0026amp; Rudnicki, W. R. Feature Selection with the Boruta Package. \u003cem\u003eJ Stat Softw\u003c/em\u003e \u003cstrong\u003e36\u003c/strong\u003e, 1\u0026ndash;13 (2010). https://doi.org/DOI 10.18637/jss.v036.i11\u003c/li\u003e\n\u003cli\u003eButler, A., Hoffman, P., Smibert, P., Papalexi, E. \u0026amp; Satija, R. Integrating single-cell transcriptomic data across different conditions, technologies, and species. \u003cem\u003eNat Biotechnol\u003c/em\u003e \u003cstrong\u003e36\u003c/strong\u003e, 411\u0026ndash;420 (2018). https://doi.org/10.1038/nbt.4096\u003c/li\u003e\n\u003cli\u003ePalla, G.\u003cem\u003e et al.\u003c/em\u003e Squidpy: a scalable framework for spatial omics analysis. \u003cem\u003eNat Methods\u003c/em\u003e \u003cstrong\u003e19\u003c/strong\u003e, 171\u0026ndash;178 (2022). https://doi.org/10.1038/s41592-021-01358-2\u003c/li\u003e\n\u003cli\u003eWolf, F. A., Angerer, P. \u0026amp; Theis, F. J. SCANPY: large-scale single-cell gene expression data analysis. \u003cem\u003eGenome Biol\u003c/em\u003e \u003cstrong\u003e19\u003c/strong\u003e, 15 (2018). https://doi.org/10.1186/s13059-017-1382-0\u003c/li\u003e\n\u003cli\u003ePedregosa, F.\u003cem\u003e et al.\u003c/em\u003e Scikit-learn: Machine Learning in Python. \u003cem\u003eJ Mach Learn Res\u003c/em\u003e \u003cstrong\u003e12\u003c/strong\u003e, 2825\u0026ndash;2830 (2011). \u003c/li\u003e\n\u003cli\u003eJjingo, D., Conley, A. B., Yi, S. V., Lunyak, V. V. \u0026amp; Jordan, I. K. On the presence and role of human gene-body DNA methylation. \u003cem\u003eOncotarget\u003c/em\u003e \u003cstrong\u003e3\u003c/strong\u003e, 462\u0026ndash;474 (2012). https://doi.org/10.18632/oncotarget.497\u003c/li\u003e\n\u003cli\u003eLoi, S.\u003cem\u003e et al.\u003c/em\u003e The journey of tumor-infiltrating lymphocytes as a biomarker in breast cancer: clinical utility in an era of checkpoint inhibition. \u003cem\u003eAnn Oncol\u003c/em\u003e \u003cstrong\u003e32\u003c/strong\u003e, 1236\u0026ndash;1244 (2021). https://doi.org/10.1016/j.annonc.2021.07.007\u003c/li\u003e\n\u003cli\u003eHarano, K.\u003cem\u003e et al.\u003c/em\u003e Rates of immune cell infiltration in patients with triple-negative breast cancer by molecular subtype. \u003cem\u003ePLoS One\u003c/em\u003e \u003cstrong\u003e13\u003c/strong\u003e, e0204513 (2018). https://doi.org/10.1371/journal.pone.0204513\u003c/li\u003e\n\u003cli\u003eTretina, K., Park, E. S., Maminska, A. \u0026amp; MacMicking, J. D. Interferon-induced guanylate-binding proteins: Guardians of host defense in health and disease. \u003cem\u003eJ Exp Med\u003c/em\u003e \u003cstrong\u003e216\u003c/strong\u003e, 482\u0026ndash;500 (2019). https://doi.org/10.1084/jem.20182031\u003c/li\u003e\n\u003cli\u003eJang, J. H., Kim, D. H. \u0026amp; Surh, Y. J. Dynamic roles of inflammasomes in inflammatory tumor microenvironment. \u003cem\u003eNPJ Precis Oncol\u003c/em\u003e \u003cstrong\u003e5\u003c/strong\u003e, 18 (2021). https://doi.org/10.1038/s41698-021-00154-7\u003c/li\u003e\n\u003cli\u003eZhou, W., Yeerkenbieke, G., Zhang, Y., Zhou, M. \u0026amp; Li, J. Guanylate binding protein 4 shapes an inflamed tumor microenvironment and identifies immuno-hot tumors. \u003cem\u003eJ Cancer Res Clin Oncol\u003c/em\u003e \u003cstrong\u003e150\u003c/strong\u003e, 90 (2024). https://doi.org/10.1007/s00432-024-05605-9\u003c/li\u003e\n\u003cli\u003eTasiheng, Y.\u003cem\u003e et al.\u003c/em\u003e DNA hypo-methylation and expression of GBP4 induces T cell exhaustion in pancreatic cancer. \u003cem\u003eCancer Immunol Immunother\u003c/em\u003e \u003cstrong\u003e73\u003c/strong\u003e, 208 (2024). https://doi.org/10.1007/s00262-024-03786-3\u003c/li\u003e\n\u003cli\u003eLiao, H.\u003cem\u003e et al.\u003c/em\u003e Whole-Genome DNA Methylation Profiling of Intrahepatic Cholangiocarcinoma Reveals Prognostic Subtypes with Distinct Biological Drivers. \u003cem\u003eCancer Res\u003c/em\u003e \u003cstrong\u003e84\u003c/strong\u003e, 1747\u0026ndash;1763 (2024). https://doi.org/10.1158/0008-5472.CAN-23-3298\u003c/li\u003e\n\u003cli\u003eMerold, V.\u003cem\u003e et al.\u003c/em\u003e Structural basis for OAS2 regulation and its antiviral function. \u003cem\u003eMol Cell\u003c/em\u003e \u003cstrong\u003e85\u003c/strong\u003e, 2176\u0026ndash;2193 e2113 (2025). https://doi.org/10.1016/j.molcel.2025.05.001\u003c/li\u003e\n\u003cli\u003eXie, F.\u003cem\u003e et al.\u003c/em\u003e ZBP1 condensate formation synergizes Z-NAs recognition and signal transduction. \u003cem\u003eCell Death Dis\u003c/em\u003e \u003cstrong\u003e15\u003c/strong\u003e, 487 (2024). https://doi.org/10.1038/s41419-024-06889-y\u003c/li\u003e\n\u003cli\u003eZhang, T.\u003cem\u003e et al.\u003c/em\u003e ADAR1 masks the cancer immunotherapeutic promise of ZBP1-driven necroptosis. \u003cem\u003eNature\u003c/em\u003e \u003cstrong\u003e606\u003c/strong\u003e, 594\u0026ndash;602 (2022). https://doi.org/10.1038/s41586-022-04753-7\u003c/li\u003e\n\u003cli\u003eBaik, J. Y.\u003cem\u003e et al.\u003c/em\u003e ZBP1 not RIPK1 mediates tumor necroptosis in breast cancer. \u003cem\u003eNat Commun\u003c/em\u003e \u003cstrong\u003e12\u003c/strong\u003e, 2666 (2021). https://doi.org/10.1038/s41467-021-23004-3\u003c/li\u003e\n\u003cli\u003ePeng, S.\u003cem\u003e et al.\u003c/em\u003e Structure and function of an effector domain in antiviral factors and tumor suppressors SAMD9 and SAMD9L. \u003cem\u003eProc Natl Acad Sci U S A\u003c/em\u003e \u003cstrong\u003e119\u003c/strong\u003e (2022). https://doi.org/10.1073/pnas.2116550119\u003c/li\u003e\n\u003cli\u003eHolland, M., Rutkowski, R. \u0026amp; T, C. L. Evolutionary Dynamics of Proinflammatory Caspases in Primates and Rodents. \u003cem\u003eMol Biol Evol\u003c/em\u003e \u003cstrong\u003e41\u003c/strong\u003e (2024). https://doi.org/10.1093/molbev/msae220\u003c/li\u003e\n\u003cli\u003eDevi, S.\u003cem\u003e et al.\u003c/em\u003e CARD-only proteins regulate in vivo inflammasome responses and ameliorate gout. \u003cem\u003eCell Rep\u003c/em\u003e \u003cstrong\u003e42\u003c/strong\u003e, 112265 (2023). https://doi.org/10.1016/j.celrep.2023.112265\u003c/li\u003e\n\u003cli\u003eGangoso, E.\u003cem\u003e et al.\u003c/em\u003e Glioblastomas acquire myeloid-affiliated transcriptional programs via epigenetic immunoediting to elicit immune evasion. \u003cem\u003eCell\u003c/em\u003e \u003cstrong\u003e184\u003c/strong\u003e, 2454\u0026ndash;2470 e2426 (2021). https://doi.org/10.1016/j.cell.2021.03.023\u003c/li\u003e\n\u003c/ol\u003e"}],"fulltextSource":"","fullText":"","funders":[],"hasAdminPriorityOnWorkflow":false,"hasManuscriptDocX":true,"hasOptedInToPreprint":true,"hasPassedJournalQc":"","hasAnyPriority":false,"hideJournal":false,"highlight":"","institution":"","isAcceptedByJournal":true,"isAuthorSuppliedPdf":false,"isDeskRejected":"","isHiddenFromSearch":false,"isInQc":false,"isInWorkflow":false,"isPdf":false,"isPdfUpToDate":true,"isWithdrawnOrRetracted":false,"journal":{"display":true,"email":"[email protected]","identity":"npj-breast-cancer","isNatureJournal":false,"hasQc":false,"allowDirectSubmit":false,"externalIdentity":"npjbcancer","sideBox":"Learn more about [npj Breast Cancer](http://www.nature.com/npjbcancer/)","snPcode":"41523","submissionUrl":"https://mts-npjbcancer.nature.com/","title":"npj Breast Cancer","twitterHandle":"","acdcEnabled":true,"dfaEnabled":true,"editorialSystem":"ejp","reportingPortfolio":"NPJ","inReviewEnabled":true,"inReviewRevisionsEnabled":true},"keywords":"","lastPublishedDoi":"10.21203/rs.3.rs-8250160/v1","lastPublishedDoiUrl":"https://doi.org/10.21203/rs.3.rs-8250160/v1","license":{"name":"CC BY 4.0","url":"https://creativecommons.org/licenses/by/4.0/"},"manuscriptAbstract":"DNA methylation deregulation is an essential feature of tumor biology, influencing cancer formation and pathogenesis. Analyses of DNA methylation in bulk cancer specimens are challenging due to the high data dimensionality, noise, and mixture of different cell types. Here we characterized methylation dynamics in cancer by investigating the variance structure of bulk tumor DNA methylation data through identification of groups of highly correlated CpGs, termed CpG cassettes. Using triple-negative breast cancer (TNBC) as a model system our approach identified co-occurring methylation patterns linked to different intrinsic tumor processes and pathways, gene inactivation, and composition of the tumor immune microenvironment. Our framework also demonstrated the presence of high variance DNA methylation seemingly not linked to tumor biology that could be excluded using chromatin accessibility filtering. Together, this work outlines a comprehensive approach to analyze bulk tumor DNA methylation data combining tumor purity adjustment, functional CpG filtering, stratification by CpG contexts, and identification of highly correlated CpG modules to enhance tumor-intrinsic DNA methylation patterns and our understanding of processes shaping epigenetic tumor evolution.","manuscriptTitle":"Genome-Wide DNA Methylation Profiling of Triple-Negative Breast Cancer Uncovers Epigenetic Biomarkers of Tumor Identity and the Immune Microenvironment","msid":"","msnumber":"","nonDraftVersions":[{"code":1,"date":"2025-12-12 10:01:18","doi":"10.21203/rs.3.rs-8250160/v1","editorialEvents":[{"type":"communityComments","content":0},{"type":"decision","content":"Revision requested","date":"2026-01-12T16:00:38+00:00","index":"","fulltext":""},{"type":"editorInvitedReview","content":"","date":"2026-01-12T10:11:36+00:00","index":"hide","fulltext":""},{"type":"editorInvitedReview","content":"","date":"2026-01-05T15:54:32+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"111974665712481754146188191531595014192","date":"2025-12-22T22:21:39+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"236805916068547500598630070460765124918","date":"2025-12-22T04:55:16+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"330152237714325726903187936841222590886","date":"2025-12-20T07:26:05+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"10372660182463826915025672690471894757","date":"2025-12-09T16:39:11+00:00","index":"hide","fulltext":""},{"type":"reviewersInvited","content":"","date":"2025-12-09T11:18:45+00:00","index":"","fulltext":""},{"type":"editorAssigned","content":"","date":"2025-12-02T17:06:14+00:00","index":"","fulltext":""},{"type":"checksComplete","content":"","date":"2025-12-02T16:33:34+00:00","index":"","fulltext":""},{"type":"submitted","content":"npj Breast Cancer","date":"2025-12-01T11:58:31+00:00","index":"","fulltext":""}],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"npj-breast-cancer","isNatureJournal":false,"hasQc":false,"allowDirectSubmit":false,"externalIdentity":"npjbcancer","sideBox":"Learn more about [npj Breast Cancer](http://www.nature.com/npjbcancer/)","snPcode":"41523","submissionUrl":"https://mts-npjbcancer.nature.com/","title":"npj Breast Cancer","twitterHandle":"","acdcEnabled":true,"dfaEnabled":true,"editorialSystem":"ejp","reportingPortfolio":"NPJ","inReviewEnabled":true,"inReviewRevisionsEnabled":true}}],"origin":"","ownerIdentity":"a572cea6-acf8-4fad-b61e-bbf148f41594","owner":[],"postedDate":"December 12th, 2025","published":true,"recentEditorialEvents":[],"rejectedJournal":[],"revision":"","amendment":"","status":"under-review","subjectAreas":[{"id":59511998,"name":"Health sciences/Biomarkers"},{"id":59511999,"name":"Biological sciences/Cancer"},{"id":59512000,"name":"Biological sciences/Computational biology and bioinformatics"},{"id":59512001,"name":"Biological sciences/Genetics"},{"id":59512002,"name":"Health sciences/Oncology"}],"tags":[],"updatedAt":"2026-04-21T21:38:32+00:00","versionOfRecord":[],"versionCreatedAt":"2025-12-12 10:01:18","video":"","vorDoi":"","vorDoiUrl":"","workflowStages":[]},"version":"v1","identity":"rs-8250160","journalConfig":"researchsquare"},"__N_SSP":true},"page":"/article/[identity]/[[...version]]","query":{"redirect":"/article/rs-8250160","identity":"rs-8250160","version":["v1"]},"buildId":"8U1c8b4HqxoKbykW_rLl7","isFallback":false,"isExperimentalCompile":false,"dynamicIds":[84888],"gssp":true,"scriptLoader":[]}