Thursday, November 30, 2023
HomeBiologyHistone methyltransferase exercise impacts metabolism in human cells independently of transcriptional regulation

Histone methyltransferase exercise impacts metabolism in human cells independently of transcriptional regulation

Quotation: Perez MF, Sarkies P (2023) Histone methyltransferase exercise impacts metabolism in human cells independently of transcriptional regulation. PLoS Biol 21(10):

Educational Editor: Jason W. Locasale, Duke College, UNITED STATES

Obtained: April 26, 2023; Accepted: September 27, 2023; Printed: October 26, 2023

Copyright: © 2023 Perez, Sarkies. That is an open entry article distributed beneath the phrases of the Artistic Commons Attribution License, which allows unrestricted use, distribution, and replica in any medium, offered the unique creator and supply are credited.

Knowledge Availability: All related information are contained throughout the paper and/or supporting info recordsdata. The information underlying the graphs proven within the figures and the code used may be present in

Funding: This work was supported by the UAS John Fell Fund (to PS) grant quantity 0011417, which supported the wage MFP till January 2023. This work was additionally funded by a Ramon y Cajal fellowship RYC2021-034496-I (to MFP), which supported the wage of MFP from January 2023. The funders had no position in examine design, information assortment and evaluation, resolution to publish, or preparation of the manuscript.

Competing pursuits: The authors have declared that no competing pursuits exist.

1-methylnicotinamide; BRCA,
breast most cancers; CCLE,
Most cancers Cell Line Encyclopedia; COAD,
colon adenocarcinoma; FDR,
false discovery charge; GRN,
gene regulatory community; GTEx,
Genotype-Tissue Expression; HERV,
human endogenous retrovirus; HMT,
histone methyltransferase; LC–MS,
liquid chromatography–mass spectrometry; LUSC,
lung squamous cell carcinoma; MEF,
mouse embryonic fibroblast; MRN,
median-ratio normalisation; NNMT,
nicotinamide N-methyltransferase; PAAD,
pancreatic adenocarcinoma; PC,
phosphatidylcholine; PE,
phosphoethanolamine; PI,
proliferative index; PRMT,
protein arginine methyltransferase; Rb,
Retinoblastoma; SAH,
S-adenosyl homocysteine; SAM,
S-adenosyl methionine; STAD,
abdomen adenocarcinoma; TCGA,
The Most cancers Genome Atlas; TR,
transcriptional regulator; TSS,
transcription begin web site


The invention of sturdy associations between the transcriptional states of genes and methylation of histones at these loci was a seminal second within the examine of gene regulation [1]. The idea of the histone code [2] proposed that data of the actual mixture of epigenetic marks in chromatin would permit a deterministic prediction of gene expression, simply because the deciphering of the genetic code had allowed exact prediction of gene merchandise. Nevertheless, 2 many years later, regardless of sturdy correlations of some histone marks with particular transcriptional states [3], these associations may be ambiguous [4,5], whereas proof of causal hyperlinks between histone marks and transcriptional activation or repression stays equivocal [6].

In that point, it has change into extensively appreciated that histone modifications may be influenced by mobile metabolism [7]. Histone methylation is influenced by the provision of S-adenosyl methionine (SAM, additionally abbreviated as AdoMet or SAMe), the common methyl donor that’s required for mobile methylation of lipids, proteins, nucleic acids, and metabolites, and which may be modulated by dietary methionine supplementation [8]. Nevertheless, the abundance of histones within the cell presents the potential for histone modifications to influence metabolism [9]. Histones have the potential to behave as a methyl sink, as histone methylation consumes methyl teams from SAM however demethylation releases formaldehyde, which can’t be simply recycled to launch methyl teams [10]. Methyl sinks play key metabolic roles, appearing to buffer the ratio of SAM to S-adenosyl homocysteine (SAH, additionally abbreviated as AdoHcy) and supporting the synthesis of essential sulphur-containing metabolites comparable to cysteine and glutathione through the SAH-dependent transsulphuration pathway [9,11].

Right here, we found sturdy detrimental relationships between the overall expression of histone methyltransferases (HMTs) and metabolic pathways beforehand characterised as methyl sinks. In each most cancers cells and wholesome tissue, we present that these relationships affected genome-wide ranges of histone posttranslational modifications however didn’t have vital penalties for transcriptional regulation. We present that HMTs have been co-expressed and negatively regulated by the Retinoblastoma (Rb) tumour suppressor in most cancers. Our information recommend the speculation that the overall exercise of HMTs has penalties for SAM homeostasis in wholesome human tissues and tumours, unbiased of the capabilities of histone methylation in transcriptional regulation.


Histone methyltransferase expression correlates to mobile metabolite ranges

We got down to examine a possible hyperlink between metabolism and histone methylation. We reasoned that results of HMT exercise on metabolism may lead to correlations between HMT ranges and mobile metabolite concentrations. To research this risk, we used a publicly obtainable metabolomics dataset consisting of 225 metabolites profiled by liquid chromatography–mass spectrometry (LC–MS) throughout 911 cell traces from the Most cancers Cell Line Encyclopedia (CCLE), representing extra 23 most cancers varieties [12]. We associated metabolite ranges to the normalised expression of HMTs in the identical cell line. We curated an inventory of 38 HMTs (S1 Desk) and examined the correlation of every HMT to all metabolites. Throughout this set 1-methylnicotinamide (1MNA) persistently emerged because the metabolite most strongly related to HMTs, with a false discovery charge (FDR) < 0.05 for 18 particular person HMTs, a geometrical imply FDR of 0.001 (S2 Desk) and a mean Pearson’s correlation of −0.090 (vary −0.261 to 0.160). Certainly, 1MNA was the metabolite with the biggest absolute correlation to the overall degree of HMTs obtained by summing the expression of the 38 particular person enzymes (Fig 1A; Pearson’s correlation = −0.274, FDR = 1.77 × 10−14).


Fig 1. Whole HMT expression is strongly anticorrelated with the exercise of NNMT in cancers.

(A) Volcano plot exhibiting Pearson’s correlation and FDR for 225 metabolites to whole HMT expression (whole RNA-seq median of ratios-normalised pseudocounts) throughout 927 most cancers cell traces from the CCLE. (B) Volcano plot exhibiting Pearson’s correlation and FDR for expression of 10,275 expressed genes to ranges of 1MNA throughout 927 CCLE most cancers cell traces. The highest and backside 2.5% of factors are proven in darker gray. HMT-encoding genes are proven as factors colored in keeping with their affiliation with transcriptional activation (inexperienced), repression (magenta), or an unclear relationship (blue). Pearson’s r for whole HMT expression is proven as a black level. (C) PCA of metabolite ranges throughout 927 most cancers cell traces from the CCLE. 1MNA is highlighted with a crimson circle. (D) NNMT and HMTs each convert SAM to SAH and so can have an effect on mobile methylation potential by appearing as a “sink.” (E) Volcano plot exhibiting Spearman’s correlation and FDR for expression of NNMT vs. 52,440 genes in a pan-cancer evaluation of 927 CCLE cell traces throughout 23 most cancers varieties. HMT-encoding genes are proven as factors as in panel 1B. (F) Volcano plot exhibiting Spearman’s correlation and FDR for expression of NNMT vs. 60,489 genes in a pan-cancer evaluation of TCGA major tumours throughout 33 most cancers varieties. HMT-encoding genes are proven as factors as in panel 1B. (G) Violin plot exhibiting Spearman’s correlation to NNMT for HMTs (black, proper) or different genes (left, gray) in 79 major ACC tumours from the TCGA. Particular person HMT-encoding genes are proven as factors as in panel 1B. (H) Spearman’s correlation vs. NNMT expression of whole expression of pooled HMTs added to the pool in a random order, and 1,000 particular person iterations are proven as black traces, with the domestically estimated smoothing (Loess match) trendline proven in crimson. (I) TCGA pan-cancer evaluation exhibiting rank percentile place of whole HMTs amongst correlations of NNMT expression to 60,489 genes and vice versa in 33 distinct most cancers varieties. Bubble measurement is inversely proportional to the log of the “relative reciprocal rating,” the sum of squares of the ranks of whole HMTs/NNMT within the reciprocal distribution. The dashed gray field signifies correlations within the strongest 2.5% of anticorrelated genes. Underlying information for all panels may be present in 1MNA, 1-methylnicotinamide; ACC, adrenocortical carcinoma; CCLE, Most cancers Cell Line Encyclopedia; FDR, false discovery charge; HMT, histone methyltransferase; NNMT, nicotinamide N-methyltransferase; SAH, S-adenosyl homocysteine; SAM, S-adenosyl methionine; TCGA, The Most cancers Genome Atlas.

HMT expression varies reciprocally with the exercise of the 1MNA/NNMT methyl sink

The sturdy relationship between 1MNA and HMT ranges indicated that HMT expression variation throughout most cancers cell traces is likely to be related to adjustments in metabolism. To research which metabolic pathways is likely to be accountable, we carried out principal element evaluation (PCA) and clustering evaluation on all 225 metabolites. Associated metabolites from identified biochemical pathways tended to cluster collectively. 1MNA was a transparent outlier in each analyses (Figs 1C and S1), indicating that 1MNA synthesis displays a discrete metabolic course of. Certainly, 1MNA is thought to be a secure metabolic end-product that has no downstream metabolites in most cancers [13] and is excreted from cells in wholesome tissues [14].

1MNA is the product of methylation of nicotinamide by the enzyme nicotinamide N-methyltransferase (NNMT). 1MNA ranges have been strongly correlated with NNMT expression as measured by RNA-seq (S2A Fig and S2 Desk) and with NNMT protein ranges as measured by quantitative mass spectrometry (S2B Fig and S3 Desk; [15]), whereas nicotinamide ranges have been anticorrelated with NNMT (S2 and S3 Tables). We conclude that 1MNA ranges mirror the exercise of the metabolic pathway that converts nicotinamide to 1MNA, catalysed by the enzyme NNMT. Furthermore, this pathway is just not tightly coupled to the exercise of different pathways of core metabolism. We subsequently determined to analyze doable explanations for the connection between HMT ranges and 1MNA pathway exercise.

HMT and NNMT are various pathways that each devour methyl teams

The response catalysed by NNMT makes use of SAM as a cofactor, transferring a methyl group from SAM to nicotinamide to kind 1MNA (Fig 1D). NNMT has been proposed to perform as a “sink” for methyl teams. Excessive NNMT exercise can cut back the SAM:SAH ratio [13,1618].

To additional examine the connection between the 1MNA synthesis pathway and HMT exercise, we investigated the connection between NNMT expression and whole HMT expression. NNMT protein ranges in 2 most cancers cell line panels, the CCLE and NCI60 [19], correlated strongly to NNMT expression (S2C and S2D Fig, respectively). NNMT expression is subsequently a dependable indicator of NNMT protein ranges and catalytic exercise. NNMT expression and whole HMT expression within the CCLE have been negatively correlated (Fig 1E). HMT protein ranges have been negatively correlated with NNMT protein ranges each for particular person HMTs (for 20 HMTs detected in >90% of samples, imply Pearson’s r to NNMT protein ranges = −0.108, vary −0.306–0.199, FDR < 0.05 detrimental correlation for 12/20 and optimistic correlation for 1/20) and collectively (S2E Fig; imply Pearson’s r with pattern imply HMT protein Z-score = −0.244, p-value = 2.62 × 10−5). HMT protein ranges are additionally negatively correlated with 1MNA ranges (S3 Desk; imply Pearson’s r with pattern imply HMT protein Z-score = −0.102, p-value = 7.94 × 10−2). Altogether, this recommended that elevated 1MNA synthesis is related to lowered HMT exercise.

One doable rationalization for the detrimental affiliation between HMT ranges and 1MNA synthesis is that 1MNA straight represses HMT transcription. Nevertheless, partial correlation evaluation indicated that the correlation between HMTs and 1MNA was weakened from −0.280 to −0.073 after controlling for NNMT expression (S2F Fig). We concluded that 1MNA itself is unlikely to control HMT expression. As an alternative, HMT ranges have been primarily related to NNMT expression, and thus with the speed of 1MNA synthesis, somewhat than 1MNA ranges themselves. These findings are in keeping with the proposal that HMT exercise and NNMT exercise are parallel pathways able to appearing as methyl sinks [9]. NNMT exercise is lowered when HMT exercise is excessive (Fig 1E)

HMT and NNMT expression are tightly coupled in major tumours

To check if the connection we uncovered in most cancers cell traces was additionally seen in tumours, we used RNA-seq information from major tumours from 33 distinct most cancers varieties discovered within the Most cancers Genome Atlas (TCGA) database to interrogate the correlation between HMTs and NNMT expression. The HMT-NNMT expression relationship was a lot stronger than that noticed within the CCLE cell traces, each in a pan-cancer evaluation (Fig 1E and 1F) and inside particular person most cancers varieties (S1 and S2 Recordsdata). The connection was strengthened by pooling the expression of HMTs collectively (Fig 1F and S1 and S2 Recordsdata). For example, among the many 177 pancreatic adenocarcinoma (PAAD) major tumour samples, the connection of NNMT expression with expression of particular person HMTs was for essentially the most half not exceptionally sturdy. Nevertheless, when their expression was pooled, HMT expression anticorrelated with NNMT higher than nearly any single gene (rank 3, prime 0.00496%; Fig 1G). Conversely, NNMT was one of the crucial negatively correlated genes to pooled HMT expression (rank 9, prime 0.0149% of genes total).

Uncooked correlation statistics (comparable to ρ) don’t essentially present a dependable comparator of the power of the connection, because the distribution of gene correlations can differ between most cancers or tissue varieties [20]. To beat this, we computed a relative reciprocal relationship rating for HMTs and NNMT (see Strategies). The relative reciprocal relationship scores have been a lot stronger within the TCGA major tumours than within the CCLE cell traces (Figs 1I and S2G). The HMT-NNMT relationship was discovered throughout most most cancers varieties and was significantly sturdy in liver, pancreatic, and bladder cancers. The one most cancers during which this relationship was not evident in any respect was acute myeloid leukemia (LAML). This can be a results of the low expression of NNMT on this most cancers sort (S2H Fig). Certainly, throughout all cancers there was a big correlation between the power of the connection between NNMT and HMT ranges and the expression of NNMT (S2I Fig). Importantly, the HMT-NNMT relationship remained strong when controlling for immune cell infiltration, as estimated by 2 distinct gene expression deconvolution instruments [21,22] (S3 Fig).

To visualise how the correlation strengthened as HMT expression was mixed, we carried out a pan-cancer evaluation during which we added successively extra HMTs and calculated the correlation between these HMTs and NNMT. The connection between HMTs and NNMT turned stronger as extra HMTs have been added, whatever the order with which the HMTs have been mixed (Fig 1H). This relationship signifies that the anticorrelation between HMTs and NNMT is distributed throughout HMTs and never on account of 1 or 2 HMTs with unusually sturdy anticorrelations. Utilizing a easy stochastic modelling method, we decided that the reinforcement of this relationship when extra HMTs are added is in keeping with co-regulation of the HMTs on the transcriptional degree (S2J Fig).

Whole HMT expression was extra strongly anticorrelated to NNMT than any particular person HMT in 16/33 most cancers varieties (S2 File). Nevertheless, 1 notable exception was melanoma (SKCM), which regardless of being one of many most cancers varieties with the strongest NNMT-whole HMT relationships exhibited a far stronger anticorrelation of NNMT to the euchromatic H3K9me3 author SETDB1 (S4 Fig and S2 File). SETDB1 is essentially the most strongly anticorrelated gene within the genome to NNMT and vice versa. As SETDB1 is a recurrently amplified, established driver of melanoma [2326], this raises the likelihood that elevated SETDB1 exercise might need essential metabolic penalties for melanoma.

Amongst mobile methyltransferases, HMTs have the strongest relationship with NNMT

To check whether or not the sturdy anticorrelation to NNMT is particular to histone lysine methyltransferases, we calculated the relative reciprocal relationship scores for protein arginine methyltransferases (PRMTs), DNA methyltransferases, and teams of RNA methyltransferases (S5 Fig; gene units in S1 Desk). The connection of NNMT to histone lysine methyltransferases was the strongest and most widespread. There was no vital relationship between NNMT and PRMTs or tRNA/rRNA methyltransferases. In breast most cancers (BRCA) and 16 different most cancers varieties NNMT had a robust detrimental relationship to mRNA methyltransferases. In lung squamous cell carcinoma (LUSC) and 4 different most cancers varieties, there was a robust detrimental relationship between NNMT and small RNA methyltransferases. Extra notably, 9 out of 33 most cancers varieties displayed a robust relationship between NNMT and methyltransferases with unknown substrates, with very sturdy relationships evident in LUSC and abdomen and colon adenocarcinomas (STAD and COAD). We propose that different mobile methyltransferases might also act as methyl sinks in parallel to lysine HMTs and NNMT. Nevertheless, this relationship is most evident for histone lysine methyltransferases throughout most cancers.

The connection between HMTs and NNMT is particular to most cancers

We investigated the connection of HMT to NNMT in RNA-seq information from wholesome tissue samples obtainable via the Genotype-Tissue Expression (GTEx) venture. A robust relationship between NNMT expression and HMT expression was noticed in just one particular person tissue sort (muscle) (Fig 2A and 2B and S3 File). Furthermore, in 10 out of the 12 TCGA most cancers varieties for which matched tumour and regular tissue samples from at the very least 30 sufferers have been obtainable, the detrimental relationship of NNMT and HMT was stronger within the most cancers samples relative to the matched regular samples (paired Wilcoxon check on reciprocal scores, p-value = 1.47 × 10−3; S6A Fig).


Fig 2. Whole HMT expression is strongly anticorrelated with the expression of PEMT in wholesome tissues.

(A) Evaluation exhibiting rank percentile place of whole HMTs amongst correlations of NNMT expression to 56,200 genes and vice versa in 48 distinct wholesome tissue varieties from the GTEx venture. Bubble measurement is inversely proportional to the log of the “relative reciprocal rating,” the sum of squares of the ranks of whole HMTs/NNMT within the reciprocal distribution (see Strategies). The dashed gray field signifies correlations within the strongest 2.5% of anticorrelated genes, with tissues labelled. (B) Volcano plot exhibiting Spearman’s correlation and FDR for expression of NNMT vs. 56,200 genes in a pan-cancer evaluation of GTEx major tumours throughout 48 tissue varieties. HMT-encoding genes are proven as factors colored in keeping with affiliation with transcriptional regulation; correlation for whole HMT expression is proven as a black level. (C) Evaluation exhibiting rank percentile place of whole HMTs amongst correlations of PEMT expression and vice versa in wholesome tissue varieties from the GTEx venture. Bubble measurement and dashed gray field as in panel 2A. (D) Volcano plot exhibiting Spearman’s correlation and FDR for expression of PEMT vs. 56,200 genes in a cross-tissue evaluation of 18 tissue varieties with a robust HMT-PEMT relationship (throughout the gray field in panel 2C). HMT-encoding genes are proven as factors as in panel 2B. (E) Violin plot exhibiting Spearman’s correlation to PEMT of HMTs (black, proper) or different genes (left, gray) in 375 affected person samples from the gastroesophageal junction. HMT-encoding genes are proven as factors as in panel 2B. (F) Spearman’s correlation vs. PEMT expression of whole expression of pooled HMTs added to the pool in a random order in a cross-tissue evaluation of tissues with a robust HMT-PEMT relationship (throughout the gray field in panel 2C); 1,000 particular person iterations are proven as black traces, with Loess match trendline in crimson. (G) Evaluation exhibiting rank percentile place of whole HMTs amongst correlations of PEMT expression to 60,489 genes and vice versa in 33 most cancers varieties from the TCGA. Bubble measurement and dashed gray field as in panel 2A. (H) PEMT sequentially methylates phosphoethanolamine to supply PC, changing 3 molecules of SAM to SAH. (I) Evaluation exhibiting rank percentile place of HMTs categorised by their substrate histone lysine residues amongst correlations of NNMT expression and vice versa in most cancers varieties from the TCGA. Bubble measurement and dashed gray field as in panel 2A. (J) Evaluation exhibiting rank percentile place of HMT units methylating distinct histone lysine residues amongst correlations of PEMT expression to 56,200 genes and vice versa in a pan-tissue evaluation of 18 tissue varieties from the GTEx with a robust HMT-PEMT relationship (throughout the gray field in panel 2C). Bubble measurement and dashed gray field as in panel 2A. (Ok) Evaluation exhibiting rank percentile place of HMT units methylating distinct histone lysine residues amongst correlations of PEMT expression and vice versa in a pan-cancer evaluation of seven most cancers varieties from the TCGA with a robust HMT-PEMT relationship (throughout the gray field in panel 2G). Bubble measurement and dashed gray field as in panel 2J. (L) Violin plot exhibiting wholesome tissue pattern PI, a measure of proliferation inferred from pattern RNA-seq gene expression information, for 48 tissue sorts of the GTEx organized by the power of the anticorrelating relationship between PEMT and whole HMTs. Be aware the x axis is inverted as a decrease relative reciprocal rating signifies a stronger relationship. (M) Violin plot exhibiting tumour PI for 31 most cancers sorts of the TCGA organized by the power of the anticorrelating relationship between NNMT and whole HMTs. Underlying information for all panels may be present in FDR, false discovery charge; GTEx, Genotype-Tissue Expression; HMT, histone methyltransferase; NNMT, nicotinamide N-methyltransferase; PC, phosphatidylcholine; PI, proliferative index; SAH, S-adenosyl homocysteine; SAM, S-adenosyl methionine; TCGA, The Most cancers Genome Atlas.

We questioned whether or not the distinction between most cancers and wholesome tissue could possibly be on account of variations within the expression ranges of NNMT and HMTs. Nevertheless, whereas NNMT has ceaselessly been reported as overexpressed in cancers [2734], we didn’t discover sturdy assist for this throughout the TCGA. NNMT was considerably up-regulated in tumour samples in comparison with matched regular samples in 3/12 most cancers varieties, with 4/12 exhibiting vital down-regulation (S6B Fig). Equally, whole HMT expression was not persistently altered in most cancers, being up-regulated in 4/12 most cancers varieties and down-regulated in 2/12 (S6C Fig). Curiously, nevertheless, we noticed that these adjustments have been inversely correlated: cancers which up-regulate HMTs have lowered NNMT expression and vice versa (S6D Fig).

In wholesome tissues, HMTs correlate with PEMT, another methyl sink

We examined whether or not whole HMT expression may present a similar relationship to a different methyltransferase or class of methyltransferases working throughout the cell. In 18/48 tissues whole HMT expression anticorrelated strongly with PEMT (reciprocally in prime 2.5% of genes; Fig 2C and S4 File). Curiously, a very excessive proportion of mind tissues (10/13) confirmed a robust, vital relationship. Throughout the 18 tissues that confirmed sturdy associations between PEMT and HMT expression (Fig 2D), the connection turned stronger as extra HMTs have been pooled (Fig 2E and 2F). PEMT was strongly negatively correlated to HMTs in 7/33 cancers (Fig 2G).

PEMT is an enzyme that provides 3 methyl teams to the phospholipid phosphoethanolamine (PE) to supply phosphatidylcholine (PC) (Fig 2H). PC makes up round 40% of the lipid content material of the plasma membrane in eukaryotic cells [35]. PEMT exercise contributes round 30% of mobile PC synthesis [36]. The abundance of PC within the membrane recommended the likelihood that PEMT may act as a sink for methyl teams, equally to NNMT; certainly, PEMT has been recommended to be the first shopper of SAM in mammals [37].

We investigated whether or not different teams of methyltransferases additionally anticorrelated to PEMT. We discovered a big detrimental relationship of PEMT to mRNA methyltransferases in 14 tissue varieties, together with 9/13 mind tissues (S7 Fig). Likewise, 7 tissues (5 from the mind) displayed a similar relationship with DNA methyltransferases (S7 Fig). Histone lysine methyltransferases had the strongest and most constant anticorrelation to PEMT.

Differential contributions of methylated residues to HMT relationship with methyl sinks

We investigated whether or not particular methylated residues on histones contribute extra strongly to the connection with mobile methyl sinks. In most cancers, the HMTs that methylate the H3K9 residue confirmed a stronger relationship to NNMT, whereas the connection to methyltransferases focusing on different residues was weaker, albeit nonetheless extremely vital (Fig 2I). In wholesome tissues PEMT confirmed a robust reciprocal relationship to H3K4 and H3K36 methyltransferases, whereas the connection with H3K9 and H3K27 methyltransferases, related to transcriptional repression, was weaker and no sturdy relationship existed for H4K20 methyltransferases (Fig 2J). The distinction between PEMT and NNMT was not on account of a distinction between most cancers and wholesome tissue as a result of cancers with a big HMT-PEMT relationship confirmed related residue specificity as wholesome tissues (Fig 2K).

We questioned whether or not the distinction in HMTs that correlated with NNMT and PEMT could possibly be as a result of cell cycle exercise of those pathways. Transcription-associated methylation on H3K4me3 and H3K36me3 is enriched in quiescent cells relative to methylation of H3K9me3, which largely happens in late S and G2 section to revive H3K9 methylation to newly synthesised histones [38]. PEMT expression can also be reported to range strongly throughout the cell cycle, peaking in G1 section and declining in S section [39,40]. We examined the correlation between the power of HMT-PEMT relationship throughout the GTEx wholesome tissues to the tissue pattern’s proliferative index (PI), a measure of proliferation inferred from RNA-seq information [41]. The HMT-PEMT relationship tended to be stronger in tissues with a decrease PI (Fig 2L; R2 = 0.116, p-value ~ 0). This may occasionally mirror the shortcoming of proliferating cells to make use of PEMT as a methyl sink all through the cell cycle. The HMT-NNMT relationship in cancers exhibited the alternative development, albeit weakly (Fig 2M; R2 = 0.0384, p-value = 1.05 × 10−84). The aggressive proliferation of most cancers cells might clarify why whole HMT expression correlates to totally different methyl sinks in cancers and wholesome tissues.

Histone methylation in chromatin is correlated to methyl sink exercise with out adjustments in transcription

We investigated whether or not the connection between HMT exercise and PEMT or NNMT impacted histone methylation ranges in chromatin genome-wide. We used ChIP-seq information from wholesome tissues and most cancers cell traces to evaluate histone methylation ranges. In wholesome tissues, we discovered that there was a world detrimental relationship of PEMT expression with H3K4me3, H3K9me3, and H3K27me3 in any respect courses of genomic areas examined. For instance, 99% of 12,355 gene promoters marked by H3K4me3 peaks confirmed an anticorrelated relationship between PEMT expression and whole H3K4me3 sign or H3K4me3 peak width, an orthogonal measure of histone methylation ranges (Figs 3A, 3B, S8A and S8B). Equally, 99% of two,870 repetitive areas modelled confirmed an anticorrelated relationship for H3K9me3 sign and PEMT expression (Fig 3A). Nevertheless, this was not the case for H3K36me3, which had a reasonably optimistic relationship with PEMT expression at promoters and gene our bodies. No detrimental relationships with NNMT expression have been noticed in wholesome tissues (S8C Fig).


Fig 3. PEMT and NNMT expression anticorrelate globally with ranges of particular histone marks genome-wide in wholesome tissues and cancers, respectively.

(A) Boxplot reveals t-values from linear blended results mannequin for pattern PEMT expression predicting ChIP-seq sign for numerous histone marks (label left) on gene our bodies, promoters or repetitive components (subpanel headers) in affected person tissue samples collected as a part of the ENCODE venture. The variety of particular person websites is famous on the plot for every boxplot; p-values derive from paired Wilcoxon assessments towards a null distribution calculated by the imply t-value at every locus for 1,000 random expressed genes. (B) Heatmap exhibiting H3K4me3 ChIP-seq sign (log2 fold change over enter) over 1,000 random genes for 4 samples from the squamous epithelium of the esophagus organized so as of PEMT expression. (C) Boxplot reveals t-values from generalised linear fashions for NNMT expression (RNA-seq) predicting ChIP-seq sign for numerous histone marks (label left) on gene our bodies and promoters in cell traces of the NCI60 most cancers cell line panel. The variety of particular person websites is famous on the plot for every boxplot; p-values derive from paired Wilcoxon assessments towards a null distribution calculated by the imply t-value at every locus for 1,000 random expressed genes. (D) Boxplot reveals t-values from generalised linear fashions for NNMT expression (RNA-seq) predicting ChIP-seq sign for numerous histone marks (label left) on totally different courses of repetitive components in cell traces of the NCI60 most cancers cell line panel. The variety of particular person websites is famous on the plot for every boxplot; p-values derive from paired Wilcoxon assessments towards a null distribution calculated by the imply t-value at every locus for 1,000 random expressed genes. Websites proven are from bin with highest ChIP sign (cf. S9D Fig). Underlying information for all panels may be present in HERVs, human endogenous retroviruses; LINEs, lengthy interspersed nuclear components; LTRs, lengthy terminal repeats; NNMT, nicotinamide N-methyltransferase; SINEs, quick interspersed nuclear components.

Altered histone modification ranges are sometimes related to adjustments in transcription of the genes on the corresponding loci. Nevertheless, we discovered that regardless of low sign of H3K4me3 (related to transcription) and H3K9me3 (related to repression) in high-PEMT samples, expression from marked genes was not affected in both case (S8D Fig). Equally, we discovered that variation throughout samples in whole H3K4me3 sign (S8E Fig) or width (S8F Fig) at marked promoters, and H3K9me3 sign on H3K9me3-marked gene our bodies, doesn’t correlate with the expression of the corresponding genes. Thus in wholesome tissues, PEMT assorted with histone methylation ranges unbiased of results on transcription.

In most cancers, H3K9me3 and H4K20me3 at each gene our bodies and promoters have been anticorrelated to cell line NNMT expression as measured by RNA-seq (Fig 3C), microarrays (S9A Fig), and proteomics (S9B Fig); nevertheless, H3K4me3 confirmed a optimistic relationship. No detrimental relationships have been noticed for PEMT expression (S9C Fig). Equally, we discovered that each H3K9me3 and H4K20me3 ranges at a number of courses of repetitive components have been negatively correlated with NNMT expression (Fig 3D). For all courses, the detrimental correlation was stronger at genomic websites with the next common ChIP-seq sign throughout samples (S9D Fig). The anticorrelation with NNMT expression was significantly sturdy at centromeric satellites (Fig 3D), unbiased of their tendency to show increased sign of heterochromatic marks than different courses of repetitive component (S9E Fig).

Whereas variation throughout samples in H4K20me3 sign at H4K20me3-marked gene our bodies was negatively correlated with NNMT ranges, expression from these genes displayed little relationship with NNMT (S9F Fig). We additionally estimated locus-specific expression of transposable components, particularly human endogenous retroviruses (HERVs). Equally, we discovered that regardless of lowered sign of this canonically repressive histone mark at HERVs in samples with excessive NNMT expression, HERV expression was not elevated (S9F Fig). Certainly, variation in whole sign of H3K9me3/H4K20me3 at marked websites was not related to the extent of transcription from both gene our bodies or HERVs (S9G Fig).

Histone methyltransferase genes are co-expressed

We wished to know the construction and origin of the variation we observe in whole HMT expression. We discovered that the expression of HMT genes was considerably extra positively correlated to one another than to random genes (S10A and S10B Fig). The strongest co-expression was evident among the many most extremely expressed HMTs, which had a spread of goal lysine residues with divergent associations with transcription (see S10A Fig sidebars). This core of 14 to 16 extremely expressed and extremely correlated HMT genes was largely secure between wholesome tissues and most cancers (S10A Fig) and throughout particular person tissue or most cancers varieties (S5 and S6 Recordsdata). We carried out gene co-expression community evaluation, exhibiting the same community structure with a small variety of modularity courses in each most cancers and wholesome tissue (S10C and S10D Fig). The power of community edges was extremely concordant between pan-cancer and pan-tissue correlation analyses (S10C and S10D Fig; Pearson’s correlation = 0.834; Jaccard index on community with |edge power| > 0.2 = 0.569). Whereas distinct modularity clusters confirmed some similarity to annotated associations with transcription (S10C and S10D Fig), the overlap was not sturdy. Collectively, this recommended the likelihood that HMTs is likely to be coordinately regulated unbiased of their transcriptional capabilities.

Histone methyltransferase genes are regulated by E2F and Rb

Searching for to know the doable foundation for co-regulation of HMTs, we turned to Caenorhabditis elegans, a mannequin organism with an easier genetic regulatory structure [42]. We noticed a optimistic correlation amongst HMTs throughout 206 numerous pure genetic backgrounds in C. elegans, with the strongest correlation in a core of 10 genes largely consisting of essentially the most extremely expressed HMTs (S11A Fig). Genes on this cluster have been extra prone to have orthologues within the human extremely expressed cluster (S11A Fig; odds ratio = 11.07, Fisher’s actual check p = 0.0108). Moreover, we noticed very sturdy detrimental correlations between whole HMT expression and expression of the NNMT orthologues anmt-1/3 and the PEMT analogue, pmt-1 (S11B Fig). This relationship was partially on account of various ranges of anmt-1/3 and pmt-1 throughout growth (S11C Fig) and was additionally evident when controlling for developmental age (S11D and S11E Fig).

We carried out a de novo motif enrichment search on the upstream areas of co-regulated C. elegans HMT genes. Probably the most strongly enriched motif, current in 9/10 genes, resembled the binding motif of the E2F orthologue EFL-2 (Fig 4A). E2F transcription elements may be certain by the Rb protein, which represses transcription of E2F targets [43]. Utilizing beforehand revealed genome-wide ChIP-seq information [44], we noticed an enrichment for binding of LIN-35, the C. elegans Rb orthologue, near the transcription begin web site (TSS) of HMT genes within the extremely expressed cluster relative to different HMTs or to random genes (Fig 4B). Certainly, in RNA-seq information from lin-35 mutants [45,46], we noticed whole HMT expression elevated by 12% to fifteen% (1-way ANOVA, p = 0.053).


Fig 4. HMTs are regulated by E2F and Retinoblastoma, with NNMT expression lowered downstream of HMTs in Rb-mutant cancers.

(A) Above: Sequence motif enriched in C. elegans HE cluster HMT promoters, relative to different HMTs promoters. Under: beforehand reported EFL-2 binding motif. (B) Binding of the C. elegans Retinoblastoma orthologue LIN-35 upstream of the TSS of the HE cluster, different HMT genes, and random genes; p-values from Wilcoxon check. (C) Enrichment for transcription issue binding, from ENCODE ChIP-seq experiments, upstream of human HMT genes. Odds ratios and p-value derived from Fisher’s actual check. (D) Boxplots present median whole HMT or NNMT expression percentile drawn from 1,000 iterations of pan-cancer sampling of tumours with wild-type RB1 or probably deleterious RB1 mutations; p-value derived from t check. (E) Whole HMT expression in small cell lung most cancers cell traces from the CCLE with wild-type RB1 or deleterious RB1 mutations; p-value derived from t check. (F) Estimated E2F1 exercise vs. whole HMT expression (each corrected for confounders) in breast most cancers major tumours from the TCGA. (G) Potential architectures of the GRN linking RB1, NNMT, and HMTs. (H) Linear mannequin t-values explaining whole HMT and NNMT expression for RB1 mutation standing as the only explanatory variable or collectively thought of with NNMT/HMT expression respectively; p-value derived from t check. Underlying information for all panels may be present in CCLE, Most cancers Cell Line Encyclopedia; GRN, gene regulatory community; HMT, histone methyltransferase; NNMT, nicotinamide N-methyltransferase; Rb, Retinoblastoma; TCGA, The Most cancers Genome Atlas; TSS, transcription begin web site.

To check whether or not E2F transcription elements additionally regulated HMT expression in people, we carried out an enrichment evaluation for transcription issue binding websites upstream of human HMT genes decided from ENCODE venture ChIP-seq experiments for 181 distinct TFs. We noticed that E2F1 was the second most importantly enriched TF (Fig 4C).

We recognized 256 major tumours from 10 most cancers varieties within the TCGA that had probably deleterious mutations within the Retinoblastoma-encoding RB1 gene in at the very least 10 samples per most cancers sort. In all most cancers varieties RB1-mutant tumours had increased imply whole HMT expression; in a pan-cancer evaluation, we noticed that this distinction was extremely vital (Fig 4D). Amongst 31 single HMTs with expression in all samples, 19 confirmed a marked up-regulation in RB1-mutant tumours, with significantly notable up-regulation of EZH2, DOT1L, and NSD2, whereas solely 5 displayed a transparent down-regulation (S7 File).

We additionally recognized most cancers cell traces from the CCLE with deleterious Rb mutations. Round half of Rb mutations in CCLE cell traces have been present in lung most cancers cell traces, significantly small cell lung most cancers traces. Whole HMT expression was considerably up-regulated in RB1-mutant lung most cancers cell traces relative to wild-type RB1 cell traces (Fig 4E; 2-way ANOVA controlling for lung most cancers subtype, p = 0.034).

To check whether or not variability in HMT expression throughout tumours was related to variable Rb/E2F exercise, even when Rb was not mutated, we inferred the exercise of 351 transcriptional regulators from expression of goal genes in RNA-seq information from hundreds of samples throughout the TCGA and GTEx [47]. Within the TCGA, inferred E2F1 exercise was considerably (FDR < 0.1) and positively correlated with pattern HMT expression in 30/33 most cancers varieties, for instance breast most cancers (Fig 4F). E2F1 was the ninth transcriptional regulator whose exercise most strongly correlated with whole HMT expression in a pan-cancer evaluation (S6 Desk; Spearman’s rho = 0.464, 2.56%ile). Curiously, this relationship was way more notable in most cancers than wholesome tissue (S7 Desk). Altogether, these outcomes recommended that Rb exercise represses HMT transcription and that this exercise is conserved in C. elegans and people.

In Rb-mutant cancers NNMT is down-regulated downstream of HMTs

The previous outcomes recommended that variation in HMT ranges throughout cancers is related to variability in Rb and E2F exercise and their results on transcription. NNMT expression was considerably lowered in RB1-mutant tumours (Fig 4D), in keeping with its anticorrelation to HMTs. A key query is how HMT exercise is likely to be coordinated with NNMT. One doable situation is that HMTs and NNMT is likely to be coordinated by a variable regulator (e.g., Rb) that has reverse results on the transcription of HMTs and NNMT. Alternatively, NNMT expression is likely to be regulated to compensate for current variability in HMT ranges, or vice versa (Fig 4G).

We examined this by modelling HMT expression and NNMT expression in response both to RB1 mutation standing alone, or in response to the mixed results of RB1 and both HMTs or NNMT. Inclusion of NNMT within the mannequin had no influence on the statistical relationship of HMT expression and RB1 mutation standing (Fig 4H). That is in keeping with a direct transcriptional regulation of HMTs by Rb/E2F and demonstrates that expression of HMTs in these tumours is unbiased of NNMT expression. Conversely, together with HMT expression when predicting NNMT expression abrogates the detrimental relationship between NNMT expression and RB1 mutation standing in tumours (Fig 4H). That is proof that in RB1-mutant cancers, NNMT is just not straight down-regulated by lack of Rb. As an alternative, lowered NNMT expression is a secondary impact of elevated HMT ranges.

Altogether, these analyses assist the speculation that variability in E2F pathway exercise drives variation in HMT expression in most cancers and that this in flip impacts the expression and exercise of the NNMT methyl sink.

Artificially lowering SAM/SAH ratio results in Nnmt down-regulation in mammalian cells

We wished to ascertain whether or not the transcription of NNMT responded to alterations within the SAM/SAH ratios. We recognized 2 RNA-seq datasets the place the SAH hydrolase Ahcy was both knocked out or inhibited pharmacologically in mouse embryonic fibroblasts (MEFs; [48]) and rat hepatic stellate cells [49]. Ahcy loss-of-function causes each SAM and SAH to build up and strongly reduces the SAM/SAH ratio [48,50].

We discovered that lowering Ahcy exercise both via Ahcy deletion or pharmacological inhibition strongly suppressed Nnmt expression in each mouse and rat cells (S12A Fig). We additionally examined how HMT expression responded to synthetic discount of the SAM/SAH ratio. In distinction to the down-regulation of Nnmt, we noticed a average up-regulation of HMT expression in response to lack of Ahcy exercise in mouse cells and no constant change in rat cells (S12B Fig).

The excessive variety of replicates within the mouse information gave us adequate energy to seek for transcriptional regulators that may management the response of Nnmt expression to synthetic perturbation of the SAM/SAH ratio. Glyr1 was the transcriptional regulator that was most importantly activated by both Ahcy deletion or Ahcy inhibition within the mouse dataset (S13A Fig and S8 Desk), implying that its exercise depends upon the SAM/SAH ratio. Prompted by this, we examined the TCGA human major tumour dataset to check whether or not GLYR1 is likely to be concerned in regulating NNMT in human cells. Of 351 transcriptional regulators, GLYR1 exercise had each the strongest detrimental correlation with NNMT expression and the strongest optimistic correlation with whole HMT expression (S6 Desk). Taking estimated GLYR1 exercise into consideration weakened the statistical relationship between NNMT and HMTs within the TCGA greater than another transcriptional regulator (S13B Fig). GLYR1 exercise was elevated in RB1 mutant cancers (S13C and S13D Fig) and taking GLYR1 exercise into consideration abolished the detrimental statistical relationship between RB1 mutation and NNMT expression (S13E Fig). These outcomes indicate a job for GLYR1 in suppressing NNMT transcription in response to a low SAM/SAH ratio.


Utilizing information from tens of hundreds of human samples, right here we demonstrated that HMT expression was strongly anticorrelated to the exercise of two pathways identified to devour extra methyl teams (often known as methyl sinks): synthesis of 1MNA by the enzyme NNMT in cancers and manufacturing of PC by PEMT in wholesome tissues. One doable interpretation of this relationship is that HMTs might also act as a methyl sink. Variation in HMT exercise thus is likely to be correlated to the extent to which various methyl sink pathways function: excessive HMT exercise related to low NNMT/PEMT exercise and vice versa. Nevertheless, we discovered no proof that this variation in HMT exercise had an impact on transcription. Under we talk about the implications of those outcomes for understanding the roles of histone posttranslational modifications.

We’ve got proven a robust anticorrelation between HMT ranges and NNMT exercise in most cancers. These outcomes match with earlier findings that adjustments in NNMT expression may modulate histone methylation [13,1618]. This was beforehand argued to be on account of a passive impact of NNMT exercise on mobile methylation potential through the SAM/SAH ratio. Nevertheless, we confirmed that this relationship corresponds to variations in expression of HMTs and NNMT. We argue that that is higher defined by the speculation that elevated histone methylation exercise ends in a lowered SAM/SAH ratio. We introduced proof that NNMT expression adjustments downstream of HMT expression in Rb-mutant cancers. One doable mechanism for this could be a transcriptional response of NNMT to the SAM/SAH ratio. Certainly, we confirmed that artificially lowering the SAM/SAH ratio by perturbing Ahcy perform led to decreased Nnmt expression in mouse and rat cells. From these experiments, the transcriptional regulator GLYR1 [51] emerged as a possible hyperlink between the SAM/SAH ratio and NNMT expression, as we confirmed GLYR1 doubtless mediates the connection between HMT ranges and NNMT expression in human major tumours. GLYR1 has been proven to be recruited, through H3K36me3, to the our bodies of transcribed genes [52] the place it most frequently promotes gene expression of targets [51]. Since we present that NNMT is repressed when GLYR1 exercise is excessive, this means it could be an oblique regulator of NNMT. Sooner or later, it will likely be attention-grabbing to check whether or not GLYR1 itself senses the SAM/SAH ratio or whether or not additional intermediate elements are concerned.

Our outcomes point out that coordinated HMT expression is managed transcriptionally by the exercise of the Rb/E2F pathway, such that E2F1 concurrently prompts a number of HMTs related to each transcriptional activation and repression. It’s doable that the coordinated change in HMT exercise has a perform in buffering the SAM/SAH ratio. Nevertheless, this could predict that HMT expression needs to be lowered when the SAM/SAH ratio is lowered. We didn’t observe this; certainly, HMT exercise was truly considerably elevated when the SAM/SAH ratio was artificially lowered by perturbing Ahcy perform. Thus, we favour the speculation that adjustments within the SAM/SAH ratio are an essential consequence of alterations in HMT exercise somewhat than their major perform. E2F exercise will increase when cells enter S-phase, which correlates to a requirement to introduce methylation marks onto newly synthesised histones [53]. One attention-grabbing risk subsequently is that the coordinated transcriptional regulation of a number of HMTs by Rb/E2F is required to keep up epigenetic landscapes via cell division.

The connection that we now have found between whole HMT exercise and the exercise of mobile methyl group sinks means that sustaining a constant exercise of methyl sink pathways is significant for mobile homeostasis. The significance of this exercise is likely to be in buffering mobile methylation potential by changing SAM to SAH to keep up the SAM/SAH ratio. Moreover, SAH is required to assist the transsulphuration pathway, which is the cell’s solely pathway to de novo synthesise cysteine and downstream metabolites (e.g., glutathione). In major tumours entry to cysteine is restricted and cells could also be compelled to depend on transsulphuration [11,54]. Nevertheless, cultured cells get pleasure from considerable cysteine equipped in ceaselessly replenished tradition medium. In assist of this notion, we noticed that the HMT/NNMT relationship is far stronger in major tumours than cultured cells. Even so, the connection remains to be evident in most cancers cell traces and we word that within the CCLE metabolomics information, the primary and third metabolites whose ranges most strongly correlate positively to whole HMT expression are cystathionine and alpha-hydroxybutyrate (Fig 1A), each attribute markers of transsulphuration [55].

A number of research have positioned the SAH hydrolase AHCY within the nucleus related to chromatin, arguing that it maintains an area atmosphere conducive to histone methylation [5659]. Nevertheless, AHCY-catalysed SAH hydrolysis is reversible and is thermodynamically favoured solely when the breakdown merchandise are quickly metabolised [55]. This speculation would subsequently predict nuclear exercise of both the methionine cycle or the transsulphuration pathway. Provided that all the enzymes required for transsulphuration and glutathione synthesis are annotated as having nuclear localisation (S9 Desk) within the Human Protein Atlas [60], the existence of a nuclear transsulphuration pathway fuelled by histone methylation, probably supplying cysteine for nuclear glutathione synthesis, is an attention-grabbing risk for future exploration.

The histone code speculation proposed that particular histone modifications have direct and instructive results on transcription [2]. Nevertheless, there are a lot of documented examples the place histone-modifying enzymes and the histone marks that they introduce have an affect on mobile states that aren’t on account of adjustments in transcription [61]. Our outcomes present one other such instance. We confirmed that a lot noticed variation in HMT expression is related to metabolic enzymes with capabilities far faraway from gene regulation. It’s conceivable that the metabolic penalties of HMT expression are unbiased of their catalytic exercise in methylation of histone residues. Nevertheless, we discovered that the degrees of many histone modifications are inversely correlated with NNMT/PEMT expression. Thus, essentially the most easy implication is that the metabolic penalties of HMT expression are on account of their catalytic exercise, which consumes SAM. Importantly, this variation doesn’t have any detectable influence on transcription, whilst a by-product. These outcomes don’t essentially contradict a job for HMTs in instructing transcriptional regulation. For instance, it could be that the variation we observe happens inside a spread that doesn’t have an effect on transcription. Alternatively, the adjustments that happen in histone methylation ranges at specific genes might require different adjustments, comparable to mixtures of histone marks or particular transcription elements, so as to result in transcriptional responses.

Taken collectively, our outcomes recommend that histone methylation impacts mobile metabolism, unbiased of the position of histone methylation marks in regulating transcription. Histone proteins advanced in archaea, the place they’ve a restricted position in transcriptional regulation and there’s little proof of posttranslational modifications comparable to methylation [62]. It’s attention-grabbing to invest whether or not the metabolic penalties of histone posttranslational modifications may predate their extra acquainted position in transcription.


RNA-seq information

RNA-Seq information have been downloaded from the GTEx information portal for GTEx v8. Knowledge have been downloaded as uncooked counts. “Harmonised” (hg38) RNA-seq information have been downloaded for TCGA initiatives utilizing the “TCGAbiolinks” package deal in “R” as uncooked counts. CCLE RNA-seq learn counts have been downloaded from the DepMap obtain portal in August 2021 (model: DepMap Public 21Q3).

Uncooked counts have been subjected to a median-ratio normalisation (MRN) previous to all analyses. The MRN was carried out utilizing the “DESeq2” package deal in “R” [63]. Normalisations have been utilized each individually for every tissue or most cancers sort cohort and throughout all samples inside every database. Normalised pseudocounts have been obtained by changing uncooked counts information to a DEseq2DataSet object utilizing the DESeqDataSetFromMatrix() perform, making use of the estimateSizeFactors() perform to the ensuing dds object, after which retrieving the normalised pseudocounts with the perform counts() with normalised = TRUE. All correlations introduced are based mostly on these MRN-normalised pseudocounts.

For TCGA most cancers sort analyses, we solely thought of samples annotated as Main Tumours, besides the place we explicitly word in any other case (e.g., adjoining regular tissue samples).

For the CCLE, we excluded samples from Main Illness varieties with fewer than 20 cell traces.

We restricted our analyses to GTEx tissues with at the very least 100 samples or TCGA most cancers varieties with at the very least 35 Main Tumour samples.

Correlation distributions

To make sure equal illustration of every tissue or most cancers sort when combining varieties throughout a database, we randomly sampled 100 (GTEx), 36 (TCGA), or 20 (CCLE) from every tissue. Combining uncooked gene expression information for tissues or most cancers varieties might introduce artifacts even when correcting for common tissue/most cancers gene expression, as high-expressing tissues/cancers should still have higher variance within the absolute worth of the residuals. To account for this, we ranked the pattern gene expression pseudocounts for every gene throughout the pattern chosen for every tissue or most cancers sort. We then mixed the ranks for the chosen samples throughout tissue varieties, utilizing the ranks rather than the uncooked residuals; this gave us 4,800 samples for GTEx, 1,188 samples for TCGA, or 460 samples for the CCLE. The Spearman’s correlation was then computed throughout these aggregated ranks. Because the end result varies barely relying on the random sampling inside every tissue, we repeated this course of 100 occasions and plotted the median correlation for every gene.

For analyses inside a single tissue sort, such a normalisation was not required and we merely correlated the uncorrected pseudocounts for our chosen gene towards all others, utilizing all samples obtainable within the cohort.

For gene units comparable to HMTs, we added collectively pseudocounts for every pattern for all the genes earlier than conducting the evaluation above.

When computing genome-wide correlations, we didn’t right for any potential confounders. The place explicitly famous that values have been corrected for confounding variables within the textual content, we corrected for the next variables in GTEx: age, intercourse, velocity of demise, ischemic time, and sequencing batch. For TCGA we corrected for age, race, intercourse, tumour stage, and sequencing centre.

Proliferative index

The PI was calculated with the “ProliferativeIndex” package deal in “R” [41]. Briefly, the whole dataset throughout all tissues or most cancers varieties was normalised by MRN and variance stabilising transformation utilizing the varianceStabilizingTransformation() perform of “DESeq2.” Following the normalisation, the PI was calculated by making use of the readDataForPI() perform with a randomly chosen gene specified within the modelIDs argument, then working calculatePI() on the ensuing object.

ENCODE ChIP-seq information processing

We recognized publicly obtainable histone methylation ChIP-seq information from grownup human affected person samples from the ENCODE venture for tissues with a robust HMT-PEMT correlation within the GTEx information (prime 2.5% in each reciprocal correlation distributions) that had at the very least 3 samples with obtainable RNA-seq information per histone mark. This gave us 19 to 23 samples (relying on the histone mark) from 5 tissues: the esophagus muscularis, the gastroesophageal sphincter, the esophagus squamous epithelium, the sigmoid colon, and the spleen. Knowledge for these samples have been obtainable for 4 distinct histone methylation marks: H3K4me3, H3K9me3, H3K27me3, and H3K36me3. For the corresponding ENCODE samples, the next processed information recordsdata have been downloaded from the ENCODE information portal: ChIP sign fold change over management (as bigwig file) and pseudoreplicated peaks (in mattress narrowPeak format). RNA-seq recordsdata have been downloaded as uncooked counts for every pattern. RNA-seq counts for all samples from all tissues have been pooled and MRN-normalised as described above to yield pseudocounts. File names and experiment DOIs are listed in S4 Desk.

Fold change over management ChIP-seq sign recordsdata downloaded from ENCODE have been transformed to log2 fold change over management. This was carried out through the use of the bigwigCompare perform of the command line package deal “deepTools” (v3.5.0) [64] to match the fold change file towards a man-made bigWig file with a flat sign of 1 throughout all chromosomes, utilizing the argument—operation log2.

To pick genomic areas during which to mannequin ChIP-seq sign by gene expression, we appeared for areas marked by a peak in at the very least 75% of samples (for H3K9me3, H3K27me3, H3K36me3) or in 100% of samples (H3K4me3, on account of higher reproducibility of peak overlaps throughout samples for this mark). We imported the pseudoreplicated peak recordsdata into “R” and used the countOverlaps() perform of the “GenomicRanges” package deal towards the coordinates of genomic areas of curiosity to find out the variety of peaks that overlapped that area (e.g., a particular promoter) in every pattern. We then excluded areas with 0 overlapping peaks in >25% of samples throughout tissues (or any samples for H3K4me3).

The coordinates of promoters and gene our bodies have been generated utilizing the “TxDb.Hsapiens.UCSC.hg38.knownGene” package deal in “R” with the genes() and promoter() capabilities (we used the default settings by which the promoter() perform returns home windows from 2,000 bp upstream to 200 bp downstream of every gene’s TSS). The coordinates of repeat areas have been obtained utilizing the hg38 “rmsk.txt” file downloaded from

In an effort to quantify sign in our chosen genomic areas for every pattern, first we exported the genomic coordinates of areas of curiosity from “R” in mattress file format. We then used the computeMatrix perform of “deepTools” with the pattern log2 fold change bigwig file and the genomic areas mattress file as inputs. We set the arguments—binSize and—regionBodyLength to be equal (often at 100); this ends in an output of a single quantity for common log2 fold change over management ChIP sign for every genomic area.

We went on to mannequin this ChIP sign at every genomic vary by expression of chosen gene. Previous to modelling, we filtered out peaks with a mean ChIP sign throughout samples that fell beneath the extent of the enter management. From any particular person evaluation, we first excluded any tissue sort that had fewer than 3 samples obtainable for that individual mixture of histone mark and tissue, as they may not be successfully modelled.

As a number of ENCODE samples from totally different tissues typically got here from the identical particular person (8 particular person donors), we used linear blended results fashions utilizing the “lmer” package deal in “R” to regulate for this lack of true pattern independence as follows. We first corrected the pattern chosen gene expression for tissue and donor of origin by becoming a linear blended results mannequin explaining log10 RNA-seq pseudocounts with tissue as a set impact and donor as a random impact, extracting the gene expression residual from the mannequin for every pattern utilizing the residuals() perform. We then went on to suit a linear blended results mannequin for every genomic area explaining the log2 fold change over management sign by corrected gene expression residual and tissue as fastened results and donor as a random impact. From this mannequin, we extracted the mannequin t-value utilizing the abstract() perform as a measure of the explanatory energy of gene expression on ChIP sign at that locus and plotted boxplots of all the t-values for every histone mark and kind of genomic area examined.

To generate a null distribution to match towards, we repeated the identical modelling process with expression of 1,000 random genes from our gold commonplace set of ubiquitously expressed genes (excluding histone methyltransferases and demethylases). For every genomic area, we then took the imply t-value for the 1,000 random genes as our null distribution that’s plotted alongside that of the chosen gene on the boxplot; p-values have been obtained from a two-tailed paired Wilcoxon check of the noticed t-value for our chosen gene for every area towards the imply t-value for that area throughout 1,000 random genes.

To measure the ChIP-seq H3K4me3 peak width, we used the H3K4me3 pseudoreplicated peak narrowPeak mattress recordsdata downloaded from ENCODE (S4 Desk). For every genomic area of curiosity, we recognized overlapping peaks utilizing the findOverlapPairs() perform of the “GenomicRanges” package deal in “R.” We then calculated the width of those peaks from the beginning and finish coordinates of peak calls from the mattress file utilizing the width() perform. Be aware we calculated the whole width of the overlapping peak and never solely the half which overlapped the area of curiosity. We then aggregated the height widths for every genomic area within the case of a number of peaks to supply a single determine for the sum whole width of known as peaks overlapping that area. We then went on to mannequin peak width as described above for ChIP-seq sign.

When modelling or correlating gene expression from promoters, we restricted the evaluation to genes with expression detected in each pattern.

NCI60 ChIP-seq information processing

To probe the connection of histone methylation ranges in chromatin to NNMT in most cancers, we used the NCI60 [65], a panel of 60 most cancers cell traces with related RNA-seq, microarray, proteomics, and ChIP-seq information.

NCI60 ChIP-seq information have been taken from [65]. We downloaded uncooked sequencing reads from the NIH’s Sequence Learn Archive, utilizing the command line instrument “SRAtoolkit” (SRA identifier numbers present in S5 Desk). We then aligned the reads to the human genome (hg38) utilizing the “bowtie2” package deal [66]. Every experiment was obtainable as 2 replicates, along with an enter management pattern. Some sequencing replicates consisted of single-reads, whereas others have been paired-end experiments. To ensure that all analyses to be comparable, we aligned solely the ahead reads from paired-end experiments. The output from “bowtie2” was saved to a.sam file and transformed to a.bam file utilizing the “samtools” package deal (model 1.16.1). bam recordsdata have been then sorted by question identify utilizing the kind perform from “samtools.”

We used the “Genrich” package deal ( to name ChIP-seq peaks on the sorted.bam recordsdata, as “Genrich” can take 2 replicate experimental recordsdata, along with a management file, and carry out an built-in peak name relative to the management. For all marks (H3K4me3, H3K9me3, and H4K20me3), we used the next settings, equivalent to “broad” peaks: -q (max FDR) 0.1, -g (max distance between vital websites) 400. Moreover, for H3K4me3 we used the default settings, equivalent to “slim” peaks. We then imported these peak name recordsdata into “R” and used the cut back() perform of the “GenomicRanges” package deal to mix the peaks known as beneath broad and slim settings for H3K4me3.

For every experimental replicate, we used the bamCompare perform of “deepTools” with—operation log2 to return a sign file in bigWig format for log2 fold change over management. We recognized genomic areas of curiosity (promoters, gene our bodies, and repeats) as described above for ENCODE, limiting to genomic websites with a known as peak in >2/3 of samples. We then ran the computeMatrix perform of “deepTools” with the pattern log2 fold change bigwig file and the genomic areas.mattress file as inputs, with equal bin measurement and physique size as described above, to return a single determine for common log2 fold change for every area. At this level, we took the common of the two replicates and carried that ahead into our modelling method.

We downloaded uncooked RNA sequencing reads utilizing “SRAtoolkit” (SRA identifier numbers in S5 Desk) and aligned them to the human genome (hg38) utilizing “bowtie2” utilizing default settings. We calculated learn counts per gene by counting reads overlapping exons utilizing the summarizeOverlaps perform of the “GenomicAlignments” package deal in “R.” Learn counts have been then MRN-normalised earlier than use in modelling. We downloaded NNMT / NNMT SWATH mass-spectrometry proteomics values and 5-microarray gene expression Z-scores for NCI60 cell traces from CellMinerCDB (

For every genomic area, we modelled the common log2 fold change over management sign utilizing a generalised linear mannequin with the glm() perform of the “stats” package deal in “R.” Within the mannequin, we included the expression of our gene/protein of curiosity (pattern RNA-seq pseudocounts, protein degree, or mRNA Z-score), along with the cell line tissue of origin and an interplay time period between tissue and expression. We extracted and plotted the mannequin t-values for expression as described above.

Mouse embryonic fibroblast and rat liver RNA-seq information processing

For MEFs, we downloaded uncooked reads with “SRAtoolkit” as described above from GEO Accession GSE126851. We used “bowtie2” to align reads to the Mus musculus genome (construct GRCm39) and used “featureCounts” from the “Subread” package deal with the GRCm39 gtf annotation downloaded from to generate gene-level learn counts. These have been then normalised with DESEq2 previous to plotting or TR exercise estimation. For rat liver samples handled with DZNep, gene degree uncooked counts have been downloaded from the Supplementary Materials of GEO Accession GSE121736 and handled as above. M. musculus and Rattus norvegicus HMT genes have been taken to be orthologues of human HMTs and are present in S1 Desk.

Bayesian iterative reweighting evaluation of multi-mapping ChIP-seq reads in NCI60

Each ChIP-seq sign and RNA-seq expression ranges from repetitive components are tough to quantify precisely on account of ambiguous multi-mapping of sequencing reads to extremely related genomic areas. These reads are sometimes discarded in information processing pipelines (as within the ENCODE pipeline). Nevertheless, for the NCI60, we used a Bayesian iterative reweighting method (“SmartMap”) to apportion multi-mapping reads to particular person genomic loci, offering extra correct estimates of ChIP-seq sign at repetitive components [67]. We carried out this evaluation for H3K9me3 and H4K20me3.

The two replicates within the NCI60 ChIP-seq information, in addition to being paired and unpaired, have totally different learn lengths. Moreover, the enter controls are unpaired reads with shorter learn lengths (150 bp). Within the case of multi-mapping reads, higher learn size was prone to have an effect on the chance of distinctive mapping and so have an effect on the validity of comparisons to the enter management. As such, we solely made use of a single replicate for every histone mark and cell line, specifically the unpaired replicate with shorter learn lengths that matched the enter management. We tailored the SmartMapPrep script from the “SmartMap” package deal to course of uncooked single-end reads, downloaded utilizing “SRAtoolkit” as described above, for therapy and enter controls, earlier than utilizing the SmartMap perform for a single reweighting iteration as default and as beneficial by the authors. The output from “SmartMap” is a bedgraph file. We transformed the bedgraph recordsdata to bigwig recordsdata with the UCSC “bedGraphToBigWig” utility. Log2 fold sign over enter management was discovered utilizing deepTools bigwigCompare with the default “log2”—operation selection.

For peak calling, we used “deepTools” bigwigCompare with the—operation subtract setting to take away the sign from the suitable enter management monitor from every monitor. We then transformed these bigWig recordsdata again to bedgraph recordsdata with the UCSC bigWigToBedGraph utility. We used MACS3 (v3.0.0) [68] for peak calling, as it could name peaks from a bedgraph file utilizing the bdgbroadcall perform.

As above we analysed particular person repetitive components that have been marked by a peak in at the very least 40 of 60 cell traces. Repetitive components have been recognized from the “rmsk.txt” file as described above, except for HERVs, which have been taken from the annotation included within the Telescope package deal (see beneath). Customized mattress recordsdata have been created with the weather to be analysed and sign was quantified throughout the whole component utilizing “deepTools” computeMatrix as described above.

Estimation of HERV expression within the NCI60

We additionally used a separate Bayesian reweighting method (“Telescope”) to estimate locus-specific expression estimates from a set of HERVs [69]. We aligned uncooked RNA-sequencing reads (downloaded with “SRAtoolkit” as above) with “bowtie2” with choices—very-sensitive-local and -k 100 (permitting as much as 100 alignments per learn), as beneficial by the “Telescope” package deal authors. The ensuing bam recordsdata have been processed in “Telescope” utilizing the telescope assign perform name. The HERVs annotation file “HERV_rmsk.hg38.v2.gtf” was downloaded from the “telescope_annotation_db” repository on GitHub ( We analysed particular person HERVs marked by ChIP peaks (recognized utilizing “SmartMap”) in 40 of 60 cell traces and which had expression detected in at the very least 30 of 60 cell traces. When utilizing HERV expression as a response variable in a linear mannequin, we used detrimental binomial generalised linear mannequin (with the glm.nb() perform)on account of typical overdispersion of the information.

Correlations inside histone methyltransferase genes

To probe co-expression of HMT genes, we first corrected the expression values for confounding variables as described above. Moreover for the GTEx pan-tissue evaluation, we corrected for donor ID as a random variable inside a linear blended results mannequin, to account for the truth that when evaluating throughout tissues a number of samples can originate from the identical particular person donor. Corrected residuals have been rank-percentile reworked inside every tissue or most cancers sort, earlier than 100 (GTEx) or 36 (TCGA) samples have been chosen from every and mixed earlier than Spearman’s correlations among the many rank-percentile reworked values have been computed throughout the grouped samples. The method of sampling was repeated 100 occasions and the median Spearman correlation from the 100 iterations was taken for plotting. We excluded very lowly expressed HMT genes from plots by filtering in keeping with a geometrical imply expression throughout all samples of at the very least 100 pseudocounts; this accounts for barely totally different numbers of HMT genes within the GTEx/TCGA plots proven in S10 Fig.

Correlations of HMT genes to random genes have been computed as above with 100 random genes from our gold commonplace set of ubiquitously expressed genes for every HMT. Values computed have been then pooled.

Community plots have been ready from pan-cancer or pan-tissue correlation matrices by first filtering out edges with correlations of magnitude lower than 0.2. Node measurement was based mostly on its diploma and edges have been weighted by the sq. of the magnitude of the correlation. Community evaluation was carried out in “Gephi” (model 0.10.1), with the next visualisation properties: ForceAtlas2 format, edge weight vary 0.1 to 2.0 and attraction 30 for TCGA and 10 for GTEx.

CeNDR C. elegans information and processing

We downloaded RNA-seq uncooked counts and TPM-normalised values for C. elegans strains from the C. elegans Pure Variety Useful resource (CeNDR; [72]) from Gene expression Omnibus, accession quantity GSE186719. TPM values have been used just for estimating pattern ages utilizing the “RAPToR” instrument (see beneath). We MRN-normalised the uncooked counts and used these counts for all different analyses. The uncooked counts have been transcript-level counts; these have been collapsed all the way down to gene-level counts previous to all analyses.

We used the “RAPToR” package deal in “R” [73] to deduce the age of the samples in keeping with the creator’s directions. We used the Cel_YA_2 reference collection from the “wormRef” package deal. To acquire age-corrected residuals for gene expression, we fitted a spline with 6 levels of freedom utilizing the clean.spline() perform of the “stats” package deal in “R” to foretell log10 pseudocounts from inferred age, taking the residuals from the spline with the residuals() perform.

C. elegans HMT genes have been chosen in keeping with their gene descriptions on WormBase (model WS287). The listing of C. elegans HMTs may be present in S1 Desk. The overwhelming majority of strains have been represented by 3 unbiased RNA-seq samples. We took the imply of the age-corrected residuals for every pressure to plot scatterplots and compute HMT correlations for the heatmap. C. elegans orthologues of human HE cluster genes have been decided utilizing OrthoList 2 [74]; a gene was annotated as an orthologue if the orthology relationship was current in at the very least 3 of the 6 databases compiled in OrthoList 2.

C. elegans lin-35 mutants and ChIP-seq information

To establish RNA-seq datasets from lin-35 mutants, we searched the Gene Expression Omnibus for lin-35 and located 2 research; [45] (GEO accession GSE62833) for lin-35 mutant or wild-type L3 larvae and [46] (GEO accession GSE155190) for L1 larvae. We downloaded the uncooked information and aligned it to the C. elegans genome (model WS276) utilizing “bowtie2” with default settings, earlier than acquiring gene degree counts utilizing summarizeOverlaps() from “GenomicRanges” as described above for NCI60 after which MRN-normalising the ensuing counts with “DESeq2” as described above. In each research whole HMT counts have been elevated in lin-35 mutants on common by 12% to fifteen%. To evaluate statistical significance, we analysed the two research collectively, performing a 2-way ANOVA for whole HMT counts with genotype and developmental stage as explanatory variables.

LIN-35 ChIP-seq information was obtained from [44], S1 Desk. We used “GenomicRanges” in “R” to find out the space from reported vital LIN-35 peaks to the TSS of HMT genes. For “different genes,” we excluded all genes which have 0 expression in any pattern; we then randomly sampled 1,000 of the remaining roughly 11,000 genes.

Rb-mutant cancers within the TCGA and CCLE

Mutation requires TCGA samples have been downloaded as MAF recordsdata utilizing the “TCGABiolinks” package deal in “R” [75]. We recognized samples with a reported mutation in RB1 which was both a missense mutation, nonsense mutation, frame-shifting insertion, in-frame deletion, or a frame-shift deletion. Ten most cancers varieties had at the very least 10 RB1-mutant samples with obtainable RNA-seq information; BLCA, BRCA, CESC, COAD, HNSC, LIHC, LUSC, LUAD, SARC, and UCEC. We assumed samples had wild-type RB1 if they’d mutations known as in different genes within the MAF recordsdata however none known as in RB1.

We rank-percentile reworked all samples in these most cancers varieties (RB1 mutant and wild sort collectively) after which sampled 10 RB1-mutant and 10 wild-type cancers from every most cancers sort, taking the median rank-percentile of every random pattern. We repeated this sampling course of 1,000 occasions and plotted the medians from these samples for every group, evaluating the medians with a t check.

To mannequin HMT or NNMT expression by RB1 mutation standing and/or counterpart expression, we rank-percentile reworked HMT and NNMT expression by most cancers sort after correction for confounding variables. We then carried out pan-cancer sampling of RB1-mutant or wild-type cancers as above, combining all samples and becoming a linear mannequin with HMT/NNMT expression as response variable and both RB1-mutation standing alone as an explanatory variable or along with counterpart expression. We then extracted the t-values from the assorted linear fashions for the statistical affiliation of RB1-mutation to both HMT or NNMT expression.

For the CCLE, we downloaded the mutation calls utilizing the depmap_mutationCalls() perform of the “depmap” package deal in “R.” We then filtered RB1 mutations by whether or not they have been known as as deleterious or not. We discovered 90/1,236 cell traces had deleterious RB1 mutations, of which 43 have been annotated as lung most cancers cell traces, 32 particularly small cell lung most cancers. As small cell lung cancers had the next expression of HMTs than different lung most cancers subtypes, we carried out a 2-way ANOVA with genotype and lung most cancers subtype as explanatory variables.

Transcriptional regulator exercise estimation

To estimate transcriptional regulator (TR) exercise in GTEx and TCGA samples, we used the “decoupleR” package deal in “R” [47] to deduce TR exercise from RNA-seq samples. “decoupleR” requires a gene regulatory community (GRN) to make use of as a foundation for TR exercise inference; we used the “dorothea” package deal in “R” [76] beforehand developed by the identical authors. The “dorothea” package deal consists of 2 totally different human GRNs; 1 common and 1 for cancers. We used the final GRN for estimating GTEx pattern TR exercise and the most cancers GRN for TCGA pattern TR exercise. In an effort to improve our confidence within the estimates, we excluded the bottom confidence TR-target interactions that have been the end result solely of in silico predictions. We then excluded any TRs that had fewer than 10 goal genes remaining by which to deduce their exercise and likewise excluded any that have been in our listing of HMTs. This left us with 351 transcriptional regulators. Earlier than working “decoupleR,” we additionally weighted targets by the arrogance of the interplay, changing the arrogance reported by “dorothea” (letters A–D after low-confidence interactions had been eradicated) into an integer worth [14] and utilizing its inverse as an interplay weight.

For estimating TR exercise in MEFs, we used the dorothea_mm GRN from the “dorothea” package deal, which relies on the human GRN. To seek out TRs whose targets have been enriched amongst differentially expressed genes within the Ahcy loss-of-function datasets, we carried out a differential expression evaluation with “DESeq2” in “R.” We used the outcomes() perform to acquire the check statistics for all genes. We then used these check statistics for all genes as enter for “decoupleR” utilizing the multilevel marketing technique to detect enriched TRs. In an effort to plot Glyr1 exercise in S13 Fig, we estimated TR actions individually for every pattern with “decoupleR” utilizing counts normalised with “DESEq2” as enter.

ENCODE transcription issue binding enrichment

We downloaded the ENCODE transcription issue binding web site profiles from the “Harmonizome” database internet portal [77]. For the genes certain by every transcription issue, we carried out a Fisher’s actual check for enrichment of HMT genes among the many certain genes. Odds ratios and p values are extracted from the Fisher’s actual check; p-values reported within the volcano plot are uncooked and uncorrected.



Please enter your comment!
Please enter your name here

Most Popular

Recent Comments