Interactions between driver genes shape the signaling pathway landscape and direct hepatocellular carcinoma therapy

Abstract Hepatocellular carcinoma (HCC) is one of the most lethal malignancies, whose initiation and development are driven by alterations in driver genes. In this study, we identified four driver genes (TP53, PTEN, CTNNB1, and KRAS) that show a high frequency of somatic mutations or copy number variations (CNVs) in patients with HCC. Four different spontaneous HCC mouse models were constructed to screen for changes in various kinase signaling pathways. The sgTrp53 + sgPten tumor upregulated mTOR and noncanonical nuclear factor‐κB signaling, which was shown to be strongly inhibited by rapamycin (an mTOR inhibitor) in vitro and in vivo. The JAK‐signal transducer and activator of transcription (STAT) signaling was activated in Ctnnb1 mut   + sgPten tumor, the proliferation of which was strongly inhibited by napabucasin (a STAT3 inhibitor). Additionally, mTOR, cytoskeleton, and AMPK signaling were upregulated while rapamycin and ezrin inhibitors exerted potent antiproliferative effects in sgPten + Kras G12D tumor. We found that JAK‐STAT, MAPK, and cytoskeleton signaling were activated in sgTrp53 + Kras G12D tumor and the combination of sorafenib and napabucasin led to the complete inhibition of tumor growth in vivo. In patients with HCC who had the same molecular classification as our mouse models, the downstream signaling pathway landscapes associated with genomic alterations were identical. Our research provides novel targeted therapeutic options for the clinical treatment of HCC, based on the presence of specific genetic alterations within the tumor.


| INTRODUC TI ON
Hepatocellular carcinoma remains a major global health challenge. [1][2][3] Several drugs have been approved for the treatment of patients with HCC to date. 4,5 However, patient response rates continue to fall below 20%, which can be partially explained by high tumor heterogeneity. 6 The highly heterogeneous pattern of genetic alterations present in tumors is believed to contribute to the high diversity of HCC. [7][8][9] Molecular classification has provided a basis for the prognosis and treatment of certain types of tumors (e.g., human epidermal growth factor receptor 2 status in breast cancer and epidermal growth factor receptor status in lung cancer). 10 Deep-sequencing studies have confirmed that TP53 and CTNNB1 are frequently mutated in HCC tumors. 11,12 In addition, RAS and mTOR signaling pathways are highly enriched in a considerable number of patients with HCC. 13 Although our understanding of the drivers of the disease has improved, this knowledge is yet to be translated into clinical practice. 14 Indeed, dominant mutational drivers in HCC have yet to be effectively targeted by drugs. 15 It has been recently shown how different combinations of driver genes affect tumor heterogeneity. 16,17 Despite this, the question of whether different combinations of mutations affecting driver genes can upregulate specific kinase cascades and provide directions for targeted therapeutic interventions, remains unanswered.
Here, we established four distinct mouse models of HCC, each with mutations affecting two distinct driver genes, and demonstrated that the cooperation between certain driver genes leads to different transcriptomic and proteomic profiles, reflecting the intertumor complexity observed in patients with HCC. Our results provide a targeted therapeutic approach for the treatment of patients with specific molecular subclasses of HCC.

| Human samples
The study was approved by the Ethics Committee of the University of Science and Technology of China (USTCEC201600004). Tissues from patients with liver cancer (n = 6) or benign disorders (hemangioma; n = 1) were collected from the First Affiliated Hospital of Anhui Medical University. Written informed consent was obtained from all patients.
The clinical characteristics of all patients are shown in Table S1.

| Cell line establishment
All the fresh cancer tissues were resected from mouse liver tumor tissues. The tumor tissues were finely chopped with scissors into small fragments and digested using the Tumor Dissociation Kit (Miltenyi Biotec). Next, cell suspensions were filtered using a 70 μm cell strainer. The culture medium was composed of DMEM (Cytiva), supplemented with 10% FBS (Gibco), 100 U/mL penicillin, 100 μg/mL streptomycin (Solarbio), sodium pyruvate, insulin, and transferrin (Procell). All the cell lines were cultured continuously from passage 10 onwards and cells that had undergone 20 or more passages were regarded as cell lines. Any potential mycoplasma contamination was removed using the Mycoplasma Removal Kit (TransGen).

| Hydrodynamic injection
For hydrodynamic liver injection, plasmid DNA suspended in 2 mL saline was injected into the tail vein of 9-week-old male B6 mice over a 5-7 s period. The amount of injected DNA was 60 μg for sgPten, 60 μg for sgp53, and 10 μg for Cas9. All mice were dosed with carbon tetrachloride as previously described. 18

| Immunohistochemistry
The mouse tumor tissue was fixed in 4% paraformaldehyde solution overnight and 4 μm tissue sections were obtained for immunohistochemical staining. A list of the Abs used is provided in Table S2.

| Plate colony formation assay
Cells (1 × 10 6 /well) were seeded into 6-well plates. The cells were incubated at 37°C and 5% CO 2 until visible colonies appeared. The colonies were fixed and stained with crystal violet staining solution. Next, the OD value was measured at 562 nm using a microplate reader (Model 680; Bio-Rad). All of the experiments were repeated three times, and the average values were adopted.

| Impedance-based label-free test of toxicity (xCELLigence for cellular proliferation)
All compounds were tested in three independent experiments in quadruplicate using the xCELLigence platform (RTCA; ACEA Biosciences). Data from each well were normalized to the time point just before compound addition using the RTCA software, which generated the normalized cell index.

| Immunoblotting
Lysates of normal liver and tumor tissue were extracted using the  Table S2.

| Cell proliferation assay
Cells (5000 /well) were seeded into 96-well plates. Ten microliters of CCK-8 was added to each well, and the cells were then incubated at 37°C for 45 min. The OD value was measured at 450 nm using a microplate reader (Model 680; Bio-Rad). All of the experiments were repeated three times, and the average values were adopted. A list of small molecule inhibitors targeting specific kinase signaling cascades is provided in Table S3.

| Whole-exome sequencing
Fresh tumor or control (adjacent to the tumor tissue) tissue samples were isolated from patients with liver cancer (n = 6) and sequenced using commercial DNA sequencing services (Guangzhou Gene Denovo Biotechnology). Whole-exome sequencing was carried out using a targeted capture approach with the SureSelect Human All Exon V6 Kit (Agilent Technologies). The captured sequences were further amplified for 150 bp paired-end sequencing on Illumina X Ten system (Illumina).
Tumor and control tissue adjacent to tumor DNA samples had an average sequencing depth of the target exonic region of 100×.The data generated in this study are publicly available through the Sequence Read Archive using the SUB12086982 reference.

| RNA sequencing
Total RNA was extracted from normal liver or tumor tissues by TRIzol reagent (Invitrogen). Samples were library-prepped and sequenced using Illumina NextSeq (GeneWiz). The data generated in this study are publicly available through the Gene Expression Omnibus using the GSE208279 reference. Gene set enrichment analyses were carried out using the clusterProfiler R package (Guangchuang Yu) version 4.0.4.

| Cancer signaling phosphoprotein profiling using an Ab array
Tissue lysates were applied to the Phospho Explorer Antibody Array, which was designed and manufactured by Full Moon Biosystems.
Data were collected and analyzed by Wayen Biotechnologies.
The array images were scanned with the SureScan Dx Microarray Scanner at 532 nm and the fluorescence intensity was quantified using GenePix Pro (Axon Instruments) version 6.0 software. The proteome array data are shown in Table S4.

| Statistical analyses
Statistical significance was determined using Prism 9.0 software (GraphPad). Two-tailed unpaired or paired Student's t-tests between two groups and two-way ANOVA across multiple groups were used to determine significance.

| Construction of four primary mouse HCC models, guided by driver gene mutation analysis
To identify driver gene alterations in HCC, we analyzed the results of exon sequencing undertaken on 358 samples derived from The Cancer Genome Atlas database. The OncoPrint plots showed that TP53 (31%) was the most frequently mutated gene following by CTNNB1 (26%) ( Figure 1A). Regardless of early or advanced liver cancer, TP53 and CTNNB1 have consistently been shown to have the highest mutation frequency. 19 More than 75% of TP53 gene mutations could induce loss-of-function effect in HCC. 11,20 In our study, we used the inactivation of the p53 protein by CRISPR/Cas9 strategy to simulate tumors.
We refer to previous work for further details of model establishment. 18 Missense mutations of CTNNB1 in HCC patients often lead to the aberrant activation of the Wnt/β-catenin pathway. 21 We also analyzed the copy numbers of driver genes. We found that 27% of all HCC samples had a PTEN deletion ( Figure 1B). The HCC samples with a PTEN deletion showed lower levels of mRNA expression than those that had diploid PTEN. 45 These data indicate that the reduction in PTEN copy number is likely to be one of the main mechanisms that contributed to the downregulation of PTEN in HCC patients. Although the mutation rate of the KRAS gene in liver cancer is only 1.4% ( Figure 1A), the frequency of KRAS with DNA copy number gain or amplification (which correlates with increased mRNA expression) in CNV was 11% ( Figure 1B). Moreover, the mutation rate of the RAS/RTK signaling pathway where the KRAS gene is located is as high as 22%-37%. 10 In conclusion, we identified four key driver genes commonly found in HCC.
A pan-cancer analysis of whole genomes found that cancer is frequently driven by multiple oncogenic drivers, while a single gene mutation often fails to cause tumor growth. 22 We further explored the combination of the driver genes mentioned above. However, because mutations in TP53 and CTNNB1 are considered to occur in a mutually exclusive manner 23  To elucidate the cooperation between the driver genes and their combined impact on intertumor heterogeneity, we established mouse models based on the hydrodynamic tail-vein delivery of genetic elements ( Figure 1C). Tumor incidence was over 70% within 8 months. (Figure 1D). In subsequent validation studies, we observed the survival curves relating to the four models ( Figure 1E). To evaluate the tumor growth within the mouse models, we analyzed liver morphological changes at different time points. Physical images and B-ultrasound results showed significant growth in the advanced stages of the tumor models ( Figure S2A Figure S2G). Taken together, these results demonstrate that the genomic and copy numbers of the four driver genes and their combination are essentially ubiquitous in HCC, which led us to construct the four mouse models.

| Cooperation between distinct driver genes activated complex downstream signals
Most often, changes in driver genes lead to the activation of complex multiple tumorigenic signaling pathways. 24 Here, we used RNA and protein microarray techniques to identify the major signaling pathways triggered by different driver gene combinations in our HCC models. Principal component analysis of transcriptomic profiles showed a clear separation between tumors and normal tissues ( Figure 2A). A heatmap ( Figure S3) and principal component analysis ( Figure 2B) of proteomic profiles showed that there were differences in protein phosphorylation levels between the models. Volcano plots were used to exhibit differential gene expression among the HCC models ( Figure S4A). Gene set enrichment analysis of RNA sequencing data revealed that our four HCC models were enriched for different pathways ( Figure 2C). Phosphoproteome array analysis of the changes in the expression of phosphoproteins was carried out to search for differences between our HCC models, specifically focusing on signaling pathways important for tumorigenesis ( Figure 2D).
The phosphorylation levels of key molecules in each signaling pathway are shown in the form of volcano plots ( Figure S4B). Pathway analysis showed differences in pathway gene expression among different HCC models and normal liver tissues ( Figure S4C).
To verify the results of RNA sequencing and the protein phosphorylation microarray, we used western blot analyses to directly detect the expression and phosphorylation levels of key molecules in the aforementioned kinase signaling cascades in tumors and normal liver tissues ( Figure 3A). The results are summarized in a dendrogram ( Figure 3B).
These results show that the combination of driver genes with specific genetic alterations differentially promotes liver tumorigenesis. To determine whether the murine HCC tumors recapitulate downstream signaling activation features that are concordant with those of human HCC tumors, protein phosphorylation levels in the cancer tissues of the six patients with HCC were measured using cancer-adjacent normal tissues and benign tumor tissues as controls ( Figure 4E). Overall, the mouse models of HCC translationally recapitulate the pathways that are found in human subjects with HCCs, thus validating our mouse models for the screening of targeted drugs for the treatment of HCC.

| Inhibitors targeting specific kinase signaling cascades can inhibit tumor cell growth in vitro
To screen for inhibitors against the most important signaling pathways in each model, murine cell lines from tumor-bearing model mice were first isolated and cultured in vitro. Sorafenib is an oral   Here, we observed the significant activation (phosphorylation) of the downstream targets of sorafenib in the sgTrp53 + sgPten model ( Figure 5A). Moreover, the viability of sgTrp53 + sgPten cells decreased significantly with increasing concentrations of sorafenib in vitro. However, the remaining three murine models did not respond well to sorafenib in a plate colony formation assay ( Figure 5B).
In order to solve the problem of sorafenib resistance, we undertook inhibition experiments for the remaining three models that did not respond to sorafenib. Targeted inhibitors of kinase cascades were screened out based on our previous results. Treatment with STAT3 inhibitor efficiently inhibited the proliferation of Ctnnb1 mut + sgPten tumor cells ( Figure 5C,D). In sgPten + Kras G12D tumor cells, P38 and AMPK inhibitors had little to no effect, while the Erk inhibitor actually facilitated tumor cell growth ( Figure 5E,F). This might be because the inhibition of Erk triggers the negative feedback loops of Erk signaling, leading to Erk inhibitor resistance. 27,28 Of the various inhibitors tested, the mTOR inhibitor almost completely suppressed sgPten + Kras G12D tumor cell proliferation ( Figure 5G,H). The

| In vivo experiments showing how a selection of screened drugs could significantly inhibit HCC tumor growth
It has been reported that combinations of targeted agents are likely to be synergistic. 29 Thus, we implanted tumor cells isolated from the primary mouse models into healthy mice and performed experiments with either a single targeted agent or used it in combination with sorafenib. Generation of subcutaneous tumors using Ctnnb1 mut + sgPten tumor cells proved unsuccessful so we did not use these cells in in vivo assays.
Although the key roles of mTOR in HCC have been reported, 30 the results from clinical trials have remained negative as most of the clinical trials to date have been undertaken on an unselected patient population. Therefore, it is critical to identify biomarkers to allow for the selection of patients with HCC who might benefit from mTOR targeted suppression. 31 Treatment with rapamycin or sorafenib alone resulted in only a slight shrinkage of the sgPten + Kras G12D tumor volume ( Figure 6A). The combination of rapamycin and sorafenib, however, resulted in a significant increase in survival of these mice ( Figure 6B). In accordance, the subcutaneous tumors were significantly reduced in size when rapamycin and sorafenib were given simultaneously ( Figure 6C,D). Our results indicate that the development of liver tumors driven by PTEN and KRAS mutation responded well to a combination of rapamycin and sorafenib.
Signal transducer and activator of transcription 3 has been reported to be critical in HCC tumorigenesis and postsurgical recurrence. 32 Although the treatment of HCC tumors with napabucasin alone was ineffective ( Figure 6E,F), the combination of sorafenib and napabucasin (both STAT3 inhibitors) completely inhibited the growth of subcutaneous tumors in mice injected with sgTrp53 + Kras G12D tumor cells. In comparison, sorafenib treatment alone had little effect on tumor growth ( Figure 6G,H).
Although the sgTrp53 + sgPten tumor cells were significantly inhibited by sorafenib in vitro, tumor growth suppression was not as efficient in vivo ( Figure 6I). The results relating to the action of sorafenib in vitro and in vivo might be inconsistent because of the complex in vivo tumor microenvironment. Based on our previous findings, we identified that the mTOR signaling pathway was aberrantly activated in sgTrp53 + sgPten tumors ( Figure 3A,B). Therefore, we chose to perform in vivo experiments with rapamycin. Our results indicate that the combination of rapamycin and sorafenib significantly increase survival and had a strong oncorepressing activity in sgTrp53 + sgPten tumors in vivo ( Figure 6J-L).
It has been reported that the activation of the PI3K/Akt/mTOR signaling pathway promoted glycolytic metabolism in HCC and was related to immunosuppression in this context. 33 We also observed similar results in our in vivo experiments. Immunohistochemistry analyses revealed a decrease in the expression of the immunosuppressive myeloid cell marker CD163 after combination treatment in sgPten + Kras G12D tumors ( Figure 6M). An increase in immune infiltration markers (CD8 and NK1.1) and a decrease in CD163 + macrophages were observed after combination treatment in sgTrp53 + Kras G12D tumors ( Figure 6N), raising the possibility that this remodeling of the immune landscape could further improve the therapeutic effect of the combination treatment. This observation could be explained by the implication of STAT3 signaling in immunoregulation. 34 It has been reported that the inhibition of STAT3 in F I G U R E 5 Inhibition of specific kinase signaling cascades can inhibit hepatocellular carcinoma (HCC) tumor cell growth. (A) Heatmap shows the phosphorylation levels of proteins pharmacologically inhibited by sorafenib in normal liver and advanced murine HCC specimens. Results standardized by minimum-maximum normalization. (B) Mouse primary HCC cells are resistant to sorafenib treatment in vitro. Plate colony formation assay results for four mouse primary HCC cell lines. Cells were grown in the absence or presence of sorafenib at the indicated concentrations for 3 days, prior to fixing and staining. (C-K) Cell lines were treated with small molecule inhibitors of specific kinase signaling cascades, as indicated. (C, F, I) Representative CCK-8 assay, (D, G, J) mean values of quantitative CCK-8 results, and (E, H, K) RTCA assay of four mouse primary HCC cell lines. Experiments were repeated three times. Ezri, ezrin inhibitor; ns, not significant; OD, optical density; PDGFRβ, platelet-derived growth factor receptor β; VEGFR, vascular endothelial growth factor receptor.
HCC can reduce regulatory T cell infiltration and inhibit tumor macrophage differentiation. 35,36 Taken together, our results indicate that the combination of sorafenib with other kinase signaling pathway inhibitors might especially benefit patients with specific driver gene combinations (Figure 7).

| DISCUSS ION
We have generated four mouse models of HCC to investigate the heterogeneity of HCC tumors and identify targeted therapies for HCC patients with specific molecular subclasses. So far, no molecular subclass has been reported as responding to a specific targeted therapy in HCC. 37 Such research tells us that the molecular characterization of HCC tumors is urgently needed. 38,39 Here, we propose that the construction of mouse models based on the mutation of different driver genes is a promising approach for the preclinical testing of HCC-specific drug candidates.
Although comprehensive transcriptomic data have provided valuable information to guide tumor prediction and treatment selection, this method has limitations (e.g., mRNA expression does not reflect the true protein expression, the regulatory processes, or posttranscriptional modifications). 40,41 Proteomics could help us to identify clinically meaningful new treatments, particularly when applied to heterogeneous and genomically complex cancers like HCC. Therefore, our study focused more on evaluating protein rather than mRNA expression. We found that different combinations of genetic alterations contributed uniquely to tumorigenic protein kinase signaling cascades and HCC progression. These alterations affected not only the mutated genes but also their downstream signaling pathways. The sgTrp53 + sgPten tumor upregulated mTOR and noncanonical nuclear factor-κB signaling, Ctnnb1 mut + sgPten tumor upregulated JAK-STAT signaling, sgPten + Kras G12D tumor upregulated mTOR, cytoskeleton, and AMPK signaling, and sgTrp53 + Kras G12D tumor upregulated JAK-STAT, MAPK, and cytoskeleton pathways. Importantly, we identified the signaling pathway that most strongly inhibited tumor growth in each model. Such information could guide the development of future precision therapies for specific groups of HCC patients.
Ezrin is a critical structural protein involved in stabilizing membrane receptor complexes. [42][43][44] Ezrin is highly expressed and reflects an unfavorable prognosis in a number of tumors. [45][46][47] Our research has shown that NSC305787 inhibits the phosphorylation of ezrin-T567, resulting in a marked reduction in cell growth ( Figure 5F-K). However, we noted a failure of the ezrin inhibitor to inhibit subcutaneous tumors in mice (data not shown). The failure could be explained by drug metabolism or solubility. Nevertheless, F I G U R E 7 Hepatocellular carcinoma (HCC) is one of the most lethal and fastest growing malignancies worldwide. The authors aimed to investigate the differences between patients with HCC based on tumor molecular classification and provide targeted therapeutic options for them. The authors elucidated the cooperation between certain driver genes that leads to different transcriptomic and proteomic profiles, reflecting the intertumor complexity observed in HCC patients and screened out targeted inhibitors of sorafenib-resistant tumors with different molecular classifications. EZR, ezrin; HCC, hepatocellular carcinoma; NF-κB, nuclear factor-κB; STAT3, signal transducer and activator of transcription 3. our research into ezrin inhibitions has revealed a novel promising target for HCC tumor development and progression.
Taken together, our models represent an approach towards HCC treatment, whereby we dissect the regulatory mechanisms linking different driver genes and target the complex downstream networks to develop potent personalized therapies.

ACK N OWLED G EM ENTS
We wish to thank all the patients, healthy donors and doctors for their participation in this study.

FU N D I N G I N FO R M ATI O N
This work was supported by the Chinese Academy of Medical Sciences (#2019-I2M-5-073).

CO N FLI C T O F I NTER E S T S TATEM ENT
The authors declare no conflict of interest.

E TH I C S S TATEM ENTS
Approval of the research protocol by an Institutional Review Board: All human tissues used in the present study were obtained under the approval of the Ethics Committee of the University of Science and Technology of China (USTCEC201600004).