Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses

Article

Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses Graphical Abstract

Authors Elizabeth Stewart, Justina McEvoy, Hong Wang, ..., James R. Downing, Michael A. Dyer, St. Jude Children’s Research Hospital – Washington University Pediatric Cancer Genome Project

Correspondence [email protected]

In Brief Using integrated analyses, Stewart et al. reveal the developmental origin of rhabdomyosarcoma (RMS) and identify commonly deregulated pathways. Comprehensive preclinical testing helps prioritize agents for future clinical trials and shows WEE1 as a highly effective therapeutic target for high-risk RMS in vivo.

Highlights d

Integrated analyses provides insight into the developmental origin of RMS

d

ARMS occurs further along the developmental program than ERMS

d

RAS/MEK/ERK/CDK4/6, G2/M, and UPR pathways are deregulated in RMS

d

Targeting WEE1 kinase in the G2/M pathway was most effective in high-risk RMS

Stewart et al., 2018, Cancer Cell 34, 1–16 September 10, 2018 ª 2018 Elsevier Inc. https://doi.org/10.1016/j.ccell.2018.07.012

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Cancer Cell

Article Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses Elizabeth Stewart,1,2,17 Justina McEvoy,3,17 Hong Wang,4,5,17 Xiang Chen,6,17 Victoria Honnell,1,2 Monica Ocarz,1 Brittney Gordon,1 Jason Dapper,1,2 Kaley Blankenship,2 Yanling Yang,4 Yuxin Li, Formal analysis4,7 Timothy I. Shaw,6,7 Ji-Hoon Cho,7 Xusheng Wang,7 Beisi Xu,6 Pankaj Gupta,6 Yiping Fan,6 Yu Liu,6 Michael Rusch, Formal analysis6 Lyra Griffiths,1 Jongrye Jeon,1 Burgess B. Freeman III,8 Michael R. Clay,9 Alberto Pappo,2 John Easton,6 Sheila Shurtleff,9 Anang Shelat,10 Xin Zhou,6 Kristy Boggs,6 Heather Mulder,6 Donald Yergeau,6 Armita Bahrami,9 Elaine R. Mardis,11,12,13 Richard K. Wilson,11,12,14 Jinghui Zhang,6 Junmin Peng,4 James R. Downing,9 Michael A. Dyer1,15,16,18,* and for the St. Jude Children’s Research Hospital – Washington University Pediatric Cancer Genome Project 1Department of Developmental Neurobiology, St. Jude Children’s Research Hospital, 262 Danny Thomas Place, MS 323, Memphis, TN 38105-3678, USA 2Department of Oncology, St. Jude Children’s Research Hospital, Memphis, TN 38105, USA 3Departments of Molecular and Cellular Biology and Pediatrics, BIO5 Institute, University of Arizona, Tucson, AZ 85721, USA 4Department of Structural Biology, St. Jude Children’s Research Hospital, Memphis, TN 38105, USA 5Integrated Biomedical Sciences, University of Tennessee Health Science Center, Memphis, TN 38105, USA 6Department of Computational Biology, St. Jude Children’s Research Hospital, Memphis, TN 38105, USA 7Proteomics Shared Resource, St. Jude Children’s Research Hospital, Memphis, TN 38105, USA 8Preclinical Pharmacokinetics Shared Resource, St. Jude Children’s Research Hospital, Memphis, TN 38105, USA 9Department of Pathology, St. Jude Children’s Research Hospital, Memphis, TN 38105, USA 10Department of Chemical Biology and Therapeutics, St. Jude Children’s Research Hospital, Memphis, TN 38105, USA 11The McDonnell Genome Institute, Washington University, St. Louis, MO 63108, USA 12Department of Genetics, Washington University, St. Louis, MO 63108, USA 13Department of Medicine, Washington University, St. Louis, MO 63108, USA 14Siteman Cancer Center, Washington University School of Medicine in St. Louis, St. Louis, MO 63108, USA 15Department of Ophthalmology, University of Tennessee Health Science Center, Memphis, TN 38105, USA 16Howard Hughes Medical Institute, Chevy Chase, MD 20815, USA 17These authors contributed equally 18Lead Contact *Correspondence: [email protected] https://doi.org/10.1016/j.ccell.2018.07.012

SUMMARY

Personalized cancer therapy targeting somatic mutations in patient tumors is increasingly being incorporated into practice. Other therapeutic vulnerabilities resulting from changes in gene expression due to tumor specific epigenetic perturbations are progressively being recognized. These genomic and epigenomic changes are ultimately manifest in the tumor proteome and phosphoproteome. We integrated transcriptomic, epigenomic, and proteomic/phosphoproteomic data to elucidate the cellular origins and therapeutic vulnerabilities of rhabdomyosarcoma (RMS). We discovered that alveolar RMS occurs further along the developmental program than embryonal RMS. We also identified deregulation of the RAS/MEK/ERK/ CDK4/6, G2/M, and unfolded protein response pathways through our integrated analysis. Comprehensive preclinical testing revealed that targeting the WEE1 kinase in the G2/M pathway is the most effective approach in vivo for high-risk RMS.

Significance This study demonstrates the value of integrating transcriptomic, epigenomic, and proteomic/phosphoproteomic data to further refine the cellular origins of developmental tumors and to identify tumor vulnerabilities that extend beyond somatic mutations in the genome. Our results demonstrate the importance of performing comprehensive preclinical testing using multiple orthotopic patient-derived xenografts to validate and prioritize vulnerabilities identified through genomic or integrated analyses. Based on these data, patients with high-risk and recurrent rhabdomyosarcoma may benefit from combination chemotherapy that includes AZD1775, irinotecan, and vincristine. To facilitate the dissemination of these data and future research, we have developed an online portal (https://pecan.stjude.org/proteinpaint/study/RHB2018) for public use. Cancer Cell 34, 1–16, September 10, 2018 ª 2018 Elsevier Inc. 1

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

INTRODUCTION Rhabdomyosarcoma (RMS) is a childhood solid tumor with features of skeletal muscle. Approximately 75% of patients with localized disease are cured with conventional multimodal therapy, but patients with recurrent or metastatic disease have overall survival rates of 17% and 30%, respectively (Pappo et al., 1999). Historically, RMS was divided into embryonal RMS (ERMS), comprising approximately 60% of all RMS cases, and alveolar RMS (ARMS), comprising 25% of cases (Newton et al., 1988). Most histopathologically classified ARMS contain a translocation between FOXO1 on chromosome 13q14 and either PAX3 on chromosome 2q35 or PAX7 on chromosome 1p36 (Raney et al., 2001). However, a subset of ARMS is translocation negative and has molecular features more similar to those of ERMS (Williamson et al., 2010). Patients with ERMS have a relatively good prognosis (Ognjanovic et al., 2009), and current clinical trials use PAX3/7-FOXO1 translocation status for risk assignment rather than histologic subgroup. Previous genomic characterization revealed that translocation-negative ERMS has a high rate of structural, copy-number, and singlenucleotide variations whereas translocation-positive ARMS has relatively few recurrent mutations in cancer consensus genes (Chen et al., 2013; Shern et al., 2014). To improve understanding of recurrent RMS and provide additional preclinical models, we have established a collection of orthotopic patient-derived xenografts (O-PDX) of pediatric solid tumors (Stewart et al., 2017). All of our O-PDXs and their associated data are freely shared through the childhood solid tumor network (CSTN) with no obligation to collaborate (www.stjude. org/CSTN/). Because the overall survival rate for patients with RMS has not appreciably improved in the past 15 years, such a streamlined approach is particularly important in elucidating efficacious RMS therapeutics (Smith et al., 2014). RESULTS Epigenetic Profiling of Rhabdomyosarcoma Previously, we generated O-PDXs for 17 RMS and characterized them along with their corresponding patient tumors (Stewart et al., 2017), 11 of which have potentially actionable mutations based on the current adult and pediatric NCI-MATCH clinical trials (Table S1). Here, we extended our analysis to include wholegenome bisulfite sequencing (WGBS) of these samples. We also performed WGBS and RNA sequencing (RNA-seq) of primary human myoblasts, myotubes, and fetal skeletal muscle (week 21 of gestation). Most patient tumors and matched O-PDXs were correlated for gene expression (r = 0.926; range, 0.76– 0.99) and DNA methylation (r = 0.874; range, 0.57–0.97), and the differences among tumor types were associated with infiltrating nontumor cells in patient tumors (Figures 1A and 1B; Table S2). We next determined the difference in correlation between DNA methylation and gene expression for ERMS and ARMS. Of the 1,782 genes that were upregulated in ERMS relative to ARMS, 427 (24%) were negatively correlated with DNA methylation, whereas 352 (20%) were positively correlated, and 279 (16%) had regions of both positive and negative correlation (Figure 1C and Table S2). Of the 1,576 genes that were downregulated in ERMS relative to ARMS, 134 (9%) were nega2 Cancer Cell 34, 1–16, September 10, 2018

tively correlated with DNA methylation, 386 (24%) were positively correlated, and 368 (23%) were both positively and negatively correlated (Figure 1C and Table S2). Many genes with negatively correlated DNA methylation and mRNA expression were identified in regions immediately downstream of promoters, known as downstream correlated regions (Hovestadt et al., 2014) (Figures 1D and 1E). In contrast, many genes with positively correlated DNA methylation and mRNA expression were found along the entire locus (Figures 1F and 1G) and contained partially methylated domains (Lister et al., 2009). Pathway analysis of the differentially expressed and methylated genes revealed that negatively correlated genes were primarily mediators of embryonic skeletal muscle development, whereas positively associated genes were predominantly mediators of nervous system development (Table S2). To determine whether epigenetic changes other than DNA methylation accompany changes in gene expression during RMS tumorigenesis, we performed chromatin immunoprecipitation and sequencing (ChIP-seq) of eight histone modifications (i.e., H3K4me1, H3K4me2, H3K4me3, H3K27me3, H3K27Ac, H3K36me3, H3K9/14Ac, and H3K9me3) and of CTCF, BRD4, and RNAPII in all O-PDXs, myoblasts, myotubes, and fetal skeletal muscle. In total, we used 19 antibodies and sequenced/ analyzed 756 ChIP-seq libraries for which data was used to perform chromatin hidden Markov modeling (chromHMM) (Ernst and Kellis, 2012). Through a combination of computational methods and manual review, we selected 18 chromHMM states for this analysis (Figure 2A). States 1 through 7 had epigenetic marks consistent with enhancers and active promoters; state 8 contained epigenetic marks of bivalent (H3K27me3/H3K4Me3) poised promoters/enhancers; and states 9 through 12 included epigenetic marks of gene bodies. State 13 was consistent with the large class of zinc-finger genes in the genome with H3K27Ac, H3K9me3, and H3K27me3 modifications (O’Geen et al., 2007). State 14 comprised polycomb repressed (H3K27me3) genomic regions, and states 15 and 16 were quiet/empty states in our analysis. State 17 consisted of heterochromatin, and state 18 represented repressed chromatin (H3K27me3 and H3K9me3). The distribution of these 18 chromHMM states across the genome was very similar for both ERMS and ARMS (Figure 2B). Differences in chromHMM states were strongly associated with the differences in gene expression between these tumors (Figure 2C), suggesting that the protein and histone modifications analyzed by ChIP-seq and chromHMM provided a more complete view of the epigenomic landscape than DNA methylation alone. We focused our epigenomic analysis on three broad categories in combination with DNA methylation: (1) genes with differential expression between ERMS and ARMS that had corresponding changes in chromHMM states (Figure 2D and Table S3), (2) genes with differential expression between the tumor types that were also associated with a superenhancer within 500 kb of the transcriptional start site (Figure 2E and Table S3), and (3) genes with differences in H3K27me3 (chromHMM states 8, 13, 14, 18) that were negatively associated with gene expression (Figure 2F and Table S3). The differences in chromHMM states accounted for 749 of 1,782 (42%) upregulated genes and 522 of 1,576 (33%) downregulated genes in ERMS relative to that of ARMS (Figure 2G). Only 36 (2%) of the upregulated

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Figure 1. Differences in DNA Methylation Relative to Gene Expression (A and B) Representative scatterplot of RNA-seq data (log2 FPKM, A) and the b value (B) from WGBS analysis for an O-PDX sample and corresponding RMS tumor. (C) Venn diagrams depicting overlap of upregulated and downregulated genes (green) in ERMS relative to ARMS tumors that are hypomethylated (blue) or hypermethylated (red). Gene number in each group is noted. (D–G) Examples of a gene that is hypermethylated downstream of the promoter in ARMS relative to ERMS (boxed region) and has reduced gene expression (MFAP2, D), that is hypomethylated in the region flanking the transcriptional start site in ARMS relative to ERMS (boxed region) with increased gene expression (NHLH1, E), that is hypermethylated in ARMS relative to ERMS with increased gene expression in ERMS (MMP16, F), and that is hypermethylated in ARMS relative to ERMS and has increased gene expression in ARMS (NOS1, G). Transcript boundaries are indicated by dashed lines, and b values for each CpG are plotted for ERMS (blue) and ARMS (red). Expression in each sample from RNA-seq is shown in the histogram on the right (FPKM). FPKM, fragments per kilobase per million reads. See also Tables S1 and S2.

genes and 90 (6%) of the downregulated genes contained the H3K27me3 modification (Table S3). Importantly, 98 genes that were upregulated in ERMS contained an ERMS-specific superenhancer, and 174 genes that were upregulated in ARMS contained an ARMS-specific superenhancer (Figure 2G and Table S3), many of which are implicated in myogenesis. To integrate ChIP-seq data with transcriptome and DNA methylation data, we added the ERMS- and ARMS-specific H3K27me3, chromHMM, and superenhancer information to Table S2. Among the 427 genes that were upregulated in ERMS relative to ARMS and hypomethylated, four had a relative decrease in H3K27me3, 47 had an ERMS-specific superenhancer, and 174 had differences in chromHMM that correlated with expression (Table S2). For example, KCNK5 was upregu-

lated in ERMS relative to ARMS with changes in chromHMM consistent with differential gene expression, was associated with an ERMS-specific superenhancer, and was hypomethylated with reduced H3K27me3 (Table S2). KCNK5 encodes a potassium channel implicated in muscle differentiation and myocyte fusion (Afzali et al., 2016) but has not been previously characterized in RMS. ERMS tumors epigenetically upregulated genes involved in extracellular matrix organization and embryonic morphogenesis, particularly of the limbs and the skeletal system (Table S3). In contrast, ARMS tumors epigenetically upregulated genes involved in muscle differentiation, suggesting that it is arrested at a later stage of muscle development (Table S3). We used core regulatory circuit (CRC) mapping (Saint-Andre et al., 2016) of our ChIP-seq and RNA-seq data to identify core Cancer Cell 34, 1–16, September 10, 2018 3

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Figure 2. Differences in Epigenetic Profiles Relative to Promoter/Enhancer Activity (A) Heatmap of the 18 chromHMM states used in this analysis. (B) Heatmap depicting proportion of the 18 chromHMM states in ERMS and ARMS for annotated regions of the genome. Proportion of individual marks or regions is directly proportional to the intensity of the blue color. (C) ChIP-seq peaks and chromHMM for MYOG. Scales are indicated on the right, and gene boundaries are marked by dashed lines. (D–F) Representative chromHMM for GAS2 that is selectively expressed in ERMS (D), NOS1 that is selectively expressed in ARMS (E), and AXL that is selectively expressed in myoblasts (F). NOS1 contains a tumor type-specific superenhancer upstream of the gene (boxed H3K27Ac region). Gene boundaries are marked by dashed lines, and colors for each state are indicated as in (A). The two highest proportion states are full-height bars, and remaining states are half height. Intensity of each bar is proportional to the percentage of each state across all stages for that gene. For half-height bars, intensity is scaled starting at 50% of maximum intensity. (G) Venn diagram for genes upregulated or downregulated (green) in ERMS relative to ARMS. Overlap with genes that have differences in chromHMM state (blue) or tumor type-specific superenhancers (red) are shown. See also Table S3.

transcription factors in developing skeletal muscle and identified MYOD1, PAX7, FOXO1, MYF6, MEF2A, and MYOG that are important in myogenesis (Table S3). Other less-studied core transcription factors included FOSL2, which was the most prevalent CRC. This combined with the CRC identified for SMAD3 and TGIF1 suggest that transforming growth factor b signaling may be important in RMS (Carthy et al., 2015) (Table S3). Similarly, hedgehog (HH), notch, and retinoic acid signaling were implicated in RMS (Table S3). Quantitative Analysis of the Rhabdomyosarcoma Proteome and Phosphoproteome We used a 10-plex isobaric labeling mass spectrometry (MS) pipeline with extensive peptide separation power and high 4 Cancer Cell 34, 1–16, September 10, 2018

mass resolution (Bai et al., 2017) to quantify the proteome and phosphoproteome of 12 O-PDXs and human myoblasts and myotubes in replicated experiments (Table S4). Three batches of ten samples were lysed, digested, and labeled with ten distinct tandem mass tags (TMTs) (Figure 3A), and a technical replicate was included in each batch as a control (Table S4). Equally pooled samples were fractionated with 5% of the pool for whole-proteome analysis (15.9 3 106 MS/MS spectra) and the remaining 95% for the phosphoproteome (8.0 3 106 MS/MS spectra) (Figure 3A). We analyzed 16,523 human proteins (12,926 genes) and 62,965 phosphorylation events on 7,177 proteins. We quantitated protein products from 10,333 genes and 12,653 phosphosites across batches, and validated several proteins and phosphoproteins by immunoblot (Figures 3B and 3C).

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

(legend on next page)

Cancer Cell 34, 1–16, September 10, 2018 5

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

The intra-batch reproducibility was 0.98 for both the whole proteome and phosphoproteome (Figures 3D and 3E), and principal component analysis revealed clear separation of each sample type (Figures 3F and 3G). Weighted gene coexpression network analysis (Langfelder and Horvath, 2008) identified six groups of proteins and phosphoproteins differentially expressed in ARMS and ERMS and developing muscle (Table S4). We performed pathway analysis on the differentially expressed proteins and phosphoproteins between myoblasts and myotubes and RMS. We identified several pathways important for myogenesis, including Wnt (Cossu and Borello, 1999), HH (Borycki et al., 1999), phosphoinositide 3-kinase (PI3K) (Jiang et al., 1998), p38/mitogen-activated protein kinase (MAPK) (Bennett and Tonks, 1997), bone morphogenetic protein (BMP) (Reshef et al., 1998), and adenylyl cyclase (AC) (Chen et al., 2005) (Figures S1–S5). MYOG is a transcription factor associated with muscle differentiation that is expressed in both ERMS and ARMS (Hostein et al., 2004) (Figures S1B and S1C). ARMS expressed higher levels of MYOG mRNA and protein than did ERMS, which showed graded expression (Figures 3H and S1C). Because MYF5 is transiently expressed by proliferating myoblasts during myogenic determination, we determined the expression variability and its transcriptional targets in individual tumors. All ARMS expressed lower levels of MYF5 relative to that of muscle, and a subset of ERMS exhibited higher levels of MYF5 mRNA and protein (Figures 3I and S1C). MYF5 and MYOG protein expression was inversely correlated in most RMS (Figure 3J), and MYOG transcriptional target genes were activated in the same pattern (Figures S1D and S1E; Table S5). We performed immunohistochemical staining of MYF5 and MYOG for 17 RMS tumors and their matched O-PDXs and found agreement between staining, gene expression, and proteomic data (Figure 3K and Table S5). Integrated Analyses Using Multiple Platforms To integrate our data across multiple platforms, we used two different approaches. First, we identified differentially expressed proteins (7,546/10,333) and mRNAs (7,317/12,362) in the tumors relative to human myoblasts and performed coexpression clustering of the genes/proteins conserved across both platforms after integration, comprising four integrated protein clusters (IPCs) (Figures 4A and 4B; Table S6). IPC1 contained genes/proteins such as WEE1, CHEK2, and RB1 that were upregulated in

tumors relative to myoblasts. IPC2 contained genes/proteins such as MYC, GLI1/2, and ATM that were upregulated in tumors, but at higher levels in ERMS than in ARMS. IPC3 contained genes/proteins such as NOS1, MYOG, and CDK4 that were upregulated specifically in ARMS. IPC4 contained genes/proteins such as AXL, MYH9, and MURC that were downregulated in tumors relative to myoblasts. Each IPC was used for pathway analysis, which were overlaid with protein-protein interaction (PPI) network data (Figures 4C–4E). In our second approach, we integrated the transcriptomic, epigenomic, proteomic, and phosphoproteomic data with order statistics. For epigenetic analysis, the fractional difference between active states (1, 2, 3, 4, 5, 6, and 10) and repressive/inactive states (14, 17, and 18) was calculated for each gene in each sample and ranked by the log2 fold change for each tumor relative to myoblasts. For proteomic and phosphoproteomic data, the total or relative phosphorylation level of each protein was calculated, and individual proteins were ranked by the log2 ratio of each tumor to myoblasts (Table S6). MYOG was one of the top ranked genes/proteins in ARMS and MYF5 was in the top in ERMS (Figures 3H–3K; Table S6). The top upregulated pathway in ARMS relative to ERMS was myogenesis (Table S6). We observed differential expression of NOS1 mRNA/protein between ARMS and ERMS (Figures 2E and 3B). NOS1 was highest ranked in ARMS (Table S6). We also identified upregulated G2/M and E2F targets in ARMS and ERMS relative to myoblasts (Table S6). Targeting the RAS-CDK4/6 Pathway in Rhabdomyosarcoma We analyzed whole-genome sequencing and whole-exome sequencing data for patient germline tissue, RMS, and matched O-PDXs. We found 18 genes with pathogenic cancer mutations in RMS, 10 of which were included in the oncomine comprehensive assay panel used for the NCI-COG pediatric MATCH trial (Tables S1 and S7). The limited druggability of NOTCH1, NOTCH3, and p53 precludes these targets from pediatric or adult NCI-MATCH trials. However, the remaining seven genes (NRAS, HRAS, NF1, CDKN2A, PIK3CA, FGFR1, and FGFR4) fall into a signal transduction cascade that regulates proliferation through the Rb/E2F pathway (Figure 5A). Although most individually mutated genes in this pathway did not rank highly in our analysis, the RAS, mammalian target of rapamycin (mTOR) complex 1, and PI3K-AKT-mTOR signaling pathways were

Figure 3. Profiling the Rhabdomyosarcoma Proteome and Phosphoproteome (A) Workflow of proteomic and phosphoproteomic profiling of RMS tumors, myotubes, and myoblasts. LC, liquid chromatography; MS, mass spectrometry. (B) Representative immunoblots validating protein and phosphoprotein expression of NOS1, phospho-Ser235/236 of RPS6, and b-actin (loading control). (C) Bar plot of normalized relative fold expression of the proteins in (B) and HMGA2. mb, myoblasts; mt, myotubes. (D and E) Scatterplots of quantitation of proteins (D) and phosphoproteins (E) across replicate samples of SJRHB000026_X1 in the same and in different batches. (F and G) Principal component analyses of myoblasts (light green), myotubes (dark green), ERMS (red), and ARMS (purple) for the whole proteome (F) and phosphoproteome (G). (H and I) Bar plots of RNA expression (FPKM) for MYOG (H) and MYF5 (I) in developing muscle and ARMS and ERMS tumors. Red line indicates expression level in human myoblasts. The ERMS-high group is the samples with highest levels of MYOG and the ERMS-low group the remaining samples. FQ21, fetal quadriceps at 21 weeks of gestation. (J) Scatterplot of the log10 protein expression intensity from proteomic analyses of MYOG and MYF5. The two outlier samples (black boxes) were excluded from linear regression (dashed line). (K) Representative immunostaining of samples for MYOG and MYF5 and scatterplot of the percentage of MYOG and MYF5 immunopositive cells for each O-PDX tumor with a best fit linear regression (dashed line). Scale bar, 20 mm. See also Figures S1–S5; Tables S4 and S5.

6 Cancer Cell 34, 1–16, September 10, 2018

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Figure 4. Integration of Transcriptomic, Epigenetic, and Proteomic/Phosphoproteomic Data (A) Scheme for integration of proteomic and transcriptomic data based on differential expression clustering and PPI network modules. (B) Box plot of normalized relative abundance of individual proteins in each IPC from the integrated analyses. The rectangle represents the interquartile range (IQR) for first to third quartile and the horizontal line within the box is the median. The vertical lines that extend above and below the rectangle represent the maximum and minimum values excluding the outliers that are more than three times the IQR. Data were normalized to myoblasts, and number of genes/proteins in each IPC are shown. Representative gene/proteins in each IPC are shown in gray. (C) Heatmap of the
log10 enrichment false discovery rate for pathway analysis of the genes/proteins in each IPC. For IPC3, no pathway reached statistical significance (p < 0.05) in that cluster. (D) Representative PPI modules for IPC1 showing E2F modules (blue circle) and G2/M checkpoint modules (dashed circle). (E) Expanded view of the G2/M checkpoint module showing representative genes/proteins. See also Table S6.

deregulated in RMS relative to myoblasts (Table S7). Moreover, we found strong deregulation of E2F targets downstream of those pathways (Figures 4C, 4D, and S6A; Table S6). Epigenetic profiling confirmed that E2F transcriptional targets such as CDK1 were activated in RMS (Figures 5B and S6). Proteomic/phosphoproteomic and immunohistochemical staining data indicated that the cyclin D/CDK4/6-RB/E2F pathway was active in O-PDXs (Figures 5C, 5D, S6B, and S6D; Table S7). Combining CDK4/6 inhibitors with MEK inhibitors is an established strategy used in NCI-MATCH trials (Vora et al., 2014) and is being considered for the NCI-COG pediatric MATCH trial (Allen et al., 2017). We tested an MEK inhibitor (trametinib) combined with CDK4/6 inhibitors (palbociclib or abemaciclib) and found

that these combinations synergistically reduced survival of RAS-mutant RMS (RD) and lung cancer (A549) cell lines but did not affect survival of an RB-deficient osteosarcoma cell line (SAOS-2) (Figure 5E and Table S7; data not shown). However, there was not a significant correlation between RAS pathway mutation and sensitivity to MEK inhibitors or CDKN2A mutations and sensitivity to CDK4/6 inhibitors as single agents (https:// braid.stjude.org/masttour/). To establish clinically relevant murine equivalent doses (MEDs), we performed pharmacokinetic (PK) studies of each drug in tumor-bearing mice. We also measured the tumor penetration of each drug (Figure 5F) to relate in vivo exposure to in vitro sensitivity using a previously described mathematical Cancer Cell 34, 1–16, September 10, 2018 7

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

(legend on next page)

8 Cancer Cell 34, 1–16, September 10, 2018

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

model (Stewart et al., 2014). We predicted that sufficient levels of palbociclib and trametinib or abemaciclib and trametinib would be achieved to produce an antitumor effect in vivo (Table S7). Preclinical phase 1 tolerability studies were performed as previously described (Stewart et al., 2014). We next established luciferase-labeled O-PDXs of two ERMS (SJRHB000026_X1 and SJRHB012_Y) and one ARMS (SJRHB013759_X1) (Table S1). SJRHB000026_X1 has oncogenic HRAS and PIK3CA mutations, whereas SJRHB012_Y has no mutations in this pathway (Stewart et al., 2017). We tested palbociclib + trametinib and abemaciclib + trametinib and compared them with a recurrent RMS standard-of-care regimen—irinotecan (IRN) + vincristine (VCR)—in a randomized preclinical phase 2 trial as described previously (Stewart et al., 2014). Mice treated with IRN + VCR had stable disease (SD), whereas mice treated with either palbociclib + trametinib or abemaciclib + trametinib exhibited progressive disease (PD; Figures 5G and 5H; Table S7). Targeting the G2/M and Unfolded Protein Response Pathways in Rhabdomyosarcoma Because targeting the RAS-CDK4/6 pathway in vivo was ineffective and few druggable mutations exist in RMS, we used our integrated dataset to identify therapeutically relevant tumor vulnerabilities. We focused on pathways elucidated from the integrated analyses that also exhibited RMS-selective drug sensitivity using primary cultures of O-PDXs and cell lines (https:// braid.stjude.org/masttour/) (Stewart et al., 2017). In total, 1,767 drug/tumor combinations were analyzed in RMS O-PDXs. The most active drugs were HDAC and proteasome inhibitors (Table S8). Those pathways were not significantly perturbed in our integrated analyses, and previous preclinical phase 2 studies with these drugs (panobinostat and bortezomib) failed to show efficacy in vivo (Stewart et al., 2017). The HSP90 inhibitor ganetespib (GSP) and the WEE1 inhibitor AZD1775 that target unfolded protein response (UPR) and G2/M checkpoint pathways, respectively, found by our integrated analyses also showed significant in vitro activity (Table S8). Many HSP90 clients are key regulators of signal transduction pathways (Karagoz and Rudiger, 2015) (Figure 6A). As hyperphosphorylation of the heat-shock transcription factor HSF1 is a hallmark of cellular stress (Akerfelt et al., 2010), we found augmented protein levels and phosphorylation of HSF1 (Figure 6B and Table S4), as well as HSF1 target gene expression, in RMS tumors (Figures 6C–6E and Table S8). RMS cell lines and primary cultures of O-PDX cells were also more sensitive

to GSP than were other pediatric solid tumors (Figure 6F). Treatment of RD cells with GSP rapidly induced expression of HSP70 and other HSF1 target genes (Figure 6G). We screened 156 drugs in combination with several concentrations of GSP (100 nM, 33 nM, 10 nM, and 3 nM). We used a bivariate response to additive interacting doses model (Twarog et al., 2016) to analyze the effects of each drug combination and found that GSP potentiated chemotherapeutic-induced death of RD and O-PDX cells (Table S8). The target of AZD1775, WEE1, is expressed at higher levels in RMS relative to other pediatric cancers, except Wilms’ tumor and retinoblastoma (https://pecan.stjude.org) (Zhou et al., 2016) (Figure 7A). WEE1 regulates the G2/M checkpoint and DNA replication during the S phase of the cell cycle (Do et al., 2013) (Figure 7B). Epigenetic, gene expression, and proteomic/ phosphoproteomic data for RMS relative to myoblasts revealed upregulation of this pathway (Figures 7C and 7D; Table S6). Treatment of RMS cell lines with AZD1775 led to G2/M phase arrest, mitotic catastrophe, and nuclear fragmentation (Figures 7E–7G). AZD1775 combined with IRN or VCR led to high rates of DNA damage (Figure 7G), suggesting that both WEE1 and HSP90 may be druggable targets in RMS. Targeting the G2/M and UPR Pathways In Vivo To establish a clinically relevant MED for GSP, we performed PK studies with tumor-bearing mice (Figure 8A); PK for AZD1775 was done previously (Stewart et al., 2017). Sustained levels of GSP and AZD1775 in tumors were above those required to potentiate cytotoxic effects in culture (Figures 8A and 8B). We performed preclinical phase 1 studies of clinically relevant GSP and AZD1775 combinations (Figure 8C) as described previously (Stewart et al., 2014, 2017). Two schedules of IRN commonly used to treat children with RMS and other solid tumors (Furman et al., 1999) were included. The MED of GSP (10 mg/kg) and AZD1775 (60 mg/kg twice daily) was well tolerated in all groups (data not shown). We next performed preclinical phase 2 trials using luciferaselabeled O-PDXs as described previously (Stewart et al., 2014, 2017). Mice were randomized into eight treatment groups— AZD1775 + VCR + IRN, AZD1775 + IRN, GSP + VCR + IRN, GSP + IRN, VCR + IRN, GSP, IRN, or placebo—and treated for four to six courses (21 days per course) (Table S9). Tumor burden was measured weekly by bioluminescence, and some mice treated with GSP + IRN + VCR or AZD1775 + IRN + VCR exhibited partial responses (PRs) or complete responses (CRs) (Stewart et al., 2017) (Figures 8D and 8E; Table S9).

Figure 5. Targeting RAS and CDK4/6 in Rhabdomyosarcoma (A) Simplified pathway map for the RAS and CDK4/6 pathways. Recurrent mutations in genes encoding proteins highlighted in purple and those in bold were targeted with small-molecule inhibitors. (B) Representative ChIP-seq tracks and chromHMM for CDK1. Scales are indicated on the right. Only one splice variant is shown for simplicity. (C) Bar plot of normalized relative CDK1 protein expression for myoblasts, myotubes, and RMS tumors. (D) Representative micrographs of ERMS and ARMS tumor sections with H&E staining or immunohistochemical staining of cyclin D2 (ERMS), cyclin D3 (ARMS), and phosphoRB1 (pRB1). Scale bars, 50 mm in the full image and 10 mm in the magnified inset. (E) Representative combination drug sensitivity study of ten concentrations of palbociclib and four concentrations of trametinib in RD cells. Survival is plotted relative to untreated cells (100% survival) and complete ablation (0% survival) with a positive-control cytotoxic compound. (F) PK plot of palbociclib in O-PDX tumors (blue) and plasma (red), n = 3 mice per treatment group. Mean and SD are plotted for each time point. (G and H) Representative plot of SJRHB013759_X1 O-PDX tumor burden, as measured by bioluminescence over time (G), and images of mice with PD in the placebo control group and the palbociclib + trametinib treatment group (H). See also Figure S6 and Table S7.

Cancer Cell 34, 1–16, September 10, 2018 9

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Figure 6. Targeting HSP90 in Rhabdomyosarcoma (A) Illustration of the six myogenic pathways deregulated in RMS and modulated by HSP90. (B) Bar plots of HSF1 and S368 phosphorylated HSF1 in myoblasts, myotubes, and RMS tumors. Data normalized to myoblasts (red line). (C and D) Bar plots of quantitative expression of 396 HSF1 target genes (C) and representative examples of HSF1 target gene expression (D) in RMS tumors normalized relative to that of myoblasts, as measured by RNA-seq. Bars denote mean and SD across multiple tumors. (E) ChromHMM of an HSF1 target gene (HSPA1A) and corresponding gene expression by RNA-seq (FPKM). (F) Heatmap of the log10 50% effective concentrations of a subset of drugs tested. Bolded terms denote primary cultures of O-PDXs, whereas nonbolded terms denote established cell lines. (G) Immunoblot of HSP70 expression in RD cells after treatment with GSP alone or IRN + VCR. Control indicates untreated cells. HSP70 intensities were quantified and normalized to b-actin before calculating fold differences between treated and untreated samples. See also Table S8.

To further explore the efficacy of combination therapy using GSP and AZD1775, we performed a blinded, randomized, and placebo-controlled preclinical phase 3 trial. To recapitulate the heterogeneity of RMS patients, we used two ERMS O-PDXs derived from recurrent patients (SJRHB000026_X1, SJRHB012_Y) and two ARMS O-PDXs, one derived from initial diagnosis (SJRHB013757_X2) and one from recurrence (SJRHB013759_X1) (Table S1). We enrolled 499 mice and randomized them to 14 treatment groups (Table S9). We included three different standard-of-care groups (VCR + IRN) with different doses and schedules of IRN (Table S9). Four treatment groups incorporated GSP and an additional four groups incorpo10 Cancer Cell 34, 1–16, September 10, 2018

rated AZD1775 with standard-of-care agents (Table S9). A subset of these data (four treatment groups, 140 mice) has been published elsewhere (Stewart et al., 2017) but is included here for direct comparison because the mice were in the cohort of 499 that were randomized. The overall survival for mice in single-agent GSP and AZD1775 groups was not different from placebo (Table S9). Although CRs and PRs occurred in the GSP + VCR + IRN for some individual O-PDXs, overall survival of the GSP + VCR + IRN group was not different than that of the VCR + IRN group for any of the four O-PDXs tested (Table S9). All four O-PDXs had improved outcomes when treated with AZD1775 + VCR + IRN relative to VCR + IRN, AZD1775 alone,

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Figure 7. Targeting WEE1 in Rhabdomyosarcoma (A) Box plot of WEE1 gene expression in pediatric cancers from https://pecan.stjude.org/home. The rectangle represents the IQR for first to third quartile, and the vertical line within the box is the median. The horizontal lines that extend from the rectangle represent the maximum and minimum values excluding the outliers (circles) that are more than three times the IQR. CPC, choroid plexus carcinoma; EPD, ependymoma; HGG, high-grade glioma; LGG, low-grade glioma; MB, medulloblastoma; AML, acute myeloid leukemia; B-ALL, B cell acute lymphoblastic leukemia; MLL, mixed lineage leukemia; T-ALL, T cell acute lymphoblastic leukemia; ACT, adrenocortical tumor; MEL, melanoma; NBL, neuroblastoma; OS, osteosarcoma; RB, retinoblastoma; WT, Wilms’ tumor. (B) Simplified drawing of the pathway regulated by WEE1 at the G2/M checkpoint. Arrows indicate activation and bars indicate repression. (C) ChromHMM for WEE1 with corresponding RNA-seq expression across samples. (D) Bar plot of WEE1 protein levels in myoblasts, myotubes, and RMS tumors normalized to myoblasts (red line). (E–G) DNA content and cell-cycle distribution for RD cells before (black) and after treatment (red) with AZD1775 (E), representative micrographs of DAPI-stained nuclei of these cells highlighting fragmented nuclei/mitotic catastrophe (F), and proportion of cells with nuclear features consistent with mitotic catastrophe after the treatment shown (G). Bars in (G) denote mean and SD of duplicate scoring. Scale bars in (F), 2 mm.

or placebo (Figures 8F–8I and Table S9). When we tested AZD1775 + IRN in dosing groups that most closely replicate the NCI-COG pediatric study (NCT02095132), we found that only one O-PDX (SJRHB012_Y) had a significantly better response (p = 0.02) compared with IRN + VCR; however, the addition of VCR to these dosing groups significantly improved survival for both ERMS O-PDXs tested (Figures 8F–8I and Table S9). Compiled responses from the four O-PDXs revealed that AZD1775 + VCR + IRN produced 70% CR + PR (26/37), whereas AZD1775 + IRN produced 39% CR + PR (14/36) (Table S9). Discernible differences in response to VCR + IRN, GSP + VCR + IRN, and AZD1775 + VCR + IRN were noted between the individual O-PDXs, signifying the importance of using a variety of O-PDXs and relevant standard-of-care regimens for comparison (Figure 8I and Table S9). DISCUSSION Although several trials have demonstrated the feasibility of precision medicine for children with cancer (Harris et al., 2016; Parsons et al., 2016), these trials pose several challenges unique to pediatric cancer. First, because the number of patients is usually small, molecularly targeted therapy is prescribed on the basis of genetic lesions rather than cancer type. If the cellular context of the mutation affects the outcome of therapy, the responses will be difficult to interpret (Schneider et al., 2017). Second, few actionable mutations in pediatric tumors are currently druggable with available therapeutics, making it diffi-

cult to enroll enough patients to achieve statistical significance (Allen et al., 2017; Harris et al., 2016; Parsons et al., 2016). These barriers may preclude the ability to discern a sizable benefit of precision medicine in children with cancer from targeted clinical trials alone. Previous genomic studies indicate that approximately 50%– 75% of patients with intermediate- and high-risk RMS have a mutation in the RAS pathway (Chen et al., 2013). To encompass the largest group of recurrent RMS patients, we used a broad definition of the RAS pathway to include mutations in NRAS, HRAS, KRAS, NF1, FGFR1/4, and PIKC3A. In some adult RAS-mutant cancer types, targeting kinases in the RAS and PI3K pathways has shown promise (Alagesan et al., 2015; Jokinen and Koivunen, 2015). Our integrated analyses revealed that pediatric RMS also exhibit robust signaling through the RAS pathway. Although combined CDK4/6 and MEK inhibitors synergistically induced RMS cell death in culture, they failed to achieve efficacy in vivo. Our PK studies suggested this was not due to failure of the drugs to penetrate the tumor. It is possible that compensatory changes in other signal transduction pathways allowed the tumors to continue to grow (Renshaw et al., 2013). To move beyond genomic profiling and identify fundamental pathways that may be deregulated in RMS, we performed integrated transcriptomic, epigenetic, and proteomic/phosphoproteomic analysis. We identified the UPR and G2/M-mitotic spindle checkpoint pathways as the most altered in RMS integrated analysis. Others have recently confirmed that the protein Cancer Cell 34, 1–16, September 10, 2018 11

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

(legend on next page)

12 Cancer Cell 34, 1–16, September 10, 2018

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

homeostasis pathway is a fundamental vulnerability of RMS (Sabnis et al., 2016). By targeting HSP90 with GSP, we found that the UPR pathway is vulnerable in RMS tumors due in part to an augmented burden on HSP90 to maintain the six myogenic signal transduction pathways. HSP90 also regulates the G2/M checkpoint and mitotic spindle assembly pathways through WEE1. Because HSP90 levels are barely sufficient to manage the cellular stress imposed by oncogenic transformation in RMS tumors, increasing cellular stress with chemotherapeutics, especially those targeting HSP90, may achieve synthetic lethality. Indeed, we observed that RMS cell lines and O-PDXs were more sensitive to GSP than other pediatric solid tumors, and the same was true for WEE1 inhibition. Moreover, combining GSP with subtherapeutic doses of chemotherapeutics with no effect on their own synergistically ablated RMS cells in culture. We found that GSP + VCR + IRN was well tolerated in preclinical phase 1 studies, partially efficacious in phase 2 studies, but not efficacious in a blinded, randomized, and placebo-controlled phase 3 study. This highlights the value of our phase 1/2/3 approach to validate therapeutic combinations before moving into clinical trials for children with RMS. Although targeting the UPR pathway was not effective, targeting additional components of this pathway, such as HSP70, may improve efficacy in future studies (Sabnis et al., 2016). A limitation of our study is bias toward aggressive refractory RMS tumors. Our integrated analyses were not comprehensive as they did not include low-risk RMS and most O-PDXs were derived from recurrent tumors. This limitation may extend to the correlation between putative actionable mutations and drug sensitivity. There was no correlation between actionable mutation and drug sensitivity for any drugs that were tested in our high-throughput screening, including AZD1775 and MEK inhibitors that are thought to be relevant for tumors with mutations in the RAS pathway (https://braid.stjude.org/ masttour/). Another limitation of our preclinical study was a lack of pharmacodynamics. Although we observed sufficient drug penetration to elicit antitumor effects, we did not directly demonstrate this in vivo with a validated pharmacodynamic assay. Therefore, HSP90 and other components of the UPR pathway may still be rational targets for RMS with future generations of drugs. Nevertheless, AZD1775 + VCR + IRN was efficacious in our preclinical study and our findings may inform future clinical trials for patients

with high-risk RMS. It is important to emphasize the need to include VCR with AZD1775 and IRN to improve the likelihood of response in RMS patients, based on our extensive phase 3 preclinical study. STAR+METHODS Detailed methods are provided in the online version of this paper and include the following: d d d d

d

KEY RESOURCES TABLE CONTACT FOR REAGENT AND RESOURCE SHARING EXPERIMENTAL MODEL AND SUBJECT DETAILS B Animals and Orthotopic Patient-Derived Xenografts METHOD DETAILS B Genomic DNA Extraction B Chromatin Immunoprecipitation and Sequencing (ChIP-Seq) B Whole Genome Sequencing and Analysis B RNA Sequencing and Analysis B DNA Methylation and Analysis B Proteomics and Phosphoproteomics B Analysis of Transcriptional Targets B Immunohistochemistry B General Drug Screening B Combination Drug Screening with Ganetespib B Response Surface Model Analysis of Select Drug Combinations B Pharmacokinetic Studies B Western Blots of Treated RD Cells B Flow Cytometry B QPCR B Microscopy B Preclinical Testing B SJRHB000026_X1 (ERMS) B SJRHB012_Y (ERMS) B SJRHB013759_X1 (ARMS): B Xenogen Imaging and Quantification QUANTIFICATION AND STATISTICAL ANALYSIS B ChIP Sequencing Analysis B Hidden Markov Modeling analysis B Core Regulatory Circuit (CRC) Analysis B DNA Methylation and Analysis

Figure 8. Preclinical Testing of AZD1775 and Ganetespib (A) PK plot of GSP in O-PDX tumors (blue) and plasma (red), n = 3 mice per treatment group. Mean and SD are plotted for each time point. (B) Dose-response curve for GSP and AZD1775 using primary cultures of SJRHB000026_X1 O-PDX. (C) Drug schedule selected for each combination treatment regimen. The daily 3 5 3 2 low-dose protracted schedule for IRN is shown, but we also tested the alternative daily 3 5 schedule. (D) Xenogen images of mice with an O-PDX tumor (SJRHB012_Y) treated with IRN + VCR or GSP + IRN + VCR, which had PD or a CR. (E) Representative plot of tumor burden, as measured by bioluminescence over time, in a preclinical phase 2 study with each of the three O-PDX tumors indicated. Each line is a different mouse. (F) Representative images of mice and their corresponding tumors with PD in the placebo control group, AZD1775 group, and GSP + IRN + VCR group. Mice with SD and a CR are shown for the IRN + VCR and AZD1775 + IRN + VCR treatment groups, respectively. (G) Line plot of tumor burden, as measured by bioluminescence over time, in a preclinical phase 3 study of the SJRHB012_Y O-PDX model. (H) Survival curve for the SJRHB012_Y O-PDX for the four indicated treatment groups. p Value is shown for treatment groups compared with IRN + VCR. (I) Histogram of the percentage of CR and PR for each of the four O-PDX models for all 14 treatment groups. Two clinically relevant doses, high (h) and low (l), of AZD1775 and GSP were tested along with several different schedules and dosing regimens of IRN used in pediatric cancer patients. *p < 0.05 compared with IRN + VCR. See also Table S9.

Cancer Cell 34, 1–16, September 10, 2018 13

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

B

Integrated Analysis Statistical Analysis for Preclinical Phase III Study DATA AND SOFTWARE AVAILABILITY B

d

SUPPLEMENTAL INFORMATION Supplemental Information includes six figures and nine tables and can be found with this article online at https://doi.org/10.1016/j.ccell.2018.07.012. ACKNOWLEDGMENTS We thank Nisha Badders for editing the manuscript. This work was supported, in part, by Cancer Center Support (CA21765) from the NCI, grants to M.A.D. from the NIH (EY014867 and EY018599 and CA168875) and ALSAC. M.A.D. was also supported by the Alex’s Lemonade Stand Foundation for Childhood Cancer, the Tully Family Foundation, and the Peterson Foundation. E.S. was supported by the National Comprehensive Cancer Network and the National Pediatric Cancer Foundation and is a St. Baldrick’s Foundation Scholar with generous support from the Invictus Fund. J.P. is supported by funding from the NIH (AG047928). The majority of this research was supported by the Howard Hughes Medical Institute. Whole-genome sequencing, WGBS, and RNA-seq were performed as part of the St. Jude Children’s Research Hospital/Washington University Pediatric Cancer Genome Project. Preclinical studies were performed with assistance from the Animal Imaging Shared Resource; electron microscopy was performed with assistance from the Cell and Tissue Imaging Shared Resource; pharmacokinetics was performed with assistance from the Preclinical Pharmacokinetics Shared Resource; MS was performed in the Proteomics Shared Resource; and histopathologic analysis was performed with assistance from the Veterinary Pathology Shared Resource at St. Jude. AUTHOR CONTRIBUTIONS Conceptualization, M.A.D., E.S., J.M., H.W., X.C., and J.P.; Methodology, E.S., J.M., H.W., X.C., and J.P.; Formal Analysis, Y.Y., Y.L., T.I.S., J.-H.C., X.W., B.X., P.G., Y.F., Y.L., M.R., and K.B.; Investigation, M.A.D., E.S., J.M., H.W., J.P., V.H., M.O., B.G., J.D., K.B., L.G., J.J., B.B.F., M.R.C., A.S., and A.B.; Writing – Original Draft, M.A.D.; Writing – Review & Editing, M.A.D., E.S.; Supervision, M.A.D. and J.P.; Project Administration, M.A.D., J.P., A.P., J.E., S.S., X.Z., H.M., D.Y., E.R.M., R.K.W., J.Z., and J.R.D.; Funding Acquisition, M.A.D. and J.P. DECLARATION OF INTERESTS The authors declare no competing interests. Received: January 11, 2018 Revised: April 9, 2018 Accepted: July 25, 2018 Published: August 23, 2018 REFERENCES Aerts, S., Lambrechts, D., Maity, S., Van Loo, P., Coessens, B., De Smet, F., Tranchevent, L.C., De Moor, B., Marynen, P., Hassan, B., et al. (2006). Gene prioritization through genomic data fusion. Nat. Biotechnol. 24, 537–544. Afzali, A.M., Ruck, T., Herrmann, A.M., Iking, J., Sommer, C., Kleinschnitz, C., Preubetae, C., Stenzel, W., Budde, T., Wiendl, H., et al. (2016). The potassium channels TASK2 and TREK1 regulate functional differentiation of murine skeletal muscle cells. Am. J. Physiol. Cell Physiol. 311, C583–C595.

Aldiri, I., Xu, B., Wang, L., Chen, X., Hiler, D., Griffiths, L., Valentine, M., Shirinifard, A., Thiagarajan, S., Sablauer, A., et al. (2017). The dynamic epigenetic landscape of the retina during development, reprogramming, and tumorigenesis. Neuron 94, 550–568.e10. Allen, C.E., Laetsch, T.W., Mody, R., Irwin, M.S., Lim, M.S., Adamson, P.C., Seibel, N.L., Parsons, D.W., Cho, Y.J., Janeway, K., et al. (2017). Target and agent prioritization for the Children’s Oncology Group-National Cancer Institute Pediatric MATCH Trial. J. Natl. Cancer Inst. 109, https://doi.org/10. 1093/jnci/djw274. Bai, B., Tan, H., Pagala, V.R., High, A.A., Ichhaporia, V.P., Hendershot, L., and Peng, J. (2017). Deep profiling of proteome and phosphoproteome by isobaric labeling, extensive liquid chromatography, and mass spectrometry. Methods Enzymol. 585, 377–395. Benjamini, Y., and Hochberg, Y. (1995). Controlling the false discovery rate: a practical and powerful approach to multiple testing. J. R. Stat. Soc. B 57, 289–300. Bennett, A.M., and Tonks, N.K. (1997). Regulation of distinct stages of skeletal muscle differentiation by mitogen-activated protein kinases. Science 278, 1288–1291. Borycki, A.G., Brunk, B., Tajbakhsh, S., Buckingham, M., Chiang, C., and Emerson, C.P., Jr. (1999). Sonic hedgehog controls epaxial muscle determination through Myf5 activation. Development 126, 4053–4063. Carthy, J.M., Sundqvist, A., Heldin, A., van Dam, H., Kletsas, D., Heldin, C.H., and Moustakas, A. (2015). Tamoxifen inhibits TGF-beta-mediated activation of myofibroblasts by blocking non-smad signaling through ERK1/2. J. Cell Physiol. 230, 3084–3092. Chen, A.E., Ginty, D.D., and Fan, C.M. (2005). Protein kinase A signalling via CREB controls myogenesis induced by Wnt proteins. Nature 433, 317–322. Chen, X., Stewart, E., Shelat, A.A., Qu, C., Bahrami, A., Hatley, M., Wu, G., Bradley, C., McEvoy, J., Pappo, A., et al. (2013). Targeting oxidative stress in embryonal rhabdomyosarcoma. Cancer Cell 24, 710–724. Cossu, G., and Borello, U. (1999). Wnt signaling and the activation of myogenesis in mammals. EMBO J. 18, 6867–6872. Do, K., Doroshow, J.H., and Kummar, S. (2013). Wee1 kinase as a target for cancer therapy. Cell Cycle 12, 3159–3164. Ernst, J., and Kellis, M. (2012). ChromHMM: automating chromatin-state discovery and characterization. Nat. Methods 9, 215–216. Furman, W.L., Stewart, C.F., Poquette, C.A., Pratt, C.B., Santana, V.M., Zamboni, W.C., Bowman, L.C., Ma, M.K., Hoffer, F.A., Meyer, W.H., et al. (1999). Direct translation of a protracted irinotecan schedule from a xenograft model to a phase I trial in children. J. Clin. Oncol. 17, 1815–1824. Goldman, J.W., Raju, R.N., Gordon, G.A., El-Hariry, I., Teofilivici, F., Vukovic, V.M., Bradley, R., Karol, M.D., Chen, Y., Guo, W., et al. (2013). A first in human, safety, pharmacokinetics, and clinical activity phase I study of once weekly administration of the Hsp90 inhibitor ganetespib (STA-9090) in patients with solid malignancies. BMC Cancer 13, 152. Gong, J., Yuan, Y., Ward, A., Kang, L., Zhang, B., Wu, Z., Peng, J., Feng, Z., Liu, J., and Xu, X.Z.S. (2017). The C. elegans taste receptor homolog LITE-1 is a photoreceptor. Cell 168, 325. Harris, M.H., DuBois, S.G., Glade Bender, J.L., Kim, A., Crompton, B.D., Parker, E., Dumont, I.P., Hong, A.L., Guo, D., et al. (2016). Multicenter feasibility study of tumor molecular profiling to inform therapeutic decisions in advanced pediatric solid tumors: the individualized cancer therapy (iCat) study. JAMA Oncol. https://doi.org/10.1001/jamaoncol.2015.5689.

Akerfelt, M., Morimoto, R.I., and Sistonen, L. (2010). Heat shock factors: integrators of cell stress, development and lifespan. Nat. Rev. Mol. Cell Biol. 11, 545–555.

Hostein, I., Andraud-Fregeville, M., Guillou, L., Terrier-Lacombe, M.J., Deminiere, C., Ranchere, D., Lussan, C., Longavenne, E., Bui, N.B., Delattre, O., and Coindre, J.M. (2004). Rhabdomyosarcoma: value of myogenin expression analysis and molecular testing in diagnosing the alveolar subtype: an analysis of 109 paraffin-embedded specimens. Cancer 101, 2817–2824.

Alagesan, B., Contino, G., Guimaraes, A.R., Corcoran, R.B., Deshpande, V., Wojtkiewicz, G.R., Hezel, A.F., Wong, K.K., Loda, M., Weissleder, R., et al. (2015). Combined MEK and PI3K inhibition in a mouse model of pancreatic cancer. Clin. Cancer Res. 21, 396–404.

Hovestadt, V., Jones, D.T., Picelli, S., Wang, W., Kool, M., Northcott, P.A., Sultan, M., Stachurski, K., Ryzhova, M., Warnatz, H.J., et al. (2014). Decoding the regulatory landscape of medulloblastoma using DNA methylation sequencing. Nature 510, 537–541.

14 Cancer Cell 34, 1–16, September 10, 2018

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Jiang, B.H., Zheng, J.Z., and Vogt, P.K. (1998). An essential role of phosphatidylinositol 3-kinase in myogenic differentiation. Proc. Natl. Acad. Sci. USA 95, 14179–14183. Jokinen, E., and Koivunen, J.P. (2015). MEK and PI3K inhibition in solid tumors: rationale and evidence to date. Ther. Adv. Med. Oncol. 7, 170–180. Karagoz, G.E., and Rudiger, S.G. (2015). Hsp90 interaction with clients. Trends Biochem. Sci. 40, 117–125. Kharchenko, P.V., Tolstorukov, M.Y., and Park, P.J. (2008). Design and analysis of ChIP-seq experiments for DNA-binding proteins. Nat Biotechnol 26, 1351–1359. Langfelder, P., and Horvath, S. (2008). WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics 9, 559. Langfelder, P., Zhang, B., and Horvath, S. (2008). Defining clusters from a hierarchical cluster tree: the Dynamic Tree Cut package for R. Bioinformatics 24, 719–720. Li, H., and Durbin, R. (2009). Fast and accurate short read alignment with Burrows-Wheeler transform. Bioinformatics 25, 1754–1760. Li, H., Handsaker, B., Wysoker, A., Fennell, T., Ruan, J., Homer, N., Marth, G., Abecasis, G., and Durbin, R; 1000 Genome Project Data Processing Subgroup (2009). The sequence alignment/map format and SAMtools. Bioinformatics 25, 2078–2079. Liang, X., Ubhayakar, S., Liederer, B.M., Dean, B., Ran-Ran Qin, A., ShahidiLatham, S., and Deng, Y. (2011). Evaluation of homogenization techniques for the preparation of mouse tissue samples to support drug discovery. Bioanalysis 3, 1923–1933. Liberzon, A., Subramanian, A., Pinchback, R., Thorvaldsdottir, H., Tamayo, P., and Mesirov, J.P. (2011). Molecular signatures database (MSigDB) 3.0. Bioinformatics 27, 1739–1740. Lister, R., Pelizzola, M., Dowen, R.H., Hawkins, R.D., Hon, G., Tonti-Filippini, J., Nery, J.R., Lee, L., Ye, Z., Ngo, Q.M., et al. (2009). Human DNA methylomes at base resolution show widespread epigenomic differences. Nature 462, 315–322. Liu, Y., Siegmund, K.D., Laird, P.W., and Berman, B.P. (2012). Bis-SNP: combined DNA methylation and SNP calling for Bisulfite-seq data. Genome Biol. 13, R61. McAlister, G.C., Nusinow, D.P., Jedrychowski, M.P., Wuhr, M., Huttlin, E.L., Erickson, B.K., Rad, R., Haas, W., and Gygi, S.P. (2014). MultiNotch MS3 enables accurate, sensitive, and multiplexed detection of differential expression across cancer cell line proteomes. Anal. Chem. 86, 7150–7158. Nesvizhskii, A.I., and Aebersold, R. (2005). Interpretation of shotgun proteomic data: the protein inference problem. Mol. Cell. Proteomics 4, 1419–1440.

domyosarcoma: a report from the Intergroup Rhabdomyosarcoma Study Group. J. Clin. Oncol. 17, 3487–3493. Parrish, K.E., Cen, L., Murray, J., Calligaris, D., Kizilbash, S., Mittapalli, R.K., Carlson, B.L., Schroeder, M.A., Sludden, J., Boddy, A.V., et al. (2015). Efficacy of PARP inhibitor rucaparib in orthotopic glioblastoma xenografts is limited by ineffective drug penetration into the central nervous system. Mol. Cancer Ther. 14, 2735–2743. Parsons, D.W., Roy, A., Yang, Y., Wang, T., Scollon, S., Bergstrom, K., Kerstein, R.A., Gutierrez, S., Petersen, A.K., Bavle, A., et al. (2016). Diagnostic yield of clinical tumor and germline whole-exome sequencing for children with solid tumors. JAMA Oncol. https://doi.org/10.1001/jamaoncol. 2015.5699. Patnaik, A., Rosen, L.S., Tolaney, S.M., Tolcher, A.W., Goldman, J.W., Gandhi, L., Papadopoulos, K.P., Beeram, M., Rasco, D.W., Hilton, J.F., et al. (2016). Efficacy and safety of abemaciclib, an inhibitor of CDK4 and CDK6, for patients with breast cancer, non-small cell lung cancer, and other solid tumors. Cancer Discov. 6, 740–753. Raney, R.B., Anderson, J.R., Barr, F.G., Donaldson, S.S., Pappo, A.S., Qualman, S.J., Wiener, E.S., Maurer, H.M., and Crist, W.M. (2001). Rhabdomyosarcoma and undifferentiated sarcoma in the first two decades of life: a selective review of intergroup rhabdomyosarcoma study group experience and rationale for Intergroup Rhabdomyosarcoma Study V. J. Pediatr. hematology/oncology 23, 215–220. Raub, T.J., Wishart, G.N., Kulanthaivel, P., Staton, B.A., Ajamie, R.T., Sawada, G.A., Gelbert, L.M., Shannon, H.E., Sanchez-Martinez, C., and De Dios, A. (2015). Brain exposure of two selective dual CDK4 and CDK6 inhibitors and the antitumor activity of CDK4 and CDK6 inhibition in combination with temozolomide in an intracranial glioblastoma xenograft. Drug Metab. Dispos. 43, 1360–1371. R Development Core Team (2011). R: A Language and Environment for Statistical Computing (Vienna, Austria: R Foundation for Statistical Computing). Renshaw, J., Taylor, K.R., Bishop, R., Valenti, M., De Haven Brandon, A., Gowan, S., Eccles, S.A., Ruddle, R.R., Johnson, L.D., Raynaud, F.I., et al. (2013). Dual blockade of the PI3K/AKT/mTOR (AZD8055) and RAS/MEK/ ERK (AZD6244) pathways synergistically inhibits rhabdomyosarcoma cell growth in vitro and in vivo. Clin. Cancer Res. 19, 5940–5951. Reshef, R., Maroto, M., and Lassar, A.B. (1998). Regulation of dorsal somitic cell fates: BMPs and Noggin control the timing and pattern of myogenic regulator expression. Genes Dev. 12, 290–303. Ritz, C., and Streibig, J.C. (2005). Bioassay analysis using R. J. Stat. Softw. 12, https://doi.org/10.18637/jss.v012.i05.

Newton, W.A., Jr., Soule, E.H., Hamoudi, A.B., Reiman, H.M., Shimada, H., Beltangady, M., and Maurer, H. (1988). Histopathology of childhood sarcomas, intergroup rhabdomyosarcoma studies I and II: clinicopathologic correlation. J. Clin. Oncol. 6, 67–75.

Robinson, J.T., Thorvaldsdo´ttir, H., Winckler, W., Guttman, M., Lander, E.S., Getz, G., and Mesirov, J.P. (2011). Integrative genomics viewer. Nat. Biotechnol 29, 24–26.

Niu, M., Cho, J.H., Kodali, K., Pagala, V., High, A.A., Wang, H., Wu, Z., Li, Y., Bi, W., Zhang, H., et al. (2017). Extensive peptide fractionation and y1 ion-based interference detection method for enabling accurate quantification by isobaric labeling and mass spectrometry. Anal. Chem. 89, 2956–2963.

Sabnis, A.J., Guerriero, C.J., Olivas, V., Sayana, A., Shue, J., Flanagan, J., Asthana, S., Paton, A.W., Paton, J.C., Gestwicki, J.E., et al. (2016). Combined chemical-genetic approach identifies cytosolic HSP70 dependence in rhabdomyosarcoma. Proc. Natl. Acad. Sci. USA 113, 9015–9020.

O’Geen, H., Squazzo, S.L., Iyengar, S., Blahnik, K., Rinn, J.L., Chang, H.Y., Green, R., and Farnham, P.J. (2007). Genome-wide analysis of KAP1 binding suggests autoregulation of KRAB-ZNFs. PLoS Genet. 3, e89.

Saint-Andre, V., Federation, A.J., Lin, C.Y., Abraham, B.J., Reddy, J., Lee, T.I., Bradner, J.E., and Young, R.A. (2016). Models of human core transcriptional regulatory circuitries. Genome Res. 26, 385–396.

Ognjanovic, S., Linabery, A.M., Charbonneau, B., and Ross, J.A. (2009). Trends in childhood rhabdomyosarcoma incidence and survival in the United States, 1975-2005. Cancer 115, 4218–4226.

Schneider, G., Schmidt-Supprian, M., Rad, R., and Saur, D. (2017). Tissuespecific tumorigenesis: context matters. Nat. Rev. Cancer 17, 239–253.

Ouellet, D., Kassir, N., Chiu, J., Mouksassi, M.S., Leonowens, C., Cox, D., DeMarini, D.J., Gardner, O., Crist, W., and Patel, K. (2016). Population pharmacokinetics and exposure-response of trametinib, a MEK inhibitor, in patients with BRAF V600 mutation-positive melanoma. Cancer Chemother. Pharmacol. 77, 807–817. Pappo, A.S., Anderson, J.R., Crist, W.M., Wharam, M.D., Breitfeld, P.P., Hawkins, D., Raney, R.B., Womer, R.B., Parham, D.M., Qualman, S.J., and Grier, H.E. (1999). Survival after relapse in children and adolescents with rhab-

Shern, J.F., Chen, L., Chmielecki, J., Wei, J.S., Patidar, R., Rosenberg, M., Ambrogio, L., Auclair, D., Wang, J., Song, Y.K., et al. (2014). Comprehensive genomic analysis of rhabdomyosarcoma reveals a landscape of alterations affecting a common genetic axis in fusion-positive and fusion-negative tumors. Cancer Discov. 4, 216–231. Shimamura, T., Perera, S.A., Foley, K.P., Sang, J., Rodig, S.J., Inoue, T., Chen, L., Li, D., Carretero, J., Li, Y.C., et al. (2012). Ganetespib (STA-9090), a nongeldanamycin HSP90 inhibitor, has potent antitumor activity in in vitro and in vivo models of non-small cell lung cancer. Clin. Cancer Res. 18, 4973–4985.

Cancer Cell 34, 1–16, September 10, 2018 15

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Smith, M.A., Altekruse, S.F., Adamson, P.C., Reaman, G.H., and Seibel, N.L. (2014). Declining childhood and adolescent cancer mortality. Cancer 120, 2497–2506.

Taus, T., Kocher, T., Pichler, P., Paschke, C., Schmidt, A., Henrich, C., and Mechtler, K. (2011). Universal and confident phosphorylation site localization using phosphoRS. J. Proteome Res. 10, 5354–5362.

Stewart, E., Federico, S.M., Chen, X., Shelat, A.A., Bradley, C., Gordon, B., Karlstrom, A., Twarog, N.R., Clay, M.R., Bahrami, A., et al. (2017). Orthotopic patient-derived xenografts of paediatric solid tumours. Nature 549, 96–100.

Twarog, N.R., Stewart, E., Hammill, C.V., and Shelat, A.A. (2016). BRAID: a unifying paradigm for the analysis of combined drug action. Sci. Rep. 6, 25523.

Stewart, E., Goshorn, R., Bradley, C., Griffiths, L.M., Benavente, C., Twarog, N.R., Miller, G.M., Caufield, W., Freeman, B.B., 3rd, Bahrami, A., et al. (2014). Targeting the DNA repair pathway in Ewing sarcoma. Cell Rep. 9, 829–841. Storey, J.D., and Tibshirani, R. (2003). Statistical significance for genomewide studies. Proc Natl Acad Sci USA 100, 9440–9445. Stuart, J.M., Segal, E., Koller, D., and Kim, S.K. (2003). A gene-coexpression network for global discovery of conserved genetic modules. Science 302, 249–255. Subramanian, A., Tamayo, P., Mootha, V.K., Mukherjee, S., Ebert, B.L., Gillette, M.A., Paulovich, A., Pomeroy, S.L., Golub, T.R., Lander, E.S., and Mesirov, J.P. (2005). Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc. Natl. Acad. Sci. USA 102, 15545–15550. Tan, H., Wu, Z., Wang, H., Bai, B., Li, Y., Wang, X., Zhai, B., Beach, T.G., and Peng, J. (2015). Refined phosphopeptide enrichment by phosphate additive and the analysis of human brain phosphoproteome. Proteomics 15, 500–507. Tan, H., Yang, K., Li, Y., Shaw, T.I., Wang, Y., Blanco, D.B., Wang, X., Cho, J.H., Wang, H., Rankin, S., et al. (2017). Integrative proteomics and phosphoproteomics profiling reveals dynamic signaling networks and bioenergetics pathways underlying T cell activation. Immunity 46, 488–503. Tate, S.C., Cai, S., Ajamie, R.T., Burke, T., Beckmann, R.P., Chan, E.M., De Dios, A., Wishart, G.N., Gelbert, L.M., and Cronier, D.M. (2014). Semi-mechanistic pharmacokinetic/pharmacodynamic modeling of the antitumor activity of LY2835219, a new cyclin-dependent kinase 4/6 inhibitor, in mice bearing human tumor xenografts. Clin. Cancer Res. 20, 3763–3774.

16 Cancer Cell 34, 1–16, September 10, 2018

Vora, S.R., Juric, D., Kim, N., Mino-Kenudson, M., Huynh, T., Costa, C., Lockerman, E.L., Pollack, S.F., Liu, M., Li, X., et al. (2014). CDK 4/6 inhibitors sensitize PIK3CA mutant breast cancer to PI3K inhibitors. Cancer Cell 26, 136–149. Wang, X., Li, Y., Wu, Z., Wang, H., Tan, H., and Peng, J. (2014). JUMP: a tagbased database search tool for peptide identification with high sensitivity and accuracy. Mol. Cell. Proteomics 13, 3663–3673. Williamson, D., Missiaglia, E., de Reynies, A., Pierron, G., Thuille, B., Palenzuela, G., Thway, K., Orbach, D., Lae, M., Freneaux, P., et al. (2010). Fusion gene-negative alveolar rhabdomyosarcoma is clinically and molecularly indistinguishable from embryonal rhabdomyosarcoma. J. Clin. Oncol. 28, 2151–2158. Xi, Y., and Li, W. (2009). BSMAP whole genome bisulfite sequence MAPping program. BMC Bioinformatics 10, 232. Yamazaki, S., Shen, Z., Jiang, Y., Smith, B.J., and Vicini, P. (2013). Application of target-mediated drug disposition model to small molecule heat shock protein 90 inhibitors. Drug Metab. Dispos. 41, 1285–1294. Zhang, B., and Horvath, S. (2005). A general framework for weighted gene coexpression network analysis. Stat Appl Genet Mol Biol 4, https://doi.org/10. 2202/1544-6115.1128. Zhang, J., Benavente, C.A., McEvoy, J., Flores-Otero, J., Ding, L., Chen, X., Ulyanov, A., Wu, G., Wilson, M., Wang, J., et al. (2012). A novel retinoblastoma therapy from genomic and epigenetic analyses. Nature 481, 329–334. Zhou, X., Edmonson, M.N., Wilkinson, M.R., Patel, A., Wu, G., Liu, Y., Li, Y., Zhang, Z., Rusch, M.C., Parker, M., et al. (2016). Exploring genomic alteration in pediatric cancer using ProteinPaint. Nat. Genet. 48, 4–6.

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

STAR+METHODS KEY RESOURCES TABLE

REAGENT or RESOURCE

SOURCE

IDENTIFIER

H3K4me1 for ChIP

Abcam

Cat#ab8895; RRID: AB_306847

H3K4me1 for ChIP

Active Motif

Cat#39297; RRID: AB_2615075

H3K4me2 for ChIP

Active Motif

Cat#39141; RRID: AB_2614985

H3K4me2 for ChIP

Abcam

Cat#ab7766; RRID: AB_2560996

H3K4me3 for ChIP

Diagenode

Cat#C15410003-50; RRID: AB_2616052

H3K4me3 for ChIP

Active Motif

Cat#39159; RRID: AB_2615077

H3K27Ac for ChIP

Abcam

Cat#ab4729; RRID: AB_2118291

H3K27Ac for ChIP

Active Motif

Cat#39133; RRID: AB_2722569

H3K9/14Ac for ChIP

Life Technologies

Cat#49-1010; RRID: AB_2533861

H3K9/14Ac for ChIP

Diagenode

Cat#C15410200; RRID: AB_2637059

H3K36me3 for ChIP

Active Motif

Cat#61101; RRID: AB_2615073

H3K36me3 for ChIP

Abcam

Cat#ab9050; RRID: AB_306966

H3K27me3 for ChIP

Millipore

Cat#07-449; RRID: AB_310624

H3K27me3 for ChIP

Active Motif

Cat#39155; RRID: AB_2561020

H3K9me3 for ChIP

Invitrogen

Cat#491008; RRID: AB_11180319

H3K9me3 for ChIP

Active Motif

Cat#39161; RRID: AB_2532132

CTCF for ChIP

Active Motif

Cat#61311; RRID: AB_2614975

BRD4 for ChIP

Bethyl Laboratories

Cat#A301-985A100; RRID: AB_2620184

RNA pol II-phospho serine II

Abcam

Cat#ab5095; RRID: AB_304749

NOS1 for immunoblot

Cell Signaling

Cat#4234; RRID: AB_2152488

HMGA2 for immunoblot

Cell Signaling

Cat#5269S; RRID: AB_10694917

Phosphor-Ser235/236 of RPS6 for immunoblot

Cell Signaling

Cat#2211; RRID: AB_331679

CDC25 for immunoblot

Cell Signaling

Cat#4688S; RRID: AB_560956

Antibodies

PhosphoCDC25C (Ser216) for immunoblot

Cell Signaling

Cat#9528; RRID: AB_2075150

MYF5 for IHC

Abcam

Cat#ab125078; RRID: AB_10975611

Myogenin for IHC

Dako

Cat#M3559; RRID: AB_2250893

P16 for IHC

Roche

Cat#705-4713

pRB for IHC

Cell Signaling

Cat#9301S; RRID: AB_330013

CCND3 for IHC

Cell Signaling

Cat#2936; RRID: AB_2070801

CCND2 for IHC

Santa Cruz

Cat#sc-593; RRID: AB_2070794

CCND1 for IHC

Cell Marque

Cat#241R-18; RRID: AB_1158233

CDK6 for IHC

Santa Cruz

Cat#sc7961; RRID: AB_627242

CDK4 for IHC

Santa Cruz

Cat#sc260; RRID: AB_631219

HSP70 for Western Blot

Cell Signaling

Cat#4872; RRID: AB_2279841

Actin for Western Blot and immunoblot

Sigma

Cat#AC-15

Anti-rabbit for Western Blot and immunoblot

Licor

Cat#926-32221; RRID: AB_621841

Anti-mouse for Western Blot and immunoblot

Licor

Cat#926-32210; RRID: AB_621842

Superscript III System

Invitrogen

Cat#18080093

SYBR Select Master Mix for CFX

Applied Biosystems

Cat#4472942

Click-iT EdU imaging kit, Invitrogen

Invitrogen

Cat#C10083

St. Jude Children’s Research Hospital Childhood Solid Tumor Network

www.stjude.org/CSTN/

Biological Samples Orthotopic Patient Derived Xenografts (O-PDX)

(Continued on next page)

Cancer Cell 34, 1–16.e1–e19, September 10, 2018 e1

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Continued REAGENT or RESOURCE

SOURCE

IDENTIFIER

Sybr Green Select Master Mix

Applied Biosystems

Cat#4472908

PhosSTOP phosphatase inhibitor cocktail

Roche

Cat#PHOSS-RO

Sep-Pak C18 cartridge

Waters

Cat#WAT054945

Chemicals, Peptides, and Recombinant Proteins

XenoLight D-Luciferin potassium salt

PerkinElmer

Cat#122799

AZD1775 (purity 98%)

AbMole

Cat#M2143

Ganetespib (purity > 98%)

MedChem Express

Cat#HY-15205, Lot 08249

Palbociclib hydrochloride (purity 99%)

MedKoo

Lot SSC40317

Abemaciclib mesylate (purity > 98%)

AbMole

Cat#M2112, Lot 2

Trametinib free base (purity >99%)

AbMole

Lot M1759

DNeasy Blood & Tissue kit

Qiagen

Cat#69504

TruChIP Chromatin Shearing kit

Covaris

Cat#520154

iDeal Chip-seq kit

Diagenode

Cat#C01010051

MinElute PCR Purification Kit

Qiagen

Cat#28606

Quant-IT dsDNA HS assay kit

ThermoFisher Scientific

Cat#Q33120

NEBNext ChIP-Seq Library Prep Reagent Set

Illumina

Cat#E6200

RNeasy Plus Mini Kit

Qiagen

Cat#74134

TruSeq Stranded Total RNA Library Prep Kit

Illumina

Cat#20020596

BCA protein assay

ThermoFisher Scientific

Cat#23225

Sequencing data

This paper

European Molecular Biology Laboratory nucleotide sequence data library: EGAS00001002967

Genomic, transcriptomic, epigenetic, and proteomic data

This paper

https://pecan.stjude.org/proteinpaint/ study/RHB2018

Human reference genome NCBI build 37, GRCh37(hg19)

Genome Reference Consortium

https://www.ncbi.nlm.nih.gov/grc/human

Uniprot mouse and human protein database (Downloaded Feb 2015)

Uniprot

http://www.uniprot.org/downloads

Drug Sensitivity Database

Stewart et al., 2017, CSTN

https://braid.stjude.org/masttour/

Critical Commercial Assays

Deposited Data

Experimental Models: Cell Lines RD cell line

ATCC

ATCC, CCL-136

A549 cell line

ATCC

ATCC, CCL-185

SAOS-2 cell line

ATCC

ATCC, HTB-85

Experimental Models: Organisms/Strains NOD/SCID/gamma female mice

Jackson Laboratories

Strain code 005557

Athymic nude female mice

Charles River Laboratories

Strain code 553

CYB5R3 forward (used for H3K4me1): CTGACCTAATGCCCACAA

This paper

N/A

CYB5R3 reverse (used for H3K4me1): GACTCTGGCACCCTGTGTC

This paper

N/A

ZNF333 forward (used for H3K27me3 and H3K9me3): TGAAGACACATCTGCGAACC

This paper

N/A

ZNF333 reverse (used for H3K27me3 and H3K9me3): TCGCGCACTCATACAGTTTC

This paper

N/A

TRABD2B forward (used for CTCF): ACTGGTGGCAGCGACAGATT

This paper

N/A

Oligonucleotides

(Continued on next page)

e2 Cancer Cell 34, 1–16.e1–e19, September 10, 2018

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Continued REAGENT or RESOURCE

SOURCE

IDENTIFIER

TRABD2B reverse (used for CTCF): TCCTGTCC CTCTCGCTGTTC

This paper

N/A

MAG1 forward (used for H3K4me1 and H3K36me3): CGCACACGTGTTTTACTCC

This paper

N/A

MAG1 reverse (used for H3K4me1 and H3K36me3): AGGGCGGTTAGGAGGCTTAG

This paper

N/A

GAPDH-TSS (used for H3K4me3, H3K27me3, H3K9me3, CTCF)

Diagenode

Cat#C01010051

APDH-TSS (used for H3K36me3)

Diagenode

Cat#C01010051

Myoglobin Exon 2 (used for H3K4me3)

Diagenode

Cat#C01010051

ACTB-2 (used for H3K4me2, H3K9/14Ac, H3K27Ac)

Activ Motif

Cat#71005

GD Chr12 (used for H3K4me2, H3K9/14Ac, H3K27Ac)

Active Motif

Cat#71001

Li and Durbin, 2009

http://bio-bwa.sourceforge.net

Picard version 1.65(1160)

Broad Institute

https://broadinstitute.github.io/picard

samtools version 0.1.18 (r982:295)

Li et al., 2009

http://www.htslib.org/

IGV version 2.3.40

Robinson et al., 2011

http://software.broadinstitute.org/ software/igv

CRCMapper 2016-10-14

Aldiri, et al., 2017

https://bitbucket.org/young_computation/ crcmapper

Bis-SNP version 0.82.2

Liu et al., 2012

https://sourceforge.net/projects/bissnp

BSMAP version 2.74

Xi and Li, 2009

https://github.com/genestack-open/bsmap

ChromHMM version 1.10

Ernst and Kellis, 2012.

http://compbio.mit.edu/ChromHMM

R version 2.14.0

R Development Core Team, 2011

https://www.r-project.org

SPP version 1.11

Kharchenko et al., 2008

http://compbio.med.harvard.edu/ Supplements/ChIP-seq

caTools version 1.17

R package

https://cran.r-project.org/web/packages/ caTools/index.html

Bitops version 1.0-6

R package

https://cran.r-project.org/web/packages/ bitops/index.html

BRAID model for drug screening

Twarog et al., 2016

https://braid.stjude.org/masttour/

Living Image Version 4.3

Caliper Life Sciences/PerkinElmer

http://www.perkinelmer.com/lab-productsand-services/resources/in-vivo-imagingsoftware-downloads.html

drc R package

Ritz and Streibig, 2005, R package

https://cran.r-project.org/web/packages/ drc/index.html

javaGSEA application

Broad Institute

http://software.broadinstitute.org/gsea/ index.jsp

JUMP peptide identification software

Wang et al., 2014, Junmin Peng laboratory at St. Jude

http://www.stjuderesearch.org/site/ lab/peng

WGCNA R package

Zang and Horvath, 2005, R package

https://cran.r-project.org/web/packages/ WGCNA/index.html

Hallmark pathway database

Liberzon et al., 2011, Broad Institute MsigDB

http://software.broadinstitute.org/gsea/ msigdb/collection_details.jsp

LIMMA R package

Storey et al., 2003, Bioconductor R Package

https://bioconductor.org/packages/ release/bioc/html/limma.html

Antibody optimization and standard protocols for ChIP sequencing

St. Jude Children’s Research Hospital Childhood Solid Tumor Network

www.stjude.org/CSTN/

SkBM-2 cell culture medium with SingleQuots supplements

Lonza

Cat#CC-3246, CC-3244

Software and Algorithms BWA version 0.5.9-r26-dev

Other

Cancer Cell 34, 1–16.e1–e19, September 10, 2018 e3

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

CONTACT FOR REAGENT AND RESOURCE SHARING Further information and requests for resources and reagents should be directed to and will be fulfilled by the Lead Contact, Michael A. Dyer ([email protected]). EXPERIMENTAL MODEL AND SUBJECT DETAILS Animals and Orthotopic Patient-Derived Xenografts This study was carried out in strict accordance with the recommendations in the Guide to Care and Use of Laboratory Animals of the National Institute of Health. The protocol was approved by the Institutional Animal Care and Use Committee at St. Jude Children’s Research Hospital. All efforts were made to minimize suffering. All mice were housed in accordance with approved IACUC protocols. Animals were housed on a 12-12 light cycle (light on 6 am and off 6 pm) and provided food and water ad libitum. Athymic nude female mice were purchased from Charles River Laboratories (strain code 553). NOD/SCID/gamma female mice were purchased from Jackson Laboratories (strain code 005557). Human O-PDXs were generated and expanded as described previously (Stewart et al., 2017). All O-PDX samples described here are freely available through the CSTN (http://www.stjude.org/CSTN/) with no obligation to collaborate. METHOD DETAILS Genomic DNA Extraction DNA was extracted from 15-20 mg of tissue using DNeasy Blood & Tissue kit following manufacturer instructions (Qiagen 69504). Chromatin Immunoprecipitation and Sequencing (ChIP-Seq) Cell Fixation O-PDX tumors in duplicate were used for ChIP-seq. Each xenograft tumor (2-3 grams) was harvested directly from the injection site in the hind leg of an immunocompromised mouse. Tissue was harvested immediately after euthanasia. Only unilateral tumors were used for biological replicates. Within 30 min after the tumor was removed, the tissue was finely minced (1-3 mm3) at room temperature in phosphate-buffered saline without calcium and magnesium (PBS) (Lonza 17-516F) using a sterile surgical blade or razor. The tissue was fixed with fresh 1% formaldehyde (Thermo scientific 28906) with gentle mixing for no more than 10 min, including the time for centrifugation. After 8 min of fixation, the tissue was centrifuged for 2 min at 1,500 RPM at 4 C. The formaldehyde was removed and the tissue was washed with freshly made 1X glycine (0.125 M) for 5 min and centrifuged again. The tissue was washed one time with cold PBS and resuspended with protease inhibitors (Diagenode C12010011) diluted 1:200 in cold PBS. The tissue was placed in a dounce tissue grinder and dounced 20-40 times on ice to break up the tissue into a single cell suspension. Cells were strained through a 40 mm mesh filter into an ice-cold 50 ml conical tube. Cells were flash frozen in 10x106 aliquots and stored at -80 C. Nuclei Extraction and Shearing Nuclei extraction and chromatin shearing was performed using the TruChIP Chromatin Shearing kit (Covaris 520154) following manufacturer instructions. Briefly, nuclei from 10x106 cells was sheared in 0.1% SDS shearing buffer for 12 min using the Covaris E-series focused ultra-sonicator (target base pair = 350 bp, duty = 5%, intensity = 4, 200 cycles per burst). Chromatin shearing efficiency was analyzed using 20 ml of chromatin treated with 0.4 mg/ml proteinase K and placed at 65 C for 20 min. Sheared chromatin was then run on a 1% agarose gel and Agilent 2200 TapeStation to check for optimal fragment size. Antibody Optimization, ChIP Procedure, and Real-Time PCR Antibody optimization and standard operating procedures for CTCF, H3K4me1, H3K4me2, H3K4me3, H3K27AC, H3K9/14AC, H3K36me3, H3K27me3, H3K9me3 are available through the Childhood Solid Tumor Network (www.stjude.org/CSTN). ChIP and library prep for RNA pol ll and BRD4 were performed by Active Motif for all O-PDX tumors and control cells. For the remaining antibodies, the ChIP procedure was performed using the iDeal Chip-seq kit (Diagenode C01010051) following the manufacturer instructions. Briefly, chromatin from 1x106 cells was incubated with antibody and protein A- coated magnetic beads overnight and rotating at 4 C. Beads were washed 4 times at 4 C then resuspended in Diagenode elution buffer with proteinase K (0.4 mg/ml final concentration) and then incubated overnight at 65 C while shaking at 900 RPM. DNA was eluted with 20 ml of RNAse/DNAse-free ChIP-grade water using MinElute PCR Purification Kit (Qiagen 28606). The DNA fragment size was analyzed using the Agilent 2200 TapeStation and quantified using Quant-IT dsDNA HS assay kit (Thermo Fisher Q33120). Briefly, 2 ml of ChIP DNA was tested for enrichment for positive and negative control target genes for each antibody using Sybr Green Select Master Mix (Applied Biosystems 4472908) and quantitative real-time PCR. Primers sequences sourced from St. Jude for each target gene are listed below. Commercially sourced primer information is available in the Key Resources Table. Cycle conditions: 5 min at 95 C; 40 cycles of: 30 s at 95 C, 30 s at 60 C, 30 s at 72 C; Add in a melting curve step.

e4 Cancer Cell 34, 1–16.e1–e19, September 10, 2018

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Positive Control Target Genes Protein

Target gene or nearest gene

Forward Primer

Reverse Primer

H3K4me1

CYB5R3

CTGACCTAATGCCCACAA

GACTCTGGCACCCTGTGTC

H3K27me3

ZNF333

TGAAGACACATCTGCGAACC

TCGCGCACTCATACAGTTTC

H3K9me3

ZNF333

TGAAGACACATCTGCGAACC

TCGCGCACTCATACAGTTTC

CTCF

TRABD2B

ACTGGTGGCAGCGACAGATT

TCCTGTCCCTCTCGCTGTTC

Negative Control Target Genes Protein

Target gene or nearest gene

Forward Primer

Reverse Primer

H3K4me1

MAG1

CGCACACGTGTTTTACTCC

AGGGCGGTTAGGAGGCTTAG

H3K36me3

MAG1

CGCACACGTGTTTTACTCC

AGGGCGGTTAGGAGGCTTAG

Library Prep and Sequencing ChIPSeq libraries were prepared using the NEBNext ChIP-Seq Library Prep Reagent Set for Illumina (NEB) per manufacturer’s instructions with the following modifications: The amount of immunoprecipitated chromatin ranged from 1 to 10 ng depending on the sample. The PCR amplification of the adaptor ligated DNA varied between 10-15 cycles, with 1 ng of chromatin input being amplified for 15 cycles, and 10 ng of input amplified for 10 cycles. A single ampure bead cleanup step was performed using a 0.8X bead ratio. Secondary size selection was not performed. Data was generated on a HiSeq 2000 using the single read 50 cycle protocol. Whole Genome Sequencing and Analysis Whole genome sequencing and library construction was performed as described previously (Chen et al., 2013) with the following modifications; 250-500 ng of genomic DNA was input for library construction using Illumina compatible adapters, and 4-6 cycles of amplification was performed with Kapa HiFi Hotstart ReadyMix (KAPA Biosystems). Identification of SNV, CNVs, SVs and Indels was performed as described previously (Chen et al., 2013). RNA Sequencing and Analysis RNA was isolated from individual TRIzol (Life Technologies) preparations via a phenol-chloroform extraction. Samples were first homogenized at 17,000 RPM for 30 s with a tissue homogenizer (Polytron, PT10-35GT). A 1:4 volume of chloroform (Sigma Aldrich) was then added to each sample and incubated at room temperature for 3 min, followed by centrifuging at 12,000 g at 4 C for 15 min. The aqueous layer was then transferred to a siliconized 1.5 mL tube, followed by the addition of 2.0 mL glycogen (Roche Applied Science) and 500 mL isopropanol (Fisher Scientific). Samples were incubated at room temperature for 10 min and centrifuged at 12,000 g at 4 C for 15 min. The samples were then washed twice with ice-cold 80% ethanol (Fisher Scientific) to remove salts, resuspended in water containing diethyl pyrocarbonate, and quantified with a NanoDrop spectrophotometer (Thermo Fisher Scientific). RNA was extracted from resected tumors, O-PDXs, primary human myoblasts and myotubes, and fetal skeletal muscle with an RNeasy Plus Mini Kit (Qiagen). Libraries were prepared from 500 ng RNA with the TruSeq Stranded Total RNA Library Prep Kit, according to manufacturer directions (Illumina). Paired-end, 100-cycle sequencing was performed with HiSeq 2000 or HiSeq 2500 sequencers, according to manufacturer directions (Illumina). RNA-Seq was performed using the TruSeq Stranded Total RNA Library Prep Kit (Illumina) with 250 ng-1 mg of total RNA as input. Gene based FPKM quantification in RNA-Seq was performed as described previously (Chen et al., 2013). DNA Methylation and Analysis Whole genome bisulfite sequencing (WGBS) data were mapped to GRCh37-lite (human) or by BSMAP2.74 with the following parameters: -z 33 -f 5 -g 3 -r 0 -m 17 -x 600 –u and raw methylation status was called by Bis-SNP (Version 0.82.2). We developed a regres
#Methylated C + 0:5 . sion tree based approach to identify regional methylation pattern of a sample using m-values m = log2 #Unmethylated C + 0:5 Adjacent segments with insignificant (p > 1E-6, Student t-test) or small difference between average b-values 
#Methylated C (< 0.1) were recursively merged. b= #Methylated C + #Unmethylated C Segmentation breakpoints from individual samples were pooled to generate an initial list of regions for testing correlation between methylation and gene expression. Correlation between regions with variable methylation (difference between highest and lowest average m-values >= 0.3) and adjacent genes with variable expression (methylation region overlapping 10 kb upstream of TSS to TES, at least 2 fold difference between highest and lowest FPKM, highest FPKM >= 1) were calculated. Gene/region pairs with FDR <= 0.5 were retained for merging, where adjacent regions for the same gene (no more than 1 kb apart) with the same orientation (positive or negative) of correlation were merged. Correlations were re-evaluated and pairs with final FDR <= 0.1 were used as differentially methylated regions that were significantly correlated with gene expression.

Cancer Cell 34, 1–16.e1–e19, September 10, 2018 e5

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Proteomics and Phosphoproteomics Quantitative Analysis of Proteome by Multiplexed TMT-LC/LC-MS/MS Isobaric labeling, such as iTRAQ and TMT, is emerging as a powerful strategy for deep multiplexed proteomics analysis with high throughput and reproducibility (Bai et al., 2017). One limitation associated with this method is that target peptide ions are often co-isolated with other co-eluted peptide ions during LC-MS/MS analysis which causes high noise level to compress quantitative ratios and decreases measurement accuracy. Fortunately, the ratio compression effect also alleviates experimental variations, and hence has only minor impact on protein differential expression analyses (Tan et al., 2017). Moreover, the ratio compression can be diminished by extensive peptides fractionation, narrow isolation window, and post-MS correction (Niu et al., 2017). Alternatively, the MS3 method can almost eliminate this ratio compression, but it requires longer duty cycles, specific MS settings, and the use of low resolution MS2 for identification, which often compromise the peptide/protein identification (McAlister et al., 2014). To balance the pros and cons associated with isobaric labeling, we implemented extensive fractionation through long gradient high resolution LC/LC-MS/MS to achieve deep proteome coverage and to reduce ion compression during quantification. In addition, we employed sample replicates to facilitate statistical inference of differentially expressed proteins, a widely adopted strategy used in proteomics analyses (Bai et al., 2017; Tan et al., 2017). RMS O-PDX Tumor Tissue, Myoblast and Myotube Cells for Proteome and Phosphoproteome Profiling 108 myoblast and myotube cells per sample were collected, washed twice with 10 ml ice cold PBS, and then snap frozen in liquid nitrogen. These processes were managed to be completed in 5 min. O-PDX engrafted mice were anesthetized and perfused with 10 ml PBS before sacrifice, center section of O-PDX tumor tissues were collected, homogenized, and snap frozen in liquid nitrogen. These steps were finished within 10 min. Protein Extraction, Digestion, Labeling and Pooling Protein extraction, digestion, labeling and pooling were performed similarly as previously described (Bai et al., 2017; Tan et al., 2017). RMS O-PDX tissues, myoblast and myotube cells were lysed in freshly prepared lysis buffer (50 mM HEPES, pH 8.5, 8 M urea, 0.5% sodium deoxycholate and phosphatase inhibitor cocktail (PhosSTOP, Roche)). Protein concentration of sample lysates were quantified by BCA protein assay (Thermo Fisher Scientific) with titrated BSA as a standard. 1 mg proteins per sample were first digested with Lys-C (Wako, 1:100 w/w) at room temperature for 2 hr, diluted 4 times with 50 mM HEPES, pH 8.5, and then further digested with trypsin (Promega, 1:100 w/w) overnight at room temperature. 1% trifluoroacetic acid was added to quench the digestion reaction, followed by desalting with Sep-Pak C18 cartridge (Waters), and the desalted peptides were dried by speedvac. Samples were then resuspended in 50 mM HEPES, pH 8.5, and were labeled with 10-plex TMT reagents following the manufacturer’s instruction. Lastly, 10 isobaric labeled samples were pooled together with equal amount, desalted again by Sep-Pak C18 cartridge and then speedvac dried. Offline Basic pH Reverse Phase Liquid Chromatography The basic pH reverse phase liquid chromatography peptides pre-fractionation were performed on Agilent 1220 LC system as previously introduced (Bai et al., 2017; Tan et al., 2017). The pooled TMT labeled sample was solubilized in buffer A (10 mM ammonium formate, pH 8) and separated on two XBridge C18 columns (3.5 mm particle size, 4.6 mm 3 25 cm, Waters) into around 180 fractions with a 220 min long gradient started from 15% to 65% buffer B (95% acetonitrile, 10 mM ammonium formate, pH 8, flow rate: 0.4 ml/min). 5% of each combined fraction was dried for whole proteome analysis and the remaining 95% was dried by speedvac for phosphoproteome analysis. Refined Phosphopeptide Enrichment by TiO2 Phosphopeptide enrichment was performed following the refined protocol as previously described (Tan et al., 2015). Briefly, peptides were added to clean TiO2 beads (GL sciences) with a peptide-to-beads weight ratio of 1:4 in binding buffer (65% acetonitrile, 2% TFA, and 0.5 mM KH2PO4) and incubated for 20 min. Enriched phosphopeptides were washed, eluted, dried, and dissolved in 5% formic acid for LC-MS/MS analysis. Long Gradient Acidic pH Reverse Phase LC-MS/MS The analysis was carried out based on our optimized platform as previously described (Bai et al., 2017; Tan et al., 2017). The dried peptide fractions were reconstituted in loading buffer (5% formic acid), loaded on a reverse phase column (75 mm 3 50 cm, 1.9 mm C18 resin (Dr. Maisch GmbH, Germany)) and interfaced with a FUSION or Q Exactive HF mass spectrometer (Thermo Fisher Scientific). Peptides were eluted up to 6 hr by a 15-65% gradient of buffer B. (buffer A: 0.2% formic acid, 5% DMSO; buffer B: buffer A plus 65% acetonitrile, flow rate: 0.25 ml/min). A butterfly portfolio heater (Phoenix S&T) was applied to heat the column at 65 C to reduce backpressure. The mass spectrometer was operated in data-dependent mode with MS1 settings of 60,000 resolution, 1 3 106 AGC target and 50 ms maximal ion time and top 20 MS/MS high resolution scans with MS2 settings of 1 m/z isolation window with offset 0.2, 60,000 resolution, 100 ms maximal ion time, 1 3 105 AGC target, HCD, 33 normalized collision energy, and 40 s dynamic exclusion (35 normalized collision energy and 20 s dynamic exclusion for phosphoproteome). Peptides Identification by JUMP, a Tag-Based Hybrid Search Engine Peptide identification was performed using our recently developed JUMP search engine with improved sensitivity and specificity (Wang et al., 2014). Commercially available database search engines can be divided into two categories: tag-based De novo sequencing (e.g. PEAKS with limited sensitivity) and pattern-based database search (e.g. SEQUEST, MASCOT). The JUMP software integrates these two methods to score putative peptides, showing significant improvements compared with these commercially available tools (Wang et al., 2014). The JUMP software has previously been utilized in numerous publications (Gong et al., 2017). Analysis was done similarly as previously described (Tan et al., 2017), MS/MS raw files were first converted into mzXML format e6 Cancer Cell 34, 1–16.e1–e19, September 10, 2018

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

and searched against a composite target/decoy database FDR estimation. The target protein database was compiled from the Uniprot mouse and human database (Human database: 88,965 protein entries; Mouse database: 52,738 protein entries, downloaded in February 2015), the decoy database was generated by reversing target protein sequences. Spectra were searched with ± 10 ppm mass tolerance for both precursor ions and product ions with fully tryptic restriction, static modification for TMT tag on N-terminus and lysine (+229.16293), dynamic modification for serine, threonine and tyrosine (+79.96633, for phosphoproteome analysis), three maximal modification sites, two maximal missed cleavages, and the assignments of a, b, and y ions. Peptide spectrum matches (PSM) were first filtered by MS mass accuracy (2 ppm, ± 4 standard deviations). PSMs of doubly charged peptides with JUMP Jscore of > 30 were applied for global mass recalibration prior to the filtering. The survived PSMs were first grouped by precursor ion charge state and then further filtered by Jscore and dJn values. Cutoffs were applied on these values and were adjusted until a protein FDR < 1% was achieved. If one peptide was shared by multiple proteins, the protein with the highest PSM will represent the peptide according to the rule of parsimony (Nesvizhskii and Aebersold, 2005). Phosphosite Assignment by the Lscore from the JUMP Software Suite To determine the reliability of phosphosite localization on peptides, we adopted the concept of the phosphoRS algorithm (Taus et al., 2011) to compute phosphosite localization scores (Lscore, 0-100%) in each PSM same as previously described (Tan et al., 2017). Phosphosites were aligned to protein sequences to generate protein level Lscores in addition to the PSM Lscores. The protein Lscore was represented by the highest PSM Lscore if multiple PSMs were identified for one specific phosphosite. Since random assignments of PSMs containing ambiguous phosphosites often causes excessively high number of unreliable phosphosites on proteins, we implemented series rules to alleviate the problem: (i) If the gap of PSM Lscores between the 1st and 2nd site > 10% for a singly phosphorylated peptide in one PSM, the top Lscore site was selected; (ii) Otherwise, we inspect the phosphosites in the corresponding proteins instead to select the sites with the highest protein Lscore. This allows low quality PSMs to borrow information inferred from the high quality PSMs; (iii) If both PSM and protein level Lscores were indistinguishable, a heuristic priority was assigned to phosphosites according to the order of occurrence: SP-motif, S, T and Y; (iv) If the PSMs did not satisfy any of the rules above, these PSMs were first sorted by JUMP Jscores, and then we selected protein phosphosites that had been determined by other PSMs of high Jscores. TMT-Based Protein and Phosphosite Quantification Using the JUMP Software Suite This analysis was carried out in the following steps similarly as previously reported (Tan et al., 2017). (i) TMT reporter ion intensities of each PSM were extracted; (ii) the raw intensities were corrected according to isotopic distribution of each labeling reagent (Tan et al., 2017); (iii) PSMs with very low reporter ion intensities were excluded (e.g. minimum intensity < 1,000 and median intensity < 5,000); (iv) sample loading bias was corrected by normalization with the trimmed median intensity of all PSMs; (v) the mean-centered intensities across samples were calculated (Tan et al., 2017); (vi) protein or phosphosite relative intensities were summarized by averaging related PSMs; (vii) protein or phosphosite absolute intensities were derived by multiplying the relative intensities by the grand-mean intensity of the top three most highly abundant PSMs. To generate a combined quantification table for multiple batches, a common sample (SJRHB10468_X1) was included in each batch as an internal standard. MS intensities from different batches were normalized according to this internal standard. Differential Expression Analyses of Proteome and Phosphoproteome Differential expression analyses of whole proteome and phosphoproteome were carried out using LIMMA R package, a software designed for the analysis of gene expression involving comparisons between many gene targets simultaneously. LIMMA borrows information across genes by fitting linear models to overcome the problem of small sample size and complex experimental design, hence it is ideal for differential expression analysis of TMT-based deep proteomic data. Briefly, (i) Linear models were fitted for expression data of each protein or phosphosite; (ii) Empirical Bayes method was used to borrow information across genes; (iii) p values were adjusted by the Benjamin Hochberg method; (iv) The adjusted p value cutoff of 0.05 was then applied; (v) Remaining proteins were further filtered by a fold change of 1.5 and 2.0 (equivalent to around 3 x SD of tumor biological replicates) in at least one group comparison for proteome and phosphoproteome respectively. Weighted Gene Co-expression Network Analysis (WGCNA) and Pathway Annotation The analysis was done by the WGCNA R package similarly as previously described (Tan et al., 2017). Only DE proteins and phosphosites were applied to define the proteome and phosphoproteome co-expression clusters (i.e. WPCs and PPCs) respectively. (i) Pearson correlation matrix were generated by calculating correlation between proteins (only positive correlations were considered), and were further raised to a power of 16 using the scale free topology principle to calculate an adjacency matrix.(ii) Co-expression clusters were then determined by the hybrid dynamic tree-cutting method (Langfelder et al., 2008) with a height cutoff (e.g. 0.2) for merging modules; (iii) A consensus trend (eigengene) was calculated based on the first principal component for each co-expression cluster. Proteins were then assigned to the co-expression cluster with the highest correlation; (iv) Each co-expression cluster was then annotated using the Hallmark pathway database downloaded from MsigDB (Liberzon et al., 2011) (Myogenic regulatory pathways that were not annotated in Hallmark database were manually extracted from KEGG database and added in) by Fisher’s exact test. Pathways with a B.H. adjusted p value less than 0.05 were selected as deregulated pathways in each co-expression cluster. Immunoblot Validation of Proteomics Xenograft tissue was lysed in RIPA lysis buffer (Sigma, catalog R0278) containing protease inhibitors (Halt Inhibitor Cocktail, catalog 87785). Proteins were resolved in 4-15% SDS-PAGE gel and transferred to PVDF membrane (Millipore). Membranes were blocked in Odyssey blocking buffer (Licor, catalog 927-4000) and then probed for the following primary antibodies: NOS1 (1:1000, Cell Cancer Cell 34, 1–16.e1–e19, September 10, 2018 e7

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Signaling, catalog 4234), HMGA2 (1:1000, Cell Signaling, catalog 5269), phosphor-Ser235/236 of RPS6 (1:1000, Cell Signaling, catalog 2211), CDC25 (1:1000, Cell Signaling, catalog 4688), Phospho CDC25C (Ser216) (1:1000, Cell Signaling, catalog 9528), and actin (1:5000, Sigma, AC-15. Secondary antibodies were anti-rabbit (Licor, catalog 926-32221); anti-mouse (Licor, catalog 926-32210). Analysis of Transcriptional Targets Transcriptional targets were curated from the literature and the list of individual targets is provide in the Supplemental Information. The FPKM values for each target in each sample were then averaged per sample and normalized to myoblasts. To be included in the analysis, the FPKM had to be above 1.0 for at least one of the normal samples. Immunohistochemistry Histopathologic features of the patient H&E-stained slides and the O-PDX derived slides were compared for similarity by a number of histologic features. These patient-xenograft pairs were evaluated via the following attributes: overall tumor cellularity, growth pattern, cytomorphology including degree of pleomorphism, and mitotic activity. Formalin-fixed paraffin embedded tissue (4 mm in thickness) was utilized for staining with appropriate positive and negative controls evaluated. Information regarding the antibodies used can be found below: All of the patient samples and the xenograft tissues were assessed for MYF5, CDK4, CDK6, CCND1, CCND2, CCND3, pRB, p16, and Myogenin by immunohistochemistry. In total, 34 hematoxylin and eosin stained slides and 306 immunohistochemical stained slides were reviewed. Antibody information can be found below and in the Key Resources Table.

Antibodies Utilized in Immunohistochemical Analysis Antibody

Validation

Stainer

Host

Dilution

Retrieval

Amplification/Blocking Reagent

Vendor/ Cat #

MYF5

RUO

Ultra

RM

1:25

CC2-92 min

None

abcam; ab125078

Myogenin

IVD

Ultra

MP

1:500

CC1-64 min

None

Dako; M3559

P16

IVD

Ultra

MM

RTU

CC1-48 min

None

Roche; 705-4713

pRB

RUO

Ultra

RP

1:100

CC1-64 min

None

Cell Signaling; 9301S

CCND3

RUO

Omnis

MM

1:400

High pH TR 30 min

Mouse Linker - 10 min

Cell signaling; 2936

CCND2

RUO

Bond-III

RP

1:100

ER1-20 min

None

Santa Cruz; sc-593

CCND1

IVD

Ultra

RM

RTU

CC1-36 min

None

Cell Marque: 241R-18

CDK6

RUO

Ultra

MM

1:25

CC1-64 min

None

Santa Cruz; sc7961

CDK4

RUO

Ultra

MM

1:50

CC1-64 min

AB Block(Roche)

Santa Cruz; sc260

IVD, In Vitro Diagnostic; RUO, Research use only; Bond-III, Leica; Omnis, Dako; Ultra, Ventana BenchMark ULTRA (Roche); MP, mouse polyclonal; RM, rabbit monoclonal; MM, mouse monoclonal; RP, rabbit polyclonal; CC1, cell conditioning solution 1; TR, thermal retrieval; RTU, Ready-to-Use.

The antibodies were evaluated for intensity of staining from non-reactive to 3+. The percentage of tumor cells with the expected localization of staining was recorded. This was restricted to nuclear staining for MYF5, CCND1, CCND2, CCND3, pRB, p16, and Myogenin. Either cytoplasmic or nuclear staining was acceptable for CDK4 and CDK6. Discordance between the patient and xenograft samples was evaluated based on expected staining localization. Discordance was noted when either there was a complete discrepancy (one sample was positive and the other was negative or had < 5% staining), or when both samples were positive but there was a significant difference (> 50%) between the samples. If staining failed completely the comparison was considered not applicable. General Drug Screening Drug screening data for solid tumor xenografts and cell lines was performed and obtained from the Childhood Solid Tumor Network (www.stjude.org/CSTN) (Stewart et al., 2017). All xenograft lines are assayed for optimal growth conditions, including plating density, DMSO-sensitivity, and positive control selection as previously described (Stewart et al., 2017) before being validated and assayed for drug sensitivities. All RMS xenograft lines are plated in vitro on gelatinized plates in SkBM-2 cell culture medium (Lonza Cat#CC3246) and supplemented with SingleQuots supplements (Lonza Cat#CC-3244). After the assay conditions have been validated, the xenograft cell lines may be screened against compounds for activity. Disaggregated xenograft cells are plated at previously determined optimal density in 384-well plates. Twenty-four hours post-plating, assay plates are drugged in triplicate with compounds of interest. For dose-response drugging, each compound is drugged against multiple wells per plate at different concentrations per well. Seventy-two hours post-drugging, RLUs are determined by the CellTiter-Glo(CTG)/EnVision system. Data is then analyzed to determine active compounds in the case of single-point drugging as well as dose-response curves.

e8 Cancer Cell 34, 1–16.e1–e19, September 10, 2018

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Combination Drug Screening with Ganetespib Raw Data Processing Raw luminescence RLU (relative light unit) values for each compound at each concentration were pooled from replicate experiments, log2 transformed, and normalized to obtain % activity using the following equation: 100* [(negctrl – compound) / (negctrl – posctrl)]. Dose-Response Analysis Based on % Activity Dose-response curves were fit using the drc (Ritz and Streibig, 2005) package in R. The four parameters sigmoidal function LL2.4 was used with the following constraints: -10 % hill slope % 0; 0 % y0 % max median normalized % activity; 0 % yFin % max median normalized % activity; 10-10 % ec50 % 10-4 (which equates to the range of drug concentrations tested in this experiment). Values from the highest concentration tested for each compound were weighted at 10% to reduce curve fitting artifacts. In the event of a failure to fit a sigmoidal dose response curve, the smooth.spline option in R was used to fit a curve that could be used to determine area under the curve. AUC Analysis of Drugs Combined with Ganetespib AUC values were obtained for each drug combined with DMSO and 4 concentrations of Ganetespib (3 nM, 10 nM, 33 nM, and 100 nM). The delta AUC value is the difference in AUC at the current concentration of Ganetespib and AUC with DMSO only. Higher value means greater inhibition of proliferation in the presence of Ganetespib. Response Surface Model Analysis of Select Drug Combinations Seven drugs were studied in combination using the BRAID model of combined drug action (Twarog et al., 2016). RD cells were plated at 1,000 cells per well in 96 well plates and read using CTG with envision after 8 days as previously described (Stewart et al., 2017). Bioreplicates were performed on separate days. Raw luminescence RLU (relative light unit) values for each compound at each concentration were pooled from replicate experiments, log10 transformed, and normalized by subtracting the value of the DMSO control wells. The BRAID kappa parameter for a drug combination that produces a 95% confidence interval was used to indicate the type of interaction: If kappa_l and kappa_h bound 0 = additivity; If kappa_l and kappa_h are both positive = synergy; If kappa_l and kappa_h are both negative = antagonism. The Index of Achievable Efficacy at 90% efficacy (IAE90) is comparable to a two-dimensional therapeutic window and enables direct comparison of the efficacy of drug combinations: higher IAE90 equates to greater combined efficacy (Twarog et al., 2016). Pharmacokinetic Studies Palbociclib The total plasma and tumor PK of palbociclib in female athymic nude mice (Charles River Laboratories, aged 12 -16 weeks) was assessed after a single oral gavage of 16 mg/kg of palbociclib free base equivalents. Palbociclib hydrochloride (MedKoo, Lot SSC40317, purity 99%) was suspended in 1% methylcellulose (MC, type 400 cPs) and 1% Tween 80 (1% MCT) in ultrapure water at a final nominal concentration of 1.6 mg/mL free base equivalents for a 10 mL/kg gavage volume. Mice were sacrificed using an IACUC-approved method at 10 min, 1, 4, 8, and 24 hr post-dose, with 3 mice per time point. Whole blood by cardiac puncture was collected into plastic microcentrifuge tubes containing 10 mL of sodium heparin, vortexed to mix anticoagulant, immediately centrifuged to plasma, and stored on dry ice for remainder of study. Mice were then perfused with calcium- and magnesium-free PBS via the right ventricle, the orthotopic xenografts excised, rinsed with PBS, weighed, divided into aliquots and placed on dry ice. At the end of the in vivo procedures, all samples were transferred from dry ice and placed at -80 C until analysis. Total plasma and tumor homogenate palbociclib concentrations were assessed using a sensitive and specific liquid chromatography, tandem mass spectrometry assay. First, tumor samples were diluted with a 5:1 volume of ultrapure water, and homogenized with a bead-based technique (Liang et al., 2011) on a FastPrep-24 system (MP Biomedicals, Santa Ana, CA). 1.4 mm ceramic spheres (MP Biomedicals, Lysing matrix D, 10 mg per mg of tumor) were added to the microcentrifuge tubes containing samples. The samples were then subjected to three 6.0 M/S vibratory cycles of 1 min each on the FastPrep-24 system. To prevent over-heating due to friction, samples were placed on wet ice for 5 min between each cycle. The homogenates were then stored at -80 C until analysis. Palbociclib hydrochloride (MedKoo, Lot SSC40317, purity 99%) stock solutions were prepared in acetonitrile and used to spike matrix calibrators and quality controls. Plasma and tumor homogenate samples, 25 mL each, were protein precipitated with 75 mL of 10 ng/mL LEE011 (Cayman Chemical Co., Lot 0467691-2, Purity R 95%) in acetonitrile as an internal standard. A 2 mL aliquot of the extracted supernatant was injected onto a Shimadzu LC-20ADXR high performance liquid chromatography system via a LEAP CTC PAL autosampler. The LC separation was performed using a Phenomenex Gemini C6 Phenyl (3.0 mm, 30 mm x 2.0 mm) column maintained at 60 C with gradient elution at a flow rate of 0.50 mL/min. The binary mobile phase consisted of water-200 mM ammonium acetate pH 6.0 (90:10 v/v) in reservoir A and acetonitrile-water-200 mM ammonium acetate pH 6.0 (90:10:10 v/v) in reservoir B. The initial mobile phase consisted of 25% B with a linear increase to 65% B in three min. The column was then rinsed for two min at 100% B and then equilibrated at the initial conditions for two min for a total run time of seven min. Under these conditions, the analyte and IS eluted at 1.36 and 1.04 min, respectively. Analyte and IS were detected with tandem mass spectrometry using a SCIEX API 5500 Q-TRAP in the positive ESI mode with the following mass were transitions monitored: Palbociclib 448.25 -> 380.10, LEE011 435.26 -> 322.20.

Cancer Cell 34, 1–16.e1–e19, September 10, 2018 e9

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

The method qualification and bioanalytical runs all passed P-PKSR’s acceptance criteria for non-GLP assay performance. A quadratic model (1/X2 weighting) fit the calibrators across the 1 to 500 ng/mL range, with a correlation coefficient (R) of R 0.9951. The lower limit of quantitation (LLOQ), defined as a peak area signal-to-noise ratio of 5 or greater verses a matrix blank with IS, was 1 ng/mL. The intra-run precision and accuracy was < 14.9% CV and 103% to 114%, respectively. The resultant palbociclib concentration-time (Ct) data were grouped by matrix and time point, and manual imputation of data below the lower limit of quantitation (BLOQ) was as follows: IF at any time point R 2/3rds of the Ct results were above the LLOQ, the BLOQ data were replaced with a value of ½ LLOQ, ELSE the entire time point’s data were treated as missing. Then, using Phoenix WinNonlin 6.4 (Certara USA, Inc., Princeton, NJ), Ct data summary statistics (arithmetic mean, standard deviation, % CV, minimum, median, maximum) were generated, and the palbociclib arithmetic mean Ct data for each matrix was subjected to noncompartmental pharmacokinetic analysis (NCA). The extravascular model (Model 202) was applied, and area under the Ct curve (AUC) values were estimated using the ‘‘linear up log down’’ trapezoidal rule. The terminal phase was defined as the three time points at the end of the Ct profile, and the elimination rate constant (Ke) was estimated using an unweighted log-linear regression of the terminal phase. The terminal elimination half-life (T1/2) was estimated as 0.693/Ke, and the AUC from time 0 to infinity (AUCinf) was estimated as the AUC to the last time point (AUClast) + predicted Clast/Ke. Other NCA parameters estimated included observed maximum concentration (Cmax), time of Cmax (Tmax), concentration at the last observed time point (Clast), time of Clast (Tlast), apparent clearance (CL/F = Dose/AUCinf), and apparent terminal volume of distribution (Vz/F). The average concentration over a dosing interval (Cavg) was estimated as AUCinf / dosing interval in hours. The apparent partition coefficient of palbociclib from the plasma to the tissue of interest (Kp,tissue) was estimated as the ratio of the AUCinf, tissue to AUCinf plasma when available. To estimate a clinically relevant mouse dosage, the resultant mouse plasma AUCinf was compared with the reported human plasma PK values at the putative single agent palbociclib approved dosage of 125 mg PO QD. All inferences were made under the assumption of time-independent, linear and dose-proportional PK in mice and humans. Palbociclib displayed moderate within-study PK variability considering the serial sacrifice design. Concentration coefficients of variation (CV) ranged from 13.5% to 58.4% for plasma, and from 25.3% to 68.5% for tumor. Most variability occurred at the first sampling time point of 10 min, which was during the absorption phase. The apparent plasma clearance (CL/F) estimate of 5.47 L/hr/kg was similar to that derived from previously published literature (4.0 to 5.9 L/hr/kg) (Parrish et al., 2015). The apparent terminal half-life (T1/2) was slightly shorter compared with literature (2.61 hr vs. approximately 3.4 hr). Palbociclib distributed well into the orthotopic RMS xenograft, achieving concentrations approximately 6–fold higher than the plasma. The T1/2 in the tumor was 1.8-fold longer than the plasma. NCA parameter estimates and the mean (SD) palbociclib Ct plots, grouped by matrix are shown below: The common FDA-approved dose for palbociclib is 125 mg PO QD for 21 days followed by 7 days of rest. The following PK and ADME information was obtained from the palbociclib FDA NDA reviews. A population PK analysis of palbociclib revealed a mean plasma CL/F of 60.2 L/hr for a typical patient, resulting in an estimated total plasma AUCtau of 2080 hr-ng/mL at steady state. The reported plasma half-life of palbociclib in patients is 22 to 29 hr. Moreover, after multiple dosing of palbociclib 125 mg PO QD in two clinical studies, the mean total plasma Cmax ranged from 94.9 to 116 ng/mL, and the mean Ctrough was 61 ng/mL. Using standard calculations, the estimated total plasma Cavg,ss in patients is about 85 ng/mL. As there is no appreciable difference in plasma protein binding between mice and humans (Fu,p 0.159 vs 0.147 respectively), total plasma concentrations can be applied to estimate a murine equivalent dose (MED). Assuming dose- and time-linear PK and minimal inter-study PK variability in mice, and based upon the total plasma AUCinf vs the AUCtau in humans, the MED is estimated at 10 mg/kg palbociclib free base equivalents daily, with this formulation under similar study conditions. Likewise, the mouse plasma Cavg supports this dosage, with the back-estimated Cavg of 76.9 ng/mL being similar to the observed Cavg,ss of 85 ng/mL in patients. Abemaciclib The total plasma and tumor PK of abemaciclib in female athymic nude mice (Charles River Laboratories, aged 8-16 weeks) was assessed after a single oral gavage of 50 mg/kg of abemaciclib free base equivalents. Abemaciclib mesylate (LY2835219, Abmole, M2112, Lot 2, purity > 98%) was suspended in 1% hydroxyethylcellulose (HEC), 0.25% Tween 80, and 0.05% simethicone in ultrapure water at a final nominal concentration of 5 mg/mL free base equivalents for a 10 mL/kg gavage volume. Mice were sacrificed using an IACUC-approved method at 10 min, 1, 4, 8, and 24 hr post-dose, with 3 mice per time point. Whole blood by cardiac puncture was collected into Sarstedt Microvette 500 ml K3 EDTA microcentrifuge tubes, vortexed to mix anticoagulant, immediately centrifuged to plasma, and stored on dry ice for remainder of study. Mice were then perfused with PBS via the aorta, the orthotopic xenografts excised, rinsed with PBS, and placed on dry ice. At the end of the in vivo procedures, all samples were transferred from dry ice and placed at -80 C until analysis. Total plasma and tumor homogenate abemaciclib concentrations were assessed using a sensitive and specific liquid chromatography, tandem mass spectrometry assay. First, tumor samples were diluted with a 5:1 volume of ultrapure water, and homogenized with a bead-based technique (Liang et al., 2011) on a FastPrep-24 system (MP Biomedicals, Santa Ana, CA). 1.4 mm ceramic spheres (MP Biomedicals, Lysing matrix D, 10 mg per mg of tumor) were added to the microcentrifuge tubes containing samples. The samples were then subjected to three 60 M/S vibratory cycles of 1 min each on the FastPrep-24 system. To prevent over-heating due to friction, samples were placed on wet ice for 5 min between each cycle. The homogenates were then stored at -80 C until analysis. Abemaciclib mesylate (LY2835219, Abmole, M2112, Lot 2, purity > 98%) stock solutions corrected for salt content were prepared in acetonitrile and used to spike matrix calibrators and quality controls. Plasma and tumor homogenate samples, 25 mL each, were protein precipitated with 100 mL of 150 ng/mL palbociclib (LC Labs, P-7788, Lot PLH-103, Purity > 99%) in acetonitrile as an internal standard. A 2 mL aliquot of the extracted supernatant was injected onto a Shimadzu LC-20ADXR high performance liquid e10 Cancer Cell 34, 1–16.e1–e19, September 10, 2018

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

chromatography system via a LEAP CTC PAL autosampler. The LC separation was performed using a Phenomenex Kinetex 2.6 mm EVO C18 (100 A˚, 50 x 2.1 mm) maintained at 50 C with gradient elution at a flow rate of 0.5 mL/min. The binary mobile phase consisted of 20 mM ammonium acetate in H2O in reservoir A and methanol: acetonitrile: 2-propanol (40:30:30 v/v) in reservoir B. The initial mobile phase was maintained at 25% B for 0.5 min followed by a linear increase to 100% B in 2.5 min. The column was then rinsed for 2 min at 100% B and then equilibrated at the initial conditions for two min for a total run time of seven min. Under these conditions, the analyte and IS eluted at 1.37 and 0.98 min, respectively. Analyte and IS were detected with tandem mass spectrometry using a SCIEX API 5500 Q-TRAP in the positive ESI mode with monitoring of the following mass transitions: abemaciclib 507.28 -> 393.20, palbociclib 448.25 -> 380.10. The experimental bioanalytical runs were all found to be acceptable for the purpose of a singlicate non-GLP, preclinical PK assessment. A linear model (1/X2 weighting) fit the calibrators across the 5.00 to 500 ng/mL range, with a correlation coefficient (R) of R 0.9967. The lower limit of quantitation (LLOQ), defined as a peak area signal-to-noise ratio of 5 or greater verses a matrix blank with IS, was 5.00 ng/mL. The intra-run precision and accuracy was < 15.6% CV and 90.7% to 112%, respectively, across the matrices. The resultant abemaciclib concentration-time (Ct) data were grouped by matrix and time point, and manual imputation of data below the lower limit of quantitation (BLOQ) was as follows: IF at any time point R 2/3rds of the Ct results were above the LLOQ, the BLOQ data were replaced with a value of ½ LLOQ, ELSE the entire time point’s data were treated as missing. Then, using Phoenix WinNonlin 6.4 (Certara USA, Inc., Princeton, NJ), Ct data summary statistics (arithmetic mean, standard deviation, % CV, minimum, median, maximum) were generated, and the abemaciclib arithmetic mean Ct data for each matrix was subjected to noncompartmental pharmacokinetic analysis (NCA). The extravascular model (Model 202) was applied, and area under the Ct curve (AUC) values were estimated using the ‘‘linear up log down’’ trapezoidal rule. The terminal phase was defined as the three time points at the end of the Ct profile, and the elimination rate constant (Ke) was estimated using an unweighted log-linear regression of the terminal phase. The terminal elimination half-life (T1/2) was estimated as 0.693/Ke, and the AUC from time 0 to infinity (AUCinf) was estimated as the AUC to the last time point (AUClast) + predicted Clast/Ke. Other NCA parameters estimated included observed maximum concentration (Cmax), time of Cmax (Tmax), concentration at the last observed time point (Clast), time of Clast (Tlast), apparent clearance (CL/F = Dose/AUCinf), and apparent terminal volume of distribution (Vz/F). The average concentration over a dosing interval (Cavg) was estimated as AUCinf / dosing interval in hours. The apparent partition coefficient of abemaciclib from the plasma to the tissue of interest (Kp,tissue) was estimated as the ratio of the AUCinf, tissue to AUCinf plasma when available. To estimate a clinically relevant mouse dosage, the resultant mouse plasma AUCinf was compared with the reported human plasma PK values at the putative single agent abemaciclib maximum tolerated dose (MTD) of 200 mg PO BID (Patnaik et al., 2016). All inferences were made under the assumption of time-independent, linear and dose-proportional PK in mice and humans. Abemaciclib concentrations showed low variability in the mouse plasma, demonstrating coefficients of variation of 4.1% to 40.1% across the sampling time points, with the highest variability noted at the 0.167 hr time point in plasma and tumor. All plasma and tumor results were above the LLOQ. Tumor penetration appeared to be relatively rapid and extensive, with a Kp,tumor value of 6.20 based on AUCinf. The T1/2 of abemaciclib in the plasma was 5.2 hr with a similar half-life observed in the tumor (5.98 hr). Notably, our abemaciclib plasma PK at 50 mg/kg PO differed from that reported by Tate (Tate et al., 2014). Our mice showed a similar plasma Ct profile up to 4 hr post-dose, but then concentrations declined rapidly resulting in a 3.8-fold lower AUC. Our mice did not exhibit the prolonged absorption and low CL/F as described by these investigators, despite the use of the mesylate salt form and the same suspension formulation. NCA parameter estimates and the mean (SD) abemaciclib Ct plots, grouped by matrix are shown below: In a Phase 1 study of single agent abemaciclib administered continuously, the MTD was 200 mg PO BID. The geometric mean human total plasma AUCtau at steady state was 3000 hr-ng/mL, with a corresponding geometric mean trough concentration of 197 ng/mL. Assuming dose proportionality and linear PK processes in mice, a murine equivalent dose (MED) providing a similar total plasma AUC would be 14.4 mg/kg. The mean in vitro fraction unbound in the plasma (Fu,p) was determined for mice and humans, yielding values of 0.054 and 0.027, respectively (Raub et al., 2015). As these values are at the two-fold threshold defining an appreciable difference for preclinical studies, it is unclear whether abemaciclib doses should be adjusted for plasma protein binding between the two species. Therefore, the rounded and recommended MED ranges from 7.5 to 15 mg/kg PO BID, based upon unbound and total plasma AUCs respectively. Doses in this range should also provide unbound steady state trough plasma concentrations similar to humans at 200 mg PO BID. More PK studies may be necessary to determine the source of the discrepancy between our plasma PK and that reported by Tate, or to further confirm our PK findings. Trametinib The total plasma and tumor PK of trametinib in female athymic nude mice (Charles River Laboratories, aged 12 -16 weeks) was assessed after a single oral gavage of 0.1 mg/kg of trametinib free base equivalents. Trametinib free base (Abmole, Lot M1759, purity > 99%) was suspended in 5% dimethyl sulfoxide (DMSO), 10% Cremophor RH40, and 10% PEG400 in 75% ultrapure water at a final nominal concentration of 0.02 mg/mL free base equivalents for a 5 mL/kg gavage volume. Mice were sacrificed using an IACUCapproved method at 30 min, 1, 4, 8, and 16 hr post-dose, with 3 mice per time point. Whole blood by cardiac puncture was collected into Sarstedt Microvette 500 ml K3 EDTA microcentrifuge tubes, vortexed to mix anticoagulant, immediately centrifuged to plasma, and stored on dry ice for remainder of study. Mice were then perfused with calcium- and magnesium-free PBS via the right ventricle, the orthotopic xenografts excised, rinsed with PBS, weighed, divided into aliquots and placed on dry ice. At the end of the in vivo procedures, all samples were transferred from dry ice and placed at -80 C until analysis. Cancer Cell 34, 1–16.e1–e19, September 10, 2018 e11

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Total plasma and tumor homogenate trametinib concentrations were assessed using a sensitive and specific liquid chromatography, tandem mass spectrometry assay. First, tumor samples were diluted with a 5:1 volume of ultrapure water, and homogenized with a bead-based technique (Liang et al., 2011) on a FastPrep-24 system (MP Biomedicals, Santa Ana, CA). 1.4 mm ceramic spheres (MP Biomedicals, Lysing matrix D, 10 mg per mg of tumor) were added to the microcentrifuge tubes containing samples. The samples were then subjected to three 6.0 M/S vibratory cycles of 1 min each on the FastPrep-24 system. To prevent over-heating due to friction, samples were placed on wet ice for 5 min between each cycle. The homogenates were then stored at -80 C until analysis. Trametinib free-base (Abmole, Lot M1759, Purity > 99%) stock solutions were prepared in acetonitrile and used to spike matrix calibrators and quality controls. Plasma and tumor homogenate samples, 25 mL each, were protein precipitated with 75 mL of 75 ng/mL refametinib (Cayman Chemical Co., Lot 0464258-2, Purity R 95%) in acetonitrile as an internal standard. A 2 mL aliquot of the extracted supernatant was injected onto a Shimadzu LC-20ADXR high performance liquid chromatography system via a LEAP CTC PAL autosampler. The LC separation was performed using a Phenomenex Gemini C6 Phenyl (3.0 mm, 30 mm x 2.0 mm) column maintained at 60 C with gradient elution at a flow rate of 0.50 mL/min. The binary mobile phase consisted of water-200 mM ammonium acetate pH 6.0 (90:10 v/v) in reservoir A and acetonitrile-water-200 mM ammonium acetate pH 6.0 (90:10:10 v/v) in reservoir B. The initial mobile phase consisted of 25% B with a linear increase to 65% B in three min. The column was then rinsed for two min at 100% B and then equilibrated at the initial conditions for two min for a total run time of seven min. Under these conditions, the analyte and IS eluted at 2.91 and 2.8 min, respectively. Analyte and IS were detected with tandem mass spectrometry using a SCIEX API 5500 Q-TRAP in the positive ESI mode with monitoring of the following mass transitions: trametinib 616.04 -> 491.10, refametinib 572.91 -> 394.00. The experimental bioanalytical runs were all found to be acceptable for the purpose of a singlicate non-GLP, preclinical PK assessment. A linear model (1/X2 weighting) fit the calibrators across the 1.00 to 500 ng/mL range, with a correlation coefficient (R) of R 0.9985. The lower limit of quantitation (LLOQ), defined as a peak area signal-to-noise ratio of 5 or greater verses a matrix blank with IS, was 1.00 ng/mL. The intra-run precision and accuracy was < 5.73% CV and 95.8% to 103%, respectively, across the matrices. The resultant trametinib concentration-time (Ct) data were grouped by matrix and time point, and manual imputation of data below the lower limit of quantitation (BLOQ) was as follows: IF at any time point R 2/3rds of the Ct results were above the LLOQ, the BLOQ data were replaced with a value of ½ LLOQ, ELSE the entire time point’s data were treated as missing. Then, using Phoenix WinNonlin 6.4 (Certara USA, Inc., Princeton, NJ), Ct data summary statistics (arithmetic mean, standard deviation, % CV, minimum, median, maximum) were generated, and the trametinib arithmetic mean Ct data for each matrix was subjected to noncompartmental pharmacokinetic analysis (NCA). The extravascular model (Model 202) was applied, and area under the Ct curve (AUC) values were estimated using the ‘‘linear up log down’’ trapezoidal rule. The terminal phase was defined as the three time points at the end of the Ct profile, and the elimination rate constant (Ke) was estimated using an unweighted log-linear regression of the terminal phase. The terminal elimination half-life (T1/2) was estimated as 0.693/Ke, and the AUC from time 0 to infinity (AUCinf) was estimated as the AUC to the last time point (AUClast) + predicted Clast/Ke. Other NCA parameters estimated included observed maximum concentration (Cmax), time of Cmax (Tmax), concentration at the last observed time point (Clast), time of Clast (Tlast), apparent clearance (CL/F = Dose/AUCinf), and apparent terminal volume of distribution (Vz/F). The average concentration over a dosing interval (Cavg) was estimated as AUCinf / dosing interval in hours. The apparent partition coefficient of trametinib from the plasma to the tissue of interest (Kp,tissue) was estimated as the ratio of the AUCinf, tissue to AUCinf plasma when available. To estimate a clinically relevant mouse dosage, the resultant mouse plasma AUCinf was compared with the reported human plasma PK values at the putative single agent trametinib maximum tolerated dose (MTD) of 2 mg PO QD. All inferences were made under the assumption of time-independent, linear and dose-proportional PK in mice and humans. Trametinib displayed low within-study PK variability with concentration coefficients of variation (CV) ranging from 5.57% to 26.4%. An adequate terminal phase in the tumor was not evident in the sampling period, and therefore tumor terminal phase parameters were not estimable. Using the ratio of the AUClast, tumor to AUClast plasma a Kp,tumor of 5.49 was estimated. NCA parameter estimates and the mean (SD) trametinib Ct plots, grouped by matrix are shown below: The common FDA-approved dose for single agent trametinib is 2 mg PO QD. A population PK analysis of trametinib (Ouellet et al., 2016) revealed a mean plasma CL/F of 4.91 L/hr for a typical female patient and 6.19 L/hr for a typical male patient, resulting in a mean AUC of 365 hr-ng/mL across genders. The reported plasma half-life of trametinib in patients is 92.5 to 115 hr. After multiple dosing of trametinib 2 mg PO QD to steady state, the mean total plasma Cmax ranged from 17.9 to 21.2 ng/mL, and the mean Ctrough ranged from 11.7 to 15.0 ng/mL. Using standard calculations, the estimated total plasma Cavg,ss in patients is about 15.1 ng/mL. The trametinib fractions unbound in the plasma (Fu,p) in mice and humans are 0.05 and 0.03, respectively. Adjusting the key PK parameters of AUC and Cavg by the Fu,p for each species, a MED of 0.1 mg/kg PO QD in the specified cosolvent solution would be clinically comparable to 2 mg PO QD in humans. AZD1775 Mouse equivalent dosing for AZD1775 was previously completed and the mouse equivalent dosage was determined to be 60 mg/kg with twice daily administration (Stewart et al., 2017). Ganetespib The total plasma and tumor PK of ganetespib in female athymic nude mice (Charles River Laboratories, aged 8-16 weeks) was assessed after a single, slow intravenous (IV) injection of 150 mg/kg via the tail vein. Ganetespib (MedChem Express, HY-15205, Lot e12 Cancer Cell 34, 1–16.e1–e19, September 10, 2018

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

08249, purity > 98%, MolWt 364.40) was dissolved in DMSO and further diluted with Kolliphor RH 40 (Sigma) in 5% dextrose for injection, USP (D5W, Baxter) for a final nominal concentration of 15 mg/mL in 10% DMSO / 18% Kolliphor RH 40 / 72% D5W and a 10 mL/kg injection volume. Mice were sacrificed using an IACUC-approved method at 10 min, 1, 4, 8, and 16 hr post-dose, with 3 mice per time point. Whole blood was collected with sodium heparin via cardiac puncture, immediately centrifuged to plasma, and stored on dry ice for remainder of study. Mice were then perfused with PBS via the aorta, the orthotopic xenografts excised, rinsed with PBS, and placed on dry ice. At the end of the in vivo procedures, all samples were transferred from dry ice and placed at -80 C until analysis. Total plasma and tumor homogenate ganetespib concentrations were assessed using a sensitive and specific liquid chromatography, tandem mass spectrometry (LC-MS/MS) assay. First, tissue samples were macerated, diluted with a 5:1 volume of ultrapure water, and homogenized with a bead-based technique (Liang et al., 2011) on a FastPrep-24 system (MP Biomedicals, Santa Ana, CA). Ceramic lysing matrix beads (MP Biomedicals, Ceramic Bead Lysing Matrix, 3 mg per mg of tissue) were added to the microcentrifuge tubes containing tumor samples. Samples were then subjected to three 60 M/S vibratory cycles of 1 min each on the FastPrep-24 system. To prevent over-heating due to friction, samples were placed on wet ice for 5 min between each cycle. The homogenates were then stored at -80 C until analysis. Ganetespib stock solutions were prepared in 100% methanol and used to spike matrix calibrators and quality controls. Plasma and tissue homogenate samples, 50 mL each, were treated with 10 mL of internal standard (SLV320, Tocris Bioscience, purity > 99%) 100 ng/mL spiking solution and then subjected to a liquid-liquid extraction procedure using methyl tert-butyl ether, vortexed and centrifuged at 4 C, with the supernatant evaporated to dryness in a CentriVap centrifugal vacuum concentrator (Labconco). The extracts were then reconstituted with 100 mL of 100% acetonitrile, and a 5 mL aliquot of the reconstituted extract was injected onto a Shimadzu LC-20ADXR high performance liquid chromatography system via a Shimadzu SIL-20ACXR autosampler. The LC separation was performed using a Phenomenex Luna C8 80A˚ LC column (4.0 mm, 30 mm x 2 mm) maintained at ambient temperature with gradient elution at a flow rate of 0.25 mL/min. The binary mobile phase consisted of 0.1% formic acid in ultra-pure water in reservoir A and 0.1% formic acid in acetonitrile in reservoir B. The initial mobile phase consisted of 5% B with an increase to 90% B over 2.5 min. After maintenance at 90% B for 4 min, the column was then rinsed for 1 min at 100% B and then equilibrated at the initial conditions for 1.5 min for a total run time of 9 min. Under these conditions, the analyte and IS eluted at 3.52 and 2.86 min, respectively. Analyte and IS were detected with tandem mass spectrometry using a SCIEX API 4000 in the positive ESI mode with monitoring of the following mass transitions: ganetespib 365.2 -> 131.1, SLV320 309.2 -> 211.1. The experimental bioanalytical runs were all found to be acceptable for the purpose of a singlicate non-GLP, preclinical PK assessment. A linear model (1/X2 weighting) fit the calibrators across the 5 to 200 ng/mL range, with a correlation coefficient (R) of R 0.99. Above the calibration range quality control samples were diluted with adequate precision and accuracy. The lower limit of quantitation (LLOQ), defined as a peak area signal-to-noise ratio of 5 or greater verses a matrix blank with IS, was 5 ng/mL for plasma and 30 ng/mL for tissues due to the dilution factor. The intra-run precision and accuracy was % 9.62% CV and 98.5% to 111%, respectively across the matrices. The resultant ganetespib concentration-time (Ct) data were grouped by matrix and time point, and manual imputation of data below the lower limit of quantitation (BLOQ) was as follows: IF at any time point R 2/3rds of the Ct results were above the LLOQ, the BLOQ data were replaced with a value of ½ LLOQ, ELSE the entire time point’s data were treated as missing. Then, using Phoenix WinNonlin 6.4 (Certara USA, Inc., Princeton, NJ), Ct data summary statistics (arithmetic mean, standard deviation, % CV, minimum, median, maximum) were generated, and the ganetespib arithmetic mean Ct data for each matrix was subjected to noncompartmental pharmacokinetic analysis (NCA). The IV bolus model (Model 201) was applied, and area under the Ct curve (AUC) values were estimated using the ‘‘linear up log down’’ trapezoidal rule. The terminal phase was defined as the three time points at the end of the Ct profile, and the elimination rate constant (Ke) was estimated using an unweighted log-linear regression of the terminal phase. The terminal elimination half-life (T1/2) was estimated as 0.693/Ke, and the AUC from time 0 to infinity (AUCinf) was estimated as the AUC to the last time point (AUClast) + predicted Clast/Ke. Other NCA parameters estimated included the back-extrapolated initial concentration (C0), observed maximum concentration (Cmax), time of Cmax (Tmax), concentration at the last observed time point (Clast), time of Clast (Tlast), clearance (CL = Dose/ AUCinf), and steady state volume of distribution (Vss). The apparent partition coefficient of ganetespib from the plasma to the tissue of interest (Kp,tissue) was estimated as the ratio of the AUCinf, tissue to AUCinf plasma when available. To estimate a clinically relevant mouse dosage, the resultant mouse plasma AUCinf was compared with the reported human plasma PK values at the putative single agent ganetespib at 200 mg/m2 which is the recommended Phase 2 dose (RP2D) (Goldman et al., 2013). All inferences were made under the assumption of time-independent, linear and dose-proportional PK in mice and humans. Ganetespib concentrations showed moderate variability in the plasma, demonstrating coefficients of variation of 2.9% to 65.6% across the sampling time points. All plasma and tumor results were above the LLOQ. Tumor penetration appeared to be rapid yet modest, with a Kp,tumor value of 0.639 based on AUCinf. The T1/2 of ganetespib in the plasma was 2.03 hr. The apparent tumor half-life was longer than that observed in plasma (8.42 hr), suggesting a high affinity of ganetespib for the orthotopic tumor. The large apparent Vss for tumor also suggests that ganetespib has high orthotopic tumor retention. Notably, our plasma and tumor PK was in line with that published by Shimamura et al. at 125 mg/kg IV in mice (Shimamura et al., 2012). NCA parameter estimates and the mean (SD) ganetespib Ct plots, grouped by matrix are shown below: In a Phase 1 study of single agent ganetespib administered IV once weekly for 3 weeks of a 4-week cycle, 200 mg/m2 was the RP2D (Goldman et al., 2013). The major adverse events with ganetespib are diarrhea and fatigue; however, retinal toxicities have been of Cancer Cell 34, 1–16.e1–e19, September 10, 2018 e13

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

concern. Cmax and AUC values appear to increase in a dose proportional manner over the dose range studied in humans (7 to 259 mg/m2). The estimated human total plasma AUCinf at 200 mg/m2 was approximately 7000 hr-ng/mL. Assuming dose proportionality and linear PK processes in mice, a murine equivalent dose (MED) providing a similar total plasma AUCinf would be 12.8 mg/kg. The plasma protein binding of ganetespib has only been reported as ‘‘high,’’ and data regarding differences in the extent between species are lacking. Therefore a MED can only be derived using total plasma concentrations, under the assumption of similar plasma protein binding between mice and humans. It bears noting that there are data suggesting that high affinity HSP90 inhibitors such as ganetespib may not demonstrate dose proportional PK, particularly in mice (Yamazaki et al., 2013). This is thought to result from a target mediated drug disposition, where the compound’s affinity for a highly abundant and specific target drives a nonlinear tissue distribution. We have generated limited plasma PK following ganetespib 50 mg/kg IV in mice indicating that this might be the case (data not shown). In this instance, a 3-fold increase in dose from 50 to 150 mg/kg results in an approximate 17-fold increase in plasma Cmax at 10 min after injection. Therefore, it is difficult for us to assume proportionality when back-extrapolating a MED based on the AUC from the 150 mg/kg dose. Taking this into consideration, we evaluated two dose levels of ganetespib in our mouse efficacy model: 10 mg/kg (likely close to MED) and 50 mg/kg. Western Blots of Treated RD Cells The RD cell line was obtained from American Type Culture Collection (ATCC, CCL-136). Cells were treated with 0.1 mM Ganetespib, 0.25 mM Irinotecan and 3 nM Vincristine, or vehicle (DMSO) and harvested 24 or 48 hr post-treatment. Cells were lysed in RIPA lysis buffer (Sigma, catalog R0278) containing protease inhibitors (Halt Inhibitor Cocktail, Thermo Fisher, catalog 87785). Proteins were resolved in 4-15% SDS-PAGE gel and transferred to PVDF membrane (Millipore). Membranes were blocked in Odyssey blocking buffer (Licor, catalog 927-4000) and then probed for the following antibodies: HSP70 (1:1000, Cell Signaling, catalog 4872) and actin (1:5000, Sigma, AC-15). Secondary antibodies were anti-rabbit (Licor, catalog 926-32221); anti-mouse (Licor, catalog 926-32210). Flow Cytometry RD Cells were treated with 1 mM AZD1775, 0.25 mM Irinotecan, 3 nM Vincristine, or vehicle (DMSO) for 24 or 48 hr. Cells were stained with propidium iodide solution (0.05 mg/mL propidium iodide, 0.1% (w/v) sodium citrate, 0.1% (v/v) Triton X-100, 2 ug/mL RNAse) and analyzed by flow cytometry for cell cycle stages on a Becton-Dickson FACS system (Rockville, MD, USA). QPCR RD Cells were treated with 0.1 mM-0.2 mM Ganetespib, 0.25 mM-0.5 mM Irinotecan, 3-6 nM Vincristine, or vehicle (DMSO) and harvested 24 or 48 hr post-treatment. RNA was extracted using Trizol and cDNA was synthesized using Superscript III System (Invitrogen, catalog 18080093). CDNA samples were analyzed using the SYBR Select Master Mix for CFX (Applied Biosystems, catalog 4472942) in duplicate and normalized to GAPDH (CCAGCAAGAGCACAAGAGGAA, GCCCCTCCCCTCTTCAAG). QPCR was performed on HSF-1 target genes (DARS: CCTGCTAGTGCACAGGCTG, CGAACCTATACTGAATTACTCC; HSP70: CAGTTAGTCT TACTGAGCTT, GCTGTAGCTATAGCTCAACT). QPCR was performed on cDNA from xenograft tissues for HSF-1 target genes (ENY2: GAACCACTTCTGAATATAGCA, CTTAGGCAACAGTGTGGA; HNRNPA2B1: GGATGAGAGCCCAGAGGTA, TCTCTGCA TCTGCTCTGGT; LTBP4: CATCCTGCAGCCCGCTTAT, GTAGGTCCCTTCTCCAGGT; USPL1: GAACTTACAATGTGTTCAGGTa, GATTACATGCTTCATGGGCA; ZNF467: CAGAGACTCAGCGAACCCTT, CACTCTTTCCTGCCCTGGTT). Microscopy RD Cells were treated with 1 mM AZD1775, 0.25 mM Irinotecan, 3 nM Vincristine, or vehicle (DMSO) for 24 or 48 hr on chamber well slides. Cells were treated with EdU for one hr and then incubated in 4% PFA to fix cells. EdU labeling was performed per manufacturer’s instructions (Click-iT EdU imaging kit, Invitrogen, catalog C10083), and DNA was stained with 0.2 mg/ml DAPI (Sigma-Aldrich). The cells were evaluated on a Zeiss Axioplan microscope. Cells with fragmented nuclei were scored as mitotic catastrophe in a blinded manner relative to the treatment. Multiple random fields were scored in duplicate. Preclinical Testing Preclinical testing was conducted using a phase I, II, III trial design as previously described (Stewart et al., 2014). Preclinical Phase I Non-tumor bearing athymic nude mice were used for tolerability testing of the following combinations: AZD1775 + VCR + IRN, and Ganetespib + VCR + IRN. Mice were given 4 cycles of chemotherapy with 21 days per cycle. Irinotecan was tested in 2 different schedules (low dose protracted schedule as well as a 5 day higher dose schedule similar to that done in pediatric patients). The chemotherapeutic drug combinations and schedules used were designed to mimic potential human clinical trials. Each treatment group contained 3 mice. Vincristine was dosed by intraperitoneal injection once a week on days 1, 8, and 15. Irinotecan was dosed by intraperitoneal injection once daily on days 1-5 and 8-12 or on days 1-5 only. AZD1775 was dosed by oral gavage twice daily on days 1-5. Ganetespib was dosed once a week by tail vein on days 1, 8, and 15. Mouse weight and complete blood counts were monitored for each treatment group. The health of the animals was monitored daily throughout therapy. All mice in all treatment groups survived 4 courses of therapy. No significant weight loss noted in any group tested.

e14 Cancer Cell 34, 1–16.e1–e19, September 10, 2018

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Preclinical Phase II RMS orthotopic xenografts were created by injecting luciferase labeled cells from SJRHB000026_X1 (ERMS), SJRHB012_Y (ERMS), and SJRHB013759_X1 (ARMS) into recipient athymic nude mice using the intramuscular injection technique previously described (Chen et al., 2013). Mice were screened weekly by Xenogen and the bioluminescence was measured. Mice were enrolled in the study after achieving a target bioluminescence signal of 106 -107 photons/sec/cm2 or greater for 2 weeks or a palpable tumor, and chemotherapy was started the following Monday. The following preclinical phase II trials were performed with 3 different xenografts (SJRHB000026_X1, SJRHB012_Y, SJRHB013759_X1): SJRHB000026_X1 (ERMS) Mice were randomized to the following treatment groups: Control, Ganetespib (GSP), Irinotecan (IRN), GSP + IRN, AZD1775 + IRN, Vincristine (VCR) + IRN, AZD1775 + VCR + IRN, Palbociclib + Trametinib, Abemaciclib + Trametinib. The following doses and schedules were used for each drug:

GSP

50 mg/kg tail vein weekly on days 1, 8, 15

AZD1775

60 mg/kg oral gavage twice daily, days 1-5

IRN

1.25 mg/kg intraperitoneal once daily, days 1-5 and 8-12

VCR

0.38 mg/kg intraperitoneal once daily, days 1, 8, 15

Palbociclib

10 mg/kg oral gavage daily on days 1-21

Abemaciclib

15 mg/kg oral gavage daily on days 1-21

Trametinib

0.15 mg/kg oral gavage daily on days 1-21

Control

vehicle control (no chemotherapy)

Mice received 4 courses of chemotherapy (3 weeks per course) and bioluminescence was monitored weekly and at the end of therapy. Disease response was classified according to bioluminescence signal. Mice with a signal of 105 photons/sec/cm2 or less (similar to background) were classified as complete response, 105-106 photons/sec/cm2 as partial response, 107 up to 108 photons/sec/cm2 (similar to enrollment signal) as stable disease, and greater than 108 photons/sec/cm2 as progressive disease. Mice with tumor burden at any time greater than 20% of body weight were also classified as progressive disease. Mice were monitored daily while receiving chemotherapy. SJRHB012_Y (ERMS) Mice were randomized to the following treatment groups: Control, GSP + IRN, VCR + IRN, AZD1775 + VCR + IRN, AZD1775 + VCR + IRN*, GSP + VCR + IRN, Palbociclib + Trametinib The following doses and schedules were used for each drug:

GSP

50 mg/kg tail vein weekly on days 1, 8, 15

AZD1775

60 mg/kg oral gavage twice daily, days 1-5

IRN

1.25 mg/kg intraperitoneal once daily, days 1-5 and 8-12

IRN*

3.125 mg/kg intraperitoneal once daily, days 1-5

VCR

0.38 mg/kg intraperitoneal once daily, days 1, 8, 15

Palbociclib

10 mg/kg oral gavage daily on days 1-21

Trametinib

0.15 mg/kg oral gavage daily on days 1-21

Control

vehicle control (no chemotherapy)

Mice received 6 courses of chemotherapy and bioluminescence was monitored weekly. Disease response was classified according to bioluminescence signal. Mice with a signal of 105 photons/sec/cm2 or less (similar to background) were classified as complete response, 105-106 photons/sec/cm2 as partial response, 107 up to 108 photons/sec/cm2 (similar to enrollment signal) as stable disease, and greater than 108 photons/sec/cm2 as progressive disease. Mice with tumor burden at any time greater than 20% of body weight were also classified as progressive disease. Mice were monitored daily while receiving chemotherapy.

Cancer Cell 34, 1–16.e1–e19, September 10, 2018 e15

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

SJRHB013759_X1 (ARMS): Mice were randomized to the following treatment groups: Control, GSP + IRN, VCR + IRN, AZD1775 + VCR + IRN, GSP + VCR + IRN, Palbociclib + Trametinib The following doses and schedules were used for each drug:

GSP

50 mg/kg tail vein weekly on days 1, 8, 15

AZD1775

60 mg/kg oral gavage twice daily, days 1-5

IRN

1.25 mg/kg intraperitoneal once daily, days 1-5 and 8-12

VCR

0.38 mg/kg intraperitoneal once daily, days 1, 8, 15

Palbociclib

10 mg/kg oral gavage daily on days 1-21

Trametinib

0.15 mg/kg oral gavage daily on days 1-21

Control

vehicle control (no chemotherapy)

Mice received 6 courses of chemotherapy and bioluminescence was monitored weekly. Disease response was classified according to bioluminescence signal. Mice with a signal of 105 photons/sec/cm2 or less (similar to background) were classified as complete response, 105-106 photons/sec/cm2 as partial response, 107 up to 108 photons/sec/cm2 (similar to enrollment signal) as stable disease, and greater than 108 photons/sec/cm2 as progressive disease. Mice with tumor burden at any time greater than 20% of body weight were also classified as progressive disease. Mice were monitored daily while receiving chemotherapy. Preclinical Phase III We performed a randomized, double-blind placebo controlled preclinical phase III trial using 499 mice enrolled from 4 different O-PDX lines. RMS orthotopic xenografts were created by injecting luciferase labeled cells from SJRHB000026_X1 (ERMS), SJRHB012_Y (ERMS), SJRHB013757_X2 (ARMS) and SJRHB013759_X1 (ARMS) into recipient athymic nude mice using the intramuscular injection technique previously described. Mice were screened weekly by Xenogen and the bioluminescence was measured. Mice were enrolled in the study over the span of 11 weeks after achieving a target bioluminescence signal of 106 -107 photons/sec/cm2 or greater for 2 weeks or a palpable tumor, and chemotherapy was started the following Monday. Individuals giving the chemotherapy as well as those doing the bioluminescence imaging were blinded to the drugs being administered as well as the specific O-PDX tumor line. The following 14 treatment groups were tested for all four O-PDXs concurrently: The following drugs and doses were used in the preclinical phase III trial directed by AUC-guided dosing obtained from pharmacokinetic studies matched closely to human doses that are clinically relevant:

Control

Vehicle/untreated

AZD1775

AZD1775 dosed at 60 mg/kg twice daily on days 1-5.

GSP

GSP dosed at 50 mg/kg weekly on days 1, 8, 15.

GSP + IRNa

GSP dosed at 50 mg/kg weekly on days 1, 8, 15, IRN dosed at 1.25 mg/kg daily on days 1-5 and 8-12

GSP + IRNb

GSP dosed at 10 mg/kg weekly on days 1, 8, 15, IRN dosed at 1.25 mg/kg daily on days 1-5 and 8-12

GSP + VCR + IRNa

GSP dosed at 50 mg/kg weekly on days 1, 8, 15, VCR dosed at 0.38 mg/kg weekly on days 1, 8, 15 and IRN dosed at 1.25 mg/kg daily on days 1-5 and 8-12

GSP + VCR + IRNb

GSP dosed at 10 mg/kg weekly on days 1, 8, 15, VCR dosed at 0.38 mg/kg weekly on days 1, 8, 15 and IRN dosed at 1.25 mg/kg daily on days 1-5 and 8-12

AZD1775 + IRN

AZD1775 dosed at 48 mg/kg twice daily on days 1-5 and IRN dosed at 1.875 mg/kg daily on days 1-5 (dosing based on recommended phase 2 dose from ADVL1312 clinical trial)

AZD1775 + VCR + IRNa

AZD1775 dosed at 48 mg/kg twice daily on days 1-5, VCR dosed at 0.38 mg/kg weekly on days 1, 8, 15 and IRN dosed at 1.875 mg/kg daily on days 1-5

AZD1775 + VCR + IRNb

AZD1775 dosed at 48 mg/kg twice daily on days 1-5, VCR dosed at 0.38 mg/kg weekly on days 1, 8, 15 and IRN dosed at 1.25 mg/kg daily on days 1-5 and 8-12

AZD1775 + VCR + IRNc

AZD1775 dosed at 60 mg/kg twice daily on days 1-5, VCR dosed at 0.38 mg/kg weekly on days 1, 8, 15 and IRN dosed at 1.25 mg/kg daily on days 1-5 and 8-12

VCR + IRNa

VCR dosed at 0.38 mg/kg weekly on days 1, 8, 15 and IRN dosed at 1.25 mg/kg daily on days 1-5 and 8-12 (equivalent to low dose protracted schedule with IRN dose of 20 mg/m2 in children)

VCR + IRNb

VCR dosed at 0.38 mg/kg weekly on days 1, 8, 15 and IRN dosed at 1.875 mg/kg daily on days 1-5 (equivalent to 5 day IRN dose of 30 mg/m2 in children).

VCR + IRNc

VCR dosed at 0.38 mg/kg weekly on days 1, 8, 15 and IRN dosed at 3.125 mg/kg daily on days 1-5 (equivalent to 5 day IRN dose of 50 mg/m2 in children)

e16 Cancer Cell 34, 1–16.e1–e19, September 10, 2018

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

Mice received 6 courses of chemotherapy (3 weeks per course) and bioluminescence was monitored on week 3 of each course and at the end of therapy. Disease response was classified according to bioluminescence signal as described above for preclinical phase II. End of study work up was done for all mice with progressive disease as previously described (Stewart et al., 2014). Chemotherapy Vehicles for Preclinical Testing The following table shows the vehicles used to mix each chemotherapy for in vivo use:

Chemotherapy

Vehicle

AZD1775

0.5% methylcellulose

IRN

Normal saline

VCR

Normal saline

Ganetespib

10% DMSO / 18% Kolliphor RH 40 / 72% D5W

Palbociclib

1% methylcellulose (Type 400 cPs) and 1% Tween 80

Abemaciclib

1% hydroxyethylcellulose (HEC), 0.25% Tween 80, and 0.05% simethicone in ultrapure water

Trametinib

1% methylcellulose (Type 400 cPs) and 1% Tween 80

Xenogen Imaging and Quantification Mice were given intraperitoneal injections of Firefly D-Luciferin (Caliper Life Sciences 3 mg/mouse). Bioluminescent images were taken five min later using the IVIS 200 imaging system. Anesthesia was administered throughout image acquisition (isoflurane 1.5% in O2 delivered at 2 L/min). The Living Image 4.3 software (Caliper Life Sciences) was used to generate a standard region of interest (ROI) encompassing the largest tumor at maximal bioluminescence signal. The identical ROI was used to determine the average radiance (photons/s/cm2/sr) for all xenografts. QUANTIFICATION AND STATISTICAL ANALYSIS ChIP Sequencing Analysis We first employ BWA (version 0.5.9-r26-dev, default parameter) to align the ChIP-Seq reads to human genome hg19 (GRCh37-lite). Picard (version 1.65(1160)) was then used for marking duplicated reads. Then only non-duplicated reads have been kept by samtools (parameter ‘‘-q 1 -F 1024’’ version 0.1.18 (r982:295)). We followed the ENCODE criterion to quality control (QC) the data then nonduplicated version of SPP (version 1.11) have been used to draw cross-correlation and calculated relative strand correlation value (RSC) under support of R (version 2.14.0) with packages caTools (version 1.17) and bitops (version 1.0-6) and estimated the fragment size. We required > 10M unique mapped reads for point-source factor (H3K4me2/3, H3K9/14Ac, H3K27Ac, CTCF, RNAPolII, BRD4) and RSC > 1. We required 20 M unique mapped reads for broad markers (H3K9me3, H3K27me3, H3K36me3). We required 10 M unique mapped reads for INPUTs and RSC < 1. We noticed H3K4me1 is a point-source factor in some stages while broad in other stages so we QC H3K4me1 as broad markers. Then upon manual inspection of the cross-correlation plot generated by SPP, the best fragment size estimated (the smallest fragment size estimated by SPP in all our cases) were used to extend each read and generate a bigwig file to view on IGV (version 2.3.40). All profiles were manually inspected for clear peaks and good signal to noise separation. Hidden Markov Modeling analysis Non-duplicated aligned reads were extended by fragment size defined above and ChromHMM (version 1.10, with ‘‘-colfields 0,1,2,5 -center‘‘ for BinarizeBed ) was used for chromatin state modeling. To choose the state number we first modeled all mouse development stages together from 7 states to 33 states and selected the model with 16 states upon manual inspection. For better visualization of the dynamics of HMM states across stages, we normalized color intensity by the max total percentage of a state covered by a gene and flanking region. To determine the best region represent a gene, we first filtered annotated isoforms by TSS within 2 kb of any H3K4me3/H3K4me2/H3K27Ac/H3K9-14Ac peaks at any development stage and then we selected the highest expressed isoform at any development stage estimated by cuffdiff or the longest isoform if no expression level estimated by cuffdiff. At last, we reduced the interval for a HMM state to half bar and the intensity to half of the normalized intensity if it didn’t rank in the top 2 HMM states for a gene. As HMM states can be assigned by multiple genes, the max total percentage across genes have been used for these normalization. Core Regulatory Circuit (CRC) Analysis The CRC mapper was run as described previously (Aldiri et al., 2017) except we performed two separate analyses to compare data with and without promoter regions. We used our H3K4me3 ChIP-seq data to remove those putative super-enhancers that overlap with promoters. This was done because some very actively expressed genes had large regions of H3K27Ac that are called as super-enhancers and it is difficult to distinguish between these putative super-enhancers and strong promoters at genes. Cancer Cell 34, 1–16.e1–e19, September 10, 2018 e17

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

DNA Methylation and Analysis Whole genome bisulfite sequencing (WGBS) data were mapped to GRCh37-lite (human) or by BSMAP2.74 with the following parameters: -z 33 -f 5 -g 3 -r 0 -m 17 -x 600 –u and raw methylation status was called by Bis-SNP (Version 0.82.2). We developed a regres
#Methylated C + 0:5 . sion tree based approach to identify regional methylation pattern of a sample using m-values m = log2 #Unmethylated C + 0:5 Adjacent segments with insignificant (p > 1E-6, Student t-test) or small difference between average b-values
 #Methylated C b= (< 0.1) were recursively merged. #Methylated C + #Unmethylated C Segmentation breakpoints from individual samples were pooled to generate an initial list of regions for testing correlation between methylation and gene expression. Correlation between regions with variable methylation (difference between highest and lowest average m-values >= 0.3) and adjacent genes with variable expression (methylation region overlapping 10 kb upstream of TSS to TES, at least 2 fold difference between highest and lowest FPKM, highest FPKM >= 1) were calculated. Gene/region pairs with FDR <= 0.5 were retained for merging, where adjacent regions for the same gene (no more than 1 kb apart) with the same orientation (positive or negative) of correlation were merged. Correlations were re-evaluated and pairs with final FDR <= 0.1 were used as differentially methylated regions that were significantly correlated with gene expression. Integrated Analysis To prioritize genes that are potentially dysregulated in tumors, we integrated multiple types of data sources based on order statistics (Zhang et al., 2012). Concretely, this method takes as input N different sets of gene rankings (with details described below), and outputs integrative gene rankings by combining all N data source, based on the formula: Z r1 Z r2 Z rN Qðr1 ; r2 ; .; rN Þ = N! . dsN dsN
1 .ds1 (Equation 1) 0

s1

sN
1

where ri is the rank ratio of data source i, N is the total number of data source, and r0 = 0. In practice, the Q value for each gene was calculated using a more efficient algorithm (Aerts et al., 2006): Qðr1 ; r2 ; .; rN Þ = N!VN

Vk =

Xk i=1

ð
1Þi
1

Vk
i i r i! N
k + 1

(Equation 2)

(Equation 3)

Since the Q statistic is not uniformly distributed (Stuart et al., 2003), we ran permutation to estimate empirical p value for each gene by following the published procedure (Zhang et al., 2012). In each permutation, missing values were randomly assigned to k genes (with k determined by real data) and gene ranks were permuted for each data source, and Q values were calculated based on the permutated data matrix. The null distribution was constructed using all genes with at least 3 data types (n = 12 k). Empirical p values were derived from the estimated null distribution, followed by multiple testing correction using the BH method (Benjamini and Hochberg, 1995). Four types of data source (i.e., whole proteome, phosphoproteome, transcriptome and epigenome) were considered for integration. For gene ranking based on RNA (from transcriptome) and protein (from whole proteome) abundance, genes were ranked by the log2 fold change of each tumor sample against control (i.e., the myoblast samples). Note that genes can be ranked either in descending (i.e. the most up-regulated/activated genes on the top) or ascending (i.e. the most down-regulated/repressed genes on the top) order and both cases were considered in our analysis. For epigenome, we calculated the fractional difference between active states (1, 2, 3, 4, 5, 6, and 10) and repressive/inactive states (14, 17, and 18) for each gene (defined as TSS to TES plus 2kb on both 5’ end and 3’ end) in each sample, and genes were ranked by the factional difference of each tumor versus control. For phosphoproteome, the phosphorylation level for each protein was summarized by the spectral count-based weighted average of the relative intensity ratio of all associated phosphopeptides (Bai et al., unpublished data), and genes were ranked by the log2 ratio. After gene rankings were defined for each data source, an integrative ranking was generated for each tumor sample by the order statistics method as described above. To summarize gene ranking for each tumor subtype (A/E-RMS), for each data source the log2 ratio values were first averaged across samples from each subtype, and the gene ranks from multiple data source were similarly integrated. To generate tumor subtype-specific gene ranking, for each data source we calculated log2 ratio by comparing the mean of the two subtypes (i.e. ERMS versus ARMS), then performed the integration. Finally, to summarize the integrative ranking on gene level into pathway level, we performed Gene set enrichment analysis (Subramanian et al., 2005) using the pre-ranked function within javaGSEA application (http://software.broadinstitute.org/gsea/index.jsp). Briefly, genes were ranked by the integrative protein ranking (for each tumor sample, subtype, or subtype specificity), and FDR values were derived by permuting gene sets from the pathway database. Pathways with FDR < 0.05 were considered as significant. Statistical Analysis for Preclinical Phase III Study Pre-study consultation for the preclinical phase III study was done with the Department of Biostatistics at St. Jude Children’s Research Hospital. The data include a total of 499 mice treated in 14 different treatment groups from 4 O-PDX models. The goal e18 Cancer Cell 34, 1–16.e1–e19, September 10, 2018

Please cite this article in press as: Stewart et al., Identification of Therapeutic Targets in Rhabdomyosarcoma through Integrated Genomic, Epigenomic, and Proteomic Analyses, Cancer Cell (2018), https://doi.org/10.1016/j.ccell.2018.07.012

of the study was to assess and compare the tumor response and tumor progression free survival among treatment groups. Mice were randomized into each treatment group with a randomization code provided by the Department of Biostatistics. Tumor response was defined by bioluminescence at the end of therapy as described above. If a mouse was taken off study at any point in time after enrollment due to tumor size greater than 20% of mouse weight, they were automatically assigned a response of progressive disease. The survival curves of time to tumor progression were generated by Kaplan-Meier method. The log-rank tests were used to compare survival curves in each sub-group. DATA AND SOFTWARE AVAILABILITY Genomic, transcriptomic, epigenetic, and proteomic data files are available in a public online portal (https://pecan.stjude.org/ proteinpaint/study/RHB2018). All sequence data are deposited in the The accession number for the sequence data reported European Molecular Biology Laboratory nucleotide sequence data library (Accession number: EGAS00001002967). Details for software availability are in the Key Resources Table.

Cancer Cell 34, 1–16.e1–e19, September 10, 2018 e19