DNA Methylation Repels Binding of HIF Transcription Factors to Maintain Tumour Immunotolerance

Background Hypoxia is pervasive in cancer and other diseases. Cells sense and adapt to hypoxia by activating hypoxia-inducible transcription factors (HIFs), but it is still an outstanding question why cell types differ in their transcriptional response to hypoxia. Results Here, we report that HIFs fail to bind CpG dinucleotides that are methylated in their consensus binding sequence, both in in vitro biochemical binding assays and in vivo studies of differentially methylated isogenic cell lines. Based on in silico structural modelling, we show that 5-methylcytosine indeed causes steric hindrance in the HIF binding pocket. A model wherein cell-type-specific methylation landscapes, as laid-down by the differential expression and binding of other transcription factors under normoxia control cell-type-specific hypoxia responses is observed. We also discover ectopic HIF binding sites in repeat regions which are normally methylated. Genetic and pharmacological DNA demethylation, but also cancer-associated DNA hypomethylation, expose these binding sites, inducing HIF-dependent expression of cryptic transcripts. In line with such cryptic transcripts being more prone to cause double-stranded RNA and viral mimicry, we observe low DNA methylation and high cryptic transcript expression in tumours with high immune checkpoint expression, but not in tumours with low immune checkpoint expression, where they would compromise tumour immunotolerance. In a low-immunogenic tumour model, DNA demethylation upregulates cryptic transcript expression in a HIF-dependent manner, causing immune activation and reducing tumour growth. Conclusions Our data elucidate the mechanism underlying cell-type specific responses to hypoxia, and suggest DNA methylation and hypoxia to underlie tumour immunotolerance.

preferentially binds at promoters and HIF2α at enhancers, DNA methylation differences 234 do not determine their binding specificities. 235

DNA methylation directly repels HIF binding in cells 236
To more firmly establish a causal link between DNA methylation and HIF binding, we 237 excluded several confounders. Firstly, since our chromatin state analysis revealed that HIF 238 preferentially binds active enhancers and promoters, which are known to carry low levels 239 of methylation 23 , we performed HIF1β ChIP-bisulfite sequencing (HIF1β ChIP-BS-seq). 240 MCF7 cells were exposed to hypoxia, HIF1β-bound DNA was immunoprecipitated and 241 bisulfite-converted prior to sequencing to uncover its methylation pattern. This revealed 242 that, while methylation levels of input DNA (not immunoprecipitated, bisulfite-converted 243 DNA) were mostly low but with some sites displaying intermediate to high methylation 244 levels, HIF1β-bound DNA was invariably very low in methylation and this at all sites ( Figure  245 2a, Figure S4a). 246 Secondly, since TFs can drive demethylation of their binding sites both passively and 247 actively, we excluded the possibility that DNA fragments bound by HIF would undergo 248 DNA demethylation upon HIF binding. Indeed, HIF1β has previously been shown to 249 actively recruit DNA demethylases 24 . However, HIF1β ChIP-BS-seq in hypoxic ESCs 250 deficient for all DNA demethylases (Tet1, Tet2 and Tet3) showed results identical to those 251 observed in wild-type MCF7 cells: HIF1β-bound DNA was unmethylated compared to 252 input DNA subjected to whole-genome BS-seq (Figure 2b, Figure S4b). 253 Additionally, other (unknown) confounders related to the binding location of HIF, such as 254 chromatin environment or sequence context, may contribute to preferential HIF binding 255 to unmethylated DNA. To exclude this possibility, we generated isogenic murine ES cell 256 lines in which a human HIF1β binding site-encoding DNA fragment was inserted that was 257 either in vitro methylated or not (Figure 2c). Following recombination, the difference in 258 methylation state between both fragments was maintained ( Figure S4c-d). HIF1β ChIP-259 qPCR revealed that methylation was sufficient to induce a 12.4-fold reduction in HIF1β 260 binding in these isogenic cell lines (Figure 2d). 261 Finally, to directly assess methylation sensitivity of HIF binding to unchromatinized DNA, 262 we employed microscale thermophoresis, and tested the binding of recombinant co-263 purified HIF1α-HIF1β and HIF2α-HIF1β heterodimers to double-stranded DNA 264 oligonucleotides containing a methylated or unmethylated RCGTG motif. Importantly, 265 HIF1α-and HIF2α-containing heterodimers both showed a 15-fold higher affinity (KD) for 266 an unmethylated than methylated RCGTG motif, thus confirming that methylation 267 directly repels binding of HIF1α-HIF1β and HIF2α-HIF1β heterodimers (  would be poised to cause severe steric clashes with these two functionally important 275 arginine residues in HIF1α or HIF2α (Figure 2h). 276

DNA demethylation enables ectopic HIF binding
Next, we investigated which parts of the genome are protected from HIF binding by DNA 278 methylation. For this, we compared HIF1β binding in hypoxic wild-type murine ESCs 279 versus ESCs deficient for DNA methyltransferases (Dnmt-TKOs), which lack DNA 280 methylation 26 , using HIF1β ChIP-seq (n=4 replicates for each; for data quality assessment 281 see Figure S5). This revealed a marked increase in the number of HIF1β binding peaks, 282 from 7,875 in wild-type to 9,806 in Dnmt-TKO ESCs (Figure 3a). Whole-genome BS-seq 283 further revealed that, while shared binding peaks were unmethylated in both cell lines, 284 Dnmt-TKO-specific HIF1β binding peaks had high methylation levels in wild-type ESCs 285

DNA methylation represses hypoxia-induced expression of retrotransposons 295
Indeed, a substantial fraction of novel Dnmt-TKO-specific HIF1β binding peaks were found 296 in repetitive genomic regions. Particularly, repeat class analysis revealed a 1.65-fold 297 increase in binding peaks near retrotransposons (2,737 of 7,875 (34.8%) shared peaks 298 versus 1,106 of 1,931 (57.3%) Dnmt-TKO-specific peaks; Figure S6a). Although HIF1β 299 binding events were frequently observed at LINEs and SINEs, only binding at long terminal 300 repeats (LTRs) was enriched over a randomisation of HIF1β binding site positions, and this 301 both for all binding events and those distal to TSS (Figure 3i). The bulk of this increase was 302 ascribable to binding at the 5'-end of endogenous retrovirus K (ERVK) LTR sequences 303 (Figure 3j), with 344 of 1,106 (31%) novel repeat-binding peaks being at ERVKs versus only 304 3% of randomly shuffled HIF1β binding sites. These were mostly at solitary LTRs ( Figure  305 S6b-e). Given that ChIP-seq analyses rely on uniquely-mapping reads, which are 306 SQuIRE, 72% (n=2,781) of all differentially expressed repeat loci exhibited increased 336 expression under hypoxia (P<10 -16 ; Figure S6j). 337

Induction of cryptic transcripts by hypoxia 338
Hypoxia-induced transcripts were, however, often not matching the annotated repeat 339 locus, but extending well beyond their annotated end, with some transcripts 340 encompassing multiple repeat elements. Also, we noticed that for many transcripts, 341 HIF1β binding did not occur in the retrotransposon promoter, while some other 342 transcripts did even not contain a retrotransposon-associated sequence. Similar 343 transcripts were also induced by aza. We therefore refer to these as 'cryptic transcripts' 344 ( Figure S7a). To more accurately quantify them, we developed a novel analysis pipeline, 345 CRyptic Elements' Differential Expression by de Novo Transcriptome Reconstruction 346 (CREDENToR). CREDENToR first performs a de novo transcriptome assembly to define 347 cryptic transcripts and then assigns uniquely-mapping reads to them to quantify their 348 expression. The cryptic transcripts detected by CREDENToR are poorly conserved, often 349 unspliced transcripts, shorter than lincRNAs but expressed at similar levels (see Methods 350 and Figure S7b-g for benchmarking). 351 CREDENToR identified that out of 1,389 differentially expressed cryptic transcripts (1% 352 FDR), 67% were upregulated by hypoxia ( Figure S6k). As expected, focussing on HIF-bound 353 cryptic transcripts revealed an even stronger enrichment, with 82% and 91% 354 (respectively, at 1% and 0.001% FDR) differentially expressed transcripts being 355 upregulated following hypoxia (Figure 4e). HIF binding was enriched at the promoter of 356 hypoxia-induced cryptic transcripts, but far less in those induced by aza ( Figure S6l). 357 Interestingly, significant fractions of cryptic transcripts contained palindromic repeats, or 358 overlapped with other transcripts in the reverse orientation, and could thus produce 359 double-stranded (ds) RNA. HIF-bound cryptic transcripts were twice as likely to generate 360 such dsRNAs (Figure 4f). Together, this suggests HIF binding to leverage cryptic TSS 361 structures within and outside the repeat genome to express dsRNA-generating cryptic 362

transcripts. 363
Cryptic transcript expression was indeed dependent on HIF, as non-HIF-bound cryptic 364 transcripts failed to show induction following hypoxia (Figure 4g). To confirm this, we 365 assessed expression in HIF1B-knockout MCF7 cells. Here, hypoxia failed to upregulate 366 cryptic transcripts, according to both CREDENToR and RepEnrich (Figure 4g, Figure S6m). 367 As expected, aza-induced overexpression was retained, while hypoxia in HIF1B-knockout 368 MCF7 cells failed to increase the effect of aza. Pharmacological activation of HIF using 369 dimethyloxalylglycine (DMOG), a broad-spectrum inhibitor of 2-oxoglutarate-dependent 370 hydroxylases 29 , affected cryptic transcripts similar to hypoxia ( Figure S6m-o). Combined,371 these data indicate that hypoxia triggers HIF binding to unmethylated repeat regions, 372 inducing HIF-dependent expression of cryptic transcripts, most of which are associated 373 with retrotransposons. 374

Hypoxia and repeat transcript expression affect tumour immunotolerance 375
Expression of repeat transcripts has been linked to tumour foreignness 30 , interferon (IFN) 376 response 31-33 and enhanced cytolytic activity 34 , all critical determinants of response to 377 checkpoint immunotherapy. Similar to our own data, such transcripts were shown to 378 increase dsRNA formation. This triggers IFN responses through viral mimicry. Cryptic 379 transcripts induced by HIF could thus contribute to an immune-activated 380

microenvironment. 381
To study this in more detail, we reanalysed expression and DNA methylation data from 382 The Cancer Genome Atlas (TCGA). We classified 5,193 tumours from 14 tumour types as 383 hypoxic or normoxic using an established hypoxia metagene expression signature 35 . We 384 remapped all RNA-seq reads to determine expression of retrotransposon subfamilies 385 using RepEnrich, and also performed de novo transcript assembly to identify on average 386 11,654 non-overlapping cryptic transcripts per tumour type using CREDENToR ( Figure  387 S8a). While TCGA tumours were not exposed to DNA demethylating agents, they did show 388 variation in DNA methylation at TSS of cryptic transcripts. Indeed, although CpGs cryptic 389 transcript promoters showed mostly high methylation levels (median = 80.7%), there was 390 considerable variability (9.2% standard deviation), and one in 10 tumours displayed 391 median levels below 67.3%. Remarkably, and in line with our in vitro data, there was a 392 significant interaction between hypoxia and DNA methylation in determining cryptic 393 transcript expression (P=0.0109), with expression being increased in hypoxic tumours 394 having lower methylation at cryptic transcripts (Figure 5a). At least 1,279 cryptic 395 transcripts showed increased expression of 10-fold or higher (FDR<0.01, Figure S8b). A 396 reanalysis of combined single-cell methylome-and-transcriptome sequencing data from 397 colorectal cancer cells 36 moreover confirmed that cryptic transcript expression and 398 promoter methylation are inversely correlated, and this more strongly in hypoxic than 399 normoxic cancer cells (P=0.032), suggesting that the observed interactions are cancer cell-400 Table 3). 401 In TCGA, this interaction was only detected in tumour types known to respond to 402 immunotherapy 37 (see Methods; P=0.0031 in responsive versus P=0.69 in non-responsive 403 tumours; Figure 5b-c). As expected, responsive tumour types exhibited an increased 404 tumour mutation burden (TMB), elevated immune checkpoint expression, more CD8 + T-405 cells and increased cytolytic activity ( Figure S8c) 38 . Importantly, responsive types also had 406 on average lower methylation at cryptic transcripts and higher cryptic transcript 407 expression than non-responsive types (P<10 -16 for both comparisons, Figure 5d). Single-408 cell RNA-seq analyses (both from 5' and 3' end) highlighted that cancer cells show the 409 highest level of cryptic transcript expression compared to stromal cells, indicating they 410 represent the main source of cryptic transcripts expression ( Figure S8d). In line with our 411 in vitro findings, DNA hypomethylation thus underlies cryptic transcript expression in 412 hypoxic tumours, an effect that was particularly striking in immunotherapy-responsive 413

tumours. 414
Overall, these observations support a model wherein hypoxia-induced cryptic transcripts 415 are tolerated in high-immunogenic tumours, as these are characterized by high immune 416 checkpoint expression, but not in low-immunogenic tumours where their expression 417 would compromise tumour immunotolerance. This suggests that low-immunogenic 418 tumours may need to maintain high DNA methylation levels in cryptic transcripts to 419 downregulate their expression and avoid the induction of tumour immunogenicity. 420

Aza compromises tumour immunotolerance in mice via HIF
To confirm that in low-immunogenic tumours DNA methylation prohibits cryptic 422 transcript expression, we identified 59 such retrotransposons that correlate in expression 423 with cytolytic activity in immunotherapy-responsive tumour types within TCGA. 424 Remarkably, all of these were upregulated in vitro, by hypoxia alone or hypoxia in 425 combination with aza (P<0.05, Figure 5e), suggesting that hypoxia and DNA 426 demethylation can indeed enhance tumour immunogenicity. To confirm this 427 experimentally, we screened several mouse tumour models for their immunogenicity. 428 The orthotopic 4T1 breast cancer model was identified as low-immunogenic. Indeed, 4T1 429 tumours exhibited a low TMB, cytolytic activity, number of CD8 + T-cells and expression of 430 immune checkpoints (Pd1, Pdl1) compared to other models ( Figure S9a). In line with 4T1 431 grafts being low-immunogenic tumours, anti-PD1 treatment failed to affect their growth 432 (-8%, P=0.397), while significantly reducing growth of high-immunogenic tumours, as 433 described previously 39, 40 ( Figure S9b). Importantly, also the expression of cryptic 434 transcripts was lower in 4T1 than in high-immunogenic tumour models ( Figure S9c). 435 Next, we verified in low-immunogenic 4T1 cells whether DNA demethylation upregulates 436 cryptic transcripts in a HIF-dependent manner. In vitro, we observed that, similar to MCF7 437 cells, both hypoxia and aza independently increased cryptic transcript expression, both 438 using CREDENToR and RepEnrich (Figure 6a; Figure S9d). Likewise, aza increased cryptic 439 transcript expression in vivo ( Figure 6b). To confirm that this upregulation was at least 440 partially hypoxia-mediated, we investigated whether tumour hypoxia enhances aza-441 induced cryptic transcript expression. We compared aza-treated 4T1 tumour-bearing 442 mice injected either with control or anti-VEGFR-2 antibody (DC101). While vehicle-443 treated 4T1 tumours were hypoxic in ~40% of the tumour, DC101 further reduced blood 444 vessel density (-35%; P<0.05) and increased hypoxic tumour areas (68%; P<0.05; Figure  445 6c-d). Importantly, this was associated with an increase in cryptic transcripts (+ 9%; 446 P=2.6x10 -16 ; Figure 6e). 447 We then explored whether this increase also compromised immunotolerance. As 448 immunogenicity of cryptic transcripts is mediated via dsRNA formation, we first in 4T1 cells (Figure 6f). In vivo, aza reduced growth of 4T1 tumours (-32%; P=3.0x10 -3 ; 451 Figure 6g), but did not reduce cell proliferation marker expression ( Figure S9e). In 452 contrast, immune activation was enhanced in tumours treated with aza, as activated T-453 cell and natural killer cell signatures were upregulated and myeloid-derived suppressor 454 cell signatures downregulated (Figure 6h). Immunofluorescence of CD8 + T-cells confirmed 455 these changes: while T-cell infiltration was unaffected, the number of activated, granzyme 456 B-positive T-cells increased 2.1-fold (P<0.05; Figure S9f). 457 To verify HIF-dependence of these immunogenic effects, we generated polyclonal 4T1 458 cells deficient for HIF1β by CRISR-Cas9 (4T1 Hif1b-KO ; Figure S9g) and compared these cells 459 to scramble-control 4T1 cells (4T1 Hif1b-scr ), while treating with aza or vehicle. In vitro, loss 460 of HIF1β abrogated hypoxia-induced dsRNA formation and HIF1β-bound cryptic transcript 461 expression both in vehicle and aza-treated cells ( Figure S9d and h), similar to what we 462 observed in MCF7 cells. Also in vivo, 4T1 Hif1b-KO showed reduced cryptic transcript 463 expression compared to 4T1 Hif1b-scr grafts, effects that were limited to HIF-bound cryptic 464 transcripts as expected (Figure 6i; Figure S9i). Of note, 4T1 Hif1b-KO tumours grew much 465 slower than 4T1 Hif1b-scr tumours, presumably because HIF1β also has direct effects on cell 466 proliferation, thus rendering it challenging to disentangle effects on immunogenicity. 467 Nevertheless, aza induced a similar and significant upregulation of cancer testis antigen 468 expression in both cell lines ( Figure S2), suggesting similar treatment efficacy. 469 Interestingly, while 4T1 Hif1b-scr grafts also showed a significantly reduced size when 470 comparing aza to vehicle (46% reduction; P=3.3x10 -6 ), 4T1 Hif1b-KO failed to show as strong 471 a reduction (only 22%; P=7.0x10 -4 , or 1.8-fold less than the scramble effect; P=0.021; 472 Figure 6j). The differential effect of aza in 4T1 Hif1b-KO versus 4T1 Hif1b-scr grafts was also highly 473 significant in an interaction analysis (P<0.0001). Moreover, while the number of activated 474 T-cells increased in 4T1 Hif1b-scr grafts following aza, 4T1 Hif1b-KO grafts failed to show such 475 increase ( Figure 6k). This differential effect of aza, depending on the Hif1b background, 476 was similarly significant in an interaction analysis (P=6.3x10 -3 ). Together, these data 477 provide a mechanistic link between HIF1β binding, DNA methylation and immune 478 repels HIF binding, we thus highlight the importance of DNA methylation profiling, 496 especially in poorly oxygenated tissues. Since tumour hypoxia has long been associated 497 with increased malignancy, poor prognosis and resistance to radio-and chemotherapy 6 , 498 DNA methylation could especially provide insights in the processes underlying 499 therapeutic resistance. For instance, Vanharanta and colleagues recently showed an 500 association between DNA methylation near CYTIP and the survival of disseminating 501 cancer cells 43 . Combined with our observations that DNA methylation directly repels HIF 502 binding, this suggests remethylation of the CYTIP promoter as a viable avenue for 503 decreasing cancer dissemination. 504 Secondly, it has been challenging to identify a guiding principle as to why specific genes 505 are induced upon hypoxia in one, but not the other cell type 9 . Our findings suggest that 506 cell-type-specific TF binding under normoxia causes differences in DNA methylation, which subsequently determine HIF binding under hypoxia and predict the cell-type-508 specific hypoxia response. Importantly, we also confirmed earlier observations that HIF1β 509 binding peaks are characterized by an active, open chromatin structure 11 . This additional 510 requirement for functional HIF1β binding peaks probably explains why each of the RCGTG 511 consensus sequences in the genome cannot serve as an equal HIF binding substrate in 512 normal cells, or upon genetic or pharmacological demethylation. Similar observations 513 were made for other TFs, such as CTCF, for which binding was similarly limited to sites 514 containing a permissive chromatin structure 14, 23 . Importantly, binding specificities for 515 HIF1α versus HIF2α are independent of DNA methylation, but appear to be influenced by 516 chromatin context. This is in line with the identical structure of DNA binding domains of 517 HIF1α and HIF2α; swapping DNA binding domains between both proteins has no influence 518 on their binding profile 44 . Instead, the transactivation domain appears to endow 519 specificity, suggesting that accessory chromatin binding partners govern the differential 520 binding of HIF1α and HIF2α 44 . 521 Thirdly, several publications by now reported how 5-aza-2'-deoxycytidine initiates cryptic 522 TSSs in the repeat genome, leading to expression of cryptic transcripts 31, 45 . Our data add 523 to these findings by demonstrating that cryptic transcript expression is at least partly HIF-524 dependent, while more importantly, hypoxia alone is also capable of inducing their 525 expression. Based on single-cell analyses, we observed this effect to be cancer cell-526 autonomous, consistent with cancer cells being hypomethylated. Our findings reinforce a 527 growing body of evidence that highlights how during evolution transposable elements 528 have copied and amplified regulatory regions throughout the genome 22, 46-50 . Most likely, 529 transposable elements hijacked the transcriptional apparatus of their host to support 530 their germline propagation 51 . In doing so, they copied the associated TF binding site and 531 seeded it at the site of insertion. Transposable elements having binding sites for TFs that 532 are active in the germline, are more likely to hijack these and transpose. Accordingly, HIF 533 is activated in early development, when DNA methylation levels are also low 50, 52 ; 534 ancestral cooption of HIF binding sites by cryptic transcripts to increase their expression 535 is thus plausible. In line with specific TFs preferentially acting on particular 536 retrotransposon subfamilies, we observe enrichment of HIF binding and activation at 537 LTRs, particularly at the LTR of ERVKs. 538 Finally, we uncover an intriguing opportunity for cancer immunotherapy. Chiapinelli et al. 539 already demonstrated that aza-induced cryptic transcripts are highly immunogenic and 540 can sensitize tumours to checkpoint immunotherapy 31 , while Sheng et al. showed that 541 also histone demethylase LSD1-ablation increases cryptic transcripts, thereby enabling 542 checkpoint blockade 53 . The mechanism underlying immunogenicity likely depends on the 543 formation of dsRNA, which via a viral mimicry-mediated process activates the immune 544 system 33,45,53 . In addition, some of these transcripts contain open-reading frames, which 545 could translate into abnormal proteins that can be antigenic 45 . Importantly, hypoxia is 546 endemic to most solid tumours, and therefore could have a more widespread impact than 547 aza. Indeed, in hypoxic tumours with high checkpoint expression, DNA methylation at TSS 548 of cryptic transcripts was reduced and consequently, cryptic transcript expression 549 increased. Since tumours with high checkpoint expression often respond to checkpoint 550 immunotherapy, and as cryptic transcripts could sensitize tumours to checkpoint 551 blockade 31, 33 , this suggests hypoxia-induced cryptic transcripts to play an important role 552 in mediating the therapeutic effects exerted by checkpoint inhibitors. In contrast, 553 immune-cold tumours characterized by low immune checkpoint expression were much 554 less tolerant to cryptic transcript expression, showing high methylation at 555 retrotransposon promoters. In light of our findings that methylation directly repels HIF 556 binding, this suggests DNA methylation to block hypoxia-induced cryptic transcript 557 expression in immune-cold tumours to maintain immunotolerance. Pharmacological 558 demethylation of immune-cold 4T1 tumours indeed increased cryptic transcription, 559 enhanced immunogenicity and reduced tumour growth in a HIF-dependent manner. By 560 showing that low-immunogenic, hypoxic tumours can be rendered immunogenic through 561 DNA methylation inhibitors, we thus highlight a novel treatment strategy for tumours 562 otherwise refractory to immunotherapies. 563

Cell line treatment conditions 598
Cell cultures were grown under atmospheric (21%) oxygen concentrations in the presence 599 of 5% CO2, or rendered hypoxic by incubating them under 0.5% oxygen (5% CO2 and 94.5% 600 N2). For ChIP-seq experiments, hypoxia was induced during 16 hours, whereas 24 hours 601 of exposure were used when assessing effects of hypoxia on gene or protein expression 602 level. Where indicated, cells were pre-treated with 5-aza-2'-deoxycytidine (aza, 1 µM) for 603 3 days by adding the required volume to fresh culture medium. Equal volumes of the 604 carrier (DMSO) were used as control. This was followed by another day of exposure to aza 605 in hypoxia or normoxia, bringing the total aza exposure time for experiments to 4 days. 606 2mM of DMOG (dimethyloxalylglycine, Sigma) was added to culture medium for 24 hours 607 where indicated. Cells were always plated at a density tailored to reach 80-95% 608 confluence at the end of the treatment. Fresh medium was added to the cells just prior 609 to hypoxia. To prove that the extent to which cells were exposed to hypoxia was similar 610 across experiments, we assessed that induction of hypoxia marker genes (BNIP3, EGLN, 611 ALDOA, CA9) but not HIF1A occurred in each experiment ( Figure S2). For experiments 612 involving exposure to aza, we assessed the expression of cancer testis antigens as a 613 positive control ( Figure S2). 614

LC-ESI-MS/MS of DNA to measure 5mC 615
DNA was extracted and processed for LC-ESI-MS/MS to determine 5mC concentrations 616 exactly as described previously 1 . 617

Western blot 618
To assess HIF1α protein stabilization, proteins were extracted from cultured cells as 619 follows: cells were placed on ice, washed twice with ice-cold PBS and lysed in protein extraction buffer (50 mM Tris HCl, 150 mM NaCl, 1% Triton X-100, 0.5% Na-deoxycholate, 621 0.1% SDS and 1´ protease inhibitor cocktail (Roche)). Protein concentrations were 622 determined using a bicinchoninic acid protein assay (BCA, Thermo Scientific) following the 623 manufacture's protocol. An estimated 60 µg of protein was loaded per well on a NuPAGE 624 Methanol-free, Thermo Scientific) directly in the medium and incubating for 8 min on a 637 flat-bed shaker at room temperature (RT). Fixed cells were incubated with 150 mM of 638 glycine for 5 min to revert the cross-links, washed twice with ice-cold PBS 0.5% Triton-639 X100, scraped and collected by centrifugation (1000 xG, 5 min at 4 °C). The pellet was 640 resuspended in 1,400 µL of RIPA buffer (50 mM Tris-HCl pH 8, 150 mM NaCl, 2 mM EDTA 641 pH 8, 1% Triton-X100, 0.5% Sodium deoxycholate, 1% SDS, 1% protease inhibitor) and 642 transferred to a new Eppendorf tube. The lysate was homogenized by passing through an 643 insulin syringe, and incubated on ice for 10 min. The chromatin was sonicated for 3 min 644 by using a Branson 250 Digital Sonifier with 0.7 s 'On' and 1.3 s 'Off' pulses at 40% power 645 amplitude, yielding predominantly fragment sizes between 100 and 500 bps. The sample 646 was kept ice-cold at all times during the sonication. Next, samples were centrifuged (10 647 min at 16,000 ´ G at 4 °C) and supernatant transferred in a new Eppendorf tube. Protein 648 concentration was assessed using a BCA. 50 µL of sheared chromatin was used as "input" bps centered on the summit was >4-fold bigger than the local background, and as absent 676 if it was <2.5-fold smaller than the local background, with local background being defined as the read depth across regions 1.5-5 kb up-and downstream of the peak. Intermediate 678 coverage was annotated as unclassified. To compare HIF1β binding peaks between 679 murine Dnmt-WT and -TKO ESCs, the HIF1β binding peak was called as present if the 680 average coverage at the 200 bps centered on the summit was >4-fold bigger than the 681 background, and as absent if it was <4-fold smaller than the background. To compare 682 efficiency between experiments, scatter plots of reads counts at peak regions of HIF1β 683 binding regions were generated per cell line in a pairwise fashion. 684

Annotation of genomic features 685
Human sequences were mapped to genome build hg19 and murine sequences to genome 686 build mm10. Putative HIF binding sites were detected genome-wide by screening the 687 whole genome for RCGTG motifs using the regular expression search tool dreg 688 (www.bioinformatics.nl/cgi-bin/emboss/help/dreg). The frequency per bp of RCGTG 689 motifs inside HIF1β binding peaks and in the rest of the genome was calculated, and 690 enrichment of RCGTG motifs at HIF1β binding peaks quantified by overlapping RCGTG 691 positions in the genome with the HIF1β binding peak positions as defined by MACS. 692 The distances of HIF1β peaks to the nearest RCGTG motif (cumulative frequency), TSS and 693 open chromatin (frequency) were calculated by overlapping each genomic feature with 694 HIF1β peak positions using BedTools in R. Protein-coding genes were annotated as per 695 Ensembl version 92. Promoter regions were annotated as being 2 kb upstream or 500 bp 696 downstream of the start site of each gene. Chromatin state annotation of MCF7 and 697 mESCs was as described 21, 56 . HIF1β binding peaks were annotated with these features 698 and overlapped with the repeat genome using BedTools. To assess enrichment of HIF 699 binding at repeats, HIF1β binding peaks were 10,000 times either randomly distributed 700 throughout the genome, or randomly distributed while matching the distal binding peak 701 distribution. Next, the frequency of repeat binding in a random distribution was 702 compared to that in the observed distribution. Peaks randomly assigned to poorly 703 mapping regions were discarded. 704

Genome distribution of 5mC: BS-seq, SeqCapEpi BS-seq and mDIP-seq
BS-seq, SeqCapEpi BS-seq and mDIP-seq were performed as described previously 1 . To 706 quantify DNA methylation inside HIF1β binding peaks, SeqCapEpi probes with >40´ 707 coverage were overlapped with HIF1β binding peaks as defined by MACS. Methylation 708 levels at the probes overlapping and non-overlapping (rest of the genome) HIF1β binding 709 peaks were calculated using Seqmonk. 710 ChIP-BS-seq was done as ChIP-seq, except that methylated adaptors (NEB) were ligated, 711 and DNA libraries were bisulfite-converted using the EZ DNA Methylation-Lightning™ kit 712 (Zymo) prior to library amplification using HiFi Uracil+ (KAPA). Reads were mapped using 713 Bismark as described 1 . 714

RNA-seq 715
To assess the impact of HIF binding at gene promoters on their expression, strand-specific 716

RNA-seq was performed in human cell lines and murine Dnmt-WT and Dnmt-TKO ESCs. 717
Briefly, total RNA was extracted using TRIzol (Invitrogen), and remaining DNA 718 contaminants in 17-20 μg of RNA were removed using Turbo DNase (Ambion) according 719 to the manufacturer's instructions. RNA was repurified using the RNeasy Mini Kit 720

CREDENToR analysis 760
The overall strategy of CREDENToR is to perform de novo assembly of all reads and based 761 on this define all cryptic transcripts. CREDENToR will consider transcripts encompassing 762 more than one repeat element as one cryptic transcript and quantify gene expression for 763 each of them. To achieve this, fastq files of RNA-seq data were first aligned to the human 764 (build GRCh38) or the murine genome (build mm10) using STAR (version 2.5) with a 765 tolerance of two mismatches. Transcriptomes were subsequently assembled using 766 StringTie 57 (version 1.3.4d), under guidance of the transcript annotation tool Ensembl 92. 767 All de novo assembled transcription annotations from the same set of tumour samples 768 (i.e., MCF7 or 4T1 cell lines, or each of the 14 tumour types downloaded from TCGA) were 769 merged using "StringTie --merge". HTSeq-counts 58 (version 0.11.2) were used to count 770 the read numbers of known and novel genes. Non-coding transcripts (transcripts not 771 overlapping annotated coding genes) in the merged transcription annotations were 772 assigned as cryptic transcripts when any of their exons overlapped with a retrotransposon 773 repeat annotation (LTR, LINE, or SINE, based on RepeatMasker annotation from UCSC). If 774 a transcript overlapped with >1 annotated repeat, the retrotransposon with the highest 775 overlap was assigned to this cryptic transcript. 776 For the analysis of MCF7 data, the assembled annotations from all experimental 777 conditions involving MCF7 cells assessed in vitro were merged before read counting. For 778 the analysis of 4T1 data, the assembled annotations from in vitro and in vivo samples were 779 merged together. Cryptic transcripts were considered to be HIF-associated if a HIF binding 780 summit was detected within the transcript promoter (i.e. 2,000 bp upstream and 500 bp 781 downstream of the transcription start site). Per set of experiments (24 samples), we 782 further required that the read number per cryptic transcript exceeds 10 in at least 1 783 sample, and that the reads per kb per million reads (RPKM) exceeds 1 in at least one 784 sample. For these cryptic transcripts, DESeq (version 1.30.0) was used to test the 785 differences between each pair of conditions. For the TCGA data, we merged assembled 786 transcription annotations for each tumour type separately. Cryptic transcripts were 787 calculated using the total cryptic transcript read count divided by the total coding gene 788 read count. In volcano plots, individual cryptic transcripts were plotted, but in violin plots, 789 where we compare effects to the cryptic transcripts obtained by RepEnrich, we summed 790 cryptic transcript counts into retrotransposon subfamilies, log-transferred counts-per-million (normalized to total read counts) and considered those as expression values. Violin 792 plots invariably show >95% of data points. P values were corrected for multiple testing 793 following Benjamini and Hochberg correction. The CREDENToR pipeline has been made 794 available on GitHub (https://github.com/Jieyi-DiLaKULeuven/CREDENToR). 795

Gene Ontology Analysis 796
Genes were associated to ontologies as annotated in BioMART (Ensembl GRCh37 release 797 84), and enrichment of ontologies was analysed using TopGo (version 1.0) in R 59 , using 798 the classic algorithm, contrasting to all protein-coding genes. 799

Structural modeling of DNA methylation 800
The crystal structure of HIF2α:HIF1β in complex with DNA containing the RCGTG core 801 sequence 5'-ACGTG-3' ( 25 , PDB code 4ZPK) was used as a template for introducing and 802 analyzing the structural consequences of methyl groups at position 5 of the cytosines 803 using the programs PyMOL (Schrodinger, LLC) and Chimera 60 . 804

Microscale thermophoresis (MST) binding assay 805
MST measurements were performed in triplicate using the NanoTemper Monolith NT.115 806 instrument. The two protein complexes (HIF1α-HIF1β and HIF2α-HIF1β) were purified as 807 PCR product sizes were verified by gel electrophoresis, and amplicons converted into 851 sequencing libraries using the NEBNext DNA library prep master mix set (E6040L, Bioke). 852 These were next sequenced to a depth exceeding 500´, and mapped and analyzed as 853 described higher. 854 Positive colonies were expanded into 10 cm dishes and subjected to ChIP as described  Tumour types were classified as responsive or non-responsive to checkpoint 875 immunotherapy following to the classification described by Turajlik and colleagues 37 , with 876 3 notable exceptions. Firstly, kidney renal papillary cell carcinoma (KIRP) was classified as 877 non-responsive as no study has yet demonstrated responsiveness of this tumour type. 878 Secondly, clear cell kidney carcinoma (KIRC) tumours were not analysed, as the HIF-879 constitutive activation in these tumours (due to VHL-loss) precludes their classification 880 into a hypoxic and normoxic subset. Finally, also microsatellite-instable COAD tumours 881 were discarded from all analyses, as these tend to show constitutive hypermethylation, 882 precluding their stratification in high and low methylation subgroups. 61 Importantly, the 883 tumour types that we define as responsive will still contain tumours that fail to respond, 884 whereas the non-responsive tumour types will also contain a minority of tumours that do 885 respond to immunotherapy. For instance, a subset of triple-negative breast tumours 886 responds to checkpoint immunotherapy. Likewise, there is evidence that a small subset 887 of LIHC tumours also responds, and that microsatellite-instable tumours also occur in 888 UCEC tumours. Cryptic transcription loads were calculated using the total cryptic 889 transcript read count divided by the total coding gene read count. 890 For the methylation stratification, the beta values of HM450K methylation microarray 891 data were downloaded from TCGA. Probes overlapping the cryptic transcript promoter 892 (i.e. 2,000 bp upstream and 500 bp downstream of the transcription start site) were 893 regarded as cryptic transcription-associated probes. For each tumour, its promoter 894 methylation level was calculated as the mean beta value of all cryptic transcription-895 associated probes. All tumours were then classified into high and low methylation groups 896 based on the median value of methylation levels. 897 To identify which of these tumour samples were hypoxic or normoxic, we performed 898 To test the interaction between hypoxia and DNA methylation, we assessed read counts 905 for cryptic transcripts in two negative-binominal generalized linear models with both 906 oxygenation (hypoxic and normoxic; encoded as 0 and 1) and methylation (low and high 907 methylation; encoded as 0 or 1), with or without an interaction term. Both models were 908 compared to each other using DESeq. A positive interaction coefficient represents a 909 cooperative enhancement of cryptic transcript expression in low-methylation, hypoxic 910 tumours. To further enrich for tumours that are prone to respond to checkpoint 911 immunotherapy, we stratified all tumour types into high PDL1 mRNA expressing and low 912 PDL1 mRNA expressing tumours, and into tumours with a high or low tumour mutation 913 burden (TMB). Stratification was done on the third decile in both cases. TMB was 914 estimated based on the number of substitutions identified by TCGA in each tumour 915 sample. All substitutions were considered, except for those also present in non-malignant 916 samples (i.e., exclusion of germ-line variants) or those clustering within and across 917 different samples (and therefore most likely representing sequencing or mapping errors). 918

Single-cell analysis 919
We used CREDENToR to map cryptic transcript expression in each individual cell from a 920 the endpoint. Tumour volumes were monitored every two to three days by a calliper, and 980 mice were culled before tumour volumes exceeded 2,000 mm 3 . When over 20% of mice 981 were culled, the experiment was terminated (all arms). In vivo experiments in 4T1, CT26 982 and MC38 treated with aza or anti-PD1 antibody were performed three times, with at 983 least 6 mice per treatment group in each experiment. 984

Neo-epitope burden 985
To assess neo-epitope burden, we mapped RNA sequencing data of isogenic 4T1, B16 and 986 CT26 tumour models, removed duplicate reads from individual samples and merged per 987 tumour model all samples into a single file. In this file, variants were called according to 988 GATK best practices, using GATK3.4. Briefly, reads were split into exon segments and to bind for 1 h. Slides were washed 3 times for 5 min in TNT, after which signal 1019 amplification was done by 30 min incubation with peroxidase-conjugated streptavidin 1020 1:100 in TNB (for all besides pimonidazole) accompanied by nuclear staining with Hoechst 1021 (H3570, Thermo Fisher) 1:500 in TNB only for the single (CD45 or CD8a) stainings, washing 1022 (3 times 5 min in TNT) and 8 min incubation using Fluorescein Tyramide (for pimonidazole 1023 NEL701A001KT, perkin Elmer) or Cy3 (NEL704A001KT, Perkin Elmer) 1:50 in amplification 1024

diluent. 1025
Slides stained for pimonidazole and Gzmb required co-staining for CD31 and CD8a 1026 respectively and were subjected to a second indirect staining for the latter epitopes. After 1027 5 min of TNT and 5 min of TBS, slides were quenched again for peroxidase activity using 1028 H2O2 and blocked using pre-immune goat or rabbit (CD31) serum, prior to a second 1029 overnight round of primary antibody binding: rat anti-CD31 (557355, BD Biosciences) or 1030 rat anti-CD8a (14-0808-82, Thermo Fisher) 1:100 in TNB. The next day, 3 times 5 min 1031 washes with TNT were followed by a 1 h incubation with biotinylated goat anti-rat 1032 Gzmb) 1:1,000 in dilution buffer (1 % BSA in TBS) was applied for 30 min at RT, followed 1065 by 3 washes of 2 min in TBST at RT. Slides were next incubated with the secondary 1066 antibody (EnVision™+/HRP goat anti-rabbit (K4003, Dako)) for 10 min at RT, and washed 1067 3 times for 2 min in TBST at RT. The OPAL 570 fluorophore (fp1488, PerkinElmer) 10 % in 1068 amplification diluent (FP1498, PerkinElmer) was applied for 10 min at RT followed by 3 1069 washes of 2 min in TBST at RT. Slides were stripped by heating in AR6 buffer just below 1070 boiling point and cooled down in double distilled water, followed by rinsing in TBST. These 1071 steps were repeated starting from blocking for the second staining with primary antibody 1072 rat anti-CD8a 1:300, secondary antibody goat anti-rat (MP-7444, Vector) and opal 690 1073 (fp1497, PerkinElmer) and the third staining with rat anti-CD45 1:1,000, secondary 1074 antibody goat anti-rat and Opal 520 (fp1487, PerkinElmer). After the third staining, slides 1075 were incubated with spectral DAPI (fp1490, PerkinElmer) 10 % in TBST for 5 min at RT, washed for 2 min in TBST at RT and mounted with ProLong Diamond Antifade Mountant 1077 (P36961 Invitrogen). Images were acquired on a Zeiss Axio Scan.Z1 using a ´20 objective 1078 and ZEN 2 software (Zeiss) with exposure times between 10-50 ms. Image processing was 1079 done using QuPath (version 0.1.2). Specifically, following visual inspection of the staining 1080 results, cells were first automatically detected using the DAPI channel (cell size 1081 constrained between 5 and 400 µm 2 ). Next, a cell classifier was generated using QuPath. 1082 Specifically, for 1 slide out of all slides, 5 sets of cells were selected: one set that was 1083 positive for CD45, one set that was negative for CD45, and three sets of CD45 + cells 1084 positive for CD8, Gzmb and CD8, or Gzmb alone. Using these 5 sets of cells, a random 1085 trees classifier was generated. Cell classification was visually verified to have occurred 1086 correctly. Next, in each tumour section, a representative region was selected, containing 1087 at least 1,000 cells. On these cells, the random trees classifier was subsequently applied. 1088 This process was reiterated for all other tumour sections stained for the same set of 1089 markers. The resulting cell identities were then exported, and processed in R. For each 1090 tumour, average cell frequencies were generated, which were summarized using 1091

Immunofluorescence analysis of 4T1 cells 1093
For the dsRNA staining on 4T1 cells, 12.000 cells were seeded on gelatin-coated glass slips 1094 in 12-well plates on day 0. After attaching for 6h, cells were treated with aza or control 1095 (DMSO) for 3 days, with renewal of the medium after 48 hours. After 72h, the medium 1096 was refreshed and cells were incubated in hypoxia or normoxia for 16 hours. After 1097 washing 3 times with PBS, cells were fixed in 1ml of ice-cold methanol for 15 min at -20 1098 °C. Cells were washed 3 times with PBS and blocked for 1 h with blocking buffer (0.1% 1099 triton X-100 with 5% goat serum in PBS). Primary antibody (1:50 dilution in blocking 1100 buffer; clone J2, Scicons) was applied overnight at 4 °C. Cells were washed 3 times for 10 1101 min with washing buffer (0.1% triton X-100 in PBS) and secondary antibody (1:500 in 1102 secondary antibody buffer (0.1% triton X-100 with 2% goat serum in PBS)). Goat-anti-1103 mouse IgG coupled to alexa 488, Life Technologies) was applied for 1 hour in the dark. 1104 Cells were washed 3 times for 10 min with washing buffer and mounted with DAPI stain on cover glasses. Slides were imaged on an infraMouse Leica DM5500 microscope. 3 1106 slides from different treatment groups were stained in triplicate (biological replicates), 1107 while 3 pictures from different slides were used for processing with Image J. More 1108 specifically, nuclei were identified using the DAPI signal, and nucleated cells were further 1109 selected based on particle size. Signal intensities for alexa fluor 488 in the selected cells 1110 were used to detect dsRNA + cells. Analyses were exclusively performed on slide regions 1111 showing a regular density and shape of nuclei, in order to avoid inclusion of acellular or 1112 necrotic areas. Mean dsRNA expression was calculated per experiment, normalised to the 1113 background signal (secondary antibody only) and expressed as mean pixel intensity 1114 relative to the control group (normoxia + DMSO). 1115

Published data sets 1116
Published data sets were obtained from GEO under the following accession numbers:

ETHICS APPROVAL 1125
All the experimental procedures on animals were approved by the Institutional Animal 1126 Care and Research Advisory Committee of the KU Leuven. 1127

AVAILABILITY OF DATA AND MATERIALS 1128
Sequencing data generated for this study have been deposited in GEO and will be made 1129 accessible upon publication. 1130

COMPETING INTERESTS 1131
The authors declare that they have no competing interests. 1132