Genomic analysis of hepatoblastoma identifies distinct molecular and prognostic subgroups

Despite being the most common liver cancer in children, hepatoblastoma (HB) is a rare neoplasm. Consequently, few pretreatment tumors have been molecularly profiled, and there are no validated prognostic or therapeutic biomarkers for HB patients. We report on the first large‐scale effort to profile pretreatment HBs at diagnosis. Our analysis of 88 clinically annotated HBs revealed three risk‐stratifying molecular subtypes that are characterized by differential activation of hepatic progenitor cell markers and metabolic pathways: high‐risk tumors were characterized by up‐regulated nuclear factor, erythroid 2–like 2 activity; high lin‐28 homolog B, high mobility group AT‐hook 2, spalt‐like transcription factor 4, and alpha‐fetoprotein expression; and high coordinated expression of oncofetal proteins and stem‐cell markers, while low‐risk tumors had low lin‐28 homolog B and lethal‐7 expression and high hepatic nuclear factor 1 alpha activity. Conclusion: Analysis of immunohistochemical assays using antibodies targeting these genes in a prospective study of 35 HBs suggested that these candidate biomarkers have the potential to improve risk stratification and guide treatment decisions for HB patients at diagnosis; our results pave the way for clinical collaborative studies to validate candidate biomarkers and test their potential to improve outcome for HB patients. (Hepatology 2017;65:104‐121).

H epatoblastoma (HB) is the most common pediatric liver tumor. It has an annual incidence rate of approximately 1.8 diagnosed cases per million in the United States, and this rate is increasing at more than 4.3% annually. (1) HBs are embryonal neoplasms that are most commonly diagnosed during the first 3 years of life. They are believed to arise from hepatic cell precursors and are characterized by heterogeneous histological patterns reminiscent of liver developmental stages. (2) Therapeutic strategies combining surgical resection and chemotherapy have improved outcomes for children with HB, but the prognosis for patients with advanced or chemotherapy-refractory disease remains poor. (1) In addition, the most effective platinum-based agents for treatment of HB often lead to serious long-term adverse effects, including ototoxicity and nephrotoxicity. (1) We describe the results of a comprehensive genomic analysis of the largest set of clinically annotated HBs reported to date. Such efforts have previously identified prognostic groups and biomarkers for other embryonal tumors (3) as well as adult hepatocellular malignancies. (4)(5)(6) For example, survival-predictive and metastasis-predictive biomarkers based on both gene and microRNA (miRNA) expression profiles have been reported for hepatocellular carcinoma (HCC), (7) the most common liver tumor in adults. (4) Interestingly, a "stem-cell" gene-expression signature, involving several oncofetal, stem-cell markers, and pluripotent stem-cell expression profiles, has been identified in a particularly aggressive type of HCC. (4,(8)(9)(10) Results from recent international HB clinical studies suggest that underlying biological differences may be responsible for the prognostic variability and wide range of responses to chemotherapy seen in HB patients. (11,12) However, there are currently no biomarkers or international consensus regarding risk stratification for HB patients. In North America, Children's Oncology Group protocols have historically stratified patients by postsurgical stage and histological type, while in Europe and Japan, HB diagnosis is often based on tumor imaging criteria prior to therapy rather than pathologic analysis. (13) Consequently, few pretreatment HB specimens are available for molecular profiling outside of the United States. (11,12) Activation of the canonical Wnt-signaling pathway occurs in the vast majority of HBs through somatic mutations at beta-catenin (CTNNB1) and or other Wnt-signaling genes, but it can also be caused by germline alterations including adenomatous polyposis coli (APC) mutations and mutations associated with related genetic syndromes. (14,15) It remains unclear whether Wnt-signaling pathway dysregulation is required for HB genesis or whether it is prognostically significant. Studies have also identified HB-specific expression signatures that are related to liver development and HB histologic subtypes (16,17) as well as genes that are differentially expressed in HB relative to HCC and normal liver. (18) Multiple studies have speculated about the biological and prognostic importance of specific genes and signaling pathways, but these are limited by the availability of HB tumor specimens, most of which have been postchemotherapy samples with incomplete clinical annotation. (19) The goal of our study was to molecularly characterize a large cohort of pretreatment, clinically and histopathologically annotated HBs of sufficient size to identify significantly predictive diagnostic and prognostic biomarkers in this disease. Conclusions from previous efforts that focused on profiles of posttreatment tumors were limited to high-risk patients and were not able to identify prognostic biomarkers. Here, pretreatment HBs were profiled by whole-exome sequencing (WES) and targeted sequencing, copy-number single-nucleotide polymorphism (SNP) arrays, and messenger RNA (mRNA) and miRNA expression arrays to identify biomarkers that differentiate between low-risk and highrisk patients. The risk-stratification potential of candidate biomarkers identified in our 88-tumor study was evaluated prospectively using protein expression profiles of 35 clinically annotated HBs. The analysis of these data provided a more in-depth view of the landscape of HB genomes and transcriptomes and identified molecular targets for HB diagnostics and therapeutics.

Materials and Methods
We molecularly profiled 88 HBs with corresponding surgical pathology reports following histological review (Supporting Table S1). Because of limited DNA and RNA availability in some cases, not all tumors were ARTICLE INFORMATION: profiled using the same assays. In total, 34 HBs were profiled by WES and 46 by SNP arrays to identify mutations and determine copy number; 50 and 57 HBs were profiled by microarrays for mRNA and miRNA expression, respectively. All tumors were profiled for genetic alterations at CTNNB1, nuclear factor, erythroid 2-like 2 (NFE2L2), and the telomerase reverse transcriptase (TERT) promoter using targeted sequencing.

HISTOLOGICAL REVIEW
Histological review of representative glass slides and/or digital images of 94 samples was performed by M.J.F. and D.L.-T., who confirmed diagnoses and histological subtypes. There were no integrase interactor 1-negative (INI1 or SMARCB1) tumors included in the study, as recommended by the International Consensus HB Classification group. (11) Subtypes included epithelial, pure welldifferentiated fetal, pleomorphic or anaplastic fetal, embryonal, small cell, mixed, or HB with HCC features. Six samples-three non-HB malignant tumors, two samples with low HB content, and one sample with poor tumor RNA viability after extraction-were excluded from the study due to incorrect initial diagnosis or poor tumor RNA content. A total of 88 histologically confirmed cases-50 males and 38 females-were selected for profiling. CTNNB1, NFE2L2, AND TERT PROMOTER MUTATION ANALYSIS CTNNB1 exon 3 and 4 mutation analysis was completed using SuperScript One-Step reverse-transcription polymerase chain reaction (RT-PCR) with Platinum Taq reagents (Invitrogen) and primers BCAT-F and BCATlong-R to generate a 507-bp PCR product. (15) Two-directional Sanger sequencing analysis of PCR products was completed with Mutation Surveyor v.4.0.4 (Softgenetics). TOPO TA cloning (Invitrogen) of CTNNB1 PCR products was performed as needed when more than one PCR product was detected by agarose gel electrophoresis or when sequence quality was poor.

GENE AND miRNA EXPRESSION
mRNA expression was profiled by Affymetrix microarrays in a total of 50 tumors, six normal pediatric liver tissues, and a pool of five normal fetal liver samples; two of the profiled tumors (TLT-033 and TLT-079) had insufficient follow-up or clinical annotation and were excluded from biomarker discovery efforts. miRNA expression was profiled by Agilent microarrays in a total of 57 tumors and four normal pediatric liver samples. RNA was isolated from approximately 25 mg of frozen tissue specimens using the Qiagen RNeasy Mini extraction system, followed by deoxyribonuclease 1 treatment. miRNA was simultaneously isolated using the mirVana miRNA isolation kit (Ambion). Samples extracted by both methods were subsequently quantitated using the Nanodrop ND-1000 spectrophotometer, and integrity was monitored by the Agilent 2100 Bioanalyzer capillary electrophoresis system (RNA 6000 Nano Kit). Samples were selected for RNA integrity and profiled using Human Genome U133P2 Affymetrix arrays and Agilent miRNA Microarray System v2.2.

DNA AND RNA EXTRACTION
DNA from frozen tumor specimens and matched normal liver tissue or peripheral blood samples was extracted using the QIAamp DNA Mini Kit (Qiagen). Samples were treated with RNAse A and eluted in Buffer AE. DNA and RNA concentration and integrity were determined with the Nanodrop ND-1000 spectrophotometer and 0.8% agarose gel electrophoresis, using the Agilent 2100 Bioanalyzer. Paired samples demonstrating intact genomic DNA were selected for WES; these included 24 and 10 tumor samples paired with normal liver and blood DNA, respectively. RNA from frozen hepatic tumor specimens and matched normal liver tissue specimens was isolated using the RNeasy Mini Kit (Qiagen) or the mirVana miRNA Isolation Kit (Ambion). Samples were treated with deoxyribonuclease 1 and eluted in nuclease-free water. Samples with RNA integrity above 6.0 were used for genomic studies.

WHOLE-EXOME SEQUENCING
DNA was sequenced on an Illumina Hiseq 2000 with 80 million to 100 million successful 2 3 100-bp reads per sample. Putative mutations identified through exome sequencing were confirmed on a second sequencing platform. See Supporting Information for details on library construction, sequencing, and variant calling. Tumor exomes were analyzed using the Mercury v.3 pipeline, (20) including variant calling and annotation with a minimum variant ratio cutoff of 0.05. (21) Somatic mutations were annotated with information from the Catalog of Somatic Mutations in Cancer database. (22) SNP ARRAY Gene copy number was estimated from Affymetrix Genome-Wide Human SNP Array 6.0 profiles of 47 paired normal and HB tumor samples. DNA was digested with NspI and StyI enzymes (New England Biolabs), ligated to the respective Affymetrix adapters using T4 DNA ligase (New England Biolabs), amplified (Clontech), purified using magnetic beads (Agencourt), labeled, fragmented, and hybridized to the arrays. Following hybridization, the arrays were washed and stained with streptavidin-phycoerythrin (Invitrogen). GISTIC-normalized segments with log 2 ratios >0.3 or <20.3 were designated as copy-number gains or losses, respectively. Analysis of genes whose expression and copy-number variation profiles were significantly correlated identified a total of 174 genes with r > 0.58 (Pearson correlation between expression and copy-number variation), a minimum that was required to achieve both P < 0.01 and a false discovery rate (FDR) <0.1 for the selection. A chromosomal view of copy-number variation changes is given in Supporting Information.

STANDARDIZED EXPRESSION
Gene expression profiles, as estimated by individual probe sets, were standardized by transforming expression profiles, measured as maximum signal log ratio, to standard deviations from mean. Namely, for each probe set, expression mean and standard deviation were computed across all tumors-when comparing or aggregating expression-or across all nontumor samples-when identifying differentially expressed genes. Then, expression estimates in each sample were transformed to the number of standard deviations from the mean expression of this probe set. (23)

RNA EXPRESSION CLUSTERING
Profiles of our HBs together with profiles of 29 samples given by Cairo et al. (16) were processed by robust multi-array average and quantile-normalized. To exclude genes with low expression and variability across HBs, we focused on 9,835 probe sets with a maximum signal log ratio (normalized relative expression) above 6 and standard deviation above 0.3. Hierarchical clustering based on these probe sets was performed using hcluster in R with Pearson correlation and average linkage. To select the number of tumor clusters, we evaluated results using 2-10 clusters and averaged Silhouette values, Dunn's Index, Gap statistics, and cluster homogeneity and cluster separation. (24)(25)(26)(27) Analysis results, given in the Supporting Information, suggested that the optimal number of clusters is between three and five. We chose to partition the dendrogram to four clusters; partitioning to five clusters split the low-risk group into two groups with nearly identical survival and did not improve prognostic prediction.

TRANSCRIPTION FACTOR ACTIVITY INFERENCE
Transcription factor (TF) activity was estimated using average standardized expression of validated targets according to TRANSFAC. (29) We inferred activity for TFs with at least two targets, based on the average of their standardized profiles. TFs whose target sets overlapped by 50% or more were clustered, and a representative with the largest number of targets was selected. (30) Hepatic nuclear factor 1 alpha (HNF1A) activity was inferred based on expression profiles of RIPPLY1, CLDN2, CD41B, SLC22A9, and IATIL; NFE2L2 based on profiles of NQO1, CYP2A6, PREPL, BSEP, NAT2, ASL, and PRIP; and Yesassociated protein 1 (YAP1) based on profiles of CYR61, BIRC5, CTGF, and JAG1.

TESTING PREDICTIVE ACCURACY IN A VALIDATION COHORT
We quantified protein expression of alpha-fetoprotein (AFP), CTNNB1, lin-28 homolog B (LIN28B), high mobility group AT-hook 2 (HMGA2), HNF1A, and NFE2L2 by immunohistochemistry in 35 clinically annotated HB tumor specimens and matched normal liver specimens using formalin-fixed, paraffin-embedded tissue blocks available in the archives of the Department of Pathology of Texas Children's Hospital, after institutional review board approval. Immunohistochemistry was performed using formalin-fixed, paraffin-embedded tissue sections and the automated Leica Bond system. Epitope retrieval was carried out on the automated Bond system using either ER1 (Leica; AR99641) (pH 6) or ER2 (Leica; AR9640) (pH 9). Sections were incubated for 15 minutes with the primary antibody. We used Bond Polymer Refine Detection (Leica; DS9800), incubation with postprimary for 8 minutes, polymer for 8 minutes, 3,3 0 -diaminobenzidine for 10 minutes, and hematoxylin for 5 minutes. Lists of antibodies, antigen retrieval methods, working concentrations, and interpretation (scoring) guidelines for each antibody are given in Supporting Table S5.

MULTIPLE TESTING CORRECTIONS
All P values were corrected for multiple testing using Bonferroni correction. See Supporting Information for the number of tests used for each statistic.

Results
We begin by describing results from our profiling effort. Molecular profiles were used to identify prognostic biomarkers and to construct a predictive function, which was then tested prospectively using protein expression from 35 additional HB patients.  Table S2 lists all mutations identified by WES. Interestingly, younger patients were likely to have tumors with fewer somatic mutations (P < 1E-04; r 5 0.62, Pearson correlation). Somatic alterations at CTNNB1 were detected in 13 tumors, and somatic mutations were discovered in genes related to regulation of oxidative stress, including recurrent point mutations at NFE2L2 (NRF2) in 3/34 HBs. NFE2L2 mutations (p.R34G 3 2; p.D29N 3 1) occurred in hot spots that have been described in a variety of cancer types, including HCC and cancers of the lung, endometrium, and urinary tract. Somatic and germline mutations were identified in other genes that have been implicated in both pediatric and adult cancers, including APC and the chromatin-remodeling genes MLL2 and ARID1A. WES profiles are available at European Nucleotide Archive under project PRJEB11805.

TARGETED SEQUENCING
Identified alterations are given in Supporting Table  S1. RT-PCR and sequencing of CTNNB1 exons 3 and 4 identified point mutations and/or in-frame deletions in the ubiquitination domain of CTNNB1 in 78 of 88 HBs (89%; Fig. 1B); TLT-022 carried both a point mutation and a deletion. Three HBs were not successfully tested for alterations at CTNNB1 by RT-PCR due to a lack of viable RNA (TLT-029, TLT-037, and TLT-094); however, long PCR and sequencing of tumor DNA in these cases identified a p.D32Y point mutation in TLT-029 and a p.S22_G34del deletion in TLT-094. Among the 10 wild-type CTNNB1 cases, one (TLT-028) harbored a germline mutation at APC (identified by WES). Targeted sequencing of NFE2L2 exon 2 confirmed the three mutations identified by WES. In total, we identified mutations at NFE2L2 in four of 88 tumors (5%). NFE2L2 mutations (p.D29N in one patient and p.R34G in three patients; Supporting Fig. S1) clustered within the neh2 domain of the protein and were detected in tumors that also carried mutations at CTNNB1, as has been observed in adult HCC. (31) NFE2L2 p.R34G has been shown to increase NFE2L2 activity. (32) Two HBs-from patients 6 and 8 years old, the two oldest profiled-had somatic point mutations at the TERT promoter.

RNA AND miRNA EXPRESSION PROFILING
Analyses of genes that were highly expressed in HB tumors relative to normal liver tissue and of TFs and their targets across HB tumors revealed significant activation of Wnt-signaling genes, activation and upregulation of oncogenes and genes previously associated with HB and HCC prognosis, and cholangiocytic lineage differentiation genes, including genes whose expression or inferred activity was prognostically predictive (Supporting Table S3). Differential expression analysis, comparing miRNA expression profiles for HB tumors and normal liver samples, identified frequent down-regulation of lethal-7 (let-7) family members. Few signaling pathways were enriched between HB and normal liver profiles; however, multiple gene sets showed enrichment in tumors and for prognosis (Supporting Table S8).

ANALYSIS OF RNA AND miRNA EXPRESSION PROFILES: CLUSTERING OF mRNA PROFILES
We hierarchically clustered mRNA expression profiles (5,823 genes, 9,835 probe sets) of our HB tumors together with HBs reported by Cairo et al. (16) -a total of 73 tumor samples-in addition to 13 profiles of normal liver samples. Averaged Silhouette values (26) and Dunn's Index (25) suggested the presence of two to four distinct HB molecular clusters together with a cluster of normal liver profiles. We chose to represent HBs in three clusters ( Fig.  2A; Supporting Information). Survival analysis using our 51 tumor samples supported the predictive value of the resulting inferred molecular subtypes (P < 0.01), suggesting that the three clusters represent low-risk (HB1), high-risk (HB2), and intermediaterisk (HB3) groups (Fig. 2B). We found significant correlation between our molecular classification and that developed by Cairo et al. using profiles of postchemotherapy HBs (16) : six of seven of their poorprognosis tumors (C2 HBs) were included in our HB2 high-risk cluster. Moreover, while subtype prediction using the 16-gene signature of Cairo et al. was not prognostically predictive (P > 0.1), other genes that were differentially expressed in their profiled tumors relative to normal samples were differentially expressed across prognostic groups (GSEA analysis; Supporting Table S8).

DIFFERENTIALLY EXPRESSED GENES
Differential expression analysis revealed genes that were consistently dysregulated in HB tumors. The eight most up-regulated genes were previously linked to Wnt signaling (Fig. 2C). Dickkopf-1, a bellwether of Wntsignaling and beta-catenin/T-cell factor pathway activation, (15,16,(33)(34)(35) was the most up-regulated, 267-fold on average. We used average standardized expression profiles of Dickkopf-1, APC down-regulated 1, axis inhibitor 2, and cyclin D1 to infer Wnt-signaling activity. (36) Five HBs showed significant increases in both Dickkopf-1 expression and inferred Wnt-signaling activity, and the remaining HBs showed exceptional increases in both (Supporting Fig. S3); we refer to these five HBs as tumors exhibiting mild Wnt-signaling activation ( Fig. 2C; Supporting Table S3). GSEA identified multiple sets of differentially expressed genes identified by Cairo et al. but not their 16-gene signature, as enriched in differentially expressed genes across both HB and normal liver profiles and across favorable and poor prognosis. The most down-regulated genes, when compared to normal liver, were involved in metabolic and oxidation-reduction pathways, including metabolism of xenobiotics by cytochrome P450, as reported (see Supporting Fig. S4). (16,17) VARIABLY EXPRESSED miRNAs A total of 19 miRNAs were dysregulated (P < 0.01 by t test comparing expression in HB and normal, after multiple testing correction) and variably expressed (top 5% in terms of variability) in HB tumors. Hierarchical clustering suggested the presence of two clusters, including one cluster composed of let-7 family members in addition to miR-99a and miR-199a-3p ( Fig.  2D; Supporting Table S4); miR-99a is often coexpressed with let-7 miRNAs. (37) Low let-7 expression correlated with high LIN28B expression (P < 1E-3, Pearson correlation). LIN28B is a known regulator of let-7 miRNAs, (38) and its expression is inversely correlated with inferred HNF1A activity. In addition, let-7 expression directly correlated with inferred NFE2L2 activity ( Fig. 2D; P < 1e-8, Pearson correlation) and survival (P < 0.01 for let-7b by survival analysis).

GSEA
We identified gene sets that were enriched with dysregulated genes in HB tumors and in poor versus favorable prognosis tumors. Genes identified downstream from CTNNB1 (39) were among the strongest enriched in both poor-prognosis tumors and tumors versus normal tissue. In addition, sets of differentially expressed genes identified by Cairo et al., but not their 16-gene signature, were differentially expressed in both HB and prognosis. These included a 188-gene set identified as down-regulated in HB and a 676-gene set identified as up-regulated in HB. Other significantly enriched gene sets in poor-prognosis tumors include sets of differentially expressed genes in HCC by Boyault et al. (40) and Chiang et al. (39) In total, we identified 351 gene sets that were enriched in HB and/or poor-prognosis HB at P < 0.01 after multiple testing correction. These are given in Supporting Table S8.

PROGNOSTIC BIOMARKERS
Unsupervised clustering identified three distinct prognostically predictive molecular clusters of HB tumors ( Fig. 2A,B), suggesting a potential for specific prognostic biomarkers and predictive classifiers based on RNA expression. In total, we profiled RNA expression in 47 tumors with known clinical staging and survival state, including 24 tumors from patients who lived 2 or more years after diagnosis with no evidence of disease (favorable prognosis) and two and eight tumors diagnosed in patients who went on to relapse and/or died of disease (poor prognosis), respectively.

NFE2L2 ACTIVATION
NFE2L2 was recurrently mutated (5% of tumors), and NFE2L2 inferred activity was predictive of prognosis at P < 1E-03, by fluoroethyl-L-tyrosine and by survival analysis. In total, nine high-risk tumors exhibited significant increases in inferred NFE2L2 activity (>2 standard deviations from mean), including three tumors with identified NFE2L2 mutations. Note that only three of the four tumors with mutations at NFE2L2 were profiled for RNA expression and that inferred NFE2L2 activity, but not its expression, was prognostically predictive.

SALL4 UP-REGULATION
There were 174 genes that showed significant correlation between expression and copy number at r > 0.58 (P < 0.01, FDR <0.1); however, copy number at the loci of these genes was not significantly predictive of prognosis, and these genes showed no significant pathway enrichment. The only prognosis-predictive copy-number change was a recurrent gain at 20q13.2, a region that includes four genes. Analysis of RNA expression of these genes identified spalt-like transcription factor 4 (SALL4) as a candidate prognostic biomarker. No favorable and seven poor-prognosis HB tumors exhibited high SALL4 expression (TLT-009, TLT-021, TLT-036, TLT-058, and TLT-088) or SALL4 copy gain (TLT-009, TLT-012, TLT-056, and TLT-058). Amplifications at SALL4 (at least one extra copy) and its up-regulation (2-fold change and 7 standard deviations from normal) were significantly predictive of poor prognosis (P < 0.01 by survival analysis).

PROGNOSTICALLY PREDICTIVE mRNA EXPRESSION
Given the size of our data set, in order to identify prognostic biomarkers, we had to limit our pool of candidates before testing their association with prognosis to 50 or fewer. We focused on the analysis of a small subset of genes with evidence for HB-specific activity, including NFE2L2 and SALL4, identified as candidates based on DNA evidence, and a panel of 20 highly expressed genes in HB relative to normal liver samples ( Fig. 2C; Supporting Table S3). These included the prognostically predictive genes LIN28B, HMGA2, and AFP. Mild Wnt-signaling activation, identified in five tumors, was predictive of favorable prognosis (5/5 patients) but not significantly (P < 0.2); however, because Wnt signaling-associated genes were significantly up-regulated in our data, we included Wnt-signaling activation in our predictive function. Low LIN28B and high HMGA2 expression predicted favorable and poor prognosis (P < 5E-03), respectively. Higher AFP expression correlated with poor prognosis: 7 patients with tumors expressing AFP at 23 (6 standard deviations away from normal) died of the disease, and only one survived disease-free for more than 2 years (P < 5E-03). LIN28B expression levels were highly correlated, with increased expression of other oncofetal genes (AFP, glypican 3 [GPC3]) and stemcell/precursor markers (epithelial cell adhesion molecule [EPCAM], delta-like 1 homolog [DLK1]), and inversely correlated with NOTCH and phosphatase and tensin homolog signaling (Supporting Fig. S5). In addition, while DLK1 and GPC3 expression were significantly increased in HB, they were not prognostically predictive (Supporting Table S3). Interestingly, while c-MYC expression was previously implicated in determining HB risk, (16) most of our samples were negative for both c-MYC mRNA and protein expression (Supporting Fig. S6).

PROGNOSIS-PREDICTIVE FUNCTION
When combined, mRNA expression profiles of candidate biomarkers produced significant HB risk prediction. High NFE2L2 activity and high expression of LIN28B, HMGA2, SALL4, and AFP were predictive of poor prognosis and marked our high-risk group; all 10 HB patients who were profiled for mRNA expression and died of disease or relapsed were markedly high for at least one of these biomarkers. Low Wntsignaling activity and low expression of LIN28B and HMGA2 were predictive of favorable prognosis and marked our low-risk group, identifying 11 of our favorable-prognosis patients and none of our poorprognosis patients (P < 0.05). Patients who were not predicted to be at low or high risk were placed in the intermediate-risk group (Fig. 3). The resulting classification function was predictive at P < 0.05 by survival analysis.

TF ACTIVITY
TF activity was inferred for a total of 88 TFs with validated targets in TRANSFAC. It was evaluated as a predictor of age at diagnosis, survival years following diagnosis, survival state (deceased, relapsed, or free of disease), and clinical stage. Inferred HNF1A activity was the most significantly (anti) correlated with tumor stage, FDR <9e-5 (Supporting Fig. S7). NOTCH1 activity and the expression of its target, phosphatase and tensin homolog, were predictive of survival. NFE2L2 activity was also highly correlated with age at diagnosis (FDR <0.05), suggesting that high NFE2L2 activity may be a feature in older children. Interestingly, inferred YAP1 activity, based on the expression of targets identified by Tao et al., (36) suggests that YAP1 is differentially activated in high-risk HB (P < 0.001; Supporting Fig. S8).

VALIDATION SET TO TEST PREDICTIVE ACCURACY
To test whether our RNA-based predictive function is also predictive on the protein level in the prospective setting, we assembled a panel of 35 HBs (validation set), including 15 tumors from patients who died of disease or had relapsed at their last follow-up visit (poor-prognosis tumors) and 20 tumors from patients who survived at least 2 years after diagnosis with no evidence of disease (favorable-prognosis tumors). Of the 35 tumors, 15 were pretreatment and 20 were obtained after chemotherapy. Each pretreatment and posttreatment tumor sample was stained with antibodies for HNF1A, SALL4, LIN28B, AFP, HMGA2, and NFE2L2; and tumor stains were evaluated by reviewing pathologists (D.L.-T., M.J.F.) as exhibiting focal or diffuse expression and having negative or weakly or strongly positive expression of each of the proteins; see Supporting Table S6 for classifications and Fig. 4 for representative examples. We note that HBs may be highly heterogeneous and that biomarker expression can vary greatly across cellular compartments (Supporting Figs. S9 and S10). HNF1A was not a prognostic biomarker candidate, and its expression was used to predict clinical staging. In addition, tumors were tested for CTNNB1 and GPC3 protein expression as a part of standard diagnostic protocols.
A total of five tumors in this validation set were identified to exhibit low HNF1A expression, all highstage tumors; this was not significant because 90% of our validation-set tumors were high-stage. CTNNB1 expression was uniformly high and provided no predictive value. SALL4, HMGA2, and LIN28B expression were high in six, four, and three tumors, respectively, all poor prognoses. AFP expression was high, and NFE2L2 was expressed in five and six posttreatment tumors, respectively, all with poor prognoses; targeted sequencing revealed no mutations at NFE2L2.
Our classifier was permitted to select features independently for pretreatment and posttreatment tumors, and it limited use of AFP and NFE2L2 expression to posttreatment tumors only; we conjecture that this was due to successful treatment of high-risk pretreatment tumors with high expression for the two proteins. In total, two patients in the favorable-prognosis set, HBTCH28 and HBTCH8, had pretreatment tumors with high AFP and NFE2L2 expression, respectively. Both patients were cured by complete resection. Interestingly, even weak expression of NFE2L2 suggested high-risk disease: 9/9 tumors (including HBTCH8) with positive NFE2L2 expression were either poor prognosis or fully resected.
In total, 13/15 (87%) poor-prognosis tumors and 0/ 20 of the favorable-prognosis tumors were flagged by high protein expression of at least one of our five proposed prognostic biomarkers (Supporting Table S6).

FIG. 3
These results mirror conclusions from RNA-profile analysis, suggesting that our antibody panel is useful for predicting risk (FDR <0.006). In total, two patients died of disease but were not predicted to have poor prognosis by our classifier. The first case, HBTCH25, was a pretreatment macrotrabecular HB that shared histological features with HCC and had elevated NFE2L2 expression but was unclassified because NFE2L2 was used to flag posttreatment tumors only. The second case was one of two patients with tumors that were negative for GPC3; both of these patients were diagnosed with mixed HB with a mesenchymal component postchemotherapy.

Discussion
Genomic and transcriptomic analyses of 88 clinically annotated HBs facilitated the identification of distinct molecular risk groups with characteristic genomic and signaling-pathway alterations. Our analysis demonstrated that HB represents a biological spectrum of disease, ranging from tumors with few aberrations that exclusively target Wnt-pathway genes to tumors with high genomic instability and somatic alterations observed in HCC.

HALLMARKS OF HB TUMORIGENESIS
Our analysis confirmed that WNT-pathway activation is nearly universal in HB, with exceptionally high levels seen in poor-prognosis cases. Mutations and deletions of CTNNB1 were found in 78/88 (89%) tumors tested. Of the 10 CTNNB1 wild-type HBs, one had a germline APC mutation and six others showed strong Wnt-pathway activation (Figs. 1 and  3). In total, 85/88 (97%) patients had evidence from Wnt-pathway activation. Its universality may be indicative of stage of maturation arrest as the highest levels of Wnt-pathway activation in our samples were accompanied by high expression of hepatic stem-cell and progenitor markers LIN28B, EpCAM, AFP, and SALL4. (41)(42)(43)(44)(45) However, in contrast to what has been reported for HCC, EpCAM and DLK1 levels in our study were elevated in a large proportion of HBs and not clearly associated with prognosis.

PROGNOSTIC BIOMARKERS
We identified recurrent hot spot mutations at NFE2L2 as well as copy-number gains at SALL4both potential therapeutic targets. (9,46) NFE2L2 mutations and pathway activation were found in 11% of our HB patients, all with high-risk tumors. These data are consistent with a previous report that identified mutations in CTNNB1 and NFE2L2 in 12/15 (80%) and 2/15 (13%) profiled HBs, respectively. (32) CTNNB1 alterations were identified in all four NFE2L2 mutated tumors, as is typical of NFE2L2-mutated adult HCC, suggesting a biological link between NFE2L2-mutated HBs and HCCs. (31) All patients with tumors carrying NFE2L2 mutations in our study died of their disease, with the exception of one patient who underwent an early complete tumor resection.
Our analysis revealed three prognosis-predictive HB molecular subtypes that were characterized by differential expression of hepatic stem/progenitor markers (LIN28B, SALL4, AFP) hepatobiliary (HNF, NOTCH1), metabolic (NFE2L2), and cancer-related pathways (P53, TERT). Tumors with favorable prognosis were characterized by canonical Wnt-pathway activation, low levels of genetic instability, and low LIN28B and SALL4 expression and NFE2L2 activity. HNF1A activity levels directly correlated with the clinical tumor stage in our patient cohort ( Fig. 3; Supporting Fig. S7). HNF1A is known to control liverspecific transcriptional programs (47) and is associated FIG. 3. Prognostic biomarkers. Our analysis of RNA expression profiles in 48 HBs revealed three risk-stratifying HB molecular subtypes that were characterized by the degree of Wnt-pathway activation; LIN28B, HMGA2, SALL4, and AFP expression; and NFE2L2 activation. Low-risk tumors were characterized by lower Wnt-pathway activity and low LIN28 and HMGA2 expression, while high-risk tumors were characterized by high HMGA2, SALL4, and AFP expression and NFE2L2 activation. The remaining tumors were classified as intermediate. Inferred HNF1 activity was predictive of tumor stage. Mutations were identified in the TERT promoter in tumor samples taken from the two oldest children in the study. These subtypes were largely in agreement with molecular classification obtained from unsupervised hierarchical clustering. For each tumor, we present survival years after diagnosis and survival status; clinical stage, age at diagnosis, and histological classification; the degree of Wnt-pathway activation (high or very high), as estimated based on target expression; expression, activation status, or the presence of a mutation in our predictive biomarkers; and the associated molecular classification obtained from unsupervised hierarchical clustering (HB1, HB2, or HB3). Tumor profiles are ordered by patient age. Stage 1 and stage 2 tumors are considered low grade, stage 3 tumors are high-grade, and stage 4 tumors are highgrade metastatic. Abbreviations: m, mesenchymal; M, mutation; SC, small cell. with hepatic progenitor cell specification (48) and favorable prognosis in hepatic adenomas. (5,49) Poor-prognosis HBs expressed high levels of AFP, which is consistent with data regarding the adverse prognosis of high or very low serum AFP. Retrospective histological review of a small number of tumors from patients with low-AFP suggested that these may have been erroneously classified as HBs. Consequently, these types of tumors are currently under active consensus review. (50) Histologically, HBs with high AFP expression demonstrated embryonal, small cell component, anaplasia, macrotrabecular patterns or "transitional" features and were more often seen in children older than 6 years. These have been associated with greater genomic instability. (32) In contrast, HBs with low hepatic progenitor marker expression and high activation of differentiation pathways demonstrated fetal histology and clustered in the low-risk group. (51)

ASSOCIATION WITH DEVELOPMENTAL STAGES
Many of the key genes linked to prognosis by expression profiling are known to be related to liver development and disease. LIN28B and HMGA2 expression and NFE2L2 activity were predictive of prognosis, and their profiles correlated with let-7 expression (Fig. 2D). LIN28B regulates stem-cell growth and metabolism and is known to inhibit let-7 miRNA maturation; together with HMGA2, LIN28B and let-7 are known to control stem-cell selfrenewal potential. (52) HBs in the high-risk group were also characterized by increased expression of other cancer stem/precursor cells, including EpCAM, DLK1, and oncofetal markers AFP, GPC3, and SALL4, and by increased activity of YAP1. (36) In HCC, embryonic stem-cell and progenitor-cell gene expression, including expression of EpCAM and SALL4, are associated with poor prognosis. (6) Similarly, EpCAM and stem cell-like genes including AFP and Wnt-pathway genes (4) have been reported to define a subtype of HCC that is associated with poor prognosis. (53) All patients with tumors demonstrating high levels of AFP were clustered in our high-risk group (Fig. 3), which is of particular interest given that high serum AFP is one of the few clinical parameters known to be associated with poor clinical prognosis in HB. (12) GENOMIC INSTABILITY SNP array data revealed copy number alterations in previously described altered regions of the HB genome, including chromosomes 1q, 2, 8q, 12, 1p, and 20, (16,32,54) as well as novel focal alterations. Gains in 1q and losses in 20q were more frequently identified in the high-risk HBs, but our study was not sufficiently powered to confirm significance of this or other potential associations, highlighting the need for large collaborative studies to obtain sufficient numbers of clinically annotated specimens.

HISTOLOGY AND BIOLOGY
Risk classes were associated with characteristic histological types. All tumors showing fetal histology with low mitotic activity (well-differentiated fetal HB) clustered in the favorable prognostic group. These features are known to be associated with an indolent clinical course, and well-differentiated fetal HBs are treated by surgery alone according to Children's Oncology Group protocols with nearly 100% survival. (51) However, some HBs in the low-risk group displayed other epithelial histology, including three with a minor smallcell component. None of the tumors in this group showed anaplasia, significant atypia, or HCC features; and none were classified as hepatocellular neoplasm not otherwise specified following the new international consensus classification. (11) Tumors in the high-risk group included embryonal and small-cell epithelial components, and only one (TLT-074) showed FIG. 4. Immunohistochemistry targeting prognostic biomarkers. Three representative tumors, two high-risk and one low-risk, from our 35-prospective tumor set were stained using antibodies to determine AFP, CTNNB1, LIN28B, NFE2L2 (NRF2), and HNF1A expression. HBTCH11: Hematoxylin and eosin stain demonstrating a predominantly embryonal HB with strong AFP and focally nuclear and cytoplasmic CTNNB1 expression. This high-risk tumor shows strong nuclear LIN28B and NFE2L2 expression and is negative for HNF1A. HBTCH8: This is another high-risk HB showing a mixed epithelial histology arranged in an almost organoid pattern with embryonal and fetal components in the periphery of nests of small undifferentiated cells; high AFP expression and nuclear CTNNB1 are present predominantly in the central nests of small cells, while LIN28B expression is present in the other epithelial patterns. It demonstrates strong NFE2L2 nuclear expression and is HNF1A-negative. This patient survived despite having a high-risk HB as the tumor was completely resected (right hepatectomy). HBTCH12: This is a low-risk HB in the group showing a predominantly fetal histology, lacking AFP and LIN28B expression, and showing focal cytoplasmic and nuclear CTNNB1; NFE2L2 expression is exclusively cytoplasmic, and HNF1A expression is strong and nuclear. Abbreviation: H&E, hematoxylin and eosin. HBs demonstrate a diverse range of histological phenotypes and clinical behaviors that may result from the proliferation of transformed stem cells or early hepatic precursors with varying degrees of differentiation. Dysregulation of various pathways may determine the ability of mutated liver stem cells to differentiate and progress to aggressive neoplasms (Fig. 5A). Consequently, we expect that, for example, mixed epithelial and mesenchymal HBs may arise from transformed precursor cells that retained the ability to differentiate to both lineages (43) ; note that the mesenchymal component of HBs does not produce GPC3. Our hypothesis is that HB phenotypes correspond to varying degrees of differentiation and activation of oncogenic pathways (Fig. 5B). High HNF1A and NOTCH induction, representing hepatic and biliary specification, (55) is seen in more differentiated tumors, suggesting that these better-differentiated HBs develop from cells at later stages of the lineage. On the other hand, tumors expressing stem-cell and liver-precursor markers such as LIN28B, EpCAM, and SALL4 are more aggressive and often show chemoresistance. In addition, activation of the NFE2L2 signaling pathway or TERT promoter mutations, which are present in tumors of older children and with biological features that overlap those of HCC, are associated with particularly poor prognosis.
Only three high-risk HBs showed significant atypia (TLT-070) or anaplasia (TLT-009, TLT-058). Interestingly, TLT-009 and TLT-058 were the oldest children included in our series, and both carried mutations at the TERT promoter. This is of particular interest as tumors with similar histological and clinical features were described and initially designated as transitionalcell liver tumors. (56) Histological review of a large series of these tumors by a group of international experts (11) indicated that these hepatocellular neoplasms, which are usually diagnosed in older children, may represent a distinct biological group of HB with HCC features that should be included under the new hepatocellular neoplasm not otherwise specified provisional category, until further characterization. In fact, TLT-009 was diagnosed as hepatocellular neoplasm not otherwise specified and recognized to share overlapping histologic and biological features between HB and HCC. (11,32) In conclusion, our analysis suggests that testing the proposed prognostic markers using a combination of targeted DNA sequencing, RNA expression profiling, and immunohistochemistry may be useful for predicting risk and response to therapy for HB patients. We propose a new HB stratification algorithm that could be used in combination with clinical parameters and histopathology. Immunohistochemical stains incorporating diagnostic antibodies for CTNNB1, AFP, GPC3, and integrase interactor 1 were previously recommended by the International Consensus HB Classification group. (11) We have described five additional predictive biomarker candidates: HNF1A, NFE2L2, SALL4, HMGA2, and LIN28B. These, together with copy-number assessment by cytogenetic or cytogenomic analysis, should be tested within a large-scale prospective clinical study. Integration of biomarkers associated with clinical outcomes into clinical stratification algorithms may be useful for HB prognostic classification and therapeutic guidance, particularly in the case of biomarkers for which emerging therapies may be available (let-7, LIN28B, NFE2L2, SALL4). Optimally designed, prospective collaborative clinical trials that include rigorous specimen banking and molecular analysis will be needed to prospectively characterize additional HB tumor specimens and further validate the biomarkers described here. (12)