Transcription Elongation Machinery Is a Druggable Dependency and Potentiates Immunotherapy in Glioblastoma Stem Cells

Glioblastoma (GBM; World Health Organization grade IV glioma) ranks among the most aggressive and lethal types of human cancer. Treatment options for patients with GBM remain ineffective and only palliative (1). GBMs display extensive intratumoral and microenvironmental heterogeneity, in which GBM stem cells (GSC) play an important role in the tumor hierarchy (2–4). GSCs, driven by both genetic and epigenetic alterations, contribute to sustained tumor growth, invasion into normal brain, evasion of immune surveillance, and therapeutic resistance, thus representing a crucial target for GBM treatment (4–7).

Recent advances in profiling GBM at the bulk and single-cell levels have deciphered different cell states associated with GSCs, which can be regulated by genetic, epigenetic, or microenvironmental factors (8). Although genetic driver alterations, such as EGFR amplification, have failed to provide effective therapies in GBM (9), examinations of epigenetic and transcriptional regulations have identified master regulatory circuitry that defines cell state plasticity, as well as immune response in malignant cells (10–12). Dysregulated transcriptional programs in cancer cells open avenues to novel cancer therapies, illustrated by targeting MYC with BRD4 inhibitors or enhancer-associated gene expression with inhibitors of transcriptional cyclin-dependent kinases (CDK; refs. 13–15). Cell cycle–related CDK4/6 has been well studied in cancers; however, the transcriptional CDKs, including CDK7, CDK9, and CDK12, which are involved in RNA polymerase (Pol) II–mediated transcriptional regulation, are not fully explored as therapeutic targets (16).

Among the epigenetic regulations, the role of chromatin interactions in cancers remains elusive. The CCCTC-binding factor (CTCF) insulator establishes topologically associated domains (TAD; ref. 17). In gliomas, IDH1 mutations induce hypermethylation to disrupt CTCF binding sites and chromosomal neighborhoods, leading to the interaction between a constitutive enhancer with platelet-derived growth factor receptor A (PDGFRA) and the subsequent upregulation of PDGFRA expression (18). Yin Yang 1 (YY1) is another structural regulator of chromatin interaction loops that controls gene expression (19). Disrupted chromosomal neighborhoods and loops in cancers represent potential targets, but the therapeutic opportunities underlying how chromatin structural proteins regulate cancer stemness states are poorly understood (20).

Although immune checkpoint inhibitors have transformed care for several solid tumors, including those with brain metastases, their clinical efficacy in GBM has been disappointing. Effective strategies of immune modulation to enhance response to immunotherapy are not available in GBM (21, 22). Considering the lack of effective treatment and inconsistent response to immunotherapy in GBM, we hypothesized that a systematic dependency analysis of chromatin regulators in GSCs could identify targetable cell states and associated master epigenetic regulators linked to therapeutic response for GBM.

To investigate core chromatin regulators in GBM, we interrogated the landscape of chromatin regulators across GSCs and matched differentiated GBM cells (DGC), as well as their normal counterparts, neural stem cells (NSC), by considering both transcription activities and cellular dependencies (Fig. 1A; ref. 23). By integrating differentially expressed genes (DEG) between 38 GSCs and 5 NSCs (24) with chromatin immunoprecipitation followed by deep sequencing (ChIP-seq) profiling of 160 chromatin regulators, 39 chromatin regulators showed active programs inferred by their target gene expression in GSCs, including transcriptional programs driven by MYC, polycomb repressive complex 2 (PRC2; including EZH2 and SUZ12), and YY1, supporting potential roles in gliomagenesis (Fig. 1B and C; Supplementary Fig. S1A; Supplementary Table S1). To identify stemness-associated regulators, we compared expression data from GSCs and matched DGCs (25), revealing 12 regulators of which programs were enriched in GSCs, including MYC, EZH2, and YY1 (Fig. 1D and E; Supplementary Table S1). Translating master regulators as targets requires an optimal therapeutic index, that is, selective dependency of a target in neoplastic cells (in this case GSCs) relative to normal counterparts (e.g., NSCs). We then interrogated CRISPR/Cas9 dropout screens from eight GSCs and two NSCs to identify GSC-specific dependencies (26). Among 160 chromatin regulators, 10 were selectively essential in GSCs compared with NSCs, including JUN, NRF1, and YY1 (Fig. 1F and G; Supplementary Fig. S1B). Integrating these analyses indicated that MYC and PRC2 complex were important for shaping gene expression profiles but failed to show differential dependencies between GSCs and NSCs (Fig. 1B–G). Surprisingly, several components related to RNA Pol-mediated transcription, including YY1, GTF2B, GTF3C2, and POLR3A, showed differential dependencies in GSCs (Fig. 1F; Supplementary Fig. S1C). Among these regulators, YY1, a master regulator important for shaping the regulatory circuitries and conferring selective dependencies in GSCs (Fig. 1H), was selected for further studies.

YY1 is a zinc finger transcription factor and interacts with other chromatin regulators, including MYC and EZH2, suggesting YY1 may function as a regulatory hub with selective dependency in GSCs (27–29). The dependency of GSCs on YY1 was further tested. Knockdown of YY1 using three nonoverlapping short hairpin RNAs (shRNA) strongly decreased cell proliferation in a panel of five patient-derived GSCs compared with nonneoplastic NSCs (Fig. 2A–F; Supplementary Fig. S1D and S1E). YY1 has been implicated in chromatin interactions, especially those connecting active regulatory elements, such as enhancer–enhancer, enhancer–promoter, and promoter–promoter interactions (19). CTCF is a structural protein that generally insulates neighborhoods and global structures, as well as long-range interactions (17). In contrast to effects of YY1 targeting, knockdown of CTCF only modestly inhibited cell proliferation in GSCs (Fig. 2A–E; Supplementary Fig. S1D and S1E). CRISPR/Cas9-based knockout validated selective dependency on YY1 over CTCF in GSCs (Supplementary Fig. S1F and S1G). Stem cell frequency and self-renewal of GSCs were diminished upon YY1 depletion, as revealed by extreme limiting dilution assays and tumor sphere formation (Fig. 2G–L). To assess the in vivo tumor initiation dependency of GSCs on YY1, two patient-derived GSCs were transduced with a control shRNA sequence (shCONT) or one of the two nonoverlapping shRNAs targeting YY1, then implanted intracranially into immunodeficient mice. YY1 knockdown prolonged survival of intracranial tumor-bearing mice relative to the shCONT (P < 0.01 by log-rank test; Fig. 2M and N), suggesting that YY1 is a master regulator controlling proliferation and self-renewal of GSCs in vivo.

As YY1 is a structural regulator in chromatin loop formation (19), we integrated RNA sequencing (RNA-seq), ChIP-seq, and Bridge-Linker–Hi-C (BL-Hi-C) analysis to dissect the YY1-regulated program (30). Gene set enrichment analysis (GSEA) of RNA-seq data identified several pathways downregulated after YY1 knockdown, including Pol II transcription, RNA processing, MYC activity, and cell cycle, which were enriched in YY1 targets from ChIP-seq data and GSC-specific gene expression signatures (Fig. 3A; Supplementary Fig. S2A–S2C; Supplementary Table S2). The expression levels of genes involved in Pol II–mediated transcription (PAPOLA and MYC), RNA splicing (SRSF1–3 and HNRNPU), and RNA N6-methyladenosine (m⁶A) modification (METTL3 and YTHDF2) were downregulated by shRNA-mediated YY1 knockdown and CRISPR/Cas9-mediated knockout (Fig. 3B; Supplementary Fig. S2D). Whole-genome CRISPR/Cas9 dropout screens also supported that YY1 target genes were enriched as GSC-specific dependencies (Supplementary Fig. S2E).

At the regulatory level, interrogation of active transcription marks, including histone 3 lysine 27 acetylation (H3K27Ac), histone 3 lysine 4 trimethylation (H3K4me3), and open chromatin profile by Assay for Transposase-Accessible Chromatin with high-throughput sequencing (ATAC-seq), showed that YY1 largely bound to the active transcriptional regions in GSCs (Fig. 3C). About 70% of YY1 binding peaks overlapped with loop anchors (Fig. 3D), and almost all the YY1 target genes (97%) showed a decreased number of chromatin interactions (Fig. 3E; Supplementary Table S3), indicating a regulatory role of 3-D genome for YY1 targets. Pathway analysis using putative YY1 target genes from ChIP-seq profiling showed the similar enrichment of YY1-regulated biological processes in RNA-seq, including RNA Pol II transcription, RNA processing and splicing, cell cycle and DNA damage (Fig. 3F). Taken together, these results suggest that YY1 acts as a checkpoint to maintain a GSC transcription state through controlling RNA Pol II transcription and RNA processing (31).

We next focused on the transcriptional regulation of RNA processing programs, especially the m⁶A program and its effect on GSCs. Using BL-Hi-C technology, the chromatin interaction frequencies decreased after YY1 knockdown, as shown at the METTL3 locus (Supplementary Fig. S2F). Two YY1 target genes, METTL3 and YTHDF2, showed decreased promoter-anchored chromatin interactions and decreased expression levels after YY1 knockdown (Fig. 3B; Supplementary Fig. S2G and S2H), supporting the essential role of YY1-mediated chromatin interactions for transcription of RNA m⁶A regulators. In line with the downregulation of METTL3, the total m6A levels were decreased after YY1 knockdown in GSCs (Fig. 3G and H).

Inhibiting METTL3 or YTHDF2 stabilizes interferon-related genes (e.g., IFNB1, STAT1, and IRF1) and activates interferon signaling in other cell types (32). m⁶A sequencing profiling of GSCs showed m⁶A modifications on a panel of interferon-related genes (Supplementary Fig. S2I). Accordingly, GSEA of RNA-seq with YY1 knockdown revealed enrichment of interferon signaling after YY1 knockdown (Fig. 3I). The qPCR analysis confirmed that YY1 knockdown increased expression levels of IFNα and IRF1 (Fig. 3J). Supporting the role of m⁶A in interferon response, GSEA of RNA-seq after YTHDF2 knockout in GSC revealed a strong enrichment of immune response, including interferon signaling and antigen presentation, upon YTHDF2 knockout (Supplementary Fig. S2J; ref. 33). In contrast, YY1 ChIP-seq profiles showed that the interferon genes were not bound by YY1, further confirming an indirect link between YY1-mediated transcription and interferon signaling (Supplementary Fig. S2K). Taken together, these results suggest that YY1 regulated transcription process suppressed interferon signaling through regulators of m⁶A modification.

Despite the strong dependency of GSCs on YY1, targeting YY1, along with other transcription factors, remains challenging, as YY1 is not a readily druggable target. To better understand the YY1-associated mechanism and therapeutic opportunity, we interrogated drug screening data from the Cancer Therapeutics Response Portal (CTRP) to identify YY1 expression–associated drug sensitivity (34, 35). We filtered drugs whose AUC values correlated with YY1 expression in brain cancer cell lines, then clustered these drugs based on their AUCs (Fig. 4A; Supplementary Table S4). The top candidates that correlated with YY1 expression levels included ferroptosis inducers (ML162, erastin, and 1S,3R-RSL-3), transcriptional CDK inhibitors (alvocidib and dinaciclib), NAMPT inhibitors (CAY10618 and STF-31), and an Aurora kinase inhibitor (alisertib). YY1 expression negatively correlated with AUCs of these drugs (Fig. 4B), suggesting that YY1-driven transcription state was sensitive to these compounds. In contrast, there was a moderate association between high YY1 expression levels and resistance to EGFR or BRAF/MEK inhibitors (Supplementary Fig. S3A and S3B).

To confirm the therapeutic opportunities targeting YY1, six compounds associated with YY1 expression, including two transcriptional CDK inhibitors, two NAMPT inhibitors, one ferroptosis inducer, and one Aurora kinase inhibitor, underwent validation studies against five GSCs using dose–response curves (Fig. 4A). GSCs dependent on YY1 were exclusively sensitive to five of these drugs, especially transcriptional CDK9 inhibitors (alvocidib and dinaciclib; Fig. 4C; Supplementary Fig. S3C–S3E). Differential sensitivity was not observed for the cell cycle–related Aurora kinase inhibitor (Supplementary Fig. S3E). Moreover, sensitivity to inhibitors targeting other types of CDKs, including the cell cycle–related CDK4/6 inhibitor palbociclib, was limited (Fig. 4D; Supplementary Fig. S3F). In vitro extreme limiting dilution assay showed decreased sphere frequencies in GSCs upon alvocidib or dinaciclib treatment (Supplementary Fig. S3G). Consistently, alvocidib treatment decreased the expression levels of stemness markers SOX2 and OLIG2 (Supplementary Fig. S3H). Considering alvocidib and dinaciclib also inhibit other CDKs, we evaluated three potent selective CDK9 inhibitors, MC180295, AZD4573, and NVP-2, against GSCs (36–38). GSCs dependent on YY1 displayed preferential sensitivity to selective CDK9 inhibitors (Fig. 4E–G). Moreover, knockdown of CDK9 using two nonoverlapping shRNAs decreased cell proliferation (Supplementary Fig. S3I). These data suggest a correlation between YY1 dependency and transcription regulation through CDK9 in GSCs.

We next tested the in vivo antitumor efficacy of one transcriptional CDK inhibitor, alvocidib. Mice bearing intracranial patient-derived GSC tumors were treated with vehicle or alvocidib. Alvocidib treatment extended the survival of mice compared with vehicle control (Fig. 4H). Both in vivo imaging and histologic analysis of tumor-bearing brains indicated the impaired tumor growth after alvocidib treatment (Fig. 4I and J). These findings suggest a pharmacologic vulnerability of targeting YY1-driven GSC state using transcriptional CDK9 inhibitors in a cell-autonomous manner both in vitro and in vivo.

Considering the regulation of RNA Pol II–mediated transcription program by YY1 and the correlation between YY1 expression and response to transcriptional CDK inhibitors, we examined the potential link between YY1 and one of the YY1-associated transcriptional CDK inhibitors, alvocidib, using RNA-seq profiling. Similar to findings upon YY1 inhibition, alvocidib-downregulated genes were enriched in RNA Pol II transcription, RNA processing, and cell cycle (Fig. 5A). Cross-examination of alvocidib and YY1 knockdown–induced gene expression changes showed that genes downregulated after alvocidib treatment were enriched in YY1 knockdown–downregulated genes (Supplementary Fig. S4A). A panel of YY1 target genes, including METTL3 and YTHDF2, was downregulated by alvocidib (Fig. 5B). Alvocidib treatment of GSCs decreased total m⁶A levels, consistent with the downregulation of METTL3 (Fig. 5C). Similarly, shRNA-mediated knockdown of CDK9 and selective CDK9 inhibitors decreased METTL3 and YTHDF2 expression levels in GSCs, suggesting a role of CDK9 in regulating YY1 targets (Fig. 5D–F; Supplementary Fig. S4B).

Transcriptional CDKs mediate serine (Ser) phosphorylation of the RNA Pol II carboxy-terminal domain (CTD; ref. 16). Immunoblotting showed that both alvocidib treatment and YY1 loss of function dramatically downregulated Ser2 phosphorylation of Pol II CTD, a key regulatory step in transcription elongation regulated by transcriptional CDK9 (Fig. 5G and H). To evaluate transcription elongation after YY1 knockdown, we performed ChIP-seq profiling using RNA Pol II antibody to study Pol II distribution across transcripts to study the elongation process. A pausing index was calculated based on the ratio of the Pol II signal in gene promoters to the gene bodies to identify Pol II pausing (39–41). YY1 knockdown increased promoter-proximal pausing of Pol II, indicating a defect in transcription elongation, as shown by the empirical cumulative distribution function plot and the decreased Pol II signal in gene bodies of METTL3 and HNRNPU (Fig. 5I–K; Supplementary Fig. S4C).

To directly detect the interaction between YY1 and CDK9, endogenous coimmunoprecipitation (co-IP) and immunoblotting were performed in GSCs, a breast cancer cell line (MDA-MB-231), and a lung cancer cell line (A549). YY1–CDK9 interaction was detected in GSCs but not in MDA-MB-231 or A549 cell lines (Fig. 5L), indicating a potential GSC-specific YY1–CDK9 interaction. CDK9 is involved in different transcription elongation complexes, including BRD4-containing elongation complex (BEC), super elongation complex (SEC), and KAP1–7SK elongation complex (KEC; Fig. 5M; refs. 40, 42, 43). To determine which CDK9-containing complexes contain YY1, we performed an in silico analysis of GBM data sets from The Cancer Genome Atlas (TCGA), revealing YY1 expression correlated with BRD4 (BEC), AFF4 (SEC), MLLT1 (SEC), and TRIM28 (KEC; Supplementary Fig. S4D). In GSCs, endogenous co-IP experiments showed that YY1 interacted with multiple complexes, including those containing BRD4, AFF4, and TRIM28, which was not detected in MDA-MB-231 or A549 (Fig. 5L). Depletion of transcription elongation complexes showed that GSCs were most dependent on these complexes compared with MDA-MB-231 (intermediate dependent) and A549 (least dependent; Supplementary Fig. S4E–S4K). Consistent with YY1 and CDK9 inhibition, depletion of transcriptional elongation complexes decreased expression levels of METTL3 and YTHDF2 (Supplementary Fig. S4L), which was also observed in treatment of the BEC inhibitor JQ1 and the SEC inhibitor KL-2 (Fig. 5N; Supplementary Fig. S4M; refs. 41, 44). These data suggest a context-specific transcription elongation complex with YY1 and CDK9 functioning in GSCs (Fig. 5M).

Similar to YY1 knockdown (Fig. 3I and J), several immune response–related gene sets, including interferon response and antigen presentation, were enriched after treatment with the transcriptional CDK inhibitor alvocidib, as revealed by GSEA (Fig. 6A). Genes involved in the major histocompatibility complex, including HLA-A, HLA-C, and HLA-H, were upregulated (Fig. 6A and B). Concordantly, treatment with alvocidib and dinaciclib upregulated IFNα and IFNβ, as well as other interferon-stimulated genes (Fig. 6C–E; Supplementary Fig. S5A). Although knockdown of individual transcription elongation complex components induced a weak upregulation of IFNα and IFNβ, potent inhibition of transcription elongation by selective CDK9 inhibitors and SEC inhibitor consistently induced a robust activation of IFNα and IFNβ (Fig. 6F–J; Supplementary Fig. S5B). To explore the potential mechanism of interferon response activation, m⁶A levels on IFNB1 after CDK inhibitor treatment were determined using m⁶A RNA immunoprecipitation followed by qPCR, showing the decreased m⁶A level on IFNB1 after alvocidib treatment (Supplementary Fig. S5C), which resulted in the increased stability of IFNB1 mRNA (Supplementary Fig. S5D; ref. 32). Rescue experiments showed that METTL3 or YTHDF2 overexpression partially rescued the activation of interferon signaling by alvocidib treatment (Supplementary Fig. S5E), supporting the role of m⁶A in mediating the activation of interferon response.

Although immune checkpoint inhibitors have transformed care for several solid tumors, including for patients with brain metastases, clinical efficacy in GBM has been disappointing (21, 22). Effective strategies to rewire the immunosuppressive microenvironment and enhance immunotherapy response are still lacking in GBM (21, 22). Considering the activation of interferon response after YY1 inhibition and transcriptional CDK inhibitor treatment, we tested whether targeting this mechanism in vivo could reprogram the tumor microenvironment. To extensively characterize the impact of transcriptional CDK inhibitor on the tumor microenvironment, we applied cytometry by time-of-flight (CyTOF) with 34 metal-conjugated antibodies to detect major immune cell populations in a syngeneic mouse glioma model with intracranial implantation of GL261 cells (Fig. 7A; refs. 45, 46). Cell clusters were computationally defined using FlowSOM clustering and ConsensusClusterPlus, as well as visualized using t-stochastic neighbor embedding (t-SNE) plots. Major immune cell populations were identified, including CD4⁺ T, CD8⁺ T, natural killer, macrophage, and myeloid-derived suppressor cells (MDSC; Fig. 7B; Supplementary Fig. S6A and S6B). Most CD4⁺ T cells were positive for FR4 staining, a marker for regulatory T (Treg) cells (Supplementary Fig. S6B). Compared with the vehicle-treated tumors, alvocidib/anti–PD-1 antibody-treated tumors showed decreased immunosuppressive Treg cells and MDSCs and increased antitumor cytotoxic CD8⁺ T cells (Fig. 7B and C; Supplementary Fig. S6C). IHC validated that alvocidib treatment reduced the infiltration of Treg cells and increased the infiltration of CD8⁺ T cells (Supplementary Fig. S6D). YY1 expression levels were positively correlated with CD4⁺ T-cell levels but negatively correlated with CD8⁺ T-cell levels in the TCGA GBM data set (Supplementary Fig. S6E). These data suggested that transcriptional CDK inhibitor treatment rewired the tumor microenvironment in GBM models.

Considering the rewiring of the immunosuppressive microenvironment by alvocidib, we hypothesized that alvocidib could enhance response to immunotherapy in vivo. Syngeneic mouse glioma models (CT2A and GL261) were treated with control, alvocidib, anti–PD-1 antibody, or the combination. Although treatments using a single agent modestly prolonged the overall survival of mice bearing intracranial tumors, the combination of alvocidib and anti–PD-1 antibody showed the longest survival compared with control (P < 0.01 by log-rank test; Fig. 7D and E). Correspondingly, single agents alone reduced tumor growth measured by in vivo bioluminescence imaging, but the combination of alvocidib and anti–PD-1 achieved a greater suppression of tumorigenesis without evidence of increased toxicity at the indicated doses (Fig. 7F). To assess memory response, mice with long-term survival after drug treatment were rechallenged with GL261 cells (47). In contrast to the short survival in naive mice, the rechallenged group showed extended survival, indicating the development of memory response for the antitumor effect (Fig. 7G). Taken together, these data suggest that the transcriptional CDK inhibitors enhance anti–PD-1 therapy by rewiring the immunosuppressive microenvironment in GBM.

To establish the clinical relevance of the YY1 program in patient samples, we interrogated YY1 expression in glioma patient databases. Compared with nontumor brain tissues, GBM and low-grade glioma tissues showed upregulated expression levels of YY1 (Fig. 7H). Expression levels of YY1 were highly correlated with tumor grade (Fig. 7H). In both TCGA glioma and Chinese Glioma Genome Atlas (CGGA) data sets, high YY1 expression levels were associated with poor prognosis (Fig. 7I and J). In a pan-cancer analysis, GBM is one of the cancer types with the most prominent YY1 upregulation (Supplementary Fig. S7A). To further dissect the YY1-regulated program in GBM, we constructed the coexpression network of YY1 followed by GSEA. In line with the role of YY1 in RNA Pol II–mediated transcription and RNA processing, YY1 expression positively correlated with gene sets involved in Pol II transcription and RNA processing but negatively correlated with immune response gene sets (Supplementary Fig. S7B and S7C), which was also observed in the CGGA GBM data set (Supplementary Fig. S7D). YY1 expression was negatively correlated with IFNB1, IFNG, and CCL27, which are involved in interferon signaling and immune regulation (Supplementary Fig. S7E).

Considering the observed dependency in GSCs on transcription components YY1 and CDK9, we proposed a transcription state in GSCs. GSCs showed heterogeneity in transcription score, which was defined by the RNA Pol II–related transcription regulation signature (Supplementary Fig. S7F). This transcription state was negatively correlated with the interferon response score in GSCs (Supplementary Fig. S7F). The presence of this transcription state with the suppressed interferon response was confirmed in vivo using an independent panel of patient-derived xenografts from patients with GBM (Supplementary Fig. S7G), as well as TCGA GBM samples and single-cell RNA-seq data sets from five primary GBMs (refs. 48, 49; Supplementary Fig. S7H and S7I). Taken together, our findings suggest that transcription and immune axes underline the GBM trajectory and define a transcription state with suppressed interferon response, in which YY1–CDK9 complex plays a role as a gatekeeper by regulating RNA Pol II–mediated transcription elongation and RNA m⁶A modification (Fig. 7K).

Through integrative analysis of gene expression programs and dependencies underlying GSCs, DGCs, and NSCs, we surveyed the chromatin regulator landscape in GBM and identified YY1 as a selective transcriptional dependency in GSCs. YY1 regulates RNA Pol II transcription machinery through interaction with transcriptional CDK9 and mediating chromatin loops, thus acting as a gatekeeper for this transcription cell state (modeled in Fig. 7K). RNA processing programs, including those involved in RNA splicing (SRSF1–3) and m⁶A modification (METTL3 and YTHDF2), were the main downstream targets regulated by YY1-mediated transcription and chromatin interactions. Our study uncovers a role of YY1 in transcription elongation by interacting with multiple transcription elongation complexes. Targeting the YY1–CDK9 complex or other transcription elongation complexes decreased expression levels of METTL3 and YTHDF2, which induced interferon responses. Treatment with YY1-associated transcriptional CDK inhibitors or other transcription elongation inhibitor augmented interferon response, providing a combination strategy with immunotherapy to treat GBM.

Profiles of chromatin interactions in tumor specimens have been interrogated to identify rational therapeutic targets, including CD276 as a target in GBM (20). As structural proteins, including CTCF, cohesin, and YY1, are challenging drug targets, methods to pharmacologically inhibit the global oncogenic chromatin interactions have not been established. Pharmacogenomic modeling using the large data sets of drug response in hundreds of genetically validated cancer cell lines provides an opportunity to successfully identify potential gene expression–correlated drug sensitivity or synthetic lethality (34). Drug response–correlated genes might be the direct target of the drug in gene–drug pairs or represent a cell state with a functional relationship to the queried drug response. Using this approach, we propose compounds targeting a specific class of transcriptional CDKs as inhibitors of the YY1 program. An apparent overlap of transcriptional CDKs and YY1 in transcription regulation may represent a targetable cell state (14, 50).

Transcriptional CDKs—CDK7, CDK9, and CDK12, among others—promote RNA Pol II–mediated transcription initiation and elongation by phosphorylating serine residues of the CTD of RNA Pol II (16, 51). The value and mechanism of transcriptional CDKs in treating cancers remain unresolved (38, 52–54). THZ1, a covalent CDK7 and CDK12/13 inhibitor, inhibits phosphorylation of RNA Pol II and abolishes superenhancer-driven oncogenic transcription (54). However, specific targeting of CDK7 shows no effect on RNA Pol II phosphorylation but inhibits CDK1 and CDK2 phosphorylation (53). Although CDK9 contributes to RNA Pol II phosphorylation and transcription elongation, CDK9 inhibition dephosphorylates BRG1 and leads to reactivation of tumor suppressor genes (38). In our study, targeting transcriptional CDKs by alvocidib decreased RNA Pol II phosphorylation and showed a strong inhibitory effect on RNA m⁶A machinery, in addition to previously reported inhibitory effects on cell cycle, DNA damage, and DNA synthesis (53). These findings highlight the different expression programs regulated by cross-talk among transcriptional CDKs (16). The dysregulated RNA processing machineries (METTL3, YTHDF2, or HNRNPC) might regulate m⁶A modification of interferon genes or yield improperly processed RNA products, which could trigger interferon response and partially explain the interferon response in our study (32, 55). Moreover, the toxicity of transcriptional CDK inhibitors, especially multi-CDK inhibitors, remains an obstacle in applying these inhibitors in the clinic, which could be overcome by developing more selective small-molecule inhibitors (56). Indeed, several recently developed CDK9 inhibitors are being tested in clinical trials of different cancer types (56).

The immunosuppressive microenvironment in GBM tumors may contribute to the unfavorable outcome from immunotherapy (22, 57). There is strong enrichment of immune-suppressive Treg cells and tumor-associated macrophages, which release inhibitory cytokines and inhibit T cell–mediated antitumor immune responses (21, 22). However, certain therapeutic strategies, such as radiation, may be able to convert immunologically “cold” microenvironments into “hot” for GBM (58). Our data indicate a crosstalk between the YY1–CDK9-driven cell state in GSCs and the immunosuppressive microenvironment in GBM. Perturbation of the transcription cell state with transcriptional CDK inhibitors may rewire the tumor microenvironment to enhance response to immunotherapy in GBM. Transcriptional CDK inhibitors could have direct effects on the microenvironment cells to simultaneously target tumor cells and microenvironment in vivo. Future studies to test a panel of transcriptional CDK inhibitors with different target profiles may elicit robust response to immunotherapies, including anti–PD-1, anti-CTLA4, or chimeric antigen receptor T-cell therapies.

In summary, our study demonstrates that YY1, YY1-associated transcriptional CDK9, and transcription elongation complexes are essential in GSC maintenance and presents a targetable transcription cell state with suppressed interferon response in GSCs. Both genetic and pharmacogenomic targeting of this transcription cell state suppress self-renewal and proliferation of GSCs and could be combined with immunotherapy to treat patients with GBM synergistically.

GBM tissues were obtained from excess surgical resection samples from patients at the Case Western Reserve University (Cleveland, OH) after review by neuropathology with written informed consent from the patients and in accordance with an institutional review board–approved protocol (090401). All specimens were reviewed by a neuropathologist. All patient studies were conducted in accordance with the Declaration of Helsinki. To minimize in vitro cell culture–based artifacts, patient-derived xenografts were propagated as a renewable source of GSCs. Short tandem repeat analyses were performed on the tumor model for authentication. GSC lines GSC387 and GSC3565 were derived by our laboratory and transferred via a material transfer agreement from Duke University (Durham, NC). GSC23 was derived from a recurrent GBM biopsy specimen from a 63-year-old male patient and was provided as a generous gift by Dr. Erik Sulman (NYU Langone Health). GSC3028 was derived from a recurrent GBM from a 65-year-old female patient. 1,919 GSCs were derived from a GBM from a 53-year-old male patient. All GSC lines were cultured in Neurobasal media (Invitrogen) supplemented with B27 without vitamin A (Invitrogen), EGF, and basic fibroblast growth factor (bFGF; 20 ng/mL each; R&D Systems); sodium pyruvate (catalog no. 11360070; Life Technologies); and GlutaMAX (catalog no. 35050061; Life Technologies) at 37°C with 20% oxygen and 5% carbon dioxide.

All mouse experiments were performed under an animal protocol approved by the Institutional Animal Care and Use Committee, University of California, San Diego. Investigators were not blinded for randomization or treatment. Power analysis has been done according to publications and prior experience. Studies were done in female animals. Intracranial xenografts were created by implanting 10,000 human-derived GSCs into the right cerebral cortex of NSG (NOD.Cg-Prkdcscid Il2rgtm1Wjl/SzJ; IMSR catalog no. JAX:005557, RRID:IMSR_JAX:005557; The Jackson Laboratory) mice at a depth of 3.5 mm. Healthy, wild-type female mice of NSG background, 4 to 6 weeks old, were randomly selected for intracranial injection. Housing conditions and animal status were supervised by a veterinarian. Animals were monitored until neurologic signs were observed, at which point they were sacrificed. Neurologic signs included hunched posture, gait changes, lethargy, and weight loss. Mouse survival was analyzed in GraphPad Prism (RRID:SCR_002798; GraphPad Software), and the statistical significance was tested by log-rank test.

Dissociation of xenografted tumors was performed using a papain dissociation system according to the manufacturer's instructions. Cells were then cultured in Neurobasal medium supplemented with 2% B27, 1% l-glutamine, 1% sodium pyruvate, 1% penicillin–streptomycin, 10 ng/mL bFGF, and 10 ng/mL EGF for at least 6 hours to recover expression of surface antigens. GSCs were isolated immediately following dissociation or following transient xenograft passage in immunocompromised NSG mice (The Jackson Laboratory) using prospective sorting followed by assays to confirm stem cell marker expression, sphere formation, and secondary tumor initiation. When experiments were performed, we isolated AC133-positive populations using CD133/1 antibody–conjugated magnetic beads. Neurobasal media were used to wash specimens, which were then acutely dissociated to remove nontumor tissue and subjected to enzymatic dissociation using the Papain dissociation system (catalog no. LK003150; Worthington Biomedical Corp). The isolated tumor cells were briefly placed in Neurobasal media with B27 supplement (catalog no. 12587010; Life Technologies) to permit recovery following enzymatic dissociation. Cells were labeled with the CD133/2(293C)–allophycocyanin antibody kit (catalog no. 130–098–826, RRID:AB_2660882; Miltenyi Biotec), and the CD133⁺ cells were sorted and analyzed by flow cytometry. The sorted CD133⁺ cells were cultured in NBM-B27 medium containing 20 ng/mL of both EGF (catalog no. 236-EG-01M; R&D Systems) and recombinant human bFGF (cat. 4114-TC-01M; R&D Systems) for a short period before treatment and analysis (7, 59). During routine culture, all GSC models were cultured in Neurobasal media (Invitrogen) supplemented with B27 without vitamin A (Invitrogen), EGF, and bFGF (20 ng/mL each; R&D Systems); sodium pyruvate (catalog no. 11360070; Life Technologies); and Glutamax (catalog no. 35050061; Life Technologies). Expression of stem cell markers (SOX2 and OLIG2) and functional assays of self-renewal (serial neurosphere passage) and tumor propagation using in vivo limiting dilution assays were used to validate GSC phenotypes.

RNA-seq data for 38 GSC and 5 NSC cell lines were downloaded from GSE119834 (24). Sequencing reads were trimmed using Trim Galore (RRID:SCR_011847, https://www.bioinformatics.babraham.ac.uk/projects/trim_galore/) and cutadapt (http://cutadapt.readthedocs.io/en/stable/guide.html). Trimmed reads were quantified using Salmon (RRID:SCR_017036) in the quasi mapping-based mode using –seqBias, –gcBias, and –validateMappings parameters (60). The R package “tximport” (RRID:SCR_016752) was used to import Salmon quant files into R package “DESeq2” (RRID:SCR_015687) for differential expression analysis to obtain the fold change ranking metric for GSC versus NSC comparison.

RNA-seq data for three pairs of GSCs and DGCs were downloaded from GSE54791 (25). RNA-seq–derived counts-per-million matrix was obtained and differential expression analysis was performed using R package “edgeR” to obtain the fold change ranking metric for GSC versus DGC comparison.

Whole-genome CRISPR loss-of-function screening data in eight GSCs and two NSCs were obtained from a published study (26). The BAGEL algorithm was used to derive a Bayes factor (BF) for each tested gene. A greater BF value indicates increased confidence that the knockout of the gene results in a decrease in fitness. The z scores of BF value difference between GSCs and NSCs were used to indicate the specificity of gene dependency in GSCs.

The public databases, including TCGA, CGGA, and the Genotype-Tissue Expression (GTEx) project, were used to interrogate YY1 expression levels in patients with glioma and normal tissues. UCSC Xena project recomputed TCGA and GTEx RNA-seq data were obtained from GEPIA (http://gepia.cancer-pku.cn/). YY1 expression and survival of patients with TCGA GBM and low-grade glioma were obtained from GlioVis (2). Chinese glioma data sets with 325 patients with glioma were obtained from the CGGA web server. Expression levels were categorized into two groups with the median value as a cutoff. The survival analysis was analyzed using the Kaplan–Meier curve and log-rank test. The R package “GSVA” was used to calculate a gene signature score with the ssgsea method (61).

Enrichment of the transcription regulator target program was performed using the RePhine algorithm (23). RePhine analyzed the enrichment score of 160 transcription regulator programs in a preranked gene list. First, ChIP-seq data for 160 transcription regulators were downloaded from the ENCODE database (62). To quantitatively determine the targets of a transcription regulator, we applied a regulatory potential (RP) score by considering the binding site's distance to a transcription start site (TSS) of a gene and the signal strength of ChIP-seq peaks (63). Therefore, the RP score indicates the regulatory strength of a transcription regulator on the specific gene. The higher the score is, the stronger the regulatory strength is. Next, a gene list ranked by their differential expression (log₂ fold change) between 44 GSCs and 5 NSCs was obtained. Third, we determined the concordant enrichment of transcription regulator target genes in the aforementioned gene list using elastic net regression. The coefficient from the elastic net model was calculated to indicate a positive or negative enrichment in the gene list. The P value was calculated by the likelihood ratio test to indicate the statistical significance. Finally, to visualize the enrichment of a transcription regulator program, we adopted canonical GSEA for the continuous variable (RP score) in the calculation and plot of enrichment score (64). The activity of one regulator in one sample is calculated by Regression Analysis with Background Integration algorithms using the normalized gene expression data set (65).

Cell proliferation experiments were performed by plating cells at a density of 1,000 cells/well in a 96-well plate with three to six replicate wells. CellTiter-Glo (Promega) was used to measure cell viability at days 0, 1, 3, 5, and 7. All data were normalized to day 0. In vitro limiting dilution assays were used to assess neurosphere formation capacity. Cells with decreasing numbers per well (50, 20, 10, 5, and 1) were plated into each well of a 96-well plate. The presence and number of neurospheres in each well were recorded 7 days after plating. Extreme limiting dilution analysis was performed using software available at http://bioinf.wehi.edu.au/software/elda to calculate stem cell frequency and statistical significance.

Lentiviral clones expressing nonoverlapping shRNAs against human YY1 (TRCN0000019895, TRCN0000019894, TRCN0000019896, and TRCN0000019898), human CTCF (TRCN0000230190, TRCN0000230191, and TRCN0000218498), human CDK9 (TRCN0000199892, TRCN0000000494), human BRD4 (TRCN0000021427 and TRCN0000318771), human AFF4 (TRCN0000015825 and TRCN0000426769), human MLLT1 (TRCN0000019290 and TRCN0000019291), human TRIM28 (TRCN0000018001 and TRCN0000017998), and a nontargeting control shRNA (shCONT; catalog no. SCH002) were obtained from Sigma-Aldrich. Nonoverlapping shRNAs were selected based on knockdown efficiency and used for the experiments. Human YY1-expressing plasmid was obtained from VectorBuilder. The pcDNA3.1 HA-YY1 was a gift from Richard Young (Whitehead Institute for Biomedical Research; plasmid #104395, RRID:Addgene_104395; Addgene). 293FT cells (ATCC catalog no. PTA-5077, RRID:CVCL_6911) were used to generate lentiviral particles through cotransfection of the packaging vectors pCMV-dR8.2 (plasmid #8455, RRID:Addgene_8455; Addgene) and pCI-VSV-G (plasmid #1733, RRID:Addgene_1733; Addgene) using the transfection method based on polyethylenimine (catalog no. 23966-1; Polysciences) or LipoD293 In Vitro DNA Transfection Reagent (catalog no. SL100668; SignaGen Laboratories). Media were changed after 16 hours. Virus was harvested 48 hours after media change. Viral supernatants were filtered through a 0.45-μm filter and used immediately or stored at −80°C for future use. GSCs were infected with shRNA lentivirus and the media were replaced after 24 hours. Knockdown efficiency was measured by qRT-PCR after 48 hours.

Cellular total RNA was isolated using TRIzol (catalog no. 15596026; Life Technologies), and reverse transcription to cDNA was performed using the High-Capacity cDNA Reverse Transcription Kit (catalog no. 4368814; Thermo Fisher Scientific). Briefly, 1 μg total RNA was used for cDNA synthesis using random primers according to the manufacturer's instructions. Quantitative real-time PCR was performed using the SYBR Green Master mix (catalog no. 4309155; Thermo Fisher Scientific) on an Applied Biosystems 7900HT Real-Time PCR System and normalized to 18S ribosomal RNA or GAPDH. The primers in qRT-PCR were designed using Primer3 (RRID:SCR_003139, http://bioinfo.ut.ee/primer3-0.4.0/).

GSC cells after shRNA-mediated knockdown or drug treatment were subject to total RNA isolation using the RNeasy Plus Mini Kit (catalog no. 74136; Qiagen). The total RNA was used to determine the quantification of m⁶A using the EpiQuik m⁶A RNA Methylation Quantification Kit (catalog no. P-9005; EpiGentek).

Co-IP was performed following the previously reported protocol (66). Briefly, cells were lysed in IP buffer [50 mmol/L Tris-HCl (pH 7.4), 125 mmol/L NaCl, 1 mmol/L EDTA, 0.1% Triton X-100] with inhibitors (50 mmol/L NaF, 1 mmol/L phenylmethylsulfonylfluoride, 1 mmol/L Na₃VO₄, 1 μg/mL aprotinin, 1 μg/mL leupeptin, 1 μg/mL Pepstatin) followed by sonication and centrifugation to clear insoluble debris. Lysates were then incubated with protein G agarose beads and IgG antibody of the same species as the IP antibody for 2 hours at 4°C to reduce nonspecific binding. The cleared lysates were then incubated with IP antibody overnight at 4°C. IP antibody used was YY1 (catalog no. 61779, RRID:AB_2793763; Active Motif). Beads were washed three times with IP buffer, boiled in 1× SDS gel loading buffer, and subjected to Western blot analysis.

Cells were collected and lysed with RIPA buffer [50 mmol/L Tris (pH 8.0), 1 mmol/L EDTA, 150 mmol/L NaCl, 1% Triton X-100, 0.1% SDS, 0.1% sodium deoxycholate] supplemented with EDTA-free proteinase inhibitor cocktail (catalog no. 4693159001; Roche) 48 hours after lentivirus infection to extract protein. Lysate was incubated on ice for 15 minutes and quantified using BCA reaction (catalog no. 23227; Thermo Fisher Scientific) after centrifugation for 10 minutes at 4°C. Then, 25 μg lysate was resolved on SDS-PAGE and transferred to an Immobilon-P polyvinylidene difluoride membrane (catalog no. IPVH00010; Millipore). After blocking with 3% BSA in TBST buffer [20 mmol/L Tris (pH 8.0), 150 mmol/L NaCl, 0.1% Tween 20], membranes were incubated with indicated antibodies at 4°C overnight. The antibodies used were YY1 (catalog no. 61779, RRID:AB_2793763; Active Motif), CDK9 (catalog nos. sc-484, RRID:AB_2275986 and sc-13130, RRID:AB_627245; Santa Cruz Biotechnology), phospho-RNA Polymerase II Ser2 (catalog no. A300–654A, RRID:AB_519341; BETHYL), β-actin (catalog no. A2228, RRID:AB_476697; Sigma-Aldrich), BRD4 (catalog no. 67374–1-Ig, RRID:AB_2882622; ProteinTech), AFF4 (catalog no. 14662–1-AP, RRID:AB_2242609; ProteinTech), and TRIM28 (catalog no. 15202–1-AP, RRID:AB_2209890; ProteinTech). Afterward, membranes were incubated with horseradish peroxidase–linked secondary (anti-mouse IgG, catalog no. 7076, RRID:AB_330924 and anti-rabbit IgG, catalog no. 7074, RRID:AB_2099233; Cell Signaling Technology) for 1 hour at room temperature. Then, the chemical signal was detected by chemiluminescent ECL buffer (catalog no. PI34080; Thermo Fisher Scientific).

YY1-associated drug response was identified by interrogation of drug screening data sets in 835 cancer cell lines. The cancer cell line information was retrieved from the Cancer Cell Line Encyclopedia (CCLE) project (35). The expression level of YY1 was downloaded from the CCLE web portal (https://portals.broadinstitute.org/ccle). The drug screening data (481 drugs and 835 cancer cell lines) were obtained from the CTRP project (34). AUC was used to represent drug response. A small AUC value indicates sensitivity, whereas a large AUC value indicates resistance. To identify lineage-specific gene–drug interactions, only 44 brain cancer–associated cell lines were considered in downstream analysis. Pearson correlation between gene expression levels and drug response AUC values was calculated to rank gene–drug correlation. A total of 28 drugs with P < 0.01 were filtered for hierarchical analysis and heatmap visualization using the “pheatmap” package.

Cells were seeded in 96-well plates at a density of 2,000 cells/well. The cells were then treated by a 2-fold serial diluted series of drug with at least three replicate wells per dose. The plates were transferred to an incubator (37°C, 5% CO₂) for 72 hours. At the endpoint of drug treatment, 50 μL CellTiter-Glo reagent (Promega) was added to each well. After a 10-minute incubation at room temperature, the luminescent signals were measured by an EnVision Multilabel Reader (PerkinElmer) to determine the cell viabilities. The cell viabilities were normalized to the average viability of DMSO control wells. The fitted drug–response curve was plotted using GraphPad Prism software (RRID:SCR_002798; GraphPad Software).

The ChIP assay was performed using the Magna ChIP A/G Chromatin Immunoprecipitation Kit (catalog no. 17–10085; MilliporeSigma) according to the manufacturer's instructions. Formaldehyde-fixed cells were lysed and sheared. The sheared chromatin was cleared and incubated overnight at 4°C with YY1 antibody (catalog no. 61779, RRID:AB_2793763; Active Motif) or RNA Pol II antibody (catalog no. 05–623, RRID:AB_309852; Millipore). Antibody–chromatin complexes were immunoprecipitated with protein G magnetic Dyna beads (Life Technologies), washed, eluted, reverse cross-linked, and treated with RNAse A followed by proteinase K. ChIP DNA was purified using Ampure XP beads (Beckmann Coulter) and then used to prepare sequencing libraries by the NEB NEBNext Ultra DNA Library Prep Kit for Illumina for sequencing.

Raw FASTQ files were trimmed using TrimGalore (https://www.bioinformatics.babraham.ac.uk/projects/trim_galore/) and aligned to the hg19 human genome using Bowtie2 (67). SamTools (RRID:SCR_002105) was used to sort and index the BAM files. Peak calling was performed using MACS2 with the narrow peak mode (68). For subsequent visualization of ChIP-seq signal, bigwig files were generated using the DeepTools (RRID:SCR_016366) “bamCoverage” script using default settings (http://deeptools.readthedocs.io/). DeepTools “computeMatrix,” “plotHeatmap,” and “plotProfile” functions were used to generate heat maps and profile plots. Signal tracks were visualized using R package “Gviz” and Integrative Genomics Viewer (IGV). Annotation of ChIP-seq peaks was performed using R package “ChIPseeker” (69). Gene Ontology enrichment of YY1 target genes was analyzed and visualized by ClueGO (RRID:SCR_005748; ref. 70). For pausing indexes, RPKM values were generated for hg19 TSS loci (from 50 bp upstream to 200 bp downstream of TSS) and for the gene body of these transcripts (from 500 bp downstream of the TSS to the transcription termination site). The ratio of the average coverage of the promoter over the average coverage of the gene body was then taken to be the pausing index (41, 71).

Total RNA was extracted using the Zymo Direct-zol RNA MiniPrep Kit (catalog no. R2052) and TRI reagent (catalog no. R2050) according to the manufacturer's instructions. Total RNA was used for library construction using the Illumina TruSeq Stranded Total RNA Library Prep Kit. The library was sequenced as a paired-end 150-bp read. Raw FASTQ reads were trimmed using Trim Galore, and transcript quantification was performed using Salmon in the quasi-mapping mode (60). Salmon “quant” files were converted using the Tximport function, and DEG analysis was performed using DESeq2 (72). GSEA (RRID:SCR_003199) was performed on a preranked gene list using the gene expression fold change as the ranking metric from differentially expressed gene analysis by GSEA desktop application (64). Visualization of GSEA results was performed by the Enrichment Map application in Cytoscape (RRID:SCR_003032, http://www.cytoscape.org).

BL-Hi-C library construction and data analysis were performed as described previously (30). Approximately 1 × 10⁶ cells were chemically cross-linked by the addition of 1/36 volume of fresh 37% formaldehyde solution to the medium and incubation for 10 minutes at room temperature with gentle shaking. Cross-linking was stopped by adding 2.5 mol/L glycine to a final concentration of 0.2 mol/L and incubating for 10 minutes at room temperature. After rinsing twice with PBS, cells were harvested in a 1.5-mL tube by scraping and centrifugation and stored at −80°C until use. Cell pellets were resuspended in 1 mL BL-Hi-C Lysis Buffer 1 [50 mmol/L HEPES-KOH (pH 7.5), 150 mmol/L NaCl, 1 mmol/L EDTA, 1% Triton X-100, 0.1% sodium deoxycholate, 0.1% SDS, 1× protease inhibitor cocktail] and incubated for 15 minutes on ice. After centrifugation at 1,500 × g for 5 minutes at 4°C, supernatant was removed, and after washing once with 1 mL BL-Hi-C Lysis Buffer 1, cell pellet was resuspended in 1 mL BL-Hi-C Lysis Buffer 2 [50 mmol/L HEPES-KOH (pH 7.5), 150 mmol/L NaCl, 1 mmol/L EDTA, 1% Triton X-100, 0.1% sodium deoxycholate, 1% SDS, 1× protease inhibitor cocktail] and rotated for 15 minutes at 4°C. Supernatant was removed after centrifugation. Cell pellet was washed once with 1 mL BL-Hi-C Lysis Buffer 1, resuspended in 50 μL 0.5% SDS, and incubated for 5 minutes at 62°C. At the end of incubation, SDS was quenched by adding 145 μL ddH₂O and 25 μL 10% Triton X-100 and incubated for 10 minutes at 37°C. Chromatin was cleaved by adding 25 μL 10× NEBuffer 2 and 100 U HaeIII (NEB) for 2 hours at 37°C followed by additional cleavage with another 30 U HaeIII for 3 hours. Cleaved chromatin was A-tailed by adding 2.5 μL of a 10-mmol/L dATP solution (Thermo Fisher Scientific) and 2.5 μL Klenow Fragment (3′→5′ exo-; NEB) with rotation at 900 rpm for 40 minutes at 37°C in ThermoMixer C. Proximity ligation was performed by adding 750 μL ddH₂O, 120 μL 10× T4 DNA ligase buffer, 100 μL 10% Triton X-100, 5 μL T4 DNA ligase (NEB), and 4 μL 200 ng/μL biotinylated Bridge Linker S2 (annealed by 5Phos-CGCGATATC/iBIOdT/TATCTGACT and 5Phos-GTCAGATAAGATATCGCGT) with rotation for 4 hours at 16°C. After centrifugation at 3,500 × g for 5 minutes at 4°C and removal of supernatant, the pellet was resuspended in 309 μL ddH₂O, 35 μL Lambda Exonuclease buffer, 3 μL Lambda Exonuclease, and 3 mL Exonuclease I and rocked at 900 rpm for 1 hour at 37°C in the ThermoMixer C. After the addition of 45 μL 10% SDS and 55 μL 20 mg/mL Proteinase K, sample was incubated overnight at 55°C, added in the next day with 65 μL 5 mol/L NaCl, and incubated for 2 hours at 68°C. DNA was recovered by extraction with 500 μL phenol/chloroform/isoamyl alcohol (25:24:1) and ethanol precipitation with 1 μL glycoblue. After centrifugation and washing with 75% ethanol, DNA pellet was resolved in 130 μL Buffer EB (10 mmol/L Tris-HCl, pH 8.0). DNA was sheared with Covaris S220 with the setting for DNA size of 300 bp. After washing twice with 2× B&W buffer [10 mmol/L Tris-HCl (pH 7.5), 1 mmol/L EDTA, 2 mol/L NaCl], 30 μL Dynabeads M-280 Streptavidin (Thermo Fisher Scientific) was blocked with 100 μL 1× I-Block buffer [2% I-Block Protein-Based Blocking Reagent (Thermo Fisher Scientific), 0.5% SDS] for 45 minutes at room temperature in a rotating wheel. Beads were washed twice with 100 μL 1× B&W buffer [5 mmol/L Tris-HCl (pH 7.5), 0.5 mmol/L EDTA, 1 mol/L NaCl] followed by resuspension in 1× B&W buffer containing 1 μg preheated Salmon Sperm DNA solution and rotation for 30 minutes at room temperature. After washing twice with 100 μL 1× B&W buffer, beads were resuspended with 130 μL 2× B&W buffer, combined with sonicated DNA, and rotated for another 45 minutes at room temperature. Beads were washed five times with 500 μL 2× SSC containing 0.5% SDS, twice with 500 μL 1× B&W buffer, and once with 100 μL Buffer EB (Qiagen). DNA on beads was end-repaired in the reaction containing 75 μL ddH₂O, 10 μL 10× T4 DNA ligase buffer, 5 μL 10 mmol/L dNTP, 5 μL T4 polynucleotide kinase, 4 μL T4 DNA polymerase, and 1 μL large (Klenow) fragment with shaking at 900 rpm for 30 minutes at 37°C in a Thermomixer C. After washing twice with 500 μL 1× TWB [5 mmol/L Tris-HCl (pH 7.5), 0.5 mmol/L EDTA, 1 mol/L NaCl, 0.05% Tween-20] for 2 minutes at 55°C, DNA on beads was A-tailed in the reaction containing 80 μL ddH₂O, 10 μL NEBuffer 2, 5 μL 10 mmol/L dATP, and 5 μL Klenow fragment (3′→5′ exo-) with shaking at 900 rpm for 30 minutes at 37°C in a ThermoMixer C. Beads were washed twice with 500 μL 1× TWB for 2 minutes at 55°C and once with 50 μL of 1× Quick Ligase buffer. DNA on beads was ligated with adaptor in the reaction containing 6.6 μL ddH₂O, 10 μL 2× Quick Ligase buffer, 2 μL Quick Ligase (NEB), and 0.4 μL 20 μmol/L Y-Adaptor (annealed by/5Phos/GATCGGAAGAGCACACGTCTGAACTCCAGTCAC and TACACTCTTTCCCTACACGACGCTCTTCCGATCT) for 45 minutes at room temperature. Beads were washed twice with 500 μL 1× TWB for 2 minutes at 55°C and once with 100 μL Buffer EB. After resuspending in 60 μL Buffer EB, bead suspension was aliquoted in 3 × 20 μL for storage at −20°C. One aliquot of bead suspension was used as a template for PCR amplification with Q5 Hot Star DNA Polymerase, universal primer (AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGAC), and index primer (CAAGCAGAAGACGGCATACGAGAT/index/GTGACTGGAGTTCAGACGTGT) for 10 cycles. PCR products between 300 and 700 bp were purified as the BL-Hi-C library using Ampure XP beads and subjected to Illumina HiSeq ×10 for paired-end 150-bp sequencing.

The ChIA-PET2 v0.9.3 software (73) was used for quality control, mapping to the hg19 genome, and identification of chromatin interactions with the following parameter setting: -A ACGCGATATCTTATC -B AGTCAGATAAGATAT -s 1 -m 1 -t 24 -k 2 -e 1 -l 15 -S 500 -M “-q 0.05.” HiC-Pro (RRID:SCR_017643) hicpro2juicebox.sh utility was used to generate the hic file for further visualization of interaction matrix (74). A/B compartment analysis at a 50-kb resolution was performed using Juicer Tools eigenvector function and visualized in IGV. The interaction matrix and chromatin interactions were visualized in WashU Epigenome Browser (https://epigenomegateway.wustl.edu/) or Juicebox v1.11.08. The overlapping between YY1 binding peaks and chromatin interaction anchors was analyzed using bedtools (RRID:SCR_006646) pairtobed function with default parameters and further analyzed in R.

For the alvocidib treatment experiment, intracranial xenografts were created by implanting 6,000 human-derived GSC1517 cells into female immunodeficient mice. After 7 days, immunodeficient mice were randomly selected into two groups and administrated with saline or alvocidib (10 mg/kg, catalog no. HY-10006; MedChemExpress) three times a week until neurologic signs were observed.

For the alvocidib and anti–PD-1 treatment experiment, intracranial xenografts were created by implanting 20,000 mouse glioma cells (CT2A) into the right cerebral cortex of 4- to 6-week-old C57BL/6 female mice (Jackson Laboratory) at a depth of 3.5 mm. After 5 days, mice were randomly selected into four groups (control, alvocidib only, anti–PD-1 only, alvocidib + anti–PD-1). During each of the following 5 weeks, mice were treated with alvocidib (10 mg/kg) or saline for 2 consecutive days, followed by anti–PD-1 (200 μg, cat. BP0273; Bio X Cell) or IgG2a isotype control (200 μg, cat. BP0089; Bio X Cell) on the third day. Animals were monitored until neurologic signs were observed, at which point they were sacrificed.

Mouse brains implanted with GSCs or mouse glioma cells that were labeled with firefly luciferase were monitored by the bioluminescence imaging. Animals were treated with D-Luciferin (120 mg/kg, cat. L-8220; Biosynth Carbosynth) intraperitoneally and anesthetized with isoflurane for the imaging analysis. The bioluminescence images were captured by an IVIS imaging system (Spectrum CT; PerkinElmer).

Mouse glioma tissues were harvested according to the following protocol. Briefly, half of the mouse brain with tumors was chopped with a scalpel in RPMI 1640 media with 50 μg/mL Liberase TL (Roche) and 0.25 mg/mL DNase I (Roche) and incubated for 35 minutes at 190 rpm at 37°C. The digestion was stopped with cold complete RPMI 1640 media. Tumor-infiltrating immune cells were separated by a 30%–37%–70% Percoll gradient centrifugation. The red blood cells were removed by ACK Lysing Buffer (Thermo Fisher Scientific). The remaining single-cell suspension was analyzed by a routine CyTOF protocol in the Flow Cytometry Core Facility in the La Jolla Institute for Immunology. Briefly, after staining with cisplatin viability staining reagent, cells were stained with antibodies conjugated to metal isotopes. Cells were also stained with DNA intercalators for detection of viable cells. The single cell signal was recorded by a Helios CyTOF Mass Cytometer (Fluidigm). The normalized FCS files were analyzed by R package “cytofWorkflow.” Cell population clusters were identified using the FlowSOM (RRID:SCR_016899) and ConsensusClusterPlus (RRID:SCR_016954) packages and visualized by the t-SNE algorithm. The antibody panel used in this work consisted of 34 antibodies, which were designed by the Flow Cytometry Core Facility in the La Jolla Institute for Immunology to identify major immune cell populations, including antibodies to CD3 (catalog no. 100302, RRID:AB_312667; BioLegend), CD4 (catalog no. 3172003B, RRID:AB_2811242; Fluidigm), CD8a (catalog no. 3153012B; Fluidigm), CD11b (catalog no. 101202, RRID:AB_312785; BioLegend), CD11c (catalog no. 3142003B, RRID:AB_2814737; Fluidigm), CD64 (catalog no. 139302, RRID:AB_10613107; BioLegend), CD25 (catalog no. 3151007B, RRID:AB_2827880; Fluidigm), FR4 (catalog no. 125102, RRID:AB_1027724; BioLegend), MERTK (catalog no. 151502, RRID:AB_2566624; BioLegend), NK1.1 (catalog no. 108702, RRID:AB_313389; BioLegend), IgM (catalog no. 406502, RRID:AB_315052; BioLegend), B220 (catalog no. 3160012B; Fluidigm), Ly6c (catalog no. 128002, RRID:AB_1134214; BioLegend) and Ly6g (catalog no. 127602, RRID:AB_1089180; BioLegend).

Statistical tests, including Student t test, one-way ANOVA, Pearson or Spearman correlation test, and Kolmogorov–Smirnov normality test, were performed using R or GraphPad Prism v8 software (RRID:SCR_002798; GraphPad Software). Data in the bar plot or curve are presented as mean ± SD or mean + SD as denoted in each analysis. For the box-and-whisker plot in the Tukey style, the middle line of box indicates the median value, and the top and bottom lines of box represent the upper and lower quartiles (Q1 and Q3), respectively. All statistical tests were two-tailed. P < 0.05 was taken to indicate statistical significance unless otherwise stated.

All raw sequencing data and selected processed data are available on the Gene Expression Omnibus through accession number GSE155839. Additional data will be provided upon request.

Y. Zheng reports grants from NIH during the conduct of the study. J.N. Rich reports grants and personal fees from Synchronicity Pharma outside the submitted work. No disclosures were reported by the other authors.

Z. Qiu: Conceptualization, software, investigation, methodology, writing–original draft, writing–review and editing. L. Zhao: Investigation, methodology, writing–review and editing. J.Z. Shen: Resources, investigation, writing–review and editing. Z. Liang: Investigation, methodology. Q. Wu: Investigation, writing–review and editing. K. Yang: Investigation, writing–review and editing. L. Min: Investigation. R.C. Gimple: Investigation. Q. Yang: Investigation. S. Bhargava: Software. C. Jin: Investigation. C. Kim: Investigation. D. Hinz: Investigation. D. Dixit: Investigation. J.A. Bernatchez: Investigation. B.C. Prager: Investigation. G. Zhang: Investigation. Z. Dong: Investigation. D. Lv: Investigation. X. Wang: Software. L.J.Y. Kim: Software. Z. Zhu: Investigation. K.A. Jones: Resources. Y. Zheng: Resources. X. Wang: Investigation. J.L. Siqueira-Neto: Resources. L. Chavez: Investigation. X.-D. Fu: Resources. C. Spruck: Resources, supervision, writing–review and editing. J.N. Rich: Conceptualization, resources, supervision, funding acquisition, writing–original draft, project administration, writing–review and editing.

We thank the Tissue Technology Shared Resource in UCSD and Histology Core Facility in Sanford Burnham Prebys Medical Discovery Institute for histology analysis, the Flow Cytometry Core Facility in La Jolla Institute for Immunology for CyTOF analysis, Erik Sulman for providing the GSC23 GBM cell model, and Mary E. Donohoe for critical reading of the manuscript. Tissue Technology Shared Resource is supported by a National Cancer Institute Cancer Center Support grant (CCSG Grant P30CA23100). The CyTOF Helios was acquired through the Shared Instrumentation Grant (SIG) Program CyTOF Mass Cytometer S10 OD018499. Several panels were prepared in part using images from Servier Medical Art by Servier (https://smart.servier.com/), which is licensed under a Creative Commons Attribution 3.0 Unported License (https://creativecommons.org/licenses/by/3.0/).We thank our funding sources: NIH grants (CA238662, CA197718, and NS103434 to J.N. Rich; CA217065 to R.C. Gimple; CA217066 to B.C. Prager; T32CA094186 to K. Yang), Young Investigator Award in Glioblastoma from ASCO Conquer Cancer Foundation (to K. Yang), RSNA Research Resident Grant (to K. Yang), and DoD Breast Cancer Research Program (W81XWH-15-1-0383 to C. Spruck).

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.