Integrated exome and RNA sequencing of TFE3-translocation renal cell carcinoma

1TFE3-translocation renal cell carcinoma (TFE3-tRCC) is a rare subtype of kidney cancer, characterized by Xp11.2 translocations resulting in TFE3 fusion with various partner genes. To date, more than 20 partner genes have been identified in fusions with TFE3, including SFPQ, ASPSCR1, NONO, PRCC, RBM10, MED15, etc. Due to the variety of partners and fusion structures, the function of chimeric TFE3 fusion proteins is diverse, which may contribute to the high degree of heterogeneity of TFE3-tRCC, both morphologically and clinically–. Although only 1–4% RCCs in adults may have TFE3 translocations,, TFE3-tRCC constitutes 15% of RCCs in patients <45 years of age and 20–50% of pediatric RCCs,. Unfortunately, owing to the limited understanding of the underlying mechanisms, the optimal therapy for TFE3-tRCC remains to be determined, prompting an urgent clinical need to molecularly characterize TFE3-tRCC.

2Before this study, our knowledge concerning the molecular features of TFE3-tRCC was limited as previous analyses contained relatively small cohorts and were mainly confined to genomic alterations–. Although TFE3-tRCC was found to have few recurrent mutations, prior reports revealed the association between the somatic copy number alterations (SCNA) and survival outcomes–. However, these genomic alterations found in tumors were not enough to explain the high heterogeneity of TFE3-tRCC. The comprehensive genomic, especially the transcriptomic, characteristics of TFE3-tRCC are still uncovered.

3In this study, we apply whole-exome sequencing (WES) on 53 TFE3-tRCCs and RNA sequencing (RNA-seq) on 63 TFE3-tRCCs to reveal their genomic and transcriptomic characteristics and discover molecular mechanisms potentially involved in tumor progression. Our analyses reveal the prognostic value of ASPSCR1-TFE3 fusion and SCNAs for TFE3-tRCC, identify five molecular subsets with distinct transcriptional signatures.

4From 2009 to 2019, a total of 68 patients with TFE3-tRCC were identified from 4581 RCC cases (Supplementary Figure 1). The baseline clinicopathologic characteristics of the cases with TFE3-tRCC are summarized in Table 1 and Supplementary Data 1. The median age at diagnosis was 32.5 years (range: 5–70 years) and nine (13.2%) patients were younger than 18 years old. The male: female ratio was 2:3. The median tumor size was 4.7 cm (range: 1.4–19.6 cm). At initial diagnosis, 16 (23.5%) and seven (10.3%) patients presented with regional lymph node metastasis and distant metastasis, respectively. Morphologically, TFE3-tRCCs presented with diverse architectural and cytologic features, including papillary, tubular, acinar, and cystic patterns and 37 (54.4%) tumors had an international society of urological pathology (ISUP) grade ≥3. For primary kidney tumors, 27 (39.7%) and 41 (60.3%) patients underwent nephron-sparing surgery and radical nephrectomy, respectively. Ten (14.7%) patients died at the end of follow-up (median 43.8 months, 95% CI: 31.5–56.1).

5RNA-seq was performed on 63 TFE3-tRCC tumors and 14 adjacent normal kidney tissues. Gene fusions were detected in 90.5% (57/63) cases, of which 94.7% (54/57) cases showed relatively more common gene fusions, including SFPQ-TFE3 (n = 15), ASPSCR1-TFE3 (n = 13), NONO-TFE3 (n = 8), MED15-TFE3 (n = 8), PRCC-TFE3 (n = 6), and RBM10-TFE3 fusions (n = 4) (Fig. 1A). Three (5.3%, 3/57) cases were identified with rare TFE3-associated gene fusions, including FUBP1-TFE3, SETD1B-TFE3, and ZC3H4-TFE3, of which FUBP1-TFE3 and SETD1B-TFE3 fusion showed previously unreported fusion structures. FUBP1-TFE3 fusion resulted in a chimeric transcript composed of exons 1–15 of FUBP1 and exons 3–10 of TFE3, and SETD1B-TFE3 fusion was composed of exons 1–4 of SETD1B and exons 4–10 of TFE3. All three rare gene fusions were confirmed by RT-PCR and Sanger sequencing (Supplementary Figure 2). Patients with SETD1B-TFE3 and ZC3H4-TFE3 fusion developed metastasis by the end of follow-up.

6The structures of the TFE3 fusion isoforms identified in our cohort are summarized in Fig. 1B and C. All fusion genes preserved the open reading frame between the 5′ terminal of partner genes and the 3′ terminal of TFE3. According to the retained exons and functional domains of the TFE3, six types of isoforms were found, including retained fragment of TFE3 2–10 exons (5.3%, 3/57), 3–10 exons (3.5%, 2/57), 4–10 exons (15.8%, 9/57), 5–10 exons (22.8%, 13/57), 6–10 exons (50.9%, 29/57), and 7–10 exons (1.7%, 1/57). All fusions retained exons 7–10 of the TFE3 gene, containing the helix–loop-helix (bHLH) and leucine zipper (LZ) domains, but only a subset (47.4%, 27/57) of fusion isoforms contained the transcription activation (AD) domain.

7The fusion isoforms of different TFE3-tRCC subtypes were highly heterogeneous (Fig. 1C and Supplementary Data 2). Among which, SFPQ-TFE3 fusion isoforms had a maximum of six fusion structures, in contrast, ASPSCR1-TFE3 and RBM10-TFE3 isoforms had relatively fixed fusion structures. Next, we analyzed the functional domains of the retained exons of the fusion partner genes (Fig. 1C). We found 42% of fusion partner genes retained all functional domains. Interestingly, fusion partners that play regulatory roles in mRNA processing and/or mRNA splicing, including NONO, SFPQ, and RBM10, retained all RNA-recognition motifs (RRM).

8We subsequently explored the association of clinicopathologic features with different TFE3 fusion subtypes. Morphologically, we found tumors harboring ASPSCR1-TFE3 fusion showed a predominated papillary pattern, and tumors with MED15-TFE3 fusion showed a dominant multicystic pattern (Supplementary Figure 3A). Radiologically, cases with ASPSCR1-TFE3 fusion presented with features typical of hypervascularity and calcification on CT scan, while MED15-TFE3 fusion tRCC often presented as cystic masses (Supplementary Figure 3B).

9Of note, compared with other TFE3-tRCC subtypes, tumors with ASPSCR1-TFE3 fusion were highly aggressive, characterized by higher ISUP nuclear grade (ISUP ≥ 3, 11/13 vs. 19/44, P = 0.017) and more frequent lymph node metastasis (6/13 vs. 7/44, P = 0.057, Supplementary Table 1). At the end of follow-up, a total of 20 patients had presented with distal metastasis, of which over one-third were ASPSCR1-TFE3 fusion tRCCs (7/20, 35%). Survival analysis showed that those with ASPSCR1-TFE3 fusion were associated with poor overall survival (OS) compared with other subtypes (median OS: 58.5 months vs. not reached, P = 0.011, Fig. 2A and B). In multivariate analysis, after adjusting for clinicopathologic features, the ASPSCR1-TFE3 fusion was still independently associated with poor outcome (P = 0.010, Supplementary Table 2 and Supplementary Figure 4).

10WES was performed on 53 TFE3-tRCCs and matched germline samples. We identified a total of 1591 somatic mutations, including 1486 nonsynonymous single-nucleotide variants (SNVs) and 105 insertions/deletions (Indels), with a median of 14 (range: 0–238) mutations per tumor (Fig. 3A). The tumor mutation burden (median 0.28, range: 0–4.76) in our TFE3-tRCC cohort was lower than that of KIRC and KIRP in TCGA cohort (Supplementary Figure 5). We extracted three prominent mutational signatures using a non-negative matrix factorization (NMF) algorithm (Supplementary Figure 6). Signature B shows the largest contribution, which is found to be similar to SBS40, a signature correlated with age in multiple types of cancer. Signature A is similar to both the SBS87 (thiopurine exposure) and SBS1 (age-related 5-methylcytosine deamination). Signature C, which highly corresponds to SBS22 (cosine similarity = 0.91), is characterized by T > A transversions at CT [A/G] and has been associated with aristolochic acid exposure. We observed Signature C in 28.8% of patients in our cohort, indicating a potential role of aristolochic acid exposure in the development of Chinese TFE3-tRCC. The frequently mutated genes (frequency of more than four samples) included DST, DNAH8, and HMHA1, whereas the mutated loci at each gene were not recurrent (Fig. 3B and Supplementary Data 3). In addition, six tumor suppressor genes previously implicated in cancer (BTK, CHD1, FN1, NFATC2, NOTCH1, and NRP1) were found in at least two samples (Fig. 3B). Of these, 75% (9/12) were clonal mutations (Supplementary Data 3). In line with previous studies, the mutational spectrum of TFE3-tRCC was quite heterogeneous. Survival analysis indicated no relation between frequent SNVs and patient survival (Supplementary Table 3).

11The most frequently observed individual arm-level events included gain of 17q (12/53, 23%) and 19p (11/53, 21%), and loss of 19p (16/53, 30%), 14q (14/53, 26%), and 1p (11/53, 21%). The most frequent focal events were gain of 19p13.2 (17/53, 32%), 1q44 (15/53, 28%), and 8q24.3 (13/53, 25%), and loss of 19p12 (15/53, 26%), 14q21.2 (13/53, 25%, Fig. 3C and Supplementary Data 4). Previous studies reported that certain copy number events (eg. 9p loss and 17q gain) were correlated with patient outcomes, therefore, we evaluated the association of SCNA with clinicopathologic features and prognosis in our TFE3-tRCC cohort (Supplementary Table 4). We found that tumors with 22q loss were correlated with ASPSCR1-TFE3 fusion (4/4 vs. 2/38, P = 0.005), higher ISUP nuclear grade (ISUP ≥ 3, 7/19 vs. 0/27, P = 0.004) and more frequent lymph node metastasis (5/7 vs. 2/39, P = 0.004). Cases with 9p loss were associated with increased lymph node metastasis (4/8 vs. 2/39, P = 0.019). Survival analysis indicated that loss of chromosome arms 1p, 2p, 6q, 8p, 9p, and 22q were predictors for poor OS (Supplementary Figure 7 and Fig. 3D). Moreover, we identified that tumors with higher SCNA burden significantly correlated with worse survival outcomes (median OS: 59.46 months vs. 111.28 months, P = 0.006). After adjustment for clinicopathologic features, 22q loss was identified as an independent predictor for poor OS (P = 0.004, Supplementary Table 2 and Supplementary Figure 4).

12Analysis of differentially expressed genes (DEG) identified a total of 3124 over-expressed and 2143 under-expressed genes in TFE3-tRCCs compared to adjacent normal tissues (Fig. 4A and Supplementary Data 5). Among them, GPNMB, HIF1A, MET, and BIRC7, which were previously reported to have higher expression in tRCC,, were also upregulated in our TFE3-tRCC cohort (Supplementary Figure 8A). Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis identified that lysosomal, autophagy, and innate immunity pathways were significantly upregulated in TFE3-tRCC (Fig. 4B). Gene set enrichment analysis (GSEA) demonstrated enrichment of transcriptional pathways involved in proliferation (E2F targets and G2M checkpoint), PI3K/ATK/mTOR, and p53 signaling in TFE3-tRCC (Fig. 4C).

13Next, we explored the tumor immune microenvironment (TIME) in TFE3-tRCC. We found a lower immune infiltration score (IIS) and T-cell infiltration score (TIS) in TFE3-tRCC in contrast to KIRC from the TCGA cohort (Fig. 4D and Supplementary Figure 8B). Using a refined RCC immune cell gene-specific signatures, we found that compared with TCGA RCC subtypes, the T helper 2 cell (Th2) signature was increased in TFE3-tRCC, while the activated dendritic cell (aDC) and plasmacytoid dendritic cell (pDC) signatures had decreased expression (Fig. 4E and Supplementary Figure 8B). Moreover, the natural killer cell (NK) signature was increased in most TFE3-tRCC compared with KIRP and KICH. Furthermore, lower levels of CD8 + T-cell, T-cell, and macrophage signatures were identified in TFE3-tRCC relative to KIRC (Fig. 4E). Consistently, immunohistochemistry (IHC) confirmed that two-thirds (74.6%, 47/63) of tumors had low CD8 + T-cells infiltrations, while only twenty percent (20.6%, 13/63) of tumors showed massive CD8 + T cells infiltrations (Fig. 4F and Supplementary Figure 9). In tumors with a low CD8 + T-cells infiltrations, multiple immunofluorescence further revealed a low infiltration of NK (CD56 + ) cells and macrophages (CD68 + /HLA-DR + ) in the tumor stroma (Fig. 4G). In line with a previous analysis, the expression of PD-L1 mRNA was higher in TFE3-tRCC than KIRC and KIRP in our analysis (Supplementary Figure 8 C). However, different from the RNA-seq results, there was a relatively low PD-L1 positivity rate (17.5%, 11/63) in our TFE3-tRCC cohort as assessed via PD-L1 IHC staining. Interestingly, we observed higher expression of PD-L1 mRNA and protein levels in tumors with MED15-TFE3 fusion (Fig. 4H and Supplementary Figure 8D). Apart from MED15-TFE3 fusion tRCC, nearly two-thirds (63.5%, 40/63) of TFE3-tRCCs in our cohort exhibited type 1 TIME, characterized by a low CD8 + T-cells infiltration and low PD-L1 expression in tumors (Fig. 4I). Taken together, our data demonstrated an immunologically ignorant microenvironment in TFE3-tRCC.

14To move forward on understanding the biology of TFE3-tRCC, we performed unsupervised clustering analysis using NMF to identify and refine transcription-based subgroups of patients with common TFE3-tRCC subtypes (n = 54). Consequently, five distinct clusters were identified based on the 1500 most variable genes in our TFE3-tRCC cohort (Fig. 5A). To further reveal the transcriptional features driving these clusters, we performed quantitative set analyses for gene expression (QuSAGE) to compare each cluster with all others utilizing Hallmark gene sets (Fig. 5B and Supplementary Data 6). Combining this with the results of DEG analysis, which also compared clusters individually to others, we summarized and refined four signatures of representative genes involving stroma, angiogenesis, proliferation and KRAS down (Fig. 5C). Tumors in cluster 5 were characterized as highly angiogenic and proliferative, which showed high expression of angiogenesis-related genes (e.g., VEGFA, POSTN and HIF1A) and had enrichment of cell-cycle transcriptional programs (E2F targets, G2M checkpoint, and p53). Compared with tumors in the other clusters, tumors in cluster 5 had higher angiogenesis scores (Fig. 5D). IHC further confirmed that these tumors had the strongest nuclear HIF1A staining and increased vessel density (Fig. 5E). Moreover, tumors in cluster 5 also exhibited higher expression of stroma signature, exemplified by significantly higher epithelial-mesenchymal transition (EMT) score, and a high degree of collagens and activated stroma-related genes (e.g., FAP, CDH2, and CDH6; Fig. 5C and D). Cluster 1 was characterized by moderate enrichment of EMT, apical junction, TGF-β, WNT catenin, and hypoxia signaling, which have been associated with stroma and angiogenesis signatures (Fig. 5B). Tumors in cluster 3 and 4 were featured by low expression of angiogenesis modules (Fig. 5B–D). In contrast to cluster 3, cluster 4 showed moderate enrichment of E2F signaling. Cluster 2 differentiated from the other clusters by enrichment of genes repressed by KRAS down signature (Fig. 5B and C).

15Remarkably, we demonstrated that this molecular classification was associated with fusion subtypes and patient outcomes (Supplementary Figure 10). All tumors with ASPSCR1-TFE3 fusion were classified into the high angiogenesis/stroma/proliferation cluster (cluster 5), which exhibited worse survival. Cluster 1 included a half of MED15-TFE3 (62.5%, 5/8) and SFPQ-TFE3 (46.7%, 7/15) fusion tRCCs and cluster 3 was mainly comprised of NONO-TFE3 (75%, 6/8) and SFPQ-TFE3 (26.7%, 4/15) fusion tRCCs. Both the clusters showed moderate EMT score and intermediate survival (Fig. 5D). Although cluster 2 and cluster 4 consisted of mixed fusion subtypes, tumors in these clusters showed decreased EMT score, and demonstrated more favorable prognosis.

16A total of eight patients received systematic therapy in our TFE3-tRCC cohort. Of note, the median progression-free survival (PFS) for patients receiving first-line tyrosine kinase inhibitors (TKIs) treatment was only 3 months (Fig. 6A), confirming that the low angiogenesis signature interfered with the efficacy of TKIs monotherapy for patients with TFE3-tRCC. Interestingly, we observed two patients in cluster 5 had sustained control of disease to TKIs plus immune checkpoint inhibitors, with one (TFE3-68) showing partial response (tumor shrinkage 63.1%) after 5 months of first-line axitinib plus pembrolizumab and the other patient (TFE3-65) demonstrating stable disease for 12 months after third-line axitinib plus pembrolizumab (Fig. 6B). These initial clinical results indicated that tumors in cluster 5 might have a higher likelihood of response to antiangiogenic plus immunotherapy.

17This study depicted comprehensive molecular analyses of 63 untreated primary TFE3-tRCCs and investigated the molecular characteristics of this rare but highly heterogeneous RCC entity. Our findings provided important new insights into key genomic and biologic features underlying TFE3-tRCC development and progression, confirmed the promotion role of copy number alterations in tumor progression, revealed the prognostic value of TFE3 fusion patterns, and described five distinct molecular clusters associated with different transcriptional signatures.

18We identified a low degree of somatic mutation rate but a widespread existence of SCNAs in TFE3-tRCC. Copy number alterations in TFE3-tRCC have been previously explored but with inconsistent results. In a previous study of 16 tRCCs, 17q gain (44%) was identified as the most common SCNA and found to be a predictor for poor prognosis. More recently, a genomic analysis of 22 tRCCs reported loss of 9p (41%) and gain of 17q (36%) were the most frequent SCNAs, but no single SCNA event was associated with worse prognosis. Such inconsistent results may result from the relatively small sample sizes and heterogeneity of TFE3-tRCC, since both studies only included a few TFE3-tRCC subtypes and were confounded with TFEB-tRCC samples. In our cohort, gain of 17q and loss of 9p were found in 23 and 11% of TFE3-tRCC, respectively. Loss of chromosome arm 9p but not gain of arm 17q was correlated with poor survival. In addition, loss of 1p, 2p, 6q, 8p, 22q, and increased SCNA burden were also predictors for poor prognosis. More importantly, we demonstrated that 22q loss was significantly associated with aggressive clinical features and an independent predictor of worse outcomes for patients with TFE3-tRCC. Increased loss of chromosome 22q that encodes NF2, CHEK2, and SMARCB1 were observed in type 2 PRCC, which may implicate in carcinogenesis and tumor progression. A recent study suggested that the presence of copy number variations might be an earlier event in TFE3-tRCC than somatic mutations. Therefore, our results again underscore the importance of SCNA in the aggressive characteristics of TFE3-tRCC and the prognostic value of copy number aberrations.

19Our results demonstrated an immune ignorant TIME in the majority of TFE3-tRCCs, characterized by low PD-L1 expression and low CD8 + T-cell infiltration in tumor stroma. Except for tumors with MED15-TFE3 fusion, a low PD-L1 positivity rate was detected in our TFE3-tRCCs. Therefore, it is possible that the fusion subtype may impact the different PD-L1 expression in TFE3-tRCC. Given the low immunogenicity of TFE3-tRCC (mostly belonging to TIME type 1), single-agent PD-1/PD-L1 blockade could probably be inadequate to improve outcomes in those with advanced disease, which is supported by several clinical investigations–. Therefore, combination therapeutic regimens, including radiotherapy, targeted therapy, or chemotherapy plus immunotherapy, could be more effective to mediate immunogenic cell death though inducing neoantigen liberation and lymphocytes recruitment into the tumor microenvironment.

20TFE3-tRCC represents a heterogeneous disease with a wide spectrum of morphologies and varied clinical phenotypes. Several studies suggested that patients with different fusion subtypes exhibited differential prognosis. In this study, we explored the association between fusion subtypes and patient outcomes. We observed that tumors with ASPSCR1-TFE3 fusion had much worse survival compared to those with other fusions. On the contrary, patients with other common TFE3 fusions, including MED15-TFE3 and PRCC-TFE3 showed a better prognosis. Our findings highlighted the importance of the detection of TFE3 fusion partners for predicting prognosis. Also, the detailed mechanisms of TFE3 fusion partners in regulating the chimeric protein are intriguing and warrant further investigation.

21Our unsupervised transcriptomic analysis identified five molecular clusters. Tumors with ASPSCR1-TFE3 fusion were exclusively classified into cluster 5, exhibiting high expression of angiogenesis/stroma/proliferation gene signatures. The distinct genomic features of ASPSCR1-TFE3 fusion tRCC identified in the current study suggested a molecular basis for the aggressive nature of this TFE3-tRCC subtype. A preclinical study suggested that ASPSCR1-TFE3 fusion induced accumulation of lactate in tumor microenvironment, leading to tumor hypoxia and angiogenesis in alveolar soft part sarcoma. In agreement with this finding, we found that tumors with ASPSCR1-TFE3 fusion exhibited strong nuclear HIF1A staining and increased vessel density. Thus, patients with ASPSCR1-TFE3 fusion tRCC may benefit from antiangiogenic based treatment. It is important to note that the high stroma and proliferation signatures may also interfere with the efficacy of antiangiogenic monotherapy. Interestingly, we observed a favorable response to anti-PD1 plus TKI combination therapy in two ASPSCR1-TFE3 fusion tRCC patients in this study. Therefore, our data may support clinical investigation of antiangiogenic therapy in combination with immune checkpoint inhibitors in this TFE3-tRCC subtype.

22In contrast to ASPSCR1-TFE3 fusion tRCC, other relatively common TFE3-tRCC fusion subtypes were characterized by decreased angiogenesis gene signatures. These findings provide a molecular explanation for the unfavorable clinical outcome to antiangiogenic monotherapy in ours and previous clinical studies,. Tumors in cluster 1 and 3 were associated with moderate stroma signatures. Tumors in cluster 2 and 4 showed an enrichment of KRAS down signature and E2F pathway, respectively. In addition, both ours and previous studies showed the activation of autophagy, mTOR and proliferation signaling in the whole TFE3-tRCC–. Therefore, targeting these specific aberrations, such as stromal disruptors, E2F, autophagy, mTOR, and proliferation inhibitors may be options for patients with advanced TFE3-tRCCs.

23The current study is not devoid of limitations. Due to the low incidence of TFE3-tRCC and its high heterogeneity, it is difficult to collect an adequate number of cases with different subtypes for further analysis. However, our study represents the largest cohort to date investigating the molecular features of TFE3-tRCC. Additionally, this was a single-center study thus selection bias is unavoidable. Although with a relatively long follow-up time, we acknowledge that the low proportion of death cases may restrict the statistical power of survival analysis in our study. Future multicenter studies with larger sizes are needed to shed further light on the molecular mechanisms of TFE3-tRCC.

24Overall, the current study expanded our knowledge of the genomic and transcriptional landscape of TFE3-tRCC, emphasized the importance of fusion partners and structures for patient outcomes, and provided a molecular basis for the heterogeneous clinical phenotype. Our study represented a step forward in understanding the TFE3-tRCC biology. We expect that our findings will provide a genetic basis for developing personalized therapies for this disease.

25A total of 4581 cases diagnosed as RCC who underwent surgery for the kidney at our center between 2009 and 2019 were reviewed by two experienced uropathologists (Ni Chen and Mengni Zhang). Among them, 1006 suspicious non-clear cell RCC cases were identified via morphological evaluation and were selected for further TFE3 IHC. As a result, 68 TFE3 positive cases were confirmed as TFE3-tRCC by break-apart fluorescence in situ hybridization (FISH) assay (Supplementary Figure 1). Untreated primary tumor tissues and adjacent normal samples were collected. RNA-seq and WES of formalin-fixed paraffin-embedded (FFPE) tissues were subsequently performed. RNA-seq was performed on 63 tumors and adjacent normal tissues. WES was performed on 53 FFPE tumor tissues and matched adjacent normal (n = 42)/blood (n = 11) samples. The study was conducted in accordance with the Declaration of Helsinki and was approved by the Ethics Committee of West China Hospital of Sichuan University. All patients or family members provided written consent for genetic analysis.

26Clinicopathological data were retrospectively collected, including age, gender, tumor size, TNM stage, morphological features, ISUP grade, and systemic treatment strategies. For patients with localized disease, regular evaluations were performed every 6–12 months postsurgery. For those with metastatic disease, regular evaluations were carried out every 4–6 weeks after receiving systematic therapy. For patients receiving systemic treatments, PFS was defined as the time from treatment initiation to disease progression or death. Tumor response was defined by Response Evaluation Criteria in Solid Tumors (RECIST) version 1.1. OS was defined as the time from surgery to the date of death.

27IHC and multiple immunofluorescence were performed as previously described. Commercially available primary anti-TFE3 (clone MRQ-37, 1:100, MXB biotechnologies, Fujian, China), anti-HIF1A (1:5000, Novus Biologicals, Colorado, USA), anti-CD31 (clone UMAB30, ready to use, ZSGB-BIO, Beijing, ChinaMXB), CD8 (clone C8/144B, ready to use, Dako, Copenhagen, DEN) and PD-L1 (clone 22C3, 1:50, Dako) were used in this study. Multiplex immunofluorescence staining was performed with a PANO 7-plex kit (0004100100, Panovue, Beijing, CHN). CD8, macrophage, and NK cells expression were quantified as positive cell density (cell number per mm2). PD-L1 expression was assessed by tumor proportion score, which was defined as the percentage of tumor cells with membranous PD-L1 staining. PD-L1 expression >1% was defined as positivity.

28FFPE tissue sections were examined by using interphase FISH to investigate the rearrangement of the TFE3 region with an LSI dual-color break-apart probe (GSP TFE3, Anbiping company, Guangzhou, China). TFE3 gene rearrangement would result in break-apart of the normal fused green-red signals, resulting in one green/one red break-apart signal pattern (male), or one green/one red break-apart signal and one remaining normal fused green-red signal pattern (female) in a normal cell. One hundred nonoverlapping tumor nuclei were evaluated. The cut-off value was 10%.

29cDNA was synthesized using PrimeScriptTM RT Master Mix (RR036A, Takara, Shiga, Japan) according to the manufacture’s instruction and then subjected to PCR reactions. Primer sequences for fusion validations are listed in Supplementary Table 5. PCR was performed with the following thermal cycling conditions: a 98 °C for 10 s; 35 cycles of 95 °C for 10 s, 52 °C for 30 s and 72 °C for 1 min; and a 72 °C for 5 min. PCR products were separated by electrophoresis in agarose gels, purified with ChargeSwitch™ PCR Clean-Up Kit (CS12000, Invitrogen, Oberhausen, GER) and then sequenced by ABI 3730XL automatic sequencer (Life Technologies).

30Total RNA was isolated from each sample (63 tumor samples and 14 paired adjacent normal samples) using the Qiagen RNeasy formalin-fixed paraffin-embedded (FFPE) Kit (73504, Qiagen, Hilden, Germany), following the protocol from the manufacturer. The purity and quantity of total RNA were measured by Nanodrop. The integrity of RNA was evaluated using the RNA Nano6000 Assay Kit on the Bioanalyzer 2100 system (Agilent Technologies, CA, USA). 1 μg RNA of per sample was used as input for the RNA sample preparations. Strand-specific RNA sequencing libraries were generated using the Whole RNA-seq Lib Prep kit for Illumina (RK20303, ABclonal, Shanghai, China). Library quality was evaluated on the Agilent Bioanalyzer 2100 system (Agilent, USA). Final libraries were sequenced at the Novogene Bioinformatics Institute (Beijing, China) on an Illumina Hiseq X10 platform by 150 bp paired-end reads.

31The raw RNA-sequencing reads were filtered by FastQC, Reads were aligned using STAR (v2.7.0 f) with default parameters to the Ensembl human genome assembly GRCh37. Gene expression levels were estimated by raw count and Transcripts Per Kilobase Million (TPM). Annotations of mRNA in the human genome were retrieved from the GENCODE (v19) database. Paired trimmed/clipped and de-duplicated RNAseq reads were used to identify gene fusion events, and the aligned output was used as input to STAR-Fusion (v1.9.1) using the developer-supplied gencode v33 CTAT library from April 6, 2020. Fusion gene supported by at least two reads were selected.

32DEGs were determined using the R package “limma” with cut-off P value < 0.05. Upregulated genes and downregulated genes were used to perform ontology and pathway enrichment analysis based on Gene Ontology and KEGG databases using R package “ClusterProfiler”.

33DEG analysis results were used in GSEA analyses. For gene set analysis, hallmark gene sets from Msigdb database were collected. The top 1500 DEGs were chosen by P value ranked by log2 fold change, and then fed into “fgsea” R packages.

34We selected 1500 genes with the highest variability across tumors, using Median Absolute Deviation (MAD) analysis. Subclasses were then computed after reducing the dimensionality of the expression data from thousands of genes to a few metagenes using consensus NMF clustering. This method computes multiple k-factor factorization decompositions of the expression matrix, and evaluates the stability of the clustering result using the cophenetic coefficient index. The robust consensus NMF clustering of 54 patient tumors using the top 1500 variable genes, and k = 5 was identified through testing k = 2 to k = 10.

35Marker genes for renal cell carcinoma immune cell types were obtained from Wang et al. The T-cell infiltration score (TIS) was defined as the mean of the standardized values for nine T-cell subtypes, and the immune infiltration score (IIS) for a sample was similarly defined as the mean of the standardized values for innate and adaptive immune scores. A previously published angiogenesis signature was utilized to measure an angiogenesis score. And epithelial-mesenchymal transition (EMT) score was calculated following the method mentioned in Gibbons et al..

36ssGSEA was used for quantifying immune infiltration and activity in tumors using markers reported by Wang et al. The ssGSEA method is an extension of the GSEA method, which works at a single-sample level rather than a sample population. Normalized RNA-Seq data were used as input without further processing (i.e., no standardization or log transformation).

37To reveal latent biological pathways of NMF clustering, we conducted QuSAGE (R qusage v2.20.0) to compare each cluster to all others, leveraging MSigDb hallmark gene sets to identify enriched pathways within each cluster.

38Representative sections of formalin-fixed paraffin-embedded (FFPE) tumor (n = 53) and matched normal tissues (n = 42) or blood samples (n = 11) were collected. Tumor sections were reviewed by a pathologist to ensure tumor sections with >70% tumors and <10% necrosis. High-quality genomic DNA was extracted by using the GeneRead DNA FFPE Kit (180134, QIAGEN, Hilden, GER) according to the manufacturer’s instructions. Germline DNA (gDNA) was extracted from white blood cells using the Blood Genomic DNA Mini Kit (CW2087, Cwbiotech, Beijing, China).

39WES was used to determine the mutational landscape of TFE3-tRCC. Exome capture was performed using the Agilent SureSelect Human All ExonV5 kit (Technologies, CA, USA) according to the manufacturer’s instructions. This was followed by paired-end sequencing using Illumina Novaseq6000 sequencer (Illumine Inc, CA, USA). Mean clean sequence data obtained was 38.30 Gb for tumor samples and 14.95 Gb for normal tissue/blood samples. Details of somatic and germline mutation analysis are described below.

40Clean reads were aligned to the reference human genome hg19 (Genome Reference Consortium GRCh37) using the BWA (v0.7.17) (Burrows-Wheeler Aligner) MEM algorithm with default parameters. BAM was coordinate sorted and PCR duplicates were removed with Sambamba (v0.6.8).

41After the initial alignment of WES data, we followed GATK (v3.8) Best Practice to process all BAMs from the same patient together for a postalignment optimization process called “co-cleaning” which includes GATK IndelRealigner and BaseQualityScoreRecalibration (BQSR). IndelRealigner performed local realignment to further improve mapping quality across all reads at loci close to indels, and BQSR detected and fixed systematic errors made by the sequencer when it estimated the quality score of each base call–.

42The GATK MuTect2 pipeline was run for paired tumor-normal somatic mutation calling with gnomAD database and a panel of normal made from all normal samples to filter common germline mutations and recurrent technical artifacts. The resulting VCFs were filtered by Mutect2 FilterMutectCalls module, variants outside of the capture kit were removed, and FilterByOrientationBias module was used to filter out false-positive calls from OxoG and FFPE. The resulting somatic SNVs and indels were further filtered according to the flowing criteria–: read depth ≥10 in both tumor and normal samples, mapping quality ≥40 and base quality ≥20, variants allele frequency (VAF) ≥ 5% and supporting reads ≥5 in tumor, VAF in tumor was ≥5 times than that of the matched normal VAF. Variants were annotated with Oncotator v1.9.9.0. To further avoid miscalling germline variants at least 19 read depth in the normal sample in dbSNP sites. Variants were excluded with minor allele frequency (MAF) > 0.01 in public databases, including 1000 G (May 2013), Exome Sequencing Project (ESP6500) and Exome Aggregation Consortium (ExAC v0.3.1). Variants were filtered which annotated as 3’UTR, 5’UTR, 3’Flank, 5’Flank, IGR, Intron, lincRNA, Silent RNA.

43The R package MutationalPatterns (v3.0.1) was used to extract the somatic motifs of these samples. In brief, the somatic motifs for each variant were retrieved from the reference sequence and converted into a matrix. NMF was used to estimate the optimal number of mutation signatures extracted from WES samples. Cosine similarity was calculated to measure the similarity between our identified signatures and COSMIC signatures v3.2 [cancer.sanger.ac.uk/cosmic/signatures].

44FACETS (v0.5.14) was used to estimate tumor cellularity and ploidy from paired tumor and normal WES data, and calculated allele-specific somatic copy number alterations. Copy Number (CN) gains were defined as alterations showing total CN > 2 and CN losses were defined as alterations showing total CN < 2. Arm-level events were defined as any gain or loss occurring in an autosome that involved at least 10% of the arm. To identify significantly focal SCNAs, we used GISTIC2 (v2.0.23), which considers both the frequency and amplitude of every SCNA, was employed with the following modified parameters “-smallmem 1 -broad 1 -brlen 0.7 -cap 1.5 -conf 0.99 -ta 0.2 -td 0.25—armpeel 1 -genegistic 1 -savegene 1 -gcm extreme -js 4 -maxseg 2000 -qvt 0.25 -rx 0”. CN gains were defined as alterations showing 0.1 < log2 CN ratio < 0.7 and CN losses were defined as alterations showing −1.1 < log2 CN ratio < −0.1. To measure SCNA burden, fraction of copy number-altered genome (FCA) was calculated by dividing the number of bases in segments with mean log2 CN ratio >0.1 or < −0.1 by the number of bases in all segments. SCNA burden high was defined as ≥ 75th of percentile of SCNA burden in the relevant cohort. SCNA burden low was defined as < 75th of percentile of SCNA burden.

45Tumor suppressor genes (TSGs) were obtained from TSGene version 2.0 (https://bioinfo.uth.edu/TSGene/) and IntOGen (https://www.intogen.org) database. Genes mutated in at least two samples were shown in barplot and TSGs were marked with bold font. KEGG enrichment analysis was done with those genes. Clonality of mutations was determined based on the cancer cell fraction (CCF) estimated by allele-specific copy number analysis.

46All analyses were conducted using R software (v3.6.0) and SPSS (v16.0). All comparisons for continuous variables were performed using the two-sided Mann–Whitney U test for two groups. For categorical variables, Pearson’s Chi-squared test with continuity correction or Fisher’s exact test was used. Survival analyses were conducted using Kaplan–Meier method and the difference was tested using log-rank. Cox proportional hazards regression (forward likelihood ratio model) was used to determine the independent predictor of OS. All clinicopathological parameters and biomarkers at P < 0.05 were then further tested on multivariate Cox regression in three patient cohorts (WES + RNAseq cohort, WES cohort and RNAseq cohort). Least absolute shrinkage and selection operator Cox regression were also performed using all variables in the multivariate analyses to identify optimal predictors of OS. A P value less than 0.05 was considered statistically significant.

47Further information on research design is available in the Nature Research Reporting Summary linked to this article.