RBFOX2 recognizes N6-methyladenosine to suppress transcription and block myeloid leukaemia differentiation

Abstract

N⁶-methyladenosine (m⁶A) methylation can be deposited on chromatin-associated RNAs (caRNAs) by the RNA methyltransferase complex (MTC) to regulate chromatin state and transcription. However, the mechanism by which MTC is recruited to distinct genomic loci remains elusive. Here we identify RBFOX2, a well-studied RNA-binding protein, as a chromatin factor that preferentially recognizes m⁶A on caRNAs. RBFOX2 can recruit RBM15, an MTC component, to facilitate methylation of promoter-associated RNAs. RBM15 also physically interacts with YTHDC1 and recruits polycomb repressive complex 2 (PRC2) to the RBFOX2-bound loci for chromatin silencing and transcription suppression. Furthermore, we found that this RBFOX2/m⁶A/RBM15/YTHDC1/PRC2 axis plays a critical role in myeloid leukaemia. Downregulation of RBFOX2 notably inhibits survival/proliferation of acute myeloid leukaemia cells and promotes their myeloid differentiation. RBFOX2 is also required for self-renewal of leukaemia stem/initiation cells and acute myeloid leukaemia maintenance. Our study presents a pathway of m⁶A MTC recruitment and m⁶A deposition on caRNAs, resulting in locus-selective chromatin regulation, which has potential therapeutic implications in leukaemia.

Fig. 1. RBFOX2 recognizes m⁶A on paRNA.

a, Top: a protein–protein interaction network with the top-ranked (n = 50) m⁶A-associated proteins on caRNAs in K562 cells. Bottom: Venn diagram of top-ranked (n = 20) m⁶A-associated proteins on caRNAs between K562 and HepG2 cells. The highlighted protein (dotted circle) is shared between K562 and HepG2 cells. b, Average profile (top) and heat map (bottom) showing RBFOX2 binding intensity at m⁶A peak centres and the flanking 2.5 kb regions in K562 cells. rep, replicate. c, Average profile of RBFOX2 binding intensity at RBFOX2 peak centres and the flanking 2.5 kb regions in K562 cells. RBFOX2 peaks^, were categorized into two groups according to whether they overlapped with m⁶A (m⁶A) or not (non-m⁶A). d, LC–MS/MS showing m⁶A enrichment in RBFOX2-bound RNA while depleted in the FT portion (n = 6, three technical replicates over two biological replicates). Data are represented as mean values ± standard deviation. Two-sided P value was calculated by Student’s t-test. e, Distribution of m⁶A or RBFOX2 peaks^, at distinct genomic regions including promoter, exonic, intronic, transcription termination sites (TTS) and intergenic regions annotated by HOMER in K562 cells. ‘(+)’ sign represents regions harbouring m⁶A (+) or bound by RBFOX2 (+), while the sign ‘(−)’ indicates the absence of m⁶A (−) or RBFOX2 (−). f, Top consensus sequences on RBFOX2-bound sites that were marked with m⁶A (caRNA m⁶A-SAC-seq) using HOMER motif discovery algorithm. g, Oligo pulldown assay showing RBFOX2 (top) and RRM domain of RBFOX2 (bottom) bound an m⁶A-containing RNA probe with higher affinity than the unmethylated control. Source data

Fig. 2. RBFOX2 recruits m⁶A MTC to gene promoter for RNA methylation.

a, Cumulative curve and box plot (inset) of m⁶A log₂(fold change (FC)) comparing RBFOX2 KD (shRBFOX2) versus control (shNS) K562 cells. m⁶A peaks (N) were categorized into two groups according to whether they overlapped with RBFOX2 peaks (+, N = 5,089) or not (−, N = 23,381) (refs. ^,). b, Cumulative curve and box plot (inset) of RNA log₂FC comparing RBFOX2 KD versus control K562 cells. caRNAs (N) were categorized into two groups according to whether they are derived from regions with RBFOX2 binding (+, N = 5,089) or not (−, N = 23,381) (refs. ^,). c, Western blots of the immunoprecipitated RBFOX2 from K562 cells and its interactions with METTL3, METTL14 and RBM15 after RNase A/T1 treatment. d, Venn diagram of overlap between RBFOX2 ChIP–seq and RBM15 eCLIP peaks in K562 cells^,. Two-sided P value was calculated by Fisher’s exact test. e, Average profile of RBFOX2 binding intensity at RBFOX2 peak centres and the flanking 2.5 kb regions in K562 cells. RBFOX2 peaks were categorized into two groups according to whether they overlap with RBM15 eCLIP peaks (+) or not (−)^,. f, Box plot of m⁶A log₂FC comparing RBFOX2 KD versus control K562 cells. m⁶A peaks (N) were categorized into three groups according to whether they overlapped with RBFOX2 or RBM15 peaks^,. From left to right, N = 22,061, 4,049 and 1,040. g, Average profile of RBM15 binding intensity at RBFOX2 peak centres and the flanking 2.5 kb regions in K562 cells. RBM15 peaks were categorized into two groups according to whether they overlap with RBFOX2 (+) or not (−) (refs. ^,). h, The integrative genomics viewer plots showing RBFOX2 ChIP–seq signals in wild-type K562 cells^,, and RBM15 ChIP–seq signals in control and RBFOX2 KD K562 cells around H2AC8/H2BC8 (left) and SOCS1 (right) gene loci. Source data

Fig. 3. RBFOX2 regulates chromatin state through the m⁶A/RBM15/YTHDC1/PRC2 axis.

a, Average profile (top) and heat map (bottom) showing RBFOX2, RBM15, YTHDC1 and the corresponding input signal at RBFOX2 peak centres and the flanking 2.5 kb regions in K562 cells. b, Average profile of RBFOX2 binding intensity at RBFOX2 peak centres and the flanking 2.5 kb regions in K562 cells. RBFOX2 peaks were categorized into two groups according to whether they overlap with YTHDC1 (+) or not (−). c, In situ PLA assay detecting the interaction (green) between RBFOX2 and YTHDC1 in K562 cells. Nuclei were stained with Hoechst 33342 (blue). Scale bar, 5 µm. d, Western blots of the immunoprecipitated RBM15, EZH2 and SUZ12, respectively, from K562 cells and their interactions with YTHDC1 after RNase A/T1 treatment. e, Average profile showing SUZ12 binding signal at their peak centres and the flanking 2.5 kb regions in YTHDC1 KD (siYTHDC1) versus control (siControl) K562 cells. SUZ12 peaks were categorized into two groups according to whether they overlap with RBFOX2 (+) or not (−). f,g, Average profile showing SUZ12 binding signal at their peak centres and the flanking 2.5 kb regions in RBFOX2 KD versus control K562 cells. SUZ12 peaks were categorized into two groups according to whether (+) or not (−) they were overlapped with RBM15 (f) or YTHDC1 (g). The depicted genome-wide data represent an integration of all samples, including two biologically independent replicates. Source data

Fig. 4. RBFOX2 depletion promotes the differentiation of leukaemia cells.

a, Colony-forming assays of K562 cells transduced with control (shNS) or two RBFOX2 shRNAs (shRBFOX2-1 and shRBFOX2-2). n = 3. b, Effects of RBFOX2 KD on K562 cell growth. n = 3. P value was calculated by two-way ANOVA, Dunnett’s multiple comparisons test. c, Flow cytometric analysis of CD61⁺ cell populations in K562 cells transduced with shNS or RBFOX2 shRNAs (shRBFOX2-1 and shRBFOX2-2). d, Relative expression of CD61 in control (shNS) and RBFOX2 KD (shRBFOX2-1 and shRBFOX2-2) K562 cells. n = 3. Two-sided P value was calculated by Student’s t-test. e, Wright-Giemsa staining of cytospin slides of K562 cells transduced with shNS or RBFOX2 shRNAs (shRBFOX2-1 and shRBFOX2-2). Arrows indicate differentiated cells. Scale bar, 40 µm. f, Flow cytometric analysis of CD61⁺ cell populations in control group (shNS + EV, control shRNA (shNS) with empty vector (EV)), RBFOX2 KD group (shRBFOX2 + EV, RBFOX2 KD with EV) and RBFOX2 rescue group (shRBFOX2 + RBFOX2, RBFOX2 KD with RBFOX2 overexpression) in K562 cells (left), and statistics from three biological replicates (right). EV or RBFOX2 overexpressed cells were labelled with green fluorescent protein (GFP), shNS or shRBFOX2 transduced cells expressed mCherry. GFP⁺mCherry⁺ cells represent positively double transduced cells. n = 2. g, Flow cytometric analysis of CD11b⁺ cell populations in control (shNS) and RBFOX2 KD (shRBFOX2-1 and shRBFOX2-2) NB4 cells. h, Relative expression of CD11b in control (shNS) and two RBFOX2 KD (shRBFOX2-1 and shRBFOX2-2) NB4 cells. n = 3. i, Wright-Giemsa staining of cytospin slides of NB4 cells transduced with shNS or RBFOX2 shRNAs (shRBFOX2-1 and shRBFOX2-2). Arrows indicate differentiated cells. Scale bar, 20 µm. n, biologically independent samples. Data are presented as mean ± standard error of the mean. For a and h, P values were calculated by one-way ANOVA, Dunnett’s multiple comparisons test. Source data

Fig. 5. RBFOX2 is aberrantly expressed in AML and RBFOX2 depletion impairs leukaemia progression in vivo.

a, Expression levels of RBFOX2 in patients with primary AML bearing chromosomal translocations and those in bone marrow (BM) HSCs collected from healthy donors (healthy BM) (GSE13159 (ref. )). n = 73 for healthy BM, n = 38 for MLL, n = 38 for t(15;17) and n = 351 for normal karyotype. Two-sided P values were calculated by Student’s t-test. b, Kaplan–Meier survival analysis in GSE1159 (ref. ) and GSE14468 (ref. ) dataset (n = 485). The patients were divided into two groups of equal size based on RBFOX2 levels. c, In vitro LDAs. Logarithmic plot showing the percentage of non-responding wells at different doses. Non-responding wells are wells not containing colony-forming cells. The estimated LSC/LIC frequency is calculated by ELDA and shown on the right. d, Bioluminescence imaging of mice transplanted with luciferase-expressing MOLM13 cells transfected with control (shNS) and RBFOX2 shRNAs (shRBFOX2-1 and shRBFOX2-2), respectively. D, day. e, Kaplan–Meier survival curves of recipient mice transplanted with control (n = 9) and two RBFOX2 KD (n = 10 for shRBFOX2-1 and n = 9 for shRBFOX2-2) MOLM13 cells. f, Flow cytometric analysis (top) and quantification (bottom) of CD11b⁺ cell populations in control and RBFOX2 KD AML-PDX cells (2017-129 (ref. )). n = 2. g, Bioluminescence imaging of mice transplanted with luciferase-expressing AML-PDX cells (2017-129 (ref. )) transduced with shNS and RBFOX2 shRNAs (shRBFOX2-1 and shRBFOX2-2), respectively. h, Kaplan–Meier survival curves of recipient mice transplanted with control (n = 8) and RBFOX2 KD (n = 7 for shRBFOX2-1 and n = 8 for shRBFOX2-2) AML-PDX cells (2017-129 (ref. )). n, biologically independent samples. Data are presented as mean ± standard error of the mean. For b, c, e and h, P values were calculated by the log-rank test. Source data

Fig. 6. RBFOX2 suppresses TGFB1 transcription to promote tumourigenesis in human AML.

a, Kaplan–Meier survival analysis in TCGA-AML dataset (n = 106) using GEPIA. The patients were divided into two groups according to gene expression level of RBM15. P value was determined by the log-rank test. b, Top: correlation of gene expression log₂FC between RBFOX2 KD versus control and RBM15 KD versus control K562 cells^,. Here genes are required to be bound by both RBFOX2 and RBM15, and genes modified with m⁶A are highlighted. Bottom: bar plot (bottom) showing the functional enrichment analysis of genes harbouring m⁶A and regulated by RBFOX2 and RBM15. PCC was employed for correlation analysis, and the corresponding P value was obtained. c, The integrative genomics viewer plots showing the binding profiles of RBFOX2 and RBM15, and m⁶A level in K562 cells around the TGFB1 gene locus. In eCLIP-seq, ‘(+)’ represents the reads from the forward strand and ‘(−)’ represents the reads from the reverse strand to the genome. d, A representative flow cytometric (left) and quantification (right) analysis of CD11b⁺ cell populations in control, RBFOX2 KD and RBFOX2 & TGFB1 double KD NB4 cells, respectively. n = 3 biologically independent samples. Data are presented as mean ± standard error of the mean. Two-sided P values were calculated by Student’s t-test. e, A proposed model depicts RBFOX2 regulation on tumourigenesis through the m⁶A/RBM15/YTHDC1/PRC2 axis. Source data

Extended Data Fig. 1. Identification of proteins associated with m⁶A on chromatin-associated RNAs.

a, A diagram showing the analysis pipeline for identification of proteins associated with m⁶A on chromatin-associated RNAs. b-c, Heatmap showing the association score of binding profile between proteins (ENCODE project) and m⁶A of chromatin-associated RNAs (this study) in K562 (b) and HepG2 (c) cells. d, A protein-protein interaction (PPI) network with the top-ranked (n = 50) m⁶A associated proteins on chromatin-associated RNAs in HepG2 cells. The depicted genome-wide data represent an integration of all samples, including two biologically independent replicates.

Extended Data Fig. 2. RBFOX2 binds m⁶A-marked promoter-associated RNAs.

a-b, Peaks overlap between RBFOX2 (ChIP-seq, ENCODE project) and m⁶A (MeRIP-seq) of caRNAs in K562 (a) and HepG2 (b) cells. Two-sided P values by Fisher’s exact test. c, Average profile of RBFOX2 binding intensity at RBFOX2 peaks (m⁶A versus non-m⁶A, HepG2). d, LC-MS/MS showing m⁶A enrichment in RBFOX2-bound RNA while depleted in the flow-through (FT) portion (HepG2) (n = 6, three technical replicates over two biologically replicates). Mean values +/− SD. Two-sided P value by Student’s t-test. e-f, Western blot showing RBFOX2 pulled down by RNA probes containing m⁶A (e, 5′CGUGG(A/m⁶A)CUGGCUU-3′) or m⁷G/m⁶A_m/G_m/A_m/Ψ (f), as well as their corresponding unmethylated control probes. g, METTL3 expression in K562 cells. h, Average profile of RBFOX2 binding at RBFOX2 peaks in K562 cells. i, RBFOX2 binding changes at RBFOX2 binding sites (N = 2082 m⁶A methylated versus N = 20846 unmethylated ones) upon METTL3 knockdown in K562 cells. j, RBFOX2 CLIP binding log₂FC upon METTL3 KD in K562 cells. The RBFOX2 bindings on chromatin (N) were separated into three groups: those with m⁶A at their promoter regions (RBFOX2 may function as an m⁶A-binding protein, N = 272), those without m⁶A at their intronic regions (RBFOX2 may function as a splicing regulator, N = 1668), and all detected genes as background (N = 3673). k, The percentage of alternatively spliced genes by RBFOX2 knockdown in K562 cells. l, m⁶A level in RBFOX2-bound caRNAs (IP) and the corresponding input caRNAs (Input) in K562 cells (N = 275 m⁶A sites, m⁶A-SAC-seq). Mean values +/− SEM. m, The percentage of m⁶A sites in different sequence contexts (motifs) in K562 cells (m⁶A-SAC-seq). n, Calibration curve for m⁶A within the AGAUG motif is generated by linear regression. o-p, The IGV visualization (o) and quantification (p) of the mutation frequency of identified m⁶A sites (m⁶A-SAC-seq). q, Oligo pulldown assay showing full length RBFOX2 bound to an m⁶A-containing RNA probe with higher affinity than the unmethylated control. r, EMSA assay measuring the dissociation constant (Kd, nM) of His-tagged MBP-RBFOX2 protein with methylated or unmethylated RNA probes. Source data

Extended Data Fig. 3. RBFOX2 recruits m⁶A MTC to gene promoter for RNA methylation installation.

a, RBFOX2 expression in K562 cells. b, H3K4me3 level changes at H3K4me3 peaks (N) by RBFOX2 KD in K562 cells. H3K4me3 peak regions were categorized into six groups based their overlap with RBFOX2 or m⁶A. From left to right, N = 9767, 11201, 14907, 3706, 3884,178. c, Correlation of log₂FC between m⁶A level and caRNA abundance by RBFOX2 KD in K562 cells. P value by PCC. d, The elongation rate changes comparing RBFOX2 KD with control K562 cells. Transcripts (N) were categorized into three groups based on their overlap with RBFOX2 or m⁶A. From left to right, N = 130274, 13232, and 5678. e, In situ PLA detecting the close proximity (green) between RBFOX2 and METTL3 or METTL14 in K562 cells. Scale bar, 5μM. f, The binding intensity of METTL3 and METTL14 at m⁶A, RBFOX2 and RBM15 marked promoter regions, respectively. g, Distribution of RBM15 (eCLIP) or RBFOX2 (ChIP) peaks at distinct genomic regions annotated by HOMER (K562, left panel). Functional enrichment analysis of genes bound by both RBFOX2 and RBM15 (K562, right panel). One-sided P value by Fisher’s Exact test. h, m⁶A methylation level on peaks (N) occupied by either RBFOX2 or RBM15 (K562). From left to right, N = 10064, 1322 and 3187. i, Venn diagram of peak overlap between RBFOX2 (ChIP-seq) and RBM15 (eCLIP) in K562 cells (left panel). Here RBFOX2 and RMB15 peaks were required to overlap with caRNA m⁶A peaks. Barplot showing the proportion of m⁶A peaks (RBM15 bound versus unbound) that overlap with RBFOX2 peaks (right panel). Two-sided P value was calculated by Fisher’s exact test. j, m⁶A log₂FC at m⁶A peak regions (N, RBM15 bound versus unbound) in RBFOX2 KD versus control K562 cells. From left to right, N = 26110 and 2360. k, RBM15 expression in K562 cells. l, Average profile of RBFOX2 binding intensity at RBFOX2 peaks (K562). The depicted genome-wide data represent an integration of all samples, including two biologically independent replicates. Source data

Extended Data Fig. 4. RBFOX2 regulates chromatin state through the m⁶A/RBM15/YTHDC1/PRC2 axis.

a, Western blots of immunoprecipitated RBFOX2 and RBM15 from K562 cells and their interactions with YTHDC1 after RNase A/T1 treatment. b, Heatmap showing the binding intensity of METTL3 and METTL14 at YTHDC1-marked promoter regions. c, Barplot showing the proportion of YTHDC1 peaks overlapping with RBFOX2 or RMB15 in K562 cells. d, Immunofluorescence images (left panel) and intensity profile (right panel) of RBFOX2 protein (red) and YTHDC1 protein (yellow) showing their colocalization in K562 cells. Scale bar, 5μM. e, Average profile showing EZH2 level at their peak centers and the flanking 2.5 kb regions in YTHDC1 KD (siYTHDC1) versus control (siControl) K562 cells. EZH2 peaks were categorized into RBFOX2 bound (+) and unbound (−) subgroups. f, Average profile showing SUZ12 binding signals at their peak centers and the flanking 2.5 kb regions in RBFOX2 KD versus control K562 cells. SUZ12 peaks were categorized into RBFOX2 bound (+) and unbound (−) subgroups. g, Cumulative distribution of gene expression changes in siControl K562 cells upon RBFOX2 knockdown versus control, as well as in YTHDC1 knockdown K562 cells upon RBFOX2 knockdown versus control. The displayed genes (N = 1221) were those that were upregulated upon RBFOX2 knockdown in siControl K562 cells and bound by RBFOX2 at their promoter regions. h, Average profile of RBFOX2 (left panel) or RBM15 (right panel) at their respective peak centers and the flanking 2.5 kb regions in Control (siControl) versus YTHDC1 knockdown (siYTHDC1) K562 cells. Source data

Extended Data Fig. 5. RBFOX2 knockdown inhibits K562 cell growth and promotes K562 cell differentiation towards megakaryocyte lineage.

a, Western blot of RBFOX2 in K562 cells transfected with shNS or RBFOX2 shRNAs (shRBFOX2-1 and shRBFOX2-2). b, Percentages of cells with CD61⁺ staining in K562 cells transduced with shNS or RBFOX2 shRNAs (shRBFOX2-1 and shRBFOX2-2). n = 10000 cells examined over 4 biologically independent experiments. c, Western blot showing RBFOX2 expression in control group (shNS+EV, control shRNA (shNS) with empty vector (EV)), RBFOX2 KD group (shRBFOX2 + EV, RBFOX2 KD with EV) and RBFOX2 rescue group (shRBFOX2 + RBFOX2, RBFOX2 KD with RBFOX2 overexpression) in K562 cells. d-e, Relative expression of RBFOX2 during hematopoietic cell differentiation into different lineages based on GSE42519 (d) and GSE24759 (e) datasets. The sample size shown in (d) and (e) are biologically independent human samples. Data are presented as mean values +/− SEM. Two-sided P values by Student’s t-test. Source data

Extended Data Fig. 6. RBFOX2 knockdown promotes the differentiation of leukemia cells.

a, Western blot assay to quantify RBFOX2 knockdown efficiency in NB4 cells. b-c, Flow cytometric analysis of CD11b⁺ cell populations (b), and relative gene expression of CD14 and CD11b with qPCR quantification (c) in control (shNS) and RBFOX2 KD (shRBFOX2-1 and shRBFOX2-2) NB4 cells after treatment with 500 nM ATRA for 72 hr as compared to DMSO-treated cells. n = 3. d, Western blot assay to quantify RBFOX2 knockdown efficiency in MOLM13 cells. e-f, Flow cytometric analysis of CD11b⁺ cell populations (e) and CD11b mRNA expression level (f) in control (shNS) and RBFOX2 KD (shRBFOX2-1 and shRBFOX2-2) MOLM13 cells. n = 3. g, Wright-Giemsa staining of cytospin slides of control (shNS) and RBFOX2 KD (shRBFOX2-1 and shRBFOX2-2) MOLM13 cells. Arrows indicate differentiated cells. Scale bar, 20 µm. h, Western blots showing the knockdown efficiency of RBFOX2 in human cord blood-derived CD34⁺ (CB-CD34⁺) cells. i, Flow cytometric analysis of the percentage of CD34⁺ cell population in control (shNS) and RBFOX2 KD (shRBFOX2-1 and shRBFOX2-2) CB-CD34⁺ cells. Day 2 and Day 5 represent the number of days post-transduction. j, Quantification of the percentage of different hematopoietic lineage progenitor colonies in control (shNS) and RBFOX2 KD (shRBFOX2-1 and shRBFOX2-2) CB-CD34⁺ cells. n = 4. P values were calculated by two-way ANOVA, Dunnett’s multiple comparisons test. k-n, Flow cytometric analysis of megakaryocyte (k), erythroid (l), granulocyte (m) and monocyte (n) differentiation of CB-CD34⁺ cells (five days post-transduction) transfected with shNS or RBFOX2 shRNAs (shRBFOX2-1 and shRBFOX2-2). n, biologically independent samples. Data are presented as mean values +/− SEM. For c and f, two-sided P values by Student’s t-test. Source data

Extended Data Fig. 7. RBFOX2 knockdown blocks LSC self-renewal, and impairs AML progression in vivo.

a, Western blot of RBFOX2 protein in normal controls (human bone marrow CD43⁺) and AML cell lines. b, Kaplan-Meier survival analysis in TCGA-AML dataset (n = 106) using GEPIA. The patients were divided into two groups of equal size based on RBFOX2 levels. n, biologically independent patients. P value was detected by the log-rank test. c-d, Relative Rbfox2 gene expression quantified by qPCR (c, n = 3) and western blot assay (d) in control (shNS) and Rbfox2 KD (shRbfox2-1 and shRbfox2-2) mouse MLL-AF9 (MF9) cells. e-f, Image (e) and quantification (f, n = 2) of effects of Rbfox2 knockdown on the colony-forming of mouse MA9 AML cells. g, Flow cytometric analysis of LSC populations in mouse MA9 leukemia cells. h, Flow cytometric analysis of CD11b⁺ cell populations in mouse MA9 leukemia cells. i, Western blot of RBFOX2 protein in control (shNS) and RBFOX2 KD (shRBFOX2-1 and shRBFOX2-2) in in vivo bone marrow samples. j, Western blot showing RBFOX2 expression in control and RBFOX2 KD (shRBFOX2-1 and shRBFOX2-2) AML-PDX cells. k, Effects of RBFOX2 knockdown (shNS [n = 4], shRBFOX2-1 [n = 4] and shRBFOX2-2 [n = 3]) on AML-PDX cell growth. P value was calculated by one-way ANOVA, Dunnett’s multiple comparisons test. l, Flow cytometric analysis of CD11b⁺ cell populations in control and RBFOX2 KD AML-PDX cell (2017-129). m, Flow cytometric analysis of CD11b⁺ (top panel) or CD11b⁺ and CD15⁺ (bottom panel) cell populations in control and RBFOX2 KD AML-PDX cell (HBT22-0148). n, The hCD45 and hCD33 double positive cells (CD45⁺CD33⁺) were utilized to determine the engraftment of human AML-PDX cells (2017-129) in recipient mice on day 34 post xeno-transplantation (shNS [n = 8], shRBFOX2-1 [n = 7] and shRBFOX2-2 [n = 9]). o, Quantification of leukemia burden in control (shNS [n = 5]) and RBFOX2 KD (shRBFOX2-1 [n = 4] and shRBFOX2-2 [n = 3]) bone marrow (BM) cells. p, Relative gene expression of CD14 and CD15 in control (shNS [n = 4)) and RBFOX2 KD (shRBFOX2-1 [n = 4] and shRBFOX2-2 [n = 3]) BM cells. n, biologically independent mice. Data are presented as mean values +/− SEM. For c, n, o and p, two-sided P values by Student’s t-test. Source data

Extended Data Fig. 8. RBFOX2 regulation of m⁶A is not significantly related to its role in regulating alternative splicing.

a, Functional enrichment analysis of genes with their promoters bound by RBFOX2 and marked with both H3K4me3 and m⁶A. Inset boxplot shows gene expression level log₂FC comparing RBFOX2 KD with control K562 cells. Genes with RBFOX2 binding and H3K4me3 modification were categorized into two subgroups according to whether they were marked with m⁶A (+) or not (−). n = 2 biologically independent samples. RBFOX2 binding are from ENCODE project, and m⁶A MeRIP-seq and H3K4me3 ChIP-seq are from this study. b, Enrichment analysis of alternatively spliced genes upon RBFOX2 knockdown in K562 cells. c-d, Enrichment analysis of RBFOX2 bound and alternatively spliced genes upon RBFOX2 knockdown in K562 cells. Alternatively spliced genes bound by RBFOX2 were categorized into two groups according to whether those genes were methylated by m⁶A (RBFOX2(+) m⁶A(+)) or not (RBFOX2(+) m⁶A(−)). e, Enrichment analysis of genes harboring both RBFOX2 binding and m⁶A methylation, but which were not alternatively spliced upon RBFOX2 knockdown in K562 cells. f, Barplot showing the ratios of alternatively spliced genes by RBFOX2 knockdown in NB4 cells. g, Enrichment analysis of alternatively spliced genes upon RBFOX2 knockdown in NB4 cells. h-i, Enrichment analysis of RBFOX2 bound and alternatively spliced genes upon RBFOX2 knockdown in NB4 cells. RBFOX2 bound and alternatively spliced genes were categorized into two groups according to whether those genes were methylated by m⁶A (RBFOX2(+) m⁶A(+)) or not (RBFOX2(+) m⁶A(−)). j, Enrichment analysis of genes harboring both RBFOX2 binding and m⁶A methylation, but which were not alternatively spliced upon RBFOX2 knockdown in NB4 cells. For b-e and g-j, one-sided P values were calculated using Fisher’s Exact test, and the size of the circle represents the level of significance, with larger circles indicating greater significance and smaller circles indicating lower significance.

Extended Data Fig. 9. RBFOX2 functions as an m⁶A regulator in AML cells.

a, Average profile (top panel) and heatmap (bottom panel) showing RBFOX2 binding intensity at the m⁶A peak centers and the flanking 2.5 kb regions in NB4 cells. b, Average profile of RBFOX2 binding intensity at the RBFOX2 peak centers and the flanking 2.5 kb regions in NB4 cells. RBFOX2 peaks were categorized into two groups according to whether they overlap with m⁶A (+) or not (−). c, Distribution of m⁶A or RBFOX2 peaks at distinct genomic regions including promoter, exonic, intronic, transcription termination sites (TTS) and intergenic regions annotated by HOMER in NB4 cells. (+) sign represents regions harboring m⁶A (+) or bound by RBFOX2 (+), while the sign (−) indicates the absence of m⁶A (−) or RBFOX2 (−). d, Western blots of the immunoprecipitated RBFOX2 from NB4 cells and its interaction with METTL3, METTL14 and RBM15 after RNase A/T1 treatment, respectively. e, Western blots of the immunoprecipitated YTHDC1 from NB4 cells and its interaction with RBM15 after RNase A/T1 treatment. f, Western blots of the immunoprecipitated EZH2 or SUZ12 from NB4 cells, and their interactions with YTHDC1 after RNase A/T1 treatment. Source data

Extended Data Fig. 10. RBFOX2 regulates leukemia differentiation through TGF-β signaling pathway.

a, Correlation of gene expression between RBOFX2 and RBM15 in Leukemia cell line (CML) (left panel) from GTEx and in AML patients (right panel) from TCGA using GEPIA. P values by PCC. b, m⁶A level, caRNA abundance and H3K4me3 level at the TGFB1 promoter (K562, the first three panels with high-throughput sequencing data, n = 2; the forth panel with RT-qPCR, n = 4 technical replicates). c-d, Nascent RNA synthesis (c) and ELISA (d) of TGFB1 in K562 cells (n = 4, two technical replicates over two biological replicates). e, Gene expression level of TGFB1 and CD61 in shNS+EV, shRBFOX2 + EV and shRBFOX2 + RBFOX2 K562 cells. n = 3. f-g, Decreased m⁶A levels (f and g), increased RNA abundances at paRNA loci (n = 3, g) and increased mRNA expression (n = 2, g) of ARHGEF2, SPTBN1and RAB11B upon RBFOX2 KD in NB4 cells, respectively. h, The RBFOX2 binding and m⁶A level around the TGFB1 gene locus in NB4 cells. i, The m⁶A level and caRNA abundance at the TGFB1 promoter in NB4 cells (RT-qPCR, n = 4 technical replicates). j, SUZ12 or RBM15 binding at the TGFB1 promoter in NB4 cells (ChIP qPCR, n = 4 technical replicates). k, Gene expression level of TGFB1 in NB4 cells. n = 3. l, TGF-β ELISA assay in NB4 cells (shNS [n = 2], shRBFOX2-1 [n = 1] and shRBFOX2-2 [n = 2]). m, Relative gene expression of human TGFB1 in the bone marrow of leukemia mice from AML-PDX cells (2017-129) (shNS [n = 4], shRBFOX2-1 [n = 4] and shRBFOX2-2 [n = 3]). n, Flow cytometric analysis of CD11b⁺ cells in NB4 cells treated with TGF-β activator (0.6 ng/mL) or not (Untreated). o-p, The expression level of RBFOX2 (protein, o) and TGFB1 (mRNA, p) in control, RBFOX2 KD, and RBFOX2 & TGFB1 double KD (DKD) NB4 cells. n represents biologically independent samples, unless otherwise stated. For b, e, g, i, j, k, m and p, mean values +/− SEM; for c and d, mean values +/− SD. For g, m and p, two-sided P values by Student’s t-test; for e and k, P values by two-way ANOVA, Dunnett’s multiple comparisons test. Source data