An unbiased approach to measure aberrant DNA methylation alterations

Domain: biology
Published: 2026-03-27
Editorial outcome: major_revision

Abstract
The ability to accurately measure aberrant DNA methylation levels is integral to the understanding of DNA methylation biology. It is well-established that in cancer, the largest, and thus, most biologically important absolute gains of DNA methylation levels occur at CpG sites with low native levels while the largest losses occur at CpG sites with high native levels. Conventional wisdom assumes that the observed association between the degree of the alterations and the native levels are largely due to the limitations of change within the DNA methylation scale. Here, we present evidence that this association is largely caused by alterations occurring as a global rate of change relative to the native level. We show that DNA methylation alterations can be accurately compared by calculating the rate of change relative to the native level. Most importantly, this approach enables the identification of more biologically significant DNA methylation alterations.

Introduction
Methylated DNA can be found in plants, animals, fungi, and bacteria, and is a key epigenetic regulator of gene expression, and thus, is essential for normal cell function1. In general, DNA methylation studies are conducted by comparing the absolute difference in DNA methylation levels of CpG sites, which are commonly measured as the percentage or proportion (Beta value) of cytosine residues that are methylated. Studies that identify DNA methylation alterations of CpG sites often do so by using two characteristics: the probability that the DNA methylation levels are statistically different (P-value) and the absolute degree of difference between DNA methylation levels (∆Meth), which can be used as an indicator of biological significance2-4. The most basic method to calculate ∆Meth is to subtract the DNA methylation levels of CpG sites between samples, where a larger ∆Meth value assumes a more biologically significant change. However, due to multiple factors, including that DNA methylation levels can only range between 0.0-1.0 and CpG sites are naturally methylated at specific levels (which are required for cellular differentiation and proper cellular function)5-7, the assumption underlying the usage of ∆Meth may not always hold true.

Modifications of the epigenome, including DNA methylation alterations at CpG dinucleotides, have been found to occur in many diseases and are one of the earliest and most frequent aberrant changes in cancer8,9. It is well-established that many of the largest gains of DNA methylation (hypermethylation) in cancer are in CpG Islands, which are largely unmethylated in normal cells and have been shown to lead to gene silencing of tumor suppressors, and thus, contributes to the basic biological understanding of the disease10-12. Early studies were perplexed about the biological importance of the loss of DNA methylation (hypomethylation), which has been found to occur in genes in confounding biological pathways for cancer progression13. However, more recent studies have shown that large global hypomethylation events at CpG sites that are normally heavily methylated induce chromosomal instability through an increased rate of DNA double-strand breaks and the reactivation of transposable elements14-16.

Conventional wisdom assumes that the observed association between the degree of ∆Meth and the native level of the CpG sites are largely due to limitations of change within the DNA methylation scale (0.0-1.0). For example, the maximum possible loss of DNA methylation of a CpG site with a native level of 0.80 is larger than that of a CpG site with a native level of 0.10. However, since DNA methylation alterations can be caused by factors that affect the whole methylome, such as age, environmental exposures and genetic alterations17-19, it is conceivable that the CpG sites of the methylome undergo a similar global rate of change relative to their native level (relative ∆Meth). For example, a 50% global loss of DNA methylation relative to the native level would change a CpG site with a native level of 0.80 to 0.40 and a CpG site with a native level of 0.10 to 0.05. Although both scenarios would cause the degree of ∆Meth to be associated with the native level, whether ∆Meth is an inherently biased approach for measuring the level of the DNA methylation alteration depends on the nature at which alterations occur: as an absolute level change or a relative level change.

If the DNA methylation alterations occur as a global rate of change relative to the native level, then it reasons that the ∆Meth calculation would be inherently biased based on the native level. Since the majority of current approaches, including the ∆ Meth calculation, assume that alterations occur as an absolute level change, the current understanding of DNA methylation would also be biased. Therefore, determining the nature at which DNA methylation alterations occur, and using an appropriate approach to measure the level of the alteration, is fundamental for the understanding of DNA methylation biology.

Here, we propose an alternative hypothesis, that the association between the degree of ∆Meth and the native levels are due to DNA methylation alterations occurring at similar global rates of change relative to the native level. Thus, we propose that the selection for large relative ∆Meth would be an unbiased approach to enrich for more biologically important alterations in comparison to those identified using ∆Meth. To test this hypothesis, we calculated and compared the features of ∆Meth and the relative ∆Meth using a comprehensive dataset of cancer tissues from The Cancer Genome Atlas (TCGA) and a large diverse dataset of normal tissues. We found multiple lines of evidence suggesting that DNA methylation alterations occur at similar global rates of change relative to the native level and that the selection for large ∆Meth skews for CpG sites based on their native level. We also found that the relative ∆Meth can more accurately identify DNA methylation signatures and can uniquely detect DNA methylation signatures caused by multiple factors including genetic alterations, smoking and age, in comparison to ∆Meth. More significantly, we found that large relative ∆Meth identified genes in pathways that are more biologically important to cancer development and progression in comparison to ∆Meth. This result was most striking in the hypomethylation direction. While we found that the CpG sites with large ∆Meth in the hypomethylation direction were enriched in pathways important for brain function, a commonly identified and confounding pathway for cancer development, we found that the CpG sites that have large relative ∆Meth were enriched in multiple pathways important for tumor progression, including pathways for cell adhesion, increased metabolism, cell signaling and immune activation. Thus, we conclude DNA methylation alterations generally occur at a global rate of change relative to the native level and that the relative ∆Meth calculation is an unbiased approach that can identify more biologically important alterations, compared to ∆Meth.

Results
Although it is commonly observed that there is an association between the degree of ∆Meth and the native level, why this association occurs remains unknown. To determine whether this association is caused by alterations occurring as a global rate of change relative to the native level, we calculated the ∆Meth and relative ∆Meth at CpG sites with varying native levels in a large and comprehensive dataset of cancers. To highlight the difference between ∆Meth and relative ∆Meth, an illustration showing examples of calculating, comparing, and the analysis of the two approaches at CpG sites of varying native DNA methylation levels is shown in Figure 1.

Global methylation changes in cell lines
To demonstrate how the association between the degree of ∆Meth and the native DNA methylation level can be caused by alterations occurring at a similar global rate of change relative to the native levels, we investigated the 5-azacytidine induced global loss of DNA methylation in multiple cell lines (3BKO and 3ABDKO) that have deficiencies in regulating DNA methylation because of a deletion of one or multiple DNA methyltransferase (DNMT) genes17. Studies have shown that 5-azacytidine causes global DNA methylation loss by incorporating into the DNA, replacing the cytosine nucleotides and inhibiting DNMTs20. To investigate the association between the degree of ∆Meth and the native levels, we subset the CpG sites of the Illumina Human Methylation 450K array datasets (3BKO; GSE51815 and 3ABDKO; GSE68344) into four even quantile intervals based on the untreated cells’ DNA methylation level (native level). We considered the DNA methylation level of the treated cells to be changed if the difference between the untreated and treated cells (∆Meth) was ≥0.01. We found in both cell lines after 5 days of 5-azacytidine treatment, in the hypomethylation direction, that the degree of ∆Meth was positively associated with the native level, and the largest absolute degree of ∆Meth occurred at CpG sites with the highest native levels (Figure 2a and c). Conversely, we found that the relative ∆Meth (∆Meth/native DNA methylated level) was very similar across CpG sites with varying native levels (Figure 2b and d). The median global loss of DNA methylation was 35-37% relative to the native level in the treated 3BKO cells and 30-32% in the treated 3ABDKO cells. These results illustrate how a similar global rate of DNA methylation change relative to the native level can lead to an association between the degree of ∆Meth and the native levels.

Global DNA methylation changes in tumors
Next, we investigated whether CpG sites with varying native levels in tumor tissues also undergo similar global rates of DNA methylation change relative to the native level. While it is ideal to calculate the native level from matched tumor adjacent normal tissues, matched tissues are seldomly available. Thus, we used the OutlierMeth package, which was designed to overcome the limited availability of normal tissues21, to calculate the native DNA methylation levels at the upper and lower 99th quantiles from Illumina Human Methylation 450K array data derived from 25 normal tissue types consisting of 1,991 samples (Supplementary Data Table 1). We used functions in the OutlierMeth package to calculate ∆Meth and relative ∆Meth. We considered the cancer DNA methylation level to be aberrant if it was greater than the upper native level or less than the lower native level. The equation (|aberrant cancer DNA methylation level – native DNA methylation level|) was used to calculate the ∆Meth in the hypomethylation and hypermethylation directions. To calculate the relative ∆Meth in the hypomethylation and hypermethylation directions, we used the equations (∆Meth/native DNA methylated level) and (∆Meth/(native DNA unmethylated level)), respectively. To investigate whether the native levels associate with the degree of the alteration, we subset the CpG sites into four even quantile intervals based on their native level (Supplemental Figure 1a-c). Because of their known significance, we also subset the CpG sites based on their genomic region (Supplemental Figure 1d-h).

We then calculated the hypomethylation and hypermethylation ∆Meth and relative ∆Meth from TCGA Illumina Human Methylation 450K array datasets, which included 31 tumor types and 8,551 samples (Supplementary Data Table 2). To ensure all tumor types carried equal weight, each tumor type was analyzed separately prior to the results being combined. Density plots were then made to assess the levels and distribution profiles of the median ∆Meth at CpG sites with varying native levels of DNA methylation. To reduce the effects caused by differences in the frequency of the aberrant events, we grouped the CpG sites based on the frequency of the alteration (number of samples with an aberrant DNA methylation level/total number of samples) before the median level was calculated (aggregated median). As expected, with the full set of CpG sites, in both the hypomethylation and hypermethylation directions, we found that the degree of ∆Meth was associated with the native DNA methylation level (Figure 3a and b). Specifically, we found that, in the hypomethylation direction, the median ∆Meth increased as the native levels increased and, in the hypermethylation direction, the median ∆Meth decreased as the native levels increased. This trend was also observed in the individual tumor types and genomic regions, with the largest hypomethylation alterations occurring in OpenSeas and the largest hypermethylation alterations occurring in Islands (Supplemental Figure 2 and Supplemental Figure 3a and b).

Next, we calculated the relative ∆Meth in the hypomethylation and hypermethylation direction. We found that the only obvious difference was at CpG sites with the highest native DNA methylation levels in the hypomethylation direction, which had much smaller relative ∆Meth (Figure 3c). Later, we will show that this result is largely due to differences in the frequency of the aberrant events. All other relative ∆Meth distribution profiles in the hypomethylation and hypermethylation direction were remarkably similar, with the median global gain and loss of DNA methylation being 22-25% relative to the native DNA methylation level. Similar relative ∆Meth distribution profiles were also found when investigating the individual tumor types and the different genomic regions (Supplemental Figure 3c and d and Supplemental Figure 4). Again, we found that the largest relative ∆Meth hypomethylation alterations occurred in OpenSeas and the largest hypermethylation alterations occurred in Islands.

We then investigated whether the overlapping profiles of the relative ∆Meth were the result of random level changes that are simply limited by the DNA methylation scale (0.0-1.0), which is a common assumption used to explain the association between the degree of ∆Meth and the native levels. After we made a dataset of randomly generated values equivalent to the size of the compiled TCGA cancer dataset (N=8,551), we found that the median relative ∆Meth was exactly 0.50 at CpG sites with varying levels of native DNA methylation (Supplemental Figure 5a and b). We found that the TCGA relative ∆Meth distribution profiles were statistically smaller (Wilcox rank sum test P<0.001) than the random dataset, providing evidence that the overlapping distribution profiles of the TCGA data are not the results of random DNA methylation level changes.

Effects of events on methylation levels
Since the relative frequency of the alteration in tumors can affect the degree of ∆Meth and relative ∆Meth, which would affect the distribution profiles (Figure 3a-d), we controlled for this variable by calculating the rate between the degree of the alteration and the relative frequency of the alteration. To ensure that all tumor types carried equal weight, each tumor type was analyzed separately before the results were combined. Histograms were made to visualize the frequency of data with respect to the relative frequency of the alteration (Figure 3).

From the TCGA dataset, as expected, we found that the aggregated median ∆Meth levels increased as the relative frequency of the alteration increased (coefficient P<0.001), confirming the positive association between the frequency of the alteration and degree of ∆Meth (Figure 3e and f). We also found that the rate of change was associated with the native level. In the hypomethylation direction, we found that the rate of change was largest at the highest native levels. And in the hypermethylation direction, the rate of change was largest at the lowest native levels. Thus, we conclude that there is an association between the degree of ∆Meth and both the relative frequency of the alteration and the native level. These two associations were also found when each tumor type was analyzed independently (Supplemental Figure 6a and b).

Next, we calculated the rate of change between the relative ∆Meth and the relative frequency of the alteration. Similar to the ∆Meth, we found that the relative ∆Meth increased as the relative frequency of the alteration increased (coefficient P<0.001) in both the hypomethylation and hypermethylation directions (Figure 3g and h). However, unlike ∆Meth, we found that the relative ∆Meth rate of change was similar, ranging between 1.7-2.1, and was not associated with the native levels. Interestingly, the rate of change of the CpG sites with the highest native levels in the hypomethylation direction is lower than other rates of change, suggesting that the small distribution profile of these CpG sites (Figure 3c) are largely due to differences in the relative frequency of the alterations (histogram of Figure 3g). The variability of the rate of change in the different tumor types also seems to be quite similar throughout the varying native DNA methylation levels (Supplemental Figure 6c and d).

Once again, we investigated whether the relative ∆Meth similarities are the result of random level changes being confined by the DNA methylation scale. Using the previously described dataset of randomly generated values, we calculated the rate of change between the relative ∆Meth and the relative frequency of the alteration. We found that the randomly generated dataset produces non-significant coefficients (P>0.05) at CpG sites with varying levels of native DNA methylation (Supplemental Figure 5c and d), again confirming that the similar rates of change between the relative ∆Meth and the relative frequency of the alteration are not due to random DNA methylation level changes.

Detection of DNA methylation signatures
If DNA methylation alterations occur as a global rate of change relative to the native level, we would expect that the relative ∆Meth calculation would be more accurate than ∆Meth in detecting DNA methylation signatures at CpG sites with native levels close to the end of the DNA methylation scale (where the ∆Meth would be the smallest). To verify this prediction, we generated DNA methylation signatures in a dataset of 1000 randomly selected cancer samples from the compiled TCGA dataset. The DNA methylation signatures were created by altering a percentage (10, 20, 30, or 40%) of the CpG sites by a 10% gain or loss of DNA methylation relative to the cancer samples’ native level. We also altered 10% of the CpG sites by a percentage (10, 20, 30, or 40%) of gain or loss of DNA methylation relative to the cancer samples’ native level. As predicted, using the OutlierMeth package, we found that the aggregated median relative ∆Meth had higher accuracy in detecting hypomethylation signatures at CpG sites with low native levels and hypermethylation signatures at CpG sites with high native levels in comparison to ∆Meth (Supplemental Figure 7 and 8).

To determine if this phenomenon can also be found in TCGA cancer datasets, we calculated the aggregated median ∆Meth and relative ∆Meth to investigate multiple factors including genetic alterations, smoking, and age. Specifically, we investigated the well-characterized CpG Island methylator phenotype (CIMP) signature in colorectal cancers, smoking signature in lung cancers and the signature associated with age in breast cancers. In brief, colorectal CIMP tumors commonly carry a BRAF mutation and are characterized by having larger levels of hypermethylation in genomic CpG Islands at CpG sites with low native levels in comparison to non-CIMP tumors22-24. As expected, we found that high CIMP (H-CIMP) tumors had larger degrees of ∆Meth and larger relative ∆Meth levels of hypermethylation at CpG sites with low native levels in comparison to the non-H-CIMP tumors (Figure 4a). However, at CpG sites with the highest native levels, only the relative ∆Meth detected a statistical difference between H-CIMP and non-H-CIMP tumors (Figure 4b). Similarly, only the relative ∆Meth could detect hypomethylation signatures at CpG sites with the lowest native levels and hypermethylation signatures at CpG sites with the highest native levels between smokers and non-smokers in lung cancer (Figure 4c and d). This was also observed in the signatures associated with age in breast cancer (Figure 4e and f).

Identifying biologically important sites
Since the evidence overwhelmingly suggests that DNA methylation alterations occur as a global rate of change relative to the native DNA methylation level, we investigated whether selecting for large relative ∆Meth would enrich for a different and more biologically important set of CpG sites in comparison to ∆Meth. We first combined all the cancer samples into a single dataset and selected the 10% most frequently altered CpG sites in both the hypomethylation and hypermethylation directions (Fisher’s Exact P≤5.25E-78). We found that the 10% most frequently altered CpG sites had native DNA methylation levels in all quantile intervals.

From the 10% most frequently altered CpG sites, we further selected the top 10% CpG sites with the largest median ∆Meth and relative ∆Meth. In the hypomethylation direction, the largest 10% ∆Meth selected 0 CpG sites with native levels between 0.0-0.25 (Figure 5a left side). Similarly, in the hypermethylation direction, 0 CpG sites with native levels between 0.75-1.0 were selected (Figure 5a right side). In both directions, a clear association between the native levels and the degree of ∆Meth was observed (Figure 5a). Conversely, the largest 10% relative ∆Meth selected CpG sites with native DNA methylation levels in all quantile intervals in both the hypomethylation and hypermethylation directions (Figure 5b). Overall, we found that only 20% of the hypomethylation CpG sites and 29% of the hypermethylation CpG sites were selected by both ∆Meth and relative ∆Meth (Figure 5c and d and Supplementary Data Table 3 and 4). Interestingly, while we found that the CpG sites of the largest 10% ∆Meth and relative ∆Meth at different quantile intervals were mostly in OpenSeas in the hypomethylation direction, in the hypermethylation direction, we found the proportions of the genomic regions varied at different quantile intervals (Supplemental Figure 9).

To investigate the biological features of the largest 10% relative ∆Meth and ∆Meth CpG sites, we measured the correlation between DNA methylation and gene expression levels derived from RNA sequencing data of 31 tumor types and 7,910 samples. Genes with statistically significant correlation (P<0.001) were identified in the hypomethylation and hypermethylation directions by ∆Meth and relative ∆Meth. We found that ∆Meth and relative ∆Meth identified a similar number of genes with statistically significant correlation and that the strength of the correlations were also similar (Supplemental Figures 10-13). Although the largest 10% relative ∆Meth selected for CpG sites with native DNA methylation levels in all quantile intervals, no genes with statistically significant correlation were selected at the highest quantile interval in the hypomethylation direction. Overall, only a minority of genes with statistically significant (P<0.001) and strong correlations (rho  |0.2|) were identified by both ∆Meth and relative ∆Meth in the hypomethylation and hypermethylation directions, 25% and 28% respectively (Supplementary Data Tables 5 and 6 and Supplemental Figure 11).

From the Gene Set Enrichment Analysis (GSEA), using the genes with statistically significant and strong correlations, we found that the genes uniquely identified by the 10% largest ∆Meth were enriched for pathways important for the midbrain (Figure 5e and Supplementary Data Table 7 and 9). Similar brain and neuron pathways were also identified when investigating the statistically significant genes uniquely identified in the four most frequent tumor types (breast, prostate, lung and colorectal) independently (Supplemental Figure 14a). Surprisingly, in the combined cancer dataset, the genes that were uniquely identified by the 10% largest relative ∆Meth showed enrichment for multiple pathways with known importance for tumor progression including increased metabolism, cell signaling, immune activation, and cell adhesion (Figure 5e and Supplementary Data Table 8)25-28. Interestingly, when investigating the relative ∆Meth genes uniquely identified in the four most frequent tumor types independently, we found that the pathways identified using hypomethylation alterations were largely tumor specific, although similar types of pathways, such as cell signaling and immune activation, were found in multiple tumor types (Supplemental Figure 14b). In the hypermethylation direction, we found that the genes identified by both ∆Meth and relative ∆Meth were enriched for pathways important for polycomb regulated stemness (Figure 5f and Supplementary Data Tables 10-12). Unlike the alterations identified in the hypomethylation direction, in the hypermethylation direction, we found a large degree of shared polycomb-regulated stemness pathways when investigating breast, prostate, lung and colorectal cancers independently (Supplemental Figure 14c and d).

Investigating and comparing ∆M-value hypomethylation alterations
Since the relative ∆Meth approach identified more biologically significant hypomethylation alterations in comparison to ∆Meth, we investigated whether other approaches that can select for hypomethylation alterations at CpG sites that have low native levels (a distinct feature of the relative ∆Meth approach) can also be used to enrich for genes in pathways important for tumor development. To accomplish this, we used the well-established ∆M-value approach29, a method that increases the resolution of DNA methylation differences as the DNA methylation levels near 0.0 and 1.0. We compared the top 10% largest hypomethylation alterations of the ∆M-value to the top 10% largest hypomethylation alterations of relative ∆Meth. As expected, we found that the largest ∆M-value selected CpG sites with native DNA methylation levels in all quantile intervals (Figure 6a and Supplemental Figure 15).

After selecting for genes with statistically significant and strong correlation between gene expression and DNA methylation levels, we found that there were more overlapping genes identified between the largest M-value and relative ∆Meth (35%) in comparison to the overlapping genes between the largest ∆Meth and relative ∆Meth (25%) (Figure 6b and Supplemental Figure 11c). However, similar to ∆Meth, from the GSEA, we found that the genes identified using the ∆M-value approach were in pathways important for the midbrain (Figure 6c and Supplementary Data Tables 13 and 15). Once again, the genes uniquely identified by relative ∆Meth were enriched for pathways with known biological importance to cancer development (Figure 6c and Supplementary Data Table 14). Thus, we conclude that the ability of relative ∆Meth to enrich for biologically significant alterations in the hypomethylation direction is not simply due to the selection of CpG sites that have low native levels.

Hypomethylation alterations in distinct genomic compartments
Since the loss of DNA methylation at different genomic compartments has been shown to have distinct biological consequences for tumor development30-33, we investigated whether large, uniquely identified, relative ∆Meth hypomethylation alterations, that had statistically significant and strong correlation with gene expression, occur in different genomic compartments compared to large ∆Meth alterations. Since the Illumina 450K DNA methylation array does not provide the genomic coverage required to accurately call different genomic compartments, we first downloaded the genomic positions of previously identified genomic compartments, partially methylated domains (PMDs), highly methylated domains (HMDs), and genomic domains that are neither PMDs or HMDs (PMDs/HMDs), which were found to be common in multiple tissue types33. Since the genomic domains can be tissue specific, we identified the array CpG sites that showed concordance between the DNA methylation levels of normal-adjacent tissues and the previously identified genomic compartments. Because of their frequency and the availability of normal-adjacent tissues, we focused our investigation on breast cancer and colorectal cancer. We defined the CpG sites of the normal-adjacent tissues to be concordant with the previously identified genomic compartments if the median DNA methylation level was > 60% in common HMDs, between 30% and 60% in PMDs, and < 30% in neither PMDs/HMDs.

After we selected the top 10% largest ∆Meth and relative ∆Meth alterations from the CpG sites that had concordance with the genomic compartments, we found that uniquely identified alterations by ∆Meth occurred most frequently in common HMDs in both breast and colorectal cancer tissues (Figure 7a and b). We also found that the unique relative ∆Meth alterations occurred more frequently in common PMDs and in neither PMDs/HMDs in comparison to ∆Meth. Next, we investigated whether the alterations in the three genomic compartments enriched for genes in different biological pathways. Similar to our previous findings, the genes uniquely identified by the top 10% largest ∆Meth in both breast and colorectal cancers were largely enriched for pathways which are confounding for tumor development (Figure 7c and e and Supplementary Data Tables 16 and 19). Interestingly, we found that the unique genes in common PMDs identified in the top 10% largest relative ∆Meth in both breast and colorectal cancers were enriched for both pathways important for cancer development and brain development, while the genes in both common HMDs and in neither PMD/HMDs were largely enriched in pathways known to be important for tumor development (Figure 7d and f and Supplementary Data Tables 17, 18 and 20-22). Interestingly, in both breast and colorectal cancers, we found that alterations in common HMDs were largely enriched for genes in immune cell pathways. In conclusion, our data supports the findings that DNA methylation alterations in different genomic compartments have different biological roles for tumor development and that these differences are more obvious when using the relative ∆Meth approach in comparison to the ∆Meth approach.

Methods
Cell lines
For the cell line DNA methylation data, which were measured using the Illumina Human Methylation 450K array, the probe-specific native DNA methylation levels were derived from the levels of the untreated cells. We compared the cell lines 3BKO and 3ABDKO before and after 5 days of 5-azacytidine treatment. We considered the DNA methylation level to be changed if the ∆Meth was ≥0.01. For probes with DNA methylation level changes in the hypomethylation direction, the ∆Meth (untreated cell DNA methylation level– treated cell DNA methylation level) and relative ∆Meth (∆Meth/native DNA methylated level) were calculated. The R function density was used to make the density plots.

Sample selection
To determine the upper and lower probe-specific native DNA methylation levels for the normal tissues, we used the OutlierMeth package (https://github.com/bdowns4/OutlierMeth) at the significance levels P0.01 using all of the available normal tissues. A detailed description of sample acquisition and usage of the OutlierMeth package can be found in our previously published manuscript21. In short, all normal tissues available in the OutlierMeth package were included in this study (25 normal tissue types N=1,991). All TCGA cancer DNA methylation data used in this study was derived using the Illumina Human Methylation 450K array and was comprised of 31 tumor types and 8,551 samples. All paired TCGA cancer RNA-Seq and Illumina Human Methylation 450K array data was used and was comprised of 31 tumor types and 7,910 samples.

Calculating ∆Meth, relative ∆Meth, and ∆M-value for tumor tissues
We considered the cancer tissue DNA methylation level to be aberrant if it was less than the lower native DNA methylation level or greater than the upper native DNA methylation level. We used the equation (|aberrant cancer DNA methylation level – native DNA methylation level|) to calculate the ∆Meth in the hypomethylation and hypermethylation direction. We used the equations (∆Meth/native DNA methylation level) and (∆Meth/(native DNA unmethylated level)) to calculate the relative ∆Meth in the hypomethylation and hypermethylation direction respectively. The OutlierMeth package function deltMeth and relMeth was used to calculate all ∆Meth and relative ∆Meth. The M-value was calculated using the equation log (Beta/1-Beta) and the hypomethylation ∆M-value was calculated using the equation (M-value – M-value ).

Distribution profiles of ∆Meth and relative ∆Meth
To ensure that each tumor type carried equal weight, we analyzed the tumor types separately before the results were combined. To limit the effect that the relative frequency of the alteration (number of samples with an aberrant DNA methylation level/total number of samples) has on the distribution profiles of ∆Meth and relative ∆Meth, the R function aggregate was used to group the CpG sites based on the relative frequency of the aberrant events and to calculate the median value for each group (aggregated median). The R function density was used to make the density plots. The R function Wilcox.test (two-sided) was used to compare the distributions of the ∆Meth and relative ∆Meth levels. All distributions of the ∆Meth and relative ∆Meth levels were found to be statistically different (Wilcoxon rank sum test P<0.001) unless described as not statistically different (NS) (Wilcoxon rank sum test P>0.001). The R function sample was used to randomly generate DNA methylation levels for the random DNA methylation level dataset. The R package pheatmap was used to generate the heatmaps36.

The rate between the ∆Meth and relative ∆Meth and the relative frequency of the alteration
The aggregated median ∆Meth and relative ∆Meth were calculated, as described above, separately for the 31 tumor types before the results were combined. The R function hist was used to generate the histograms showing the frequency of the data based on the relative frequency of the alteration. The R functions lm was used to calculate the slopes and the function coef was used to calculate the coefficient P-values. All slopes were found to be statistically different in comparison to each other (z-statistic P<0.001) unless described as not statistically different (NS) (z-statistic P>0.001). To calculate the z-statistic we used the following equation z = (Coeff -Coeff )/√(Se 2-Se 2) 37. The R function pnorm was used to calculate the P-value using the z-statistic.

Cancer signatures
The R function sample was used to randomly select 1000 cancer samples from the combined TCGA dataset. We also used the R function sample to randomly select the CpG sites that were altered for the generated DNA methylation signatures. In total, we generated 8 test conditions by altering a range of CpG sites (10, 20, 30, or 40%) by a percentage (10, 20, 30, or 40%) of gain or loss of DNA methylation relative to the cancer samples’ native level. We generated n=10 DNA methylation signatures for each of the test conditions. The Wilcoxon rank sum test was used to determine whether the aggregated median ∆Meth and median relative ∆Meth values were statistically different (P<0.001) compared to the native levels of the selected 1000 TCGA cancer samples.

Colorectal tumors (COADREAD) from TCGA were considered to be high-CIMP (H-CIMP) if they carried a BRAF mutation or were determined H-CIMP in the TCGA manuscript24; (n=126). Non-H-CIMP tumors did not have a BRAF mutation and were not previously determined as H-CIMP; (n=503). Lung cancer patients (LUAD) from TCGA were considered a smoker if they had >5 pack years. The sample size for LUAD smokers n=343 and non-smokers n=166. Breast cancer patients from TCGA (BRCA) were divided into two groups, ≤50; (n=331) and >50 age; (n=766) at primary diagnosis.

The top 10% most altered CpG sites identified by ∆Meth and relative ∆Meth
The R function quantile was used to select the top 10% most frequently altered DNA methylation CpG sites and the CpG sites with the 10% largest ∆Meth and relative ∆Meth. The R function fisher.test (two-sided) was used to calculate the P-values for the altered CpG sites between the compiled TCGA cancer dataset (n=8,551) and the normal tissue dataset (n=1,991). The package eulerr was used to make the Venn diagrams. The R function cor was used to calculate the rho correlation values and spearman coefficient P-values. The Gene Set Enrichment Analysis (GSEA) was performed on the GSEA website (https://www.gsea-msigdb.org). The per-gene methylation information was derived from the 450k annotation file downloaded from the Illumina website. The common PMDs, HMDs and neither PMDs/HMDs genomic locations33 were downloaded from the UCSC Genome Browser website. The CpG sites of the 450K array from the normal-adjacent tissue was defined as concordant with the genomic compartment if the median DNA methylation level was > 60% in common HMDs, between 30% and 60% in PMDs, and < 30% in neither PMDs/HMDs.

Statistics & reproducibility
No statistical method was used to predetermine sample size. No data were excluded from the analyses. The R function sample was used to randomly select the CpG sites that were altered for the generated DNA methylation signatures. The experiments analyzing data from cancer tissues were not randomized. The investigators were not blinded to allocation during experiments and outcome assessment.

Data availability
3BKO and 3ABDKO cell line data was downloaded from the Gene Expression Omnibus (GEO) database (GSE51815 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc= GSE51815) and GSE68344 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc= GSE68344)). The TCGA tumor DNA methylation data, mutational data, RNA sequencing data and clinical data were downloaded from the Firebrowse site (http://firebrowse.org).

Code availability
The OutlierMeth package is on GitHub (https://github.com/bdowns4/OutlierMeth)38 and the custom code used for the analyses in this study is available upon reasonable request. The R version 4.2.3 was used for the analysis in this manuscript.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "In this manuscript, Bradley and colleagues propose the use of relative delta beta values as an alternative to absolute difference in delta beta values as a metric to measure differences between CpG methylation levels in aberrant vs. native conditions. The main working hypothesis supporting the use of relative delta beta values is that DNA methylation alterations occur as global rate of change rather than an absolute change relative to the native level and this latter aspect could potentially lead to measurement biases when using delta beta values in sites characterized by extremely low or high native beta values.\nWhile the adoption of a relative measure could certainly be helpful, especially for signature and biomarkers identification, it should be used with caution, especially when interpreting the biology of the phenomenon under study. In particular, the authors, by comparing results obtained using relative delta beta values with those obtained with delta beta, report to have provided evidence that methylation occur as a global rate of change relative to native levels, but do not provide a biological model or experimental evidence explaining how and in which conditions this can happen: contrary to gene expression, CpG methylation is a binary phenomenon in a cell (i.e. a given site can be methylated or un-methylated), thus, under this terms, working with relative rate of change would imply a model where methylation in a particular site can be found altered in a cell depending on the methylation state in other cells, that could certainly be true in some cases (e.g. under selective pressure) but not in general.\nThe manuscript needs to be revised to improve its readability as many key points are difficult to interpret.\nAlso, the methods section needs to be improved as it lacks the necessary level of detail in some parts.\nPlease find below a detailed list of comments.\nMajor points:\n1. Lines 67-68: “relative ∆Meth would be an unbiased approach to enrich for more biologically important alterations in comparison to those identified using ∆Meth”: this is questionable, one would expect that at the cell population and phenotypic level a 2 fold relative change in a locus from 0.4 to 0.8 would determine an overall effect remarkably different the one linked to a 2 fold relative change from 0,04 to 0,08.\n2. Lines 101-103: “we subset the CpG sites of the Illumina Human Methylation 450K array datasets (3BKO; GSE51815 and 3ABDKO; GSE68344) into four even quantiles based on the untreated cells’ DNA methylation level”: probes on the Infinium array are not independent and not uniformly distributed along the genome, hence beta values are not independent each other. Did the authors check for unbalance with respect to genomic regions distribution in the defined quantiles intervals?\n3. Line 166: The “relative frequency of alteration” should be introduced in the results and its exact calculation should be detailed in the methods.\n4. Lines 248-250: The text reports that a similar number of genes with statistically significant association for the absolute and relative delta Meth has been found, while in Figure S9 Panel C there is no plot for the association between relative delta methylation and gene expression levels for native levels 0.75-1.0” in the hypomethylation direction. Please comment this aspect.\n5. Lines 358-361: “The R functions lm was used to calculate the slopes and the function coef was used to calculate the coefficient P-values. All slopes were found to be statistically different in comparison to each other (z-statistic P<0.001) unless described as not statistically different (NS) (z-statistic P>0.001)”. Please provide details on the statical test used to compute significance of linear models’ coefficients and to compare the slopes.\n6. Lines 365-366: “We generated n=10 DNA methylation signatures for each test condition.”. In the result section and in Supplementary figures S7-S8 only 8 signatures are reported (4 gain and 4 loss) for each condition. Could the authors please make clear where and how these 10 signatures are used?\n7. Supplementary figures S7-S8 legend: Please explicitly state the comparison for the Wilcoxon test.\n8. Code availability: Code availability: The presented work its entirely based on in silico analyses and the code used for the analyses is not of marginal importance. It would be important to provide direct access to the code that was used to generate the results. Is there any particular reason why the authors have decided to release the code only upon reasonable request?\nMinor points:\n1. Lines 103, 128: The usage of the term quantile to indicate intervals is confusing, I suggest using the term quantile intervals.\n2. Figure 3 Legend: “Histograms show the frequency of data.” Please provide a more detailed description for the reported frequencies.\n3. Supplementary Figure S6: The association between Delta Meth and relative frequency of alteration is not present in all tumors, in particular for the native levels 0.75-1.0 there is one tumor for which the association has the opposite direction. Is there any particular reason for this?\n4. Lines 335-336: “We considered the cancer tissue DNA methylation level to be aberrant if it fell outside of the native DNA methylation level.”. Please provide a more specific definition.\n5. Lines 336-337: “We used the equation (aberrant cancer DNA methylation level – native DNA unmethylated level) and (aberrant cancer DNA methylation level – native DNA methylation level) to calculate the ∆Meth in the hypomethylation and hypermethylation direction respectively”. Shouldn’t it be the opposite for the hypomethylation direction (i.e. native DNA unmethylated level - aberrant cancer DNA methylation level) ?\n6. Lines346-347: Please provide more details on the aggregation of sites based on the relative frequency of aberrant events and, in particular, on how the grouping variable is defined.\n7. Supplementary figures S7-S8 legend: “in the hypomethylation direction a measured”. Typo here.", "reviewer_2": "Downs et al. propose a new metric for quantifying changes in DNA methylation. Instead of the standard ΔMeth (tumour – normal, or treatment – control), they introduce relative ΔMeth, which divides ΔMeth by the baseline methylated level when the change is a loss (hypomethylation) or by the baseline unmethylated level when the change is a gain (hypermethylation).\nUsing relative ΔMeth, the authors show that methylation changes appear more consistent across CpGs with different starting values, whereas absolute ΔMeth is strongly influenced by the native level. They argue that relative ΔMeth is a less biased way to detect biologically meaningful alterations. Gene-set enrichment analysis (GSEA) of the CpGs with the highest relative ΔMeth in TCGA tumours highlights cancer-specific pathways, whereas the same analysis based on absolute ΔMeth points to broader, less disease-specific pathways.\nThe core idea of treating the mean relative ΔMeth as background and identifying loci whose relative ΔMeth exceeds that background is potentially important.\nHowever, we believe several key issues need to be addressed (details below).\nA central assertion—“DNA-methylation alterations occur as a global rate of change relative to the native level”—is problematic. Although the net gain or loss is, on average, proportional to the starting level, the apparent consistency arises from averaging highly heterogeneous effects with distinct biological implications. Numerous studies have shown that DNA-methylation changes cluster in specific compartments: e.g. in cancer, focal island hypermethylation and large-scale hypomethylation tend to occur in late-replicating, lamina-associated regions (also known as PMDs; see doi: 10.1038/nsmb.2518, 10.1038/nrm.2015.31, 10.1038/ng.969). Hypomethylation of PMDs, often age-associated, also contributes to tumorigenesis (doi: 10.1038/s41588-018-0073-4). Averaging all Open-Sea probes therefore yields a figure that does not represent either highly methylated domains (HMDs) or PMDs.\nIndeed, the first experiment presented, 5-AZA treatment, may be the only one that produces a true genome-wide loss: under 5-AZA, virtually every Open-Sea CpG loses methylation by a similar relative ΔMeth, whereas in cancer only the Open-Sea CpGs within PMDs lose methylation.\nSpecific recommendations\na) Acknowledge prior evidence that methylation loss is compartment-specific, not uniform.\nb) Re-evaluate ΔMeth and relative ΔMeth separately within late-replicating / lamina-associated domains (or PMDs / hypomethylated blocks) to test whether the new metric captures biologically meaningful changes there.\nThe authors claim that relative ΔMeth identifies more relevant CpGs, based on GSEA of the top 10 % most variable sites. This apparent improvement likely reflects the inclusion of low-baseline hypomethylated CpGs and high-baseline hypermethylated CpGs. Demonstrating functional importance requires further evidence, e.g. by correlating methylation and expression specifically for this subset.\nAn established alternative to ΔMeth is ΔM, where M = logit(β). Is relative ΔMeth superior to ΔM for highlighting biologically meaningful sites? A direct comparison would strengthen the paper.\na) Baseline definition. The “native” baselines appear to pool all normal tissues. Because DNA methylation is strongly tissue-specific, this will preferentially select the small subset of CpGs that are constitutively high or low in everytissue. Using matched normal tissue would be more appropriate.\nb) CpG-to-gene mapping. For DNAm–expression correlations, how is per-gene methylation derived? Promoter methylation has very different effects from gene-body methylation; the analysis should stratify by genomic context."}, "author_response": "We would like to thank Reviewer 1 for all of their comments and suggestions. In our resubmitted\nmanuscript, we made many changes throughout our article to improve the readability and have included\nnew Supplemental Figure S11-13 to more clearly explain our findings. We believe that these modifications\nhave improved the clarity of our manuscript significantly.\nPlease find below a detailed list of comments.\nMajor points:\n1. Lines 67-68: “relative ∆Meth would be an unbiased approach to enrich for more biologically important\nalterations in comparison to those identified using ∆Meth”: this is questionable, one would expect that at the\ncell population and phenotypic level a 2 fold relative change in a locus from 0.4 to 0.8 would determine an\noverall effect remarkably different the one linked to a 2 fold relative change from 0,04 to 0,08.\nWe agree with the reviewer that there is likely a point at which the absolute level change would be\nbiologically informative. Since the degree of DNA methylation level change does not provide information\nregarding whether, and at which point, gene expression level undergoes a biologically significant change,\nwe focused on the genes that showed a statistically significant correlation between expression levels and\nDNA methylation levels.\nIn an attempt to understand whether the scenario described above is having a large effect on the ability\nof relative ∆Meth to enrich for important alterations, we included new Supplemental Figures S11-13,\nshowing the number of genes and the levels of correlation between gene expression and DNA methylation\nat CpG sites at different quantile intervals.\nOverall, we found that the relative ∆Meth was identifying a similar number of genes that had similar\nstrengths of correlation between expression levels and DNA methylation levels in comparison to ∆Meth.\nThus, we conclude that the scenario described above is rare and is having a relatively small effect on the\nability of relative ∆Meth to enrich for important alterations.\n2. Lines 101-103: “we subset the CpG sites of the Illumina Human Methylation 450K array datasets (3BKO;\nGSE51815 and 3ABDKO; GSE68344) into four even quantiles based on the untreated cells’ DNA methylation\nlevel”: probes on the Infinium array are not independent and not uniformly distributed along the genome, hence\nbeta values are not independent each other. Did the authors check for unbalance with respect to genomic regions\ndistribution in the defined quantiles intervals?\nWe would like to thank the reviewer for raising this important point. We have investigated the\ndistributions of CpG native levels across the genomes of both cell lines. While we found that the probes\nwere not uniformly distributed across the chromosomes, we did not find noticeable levels of unbalance\nwith respect to the genomic regions and the native levels in either 3BKO or 3ABDKO cell lines. Plots of\nthe CpG native levels across the whole genome and chromosome 1 of cell lines 3BKO and 3ABDKO are\nshown below.\n3. Line 166: The “relative frequency of alteration” should be introduced in the results and its exact calculation\nshould be detailed in the methods.\nWe apologize for the lack of clarity and we have included the calculation for the relative frequency of\nalterations in the Results and Methods sections. (Lines 134-135 and 413-414).\n4. Lines 248-250: The text reports that a similar number of genes with statistically significant association for the\nabsolute and relative delta Meth has been found, while in Figure S9 Panel C there is no plot for the association\nbetween relative delta methylation and gene expression levels for native levels 0.75-1.0” in the hypomethylation\ndirection. Please comment this aspect\nWe would like to apologize for the lack of clarity regarding this point. We have included a new\nSupplemental Figure S11 to show more clearly the number of genes with statistically significant and\nstrong correlation between expression level and DNA methylation level.\nThe reviewer is correct in that no genes were selected from CpG sites with native levels 0.75-1.0 in Figure\nS9 Panel C (now called Figure S10 Panel C). Although the relative ∆Meth selected CpG sites in all\nquantile intervals, very few (0.1% of the selected CpG sites; Figure 5B) were from the top quantile\ninterval in the hypomethylation direction. We found that none of these CpG sites correlated with the\nexpression level of their corresponding gene, thus, there is no data in the plot for Figure S10C. We have\nincluded a comment regarding this point in the Results section. (Lines 250-252).\n5. Lines 358-361: “The R functions lm was used to calculate the slopes and the function coef was used to\ncalculate the coefficient P-values. All slopes were found to be statistically different in comparison to each other\n(z-statistic P<0.001) unless described as not statistically different (NS) (z-statistic P>0.001)”. Please provide\ndetails on the statical test used to compute significance of linear models’ coefficients and to compare the\nslopes.\nWe apologize for the lack of details and we have provided the calculation used for this statistical test in\nthe Methods section. (Lines 430-432).\n6. Lines 365-366: “We generated n=10 DNA methylation signatures for each test\ncondition.”. In the result section and in Supplementary figures S7-S8 only 8 signatures are reported (4 gain and\n4 loss) for each condition. Could the authors please make clear where and how these 10 signatures are used?\nWe would like to apologize for the lack of clarity regarding this experiment. We generated 10 DNA\nmethylation signatures for the four levels of alterations and for the four percentages of CpG sites that\nwere altered. Thus, each condition (box) in the plots of Supplementary Figures S7-S8 has 10 dots. We\nhave provided more details regarding this experiment in the Methods section. (Lines 436-439).\n7. Supplementary figures S7-S8 legend: Please explicitly state the comparison for the Wilcoxon test.\nWe apologize for this oversight. We have stated the comparison for the Wilcoxon test in the legend of\nSupplementary Figures S7 and S8.\n8. Code availability: Code availability: The presented work its entirely based on in silico analyses and the code\nused for the analyses is not of marginal importance. It would be important to provide direct access to the code\nthat was used to generate the results. Is there any particular reason why the authors have decided to release the\ncode only upon reasonable request?\nWe would like to thank the reviewer for raising this point. While this study is based on in silico analyses,\nbesides the relative ∆Meth calculation, which is based on a previously published package and the code is\nprovided, the majority of this study uses relatively simple code, such as calculating median values and\nsorting and selecting data. Although unlikely, for transparency purposes, we decided to offer the code if\nrequested.\nMinor points:\n1. Lines 103, 128: The usage of the term quantile to indicate intervals is confusing, I suggest using the term\nquantile intervals.\nWe apologize for the confusion. We have changed the term “quantiles” to “quantile intervals” throughout\nthe manuscript.\n2. Figure 3 Legend: “Histograms show the frequency of data.” Please provide a more detailed description for\nthe reported frequencies.\nWe have provided more details regarding the histograms in the legend.\n3. Supplementary Figure S6: The association between Delta Meth and relative frequency of alteration is not\npresent in all tumors, in particular for the native levels 0.75-1.0 there is one tumor for which the association has\nthe opposite direction. Is there any particular reason for this?\nWe agree that the associations that occur in the opposite direction in Supplementary Figure S6 is an\ninteresting feature. When we investigated this anomaly, we found that this cancer type (CHOL;\ncholangiocarcinoma) had relatively few alterations at the upper quantile interval in comparison to the\nother cancer types and other quantile intervals. We have concluded that this is a unique feature of CHOL\nand not an error in the analysis.\n4. Lines 335-336: “We considered the cancer tissue DNA methylation level to be aberrant if it fell outside of the\nnative DNA methylation level.”. Please provide a more specific definition.\nWe would like to apologize for the lack of clarity regarding this point, and we have provided a more\nspecific definition in the Methods section. (Lines 402-403).\n5. Lines 336-337: “We used the equation (aberrant cancer DNA methylation level – native DNA unmethylated\nlevel) and (aberrant cancer DNA methylation level – native DNA methylation level) to calculate the ∆Meth in\nthe hypomethylation and hypermethylation direction respectively”. Shouldn’t it be the opposite for the\nhypomethylation direction (i.e. native DNA unmethylated level - aberrant cancer DNA methylation level) ?\nCorrect, this is a typo. We have corrected this error in the manuscript.\n6. Lines346-347: Please provide more details on the aggregation of sites based on the relative frequency of\naberrant events and, in particular, on how the grouping variable is defined.\nWe apologize for the lack of clarity and have included a more detailed description describing the\naggregation calculation in the Methods section. (Lines 413-417).\n7. Supplementary figures S7-S8 legend: “in the hypomethylation direction a measured”. Typo here.\nWe would like to thank the reviewer for pointing out this typo and we have corrected this error in the\nmanuscript.\nReviewer #2 (Remarks to the Author):\nDowns et al. propose a new metric for quantifying changes in DNA methylation. Instead of the standard ΔMeth\n(tumour – normal, or treatment – control), they introduce relative ΔMeth, which divides ΔMeth by the baseline\nmethylated level when the change is a loss (hypomethylation) or by the baseline unmethylated level when the\nchange is a gain (hypermethylation).\nUsing relative ΔMeth, the authors show that methylation changes appear more consistent across CpGs with\ndifferent starting values, whereas absolute ΔMeth is strongly influenced by the native level. They argue that\nrelative ΔMeth is a less biased way to detect biologically meaningful alterations. Gene-set enrichment analysis\n(GSEA) of the CpGs with the highest relative ΔMeth in TCGA tumours highlights cancer-specific pathways,\nwhereas the same analysis based on absolute ΔMeth points to broader, less disease-specific pathways.\nThe core idea of treating the mean relative ΔMeth as background and identifying loci whose relative ΔMeth\nexceeds that background is potentially important.\nHowever, we believe several key issues need to be addressed (details below).\nWe would like to thank Reviewer 2 for their insightful suggestions to better understand the relative ∆Meth\napproach. To address these comments and suggestions, we have added 3 new paragraphs in the Results\nsection, 2 new paragraphs in the Discussion section, new Figure 6 and 7, new Supplemental Figure S9\nand S15, and new Supplemental Tables S13-21. Overall, we believe that these new results have dramatically\nimproved the scientific significance of our manuscript.\nA central assertion—“DNA-methylation alterations occur as a global rate of change relative to the native\nlevel”—is problematic. Although the net gain or loss is, on average, proportional to the starting level, the\napparent consistency arises from averaging highly heterogeneous effects with distinct biological implications.\nNumerous studies have shown that DNA-methylation changes cluster in specific compartments: e.g. in cancer,\nfocal island hypermethylation and large-scale hypomethylation tend to occur in late-replicating, lamina-\nassociated regions (also known as PMDs; see doi: 10.1038/nsmb.2518, 10.1038/nrm.2015.31, 10.1038/ng.969).\nHypomethylation of PMDs, often age-associated, also contributes to tumorigenesis (doi: 10.1038/s41588-018-\n0073-4). Averaging all Open-Sea probes therefore yields a figure that does not represent either highly\nmethylated domains (HMDs) or PMDs.\nIndeed, the first experiment presented, 5-AZA treatment, may be the only one that produces a true genome-wide\nloss: under 5-AZA, virtually every Open-Sea CpG loses methylation by a similar relative ΔMeth, whereas in\ncancer only the Open-Sea CpGs within PMDs lose methylation.\nSpecific recommendations\na) Acknowledge prior evidence that methylation loss is compartment-specific, not uniform.\nb) Re-evaluate ΔMeth and relative ΔMeth separately within late-replicating / lamina-associated domains (or\nPMDs / hypomethylated blocks) to test whether the new metric captures biologically meaningful changes there.\nWe would like to thank the reviewer for raising this insightful point, and we have acknowledged\ninformation regarding the loss of DNA methylation in PMDs and HMDs in our manuscript.\nIn our revised manuscript, we further investigated large ∆Meth and relative ∆Meth hypomethylation\nalterations that occur in common PMDs, common HMDs and in neither PMDs/HMDs. We found that\nboth the CpG sites selected using large ∆Meth and relative ∆Meth were mostly in PMDs (Figure 7A).\nInterestingly, we found that relative ∆Meth selected for more CpG sites in common HMDs and in neither\nPMDs/HMDs in comparison to ∆Meth (Figure 7A).\nSimilar to our previous findings, we found that the large ∆Meth in common PMDs, common HMDs, and\nneither PMDs/HMDs enriched for genes in brain pathways (Figure 7B and Supplemental Table S16-18).\nConversely, we found that CpG sites with large relative ∆Meth in common HMDs and neither\nPMDs/HMDs domains were enriched for genes in pathways with known importance for tumor\ndevelopment in comparison to those in common PMDs (Figure 7C and Supplemental Table S16-18).\nFrom these results, we confirmed the previous findings that alterations at different genomic\ncompartments have different biological functions for tumor development. Moreover, we found that these\nbiological differences are more obvious when using the relative ∆Meth approach in comparison to the\n∆Meth approach.\nDetails describing these results can be found in the Results section lines 294-314 and Discussion section\nline 355-367).\nThe authors claim that relative ΔMeth identifies more relevant CpGs, based on GSEA of the top 10 % most\nvariable sites. This apparent improvement likely reflects the inclusion of low-baseline hypomethylated CpGs\nand high-baseline hypermethylated CpGs. Demonstrating functional importance requires further evidence, e.g.\nby correlating methylation and expression specifically for this subset.\nWe have provided new correlation plots showing DNA methylation levels and expression levels for the top\ntwo CpG sites that had the strongest correlation at different native DNA methylation quantile intervals\n(Supplementary Figure S12 and S13).\nWe found that in the hypomethylation direction, CpG sites with large relative ∆Meth in the lowest\nquantile interval show a more even spread of DNA methylation levels across the beta scale in comparison\nto the CpG sites with large ∆Meth, which show large clusters of DNA methylation levels at large beta\nvalues.\nWhen we investigated the CpG sites with large relative ∆Meth in the highest quantile interval in the\nhypermethylation direction, we found very little expression change between the DNA methylation levels\n0.0-0.6. In our experience, this feature is uncommon in ∆Meth-selected hypermethylation alterations and\nsuggests that relative ∆Meth may be revealing new and interesting DNA methylation biology.\nAn established alternative to ΔMeth is ΔM, where M = logit(β). Is relative ΔMeth superior to ΔM for\nhighlighting biologically meaningful sites? A direct comparison would strengthen the paper.\nWe would like to thank the reviewer for the suggestion to investigate whether the M = logit(β) approach\nwould also select for biologically meaningful DNA methylation alterations. The ∆M-value is a good\napproach to compare since it increases the resolution of DNA methylation loss when the levels are near\n0.0, which can lead to the selection of CpG sites with low native levels.\nAs expected, we found that large ∆M-value hypomethylation alterations occurred at CpG sites in all\nnative level quantile intervals (Figure 6A and Supplementary Figure S15). However, when we performed\nGSEA on the 10% largest ∆M-values that occurred in the hypomethylation direction, we found that the\nselected genes were enriched in pathways for the midbrain (Figure 6C and Supplementary Table S13 and\nS15). Again, we found that the unique relative ∆Meth selected genes were enriched in pathways known to\nbe important for tumor development (Figure 6C and Supplementary Table S14). Thus, we conclude that\nthe ability of relative ∆Meth to enrich for biologically significant alterations in the hypomethylation\ndirection is not simply due to the selection of CpG sites with low native levels.\nDetails describing these results can be found in the Results section lines 274-293 and Discussion section\nline 343-354).\na) Baseline definition. The “native” baselines appear to pool all normal tissues. Because DNA methylation is\nstrongly tissue-specific, this will preferentially select the small subset of CpGs that are constitutively high or\nlow in everytissue. Using matched normal tissue would be more appropriate.\nWe agree with the reviewer that using matched tissues would provide more biologically important\ninformation regarding DNA methylation alterations in comparison to using native levels derived from a\ndiverse dataset of normal tissues. Although it is uncommon to have matched tissues, there are a number\nof matched normal-tumor tissues in the TCGA database. To investigate whether using matched tissues\nwould cause a large difference in the relative ∆Meth distribution profiles, we calculated the native levels\nand relative ∆Meth using only matched breast tissues (n=97). We found that matched breast tissues had\nsimilar distribution profiles in comparison to the breast cancer relative ∆Meth distribution profiles made\nusing the native levels calculated from the pooled normal tissues, see plots below.\nWhile the overall conclusions regarding relative ∆Meth may hold true, we agree with the reviewer that\nmore biological information can be determined by the investigation of match tumor-normal tissues. In\nour next study, we plan to investigate relative ∆Meth alterations in matched tissues.\nb) CpG-to-gene mapping. For DNAm–expression correlations, how is per-gene methylation derived? Promoter\nmethylation has very different effects from gene-body methylation; the analysis should stratify by genomic\ncontext.\nWe derived the per-gene methylation information from the array annotation file from the Illumina\nwebsite. This information is now included in the Methods section. (Lines 455-456).\nWe agree that promoter methylation has different effects from gene-body methylation. To this point, we\nhave included the genomic context in relation to the CpG island of the top 10% most frequently altered\nCpG sites in Supplementary Tables S3 and S4 and included a new Supplemental Figure S9. Although the\nbiological functions of different genomic regions are important, we feel it is slightly outside the focus of\nthis manuscript and deserves its own investigation.\nIn our next study, we plan to investigate whether the relative ∆Meth approach can help illuminate new\nDNA methylation biological features in tumors. To accomplish this, we will investigate the different tumor\ntypes independently. We will also stratify the CpG sites based on their genomic regions within the genes\nand will continue to investigate PMD and HMD genomic compartments."}, {"round": 2, "reviewer_comments": {"reviewer_1": "The authors have taken action to address all the concerns raised by this reviewer.\nThe text is now clearer, and the previously missing details about various aspects have been sufficiently described.", "reviewer_2": "The authors have partially addressed my concerns from the original review. In particular, I appreciate the work done in response to the issue of methylation patterns (both normal/baseline and 'aberrant') being very different between different genomic regions. My remaining concern is that the new approach of studying changes in genomic regions defined as HMDs vs PMDs for all cancers grouped together is difficult to interpret since genomic methylation patterns are so tissue specific. (In fact a list of HMDs and PMDs for a given sample is a great predictor of which tissue type the sample comes from). In my opinion, this analysis would have to be done in a tissue-specific manner (or at least a subset of individual tissues) to be meaningful. (I.e. define HMDs, PMDs for that tissue, and then relate it to methylation alterations in the same tissue)."}, "author_response": "RESPONSE TO REVIEWERS’ COMMENTS\nReviewer #2 (Remarks to the Author):\nThe authors have partially addressed my concerns from the original review. In particular, I appreciate the work\ndone in response to the issue of methylation patterns (both normal/baseline and 'aberrant') being very different\nbetween different genomic regions. My remaining concern is that the new approach of studying changes in\ngenomic regions defined as HMDs vs PMDs for all cancers grouped together is difficult to interpret since\ngenomic methylation patterns are so tissue specific. (In fact a list of HMDs and PMDs for a given sample is a\ngreat predictor of which tissue type the sample comes from). In my opinion, this analysis would have to be done\nin a tissue-specific manner (or at least a subset of individual tissues) to be meaningful. (I.e. define HMDs,\nPMDs for that tissue, and then relate it to methylation alterations in the same tissue).\nWe would like to thank the reviewer for their thoughtful suggestion to reanalyze the genomic\ncompartment analysis in a tissue-specific manner. In our revised manuscript, we investigated large\n∆Meth and relative ∆Meth hypomethylation alterations that occur in common PMDs, common HMDs\nand in neither PMDs/HMDs in breast and colorectal cancers.\nTo accomplish this, we first identified the CpG sites of the array that were concordant with the previously\nidentified common genomic compartments. We defined the CpG sites to be concordant if the median\nDNA methylation level of the normal-adjacent breast and colorectal tissues were > 60% in common\nHMDs, between 30% and 60% in common PMDs, and < 30% in neither PMDs/HMDs.\nUsing the concordant CpG sites, we then identified the top 10% largest ∆Meth and relative ∆Meth\nalterations in the cancers. We found that in both tissue types, the CpG sites selected using large ∆Meth\nwere mostly in HMDs while the large relative ∆Meth were mostly in PMDs (Figure 7A and B). Similar to\nour previous findings, we found that the large ∆Meth in common HMDs enriched for genes in pathways\nconfounding for tumor development (Figure 7C and E and Supplemental Tables S16 and S19).\nConversely, we found that CpG sites with large relative ∆Meth in all genomic compartments were\nenriched for genes in pathways with known importance for tumor development (Figure 7D and F and\nSupplemental Tables S17, S18, and S20-22). Furthermore, we also found that in both breast and\ncolorectal cancers, alterations in HMDs were enriched in more immune pathways than alterations in\nPMDs and neither PMDs/HMDs.\nFrom these results, we confirmed the previous findings that alterations at different genomic\ncompartments have different biological functions for tumor development. Moreover, we found that these\nbiological differences are more obvious when using the relative ∆Meth approach in comparison to the\n∆Meth approach."}, {"round": 3, "reviewer_comments": {"reviewer_2": "I appreciate this additional breast and colon tissue analysis which addresses my earlier concern.\nOne final comment to the authors: For future work (i.e. in my opinion not necessary for this manuscript) it could be interesting to understand what happens specifically in regions that are discordant between tissues (e.g. regions that are HMD in breast and PMD in colon)."}, "author_response": ""}]

Decision letter