In Vivo Phosphoproteome Analysis Reveals Kinome Reprogramming in Hepatocellular Carcinoma*

Aberrant kinases contribute to cancer survival and proliferation. Here, we quantitatively characterized phosphoproteomic changes in an HBx-transgenic mouse model of hepatocellular carcinoma (HCC) using high-resolution mass spectrometry, profiled 22,539 phosphorylation sites on 5431 proteins. Using a strategy to interpret kinase- substrate relations in HCC and to uncover predominant kinases in tumors, our results, revealed elevated kinase activities of Src family kinases (SFKs), PKCs, MAPKs, and ROCK2 in HCC, representatives of which were further validated in cell models and clinical HBV-positive HCC samples. Inhibitor combinations targeting Src and PKCs or ROCK2 both synergized significantly to inhibit cell growth. In addition, we demonstrated that phosphorylation at Src Ser17 directly affects its kinase activity. Our phosphoproteome data facilitated the construction of a detailed molecular landscape in HCC and should serve as a resource for the cancer community. Our strategy is generally applicable to targeted therapeutics, also highlights potential mechanisms of kinase regulation.

Aberrant kinases contribute to cancer survival and proliferation. Here, we quantitatively characterized phosphoproteomic changes in an HBx-transgenic mouse model of hepatocellular carcinoma (HCC) using high-resolution mass spectrometry, profiled 22,539 phosphorylation sites on 5431 proteins. Using a strategy to interpret kinasesubstrate relations in HCC and to uncover predominant kinases in tumors, our results, revealed elevated kinase activities of Src family kinases (SFKs), PKCs, MAPKs, and ROCK2 in HCC, representatives of which were further validated in cell models and clinical HBV-positive HCC samples. Inhibitor combinations targeting Src and PKCs or ROCK2 both synergized significantly to inhibit cell growth. In addition, we demonstrated that phosphorylation at Src Ser17 directly affects its kinase activity. Our phosphoproteome data facilitated the construction of a detailed molecular landscape in HCC and should serve as a resource for the cancer community. Our strategy is generally applicable to targeted therapeutics, also highlights potential mechanisms of kinase regulation. Aberrant kinase initiated dysregulation of phosphorylation signaling commonly contributes to tumor cell proliferation, survival, and migration (1,2). Targeting kinases with small molecular kinase inhibitors has been demonstrated to be an efficient strategy to treat cancer. Because of intrinsic and acquired resistance via the reactivation of receptor tyrosine kinases (RTKs) 1 and by-passing signaling pathways, or other mechanisms (3)(4)(5), the emergence of agent-resistance stands as a formidable challenge for targeted cancer therapeutics. Further, for the heterogeneous tumors, resistance can also arise via the positive selection of agent-resistant subpopulations (6,7). Combination strategies employing drugs to target tumors provide an opportunity to optimize antitumor efficacy and to prevent resistance. However, the selection of appropriate inhibitor combinations for therapy is challenging in terms of scope and cost. Thus, a clarion call has emerged for the development of combinatorial approaches for targeted cancer therapies (4).
Hepatocellular carcinoma (HCC) is a highly lethal cancer, ranking as the fifth most prevalent tumor and the second leading cause of cancer-related deaths worldwide (8,9). Activation of several key carcinogenic signaling pathways, involving VEGF, EGF, PDGF, FGF, RAF, and IGF was shown to be responsible for the initiation and progression of HCC, and sorafenib has been shown to inhibit RAF, VEGF, and PDGF. Most patients are diagnosed at advanced-stages HCC, and only treatment with sorafenib significantly prolongs patients survival (10). However, very few groups of patients with advanced-stage HCC treated with sorafenib have achieved even partial responses (11). There is little evidence for a specific pathway addiction in HCC, and many patients are likely to develop acquired resistance to sorafenib (12). After the approval of sorafenib, up to seven large phase III clinical trials reported negative findings (10). Thus, there is an urgent need to identify oncogenic pathways addiction and escape pathways in HCC, and novel therapeutic strategies are required to improve the survival of patients with HCC.
To meet these challenges, we present an integrative approach combining phosphoproteomic and computational strategies to characterize phosphorylation dysregulation in tumors, to identify predominant kinases for targeted cancer therapies and to unearth novel mechanisms regulating kinase activities ( Fig. 1A and supplemental Fig. S1). As HBx plays crucial role in carcinogenesis of HBV-related HCC (13), we applied this strategy to an HBx transgenic mouse model of HCC and as well as the corresponding wild-type (WT) mice (14), and constructed the first large-scale phosphoproteome analysis profile of HCC in vivo, which revealed a total of 22,539 phosphorylation sites on 5,431 proteins. Further analysis of HCC and WT mice revealed system-wide reprogramming of kinome activity, presenting a landscape depicting the phosphorylation dysregulation network in HCC, and ultimately enabled the identification of predominant kinase families in HCC, including Src family kinases (SFKs), PKCs, and Rhoassociated protein kinase 2 (ROCK2). We confirmed the activities of specific kinases in cell line models and clinical samples, and investigated inhibitor combinations targeting these kinases, achieving significant synergistic inhibition of HCC cell growth. In addition, we found Src Ser17 is involved in HCC cell migration and mediates ROCK2 kinase activity.

EXPERIMENTAL PROCEDURES
Animal Models and Clinical Samples-HBx transgenic mouse liver tissues (HCC samples) and corresponding wild-type mouse liver tissues (WT samples) were provided by Wang et al. (14). The mouse tissues were immediately frozen in liquid nitrogen after resection. The mouse tissues were transferred with plenty dry-ice, kept frozen at Ϫ80°C until used.
A total of sixty-two HCC and pairwise adjacent nontumor liver tissue samples were collected from Shanghai Zhongshan Hospital and Beijing Cancer Hospital after surgical resection. Patient information was supplemented as supplemental Table S7. Written consents approving the use of their tissues for research purposes after the operation were obtained from each patient. The study was approved by the institute research ethics committee at the Beijing Institute of Radiation Medicine. The tissues were diagnosed by senior pathologists and the biospecimens were collected from patients who were undergoing surgical resection and without prior chemo/radiotherapy for their disease. The tissue samples were immediately frozen in liquid nitrogen for storage after surgical resection. The specimens were shipped with plenty dry-ice and kept frozen at Ϫ80°C until used.
The human HCC and counterpart tissues (used only for Western blot analysis) were minced and lysed in T-PER TM Tissue Protein Extraction Reagent (Thermo Fisher Scientific) containing protease and phosphatase inhibitors (Thermo Fisher Scientific) on ice, followed by 1 min sonication (3 s on and 3 s off, amplitude 25%) to reduce lysate viscosity. The lysate was centrifuged at 10,000 ϫ g for 20 min, and the supernatant was collected for Western blotting assay.
Phosphoproteomics Sample Preparation and MS Analysis-Protein Extraction and Peptide Digestion-HBx transgenic mouse liver tissues (HCC group) and corresponding wild-type mouse liver tissues (WT group) were collected, ground with liquid nitrogen, and pooled from 5 individuals in each group (to minimize the biological variance in each group). Samples were stored at Ϫ80°C before analysis. Thirty milligrams of pooled tissues (ϳ 2.4 mg proteins) were lysed with lysis solution (4% SDS, 100 mM Tris/HCl pH 7.6, 0.1 M DTT, 5 mM sodium fluoride, 50 mM ␤-glycerophosphate, 1 mM sodium orthovanadate, supplemented with HALT protease and phosphatase in-hibitor mixture (Thermo Fisher Scientific, 78441)). Lysates were incubated at 95°C for 4 min, and then DNA was sheared by sonication to reduce lysate viscosity. Samples were centrifuged at 16,000 ϫ g for 10 min, and insoluble materials were removed. The protein in lysates was quantified relative to a bovine serum albumin (BSA) control.
Two milligrams of protein lysate from mouse tissues was digested with the endoproteinase Lys-C followed by trypsin according to the MED-FASP method (15). The trypsin digestion was halted by acidification with trifluoroacetic acid (TFA) to 1% (v/v) on ice. The tryptic digests were then centrifuged to remove insoluble materials, and the supernatant was transferred into clean tubes, vacuum dried and stored at Ϫ80°C for further analysis.
TiO 2 Phosphopeptide Enrichment-Phosphopeptides were enriched using titanium dioxide (TiO 2 ) (16,17). Briefly, a homemade micro-column (200 l) with a small plug of C 8 material was prepared as previously described (16). TiO 2 beads in tubes and TiO 2 microcolumns packed on centrifugation adaptors (GL Sciences, 5010 -21514) were individually equilibrated with 100 l of loading buffer (70% ACN, 5% TFA, 20% lactic acid (Sigma-Aldrich, L6661)). After centrifugation, supernatants containing the TiO 2 beads were removed. Peptides from all 6 subfractions in each experiment were diluted in 100 l of loading buffer and mixed thoroughly with 5 mg of TiO 2 beads, followed by incubation for 30 min. The mixtures of TiO 2 beads and peptides were then centrifuged at 800 ϫ g for 5 min, the supernatants were transferred into new clear tubes for a second round of phosphopeptides enrichment, and enriched TiO 2 beads were suspended in 100 l of loading buffer and transferred into micro-columns, then packed into centrifugation adaptors. Each column with enriched TiO 2 beads was washed four times with 120 l of loading buffer (150 ϫ g, 10 min), followed by 100 l of washing buffer 1 (30% ACN, 0.5% TFA) (250 ϫ g, 8 min), and two washes with washing buffer 2 (0.4% TFA, 80% ACN) (250 ϫ g, 5 min).
Phosphopeptides were eluted with 150 l of elution buffer (15% NH 4 OH (Wako, 013-23355), 40% ACN), and the first and second phosphopeptides enrichment of each subfraction were combined, dried in a speedvac and stored at Ϫ80°C.
LC-MS/MS Analysis-Enriched phosphopeptides were resuspended in 10 l of 10% FA and subjected to online reverse phase nanoLC-MS/MS analysis with 50% of sample loading using an Easy nLC 1000(Thermo Fisher Scientific, San Jose, CA, USA), coupled to a Q-Exactive mass spec-trometer (Thermo Fisher Scientific). Peptide samples were concentrated on a 2 cm trap column (100 m diameter) and separated on a 12 cm capillary column (75 m diameter), both packed in-house with 1.9 m C18 reverse-phase fused silica (Michrom Bioresources, Inc., Auburn, CA). The samples were separated at a flow rate of 0.6 l/min with a 71 min linear gradient from 5% to 30% mobile phase B (phase A: 0.1% formic acid in water, phase B: 0.1% FA in ACN), followed by a quick ramp from 30% mobile phase B to 95% mobile phase B within 1 min, where samples were held for 6 min before a quick ramp down; then, the column was re-equilibrated. Eluted peptides were analyzed with a Q-Exactive mass spectrometer (Thermo Fisher Scientific). The MS survey scan was analyzed over a mass range of 300 -1400 Da with a resolution of 70000 at m/z 200. The isolation width was 3 m/z for precursor ion selection. The automatic gain control (AGC) was set to 3e6, and the maximum injection time (MIT) was 60 ms. The MS2 was analyzed using data-dependent mode searching for the 20 most intense ions fragmented in the HCD. For each scan with a resolution of 17500 at m/z 200, the AGC was set at 5e4 and the MIT was 80 ms. The dynamic exclusion was set at 18 s to suppress the repeated detection of the same fragment ion peaks. The relative collision energy for MS 2 was set at 27% for HCD.
Experimental Design, Statistical Rationale and Data Deposition-Liver tissues from 5 individual mice in the HCC and WT groups were each pooled to minimize the biological variance in each group. We then performed three technical replicates on each pool to detect the difference in phosphoproteome between HCC and WT groups. Additionally, in order to minimize bias in phosphopeptide preparations, the sample preparation of HCC and WT groups were performed in parallel during each replicate, and the three independent technical replicates were performed at different times within two months and each was run through taking samples from pools, lysates extraction of tissues, FASP digestion of proteins, high-pH HPLC fractionation of peptides, enrichment of phosphopeptides and LC-MS/MS analysis workflow. After data filter for class I phosphosites, the quantitative values in two groups required at two valid values of three replicates in at least one group before missing data imputation and statistical comparison. To determine significantly regulated phosphosites, a two-sample t test was applied and a p value of 0.01 was set as the threshold to determine significance in two groups. Pathway enrichment analysis was performed with GengGO/MetaCore software (https://portal.genego.com), false discovery rate (FDR) cut-off 0.001, we further filtered the pathway enrichment result by requiring at least three proteins in a group and coverage at least 5%, and 0.001 cut-off of p value. For the GO analysis of proteins regulated by phosphorylation, the proteins that detected in this study but not regulated by phosphorylation were used as reference set. The Benjamini-Hochberg method was used to adjust the p values for multiple test correction and a cut-off less than 0.05 was used for GO category analysis.
All raw files generated from this study have been uploaded to iPROX with the identifier IPX0000776001 (http://www.iprox.org/page/ SDV015.html?subprojectIdϭIPX0000776001).
MS data Processing and Data Analysis-Raw MS data were analyzed using MaxQuant software (18) version 1.5.0.30 and searched against Mus musculus proteins from the UniProt sequence database (v. 140903) containing 51,552 entries including 24,686 reviewed protein entries (Swiss-Prot IDs) and 26,866 unreviewed protein entries (TrEMBL IDs). The database search was performed with the following parameters: oxidation of methionine, protein N-terminal acetylation, and serine, threonine, and tyrosine phosphorylation were searched as variable modifications, and cysteine carbamidomethylation was used as a fixed modification. The enzyme specificity was set to trypsin with maximum of 2 missed cleavages. The precursor mass tolerance was set to 20 ppm for the first search and then 4.5 ppm for the main search. The false discovery rate (FDR) for peptide, protein, and site identification was filtered to 1%, and the modified peptides required a minimum peptide length of 7 amino acids and a minimum Androm-eda score of 40. To match identifications across different technical replicates, the "match between runs" feature was enabled with a retention time window of 0.7 min, information contained in phosphopeptide and recalibrated retention time is used to transfer identifications with a sufficiently low FDR (in the range of 1%) between similar LC-MS runs (19). The run-specific normalization factors N j were introduced by a "delaying normalization" algorithm, which were determined in a nonlinear optimization model that minimized overall changes for all phosphopeptides across all samples. Then the phosphopeptide ion signals were summarized over fractions in each sample with the run-specific normalization factors N j (19).
Data analysis was performed using R software and Perseus (Max Planck Institute of Biochemistry, Munich). The "Phospho(STY).txt" file generated by MaxQuant was used for quantification data analysis. The "Reverse" hits were filtered and a phosphosite localization probability of at least 0.75 was used as threshold for phosphosite localization (15,644 phosphosites), and a delta score of 8 was used as threshold for phosphosite filter (15,603 phosphosites). The intensities of phosphosites were log2 transformed and triple replicates of each group (HCC or WT) were grouped, the phosphosites required be quantified with Ն 2 observations of 3 replicates in at least one group for following normalization, we further excluded the phosphosites that were quantified with only one measurement in a group but at least two in another group whereas the former were larger than mean of the latter. The remained phosphosites were processed with imputation of missing values, and quantitative analysis (11,008 phosphosites).
The normalization of the replicates was performed by the "normalizeQuantiles" function in the R package "LIMMA" which normalizes the intensities of phosphosites by applying quantile normalization for each replicates, this resulted in identical distribution of phosphosite intensity among the replicates of two groups. For the imputation of missing values, two situations were considered: (20) (1) the phosphopeptide of the corresponding phosphosite was truly presented at an abundance the mass spectrometer should be able to detect whereas it was not detected. (2) the phosphopeptide of the corresponding phosphosite was presented at a low abundance beyond the capability of the mass spectrometer or the phosphopeptide was not presented in the sample. Function "model.Selector" from the R package "imputeLCMD" was used as model selector for hybrid missing data imputation, which was applied in each group (HCC or WT) individually to determine the missing values should be imputed as a missing at random (MAR) or missing not at random (MNAR) value. The "model.Selector" function returned an estimated censoring threshold and flag indicators of rows (1 or 0) in each group, i.e. the rows in each group with flag indicator "1" would be imputed as MAR (case (1)), and the rows in each group with flag indicator "0" would be imputed as MNAR (case (2)). For the case (1), function "impute.MAR" with method parameter "MLE" (MLE-based imputation of missing data) from the "imputeLCMD" was used to impute the missing values. For the case (2), the left-censored missing data were imputed by normally distributed values (width 0.3, down-shift 3.5). To identify significantly regulated phosphosites between the HCC and WT groups, a two-sample t test was performed and a p value of 0.01 was set as the threshold to determine significance in two groups.
Functional analysis for gene ontology (GO) enrichment was performed using PantherDB bioinformatics tool (21). For the GO analysis of proteins regulated by phosphorylation, the proteins that detected in this study but not regulated by phosphorylation were used as reference set. The Benjamini-Hochberg method was used to adjust the p values for multiple test correction and a cut-off less than 0.05 was used for GO category analysis.
Pathway enrichment analysis was performed with GengGO/Meta-Core software (https://portal.genego.com) with false discovery rate (FDR) cut-off 0.001, we further filtered the pathway enrichment result by requiring at least three proteins in a group and coverage at least 5%, and 0.001 cut-off of p value. Protein-protein interaction network analysis was performed using DAPPLE, a previously described up to date interaction network of quality controlled and network-building algorithm. Briefly, the DAPPLE uses a refined database of qualitycontrolled PPIs (the inWeb database), (22,23) combining reported protein interaction from MINT, BIND, IntAct, KEGG annotated proteinprotein interactions (PPrel), KEGG Enzymes involved in neighboring steps (ECrel), Reactome and others. All human interactions were pooled and interactions in orthologous protein pairs passing a strict threshold for orthology were included. Each interaction was assigned a probabilistic score based on the neighborhood of the interaction, the scale of the experiment in which the interaction was reported and the number of different publications in which the interaction had been cited (24). Then the DAPPLE builds interaction networks among the proteins of interest, assessing the statistical significance of a number of network connectivity parameters as well as the connectivity among individual proteins to other seeds by performing a within-degree node-label permutation (22)(23)(24)(25). The p value was adjusted for multiple testing using a Bonferroni correc-tion method. The direct connections were considered, and DAPPLE evaluated the statistical significance of the protein interactions network by comparing it to 1000 randomized networks that are generated by shuffling protein names in an interaction-degree-preserving manner. The protein-protein interactions generated by NetworkAnalyst (26) was provided as supplement to results of DAPPLE. the following parameters were used: PPIs database from STRING was used, only experimental evidence was required, a confidence cut-off of 900 was applied (medium (400)-high (1000)). The "degree" centrality and "betweenness" centrality of nodes in the network were extracted by NetworkAnalyst and supplemented in "supplemental Data S2". The extracted protein-protein interactions were visualized in Cytoscape.
Kinase Substrate Relation Analysis-Overrepresented motifs were determined by analyzing the amino acid sequence motifs of regulated phosphosites in HCC with Motif-X (27) against the IPI mouse database, a minimum number of occurrences of 20 and a minimum 3.0-fold change enrichment compared with the background, the p value threshold was set to 0.000001 to maintain a low false positive rate in standard protein motif analyses.
To investigate kinase activity in HCC, the NetworKIN algorithm (version 3.0) (28 -30), which can predict substrates for 222 kinases, was used to interrogate the kinase-substrate relations in HCC. The NetworKIN algorithm first uses neural networks and consensus sequence motifs to assign phosphosites to one or more kinase families, then the protein-association networks from STRING were used to improve the specificity from the motif-based predictions. The NetworKIN was updated and tackled bias in network structure caused by over-studied proteins (30). The mouse sequences were matched to the homologous human sequences with BLAST, and the predictions of upregulated/downregulated and the nonregulated phosphosites were generated in the human background. To reveal kinase activity in HCC, individually prediction of possible upstream kinases in upregulated phosphosites, downregulated phosphosites, and nonregulated phosphosites in HCC were performed. Only kinases predicted both in upregulated/downregulated data and nonregulated data were used for following KS test interrogation. For each kinase, the distribution for prediction scores of upregulated or downregulated phosphosites was respectively compared with the distribution of prediction scores for nonregulated phosphosites using the nonparametric Kolmogorov-Smirnov test. We then filtered the kinases with a two-sided Kolmogorov-Smirnov test (p Ͻ 0.05). Then using one-side Kolmogorov-Smirnov test, we determined whether kinase activity was amplified or depressed in HCC (p Ͻ 0.05). We further focused on the kinases that were observed peptides detection in HCC. The kinases demonstrating significantly altered activity were mapped to the human kinome tree, which was adapted with Adobe Illustrator.
To construct the phosphorylation network in HCC, predicted active kinases (only kinases with experimentally observed peptides were used) and their substrates were employed, and the top 3% of substrates for each kinase were extracted to construct the network topology. Subcellular localization annotation of substrates was extracted from UniProt and the Human Protein Reference Database (HPRD). The phosphorylation network was visualized in Cytoscape and adapted in Adobe Illustrator.
Cell Culture and Expression Vectors-Human liver LO2 and hepatoma HepG2 cells were maintained in DMEM medium supplemented with 10% FBS, 1% 100 U/ml penicillin and 100g/ml streptomycin at 37°C in a humidified atmosphere containing 5% CO 2 . The full-length HBx-expressing vector pHAGEHBx was constructed in our laboratory. The pcDNA3.1-Src wild type and pcDNA3.1-S17A mutant plasmids were purchased from Generay, Shanghai, China.
Transient Transfection and Generation of Stable Cell Lines-HepG2 and L02 cells were transfected using Lipofectamine 2000 (Invitrogen, Carlsbad, CA) with plasmids according to manufacturer's protocol. For stable transfection, transfected HepG2 and L02 cells were maintained in DMEM containing 500 g/ml G418 (pcDNA3.1-Src or pcDNA3.1-Src-S17A) or 2g/ml puromycin (pHAGE-HBx) for 48 h post-transfection. After 3ϳ4 weeks of selection, individual colonies were isolated and expanded. The overexpression of target genes in these clones was confirmed by Western blotting.
CCK8 Proliferation Assay-Approximately 5000 cells were seeded in 96-well plates. The cells were cultured in 10% DMEM for 24 h, 48 h and 72 h. Subsequently, the absorbance was measured at 450 nm. This assay was performed once per day for four consecutive days. To evaluate the proliferation of cell lines in the presence of saracatinib(AZD0530), GSK429286A and sotrastaurin, ϳ100l of L02, L02-HBx, HepG2, and HepG2-HBx cells were seeded in 96-well plates. After overnight incubation at 37°C, saracatinib(AZD0530), GSK429286A and sotrastaurin were added to wells individually or in combination. Cells were further incubated for 24 h, 48 h, 72 h, and 10l CCK8 was added to each well. Absorbance was measured at 450 nm for 1 h after the cells were incubated with CCK8.
Cell Migration Analysis-For the transwell migration assay, 2 ϫ 10e4 serum starved (0.5% serum) cells were trypsinized and re-plated onto an FN-coated upper chamber membrane (8 m pore filter, Corning Costar, Cat:3422) of the transwell. The lower chamber of the transwell was filled with 2.5% serum containing DMEM. After incubation for 24 h, the filters were removed and the cells on the membrane were fixed with methanol. The migrated cells on the underside of the membrane were stained with 0.5% Crystal Violet. The dye was washed with water, and cells were examined by microscope.
Immunoprecipitation and Western Blotting-Whole cell extracts (500g) were incubated overnight at 4°Cwith the monoclonal antibody anti-HBx(5g), and then immunocomplexes were incubated with 50l (50%) of protein ASepharose beads for 8 h at 4°C. Beads were washed three times with immunoprecipitation buffer, and then bound proteins were eluted with 50l of SDSPAGE sample buffer and boiled for 5min. Cell lysates were run on an SDS-PAGE gel, and immunoblotting was performed according to standard techniques.

Phosphoproteome Profiles of HCC Tumor and Wild Type
Tissues-To characterize the signaling dysregulation in HCC, we applied quantitative phosphoproteomic strategies via a robust label-free quantification method. We isolated liver tis-sue samples from male HBx transgenic HCC mice and corresponding wild-type mice (supplemental Fig. S2), and pooled five liver tissue samples from each group to account for biological variation. After protein extraction, we digested the liver proteins, fractionated peptides using high-pH high-performance liquid chromatography (HPLC), and subsequently enriched for phosphopeptides via TiO 2 chromatography (16,17). We then analyzed the phosphopeptides by performing nanoflow liquid chromatography-tandem mass spectrometry (nLC-MS/MS). In total, we achieved deep coverage of the phosphoproteome by analyzing 48 MS runs from three replicates, resulting in the identification of 19,627 unique phosphopeptides corresponding to 22,539 unique phosphosites on 5,431 phosphoproteins (supplemental Fig. S1, supplemental Table S1-S4). Phosphopeptide identification confidence was high, with 95% of all peptides identified with an average absolute peptide mass error of 1.95 ppm (supplemental Fig.  S3A). Evaluation of this high-quality MS data is shown in supplemental Fig. S3.
Single-, double-, triple-, and higher phosphorylated peptides represented 83.3%, 13.9%, 2.4%, and 0.4% of total phosphopeptides, respectively ( Fig. 1B and supplemental Fig.  S3D). Depending on the phosphorylation site localization probability calculated by MaxQuant, we divided the observed phosphosites into four classes based on their probability score: 15,644 phosphosites were assigned to class I with a specific residue with high confidence (Ն0.75, supplemental Table S2), 3541 were class II sites (0.5-0.75), 3021 were class III sites (0.25-0.5), and 333 were class IV sites (Ͻ0.25). A total of 69% of all phosphosites were classified as class I sites with accurate quantification and localization (with a median localization confidence of Ͼ99.9%), which were only considered in further analysis (Fig. 1C). The distribution measurements of phosphorylated amino acids in class I were classified: phosphoserine (85.5%), phosphothreonine (13.5%), and phosphotyrosine (1.0%) (Fig. 1C), which is like previous estimations and recently reported phosphorylation amino acid distribution (31,32). Comparison of the data set in this study with the PhosphositePlus database (version 20170217) revealed that 5,552 of the identified 22,539 phosphorylation sites are identified for the first time whereas 2284 of the 15,644 class I phosphosites are previously unknown, implying a considerable phosphoproteome depth of this study. We further explored the phosphorylation sites detected in this study by comparing the data set with previous reports. The Phos-phophositePlus database provides a "disease-associated sites" data set, which included information curated from the literature about phosphorylation sites shown to correlate with specific disease states, especially with various cancers. There are 762 unique phosphorylation sites annotated to correlate with diseases in the data set, among which 108 phosphorylation sites (14% of the data set) were detected in our data set (and 91 of Class I phosphorylation sites). Considering the data set includes phosphosites related to various diseases and our data set was derived from a single experiments of HCC mouse tissues, that suggests a substantial sensitivity of our study.
Phosphorylation Regulation in HCC-For the quantitative analysis, the 15,644 class I phosphosites were further filtered with a delta score cut-off of 8, and only those that were quantified with Ն 2 observations of 3 replicates in at least one group, were kept for further normalization, imputation of missing values, and quantitative analysis (11,008 phosphosites). We normalized the MS signal intensities of phosphosites to ensure quantification comparable across triplicate experiments for the HCC group and the WT group. Using a twosample t test to the 11,008 highly confident phosphosites, we identified 1451 phosphosites from 1099 proteins that were statistically significantly regulated in HCC, specifically 1088 phosphosites that were upregulated and 363 phosphosites that were downregulated, and the difference between the means of HCC and WT group of each phosphosite was reported and plotted with the corresponding p value in Fig. 1D (supplemental Table S3), representing 9.3% of the phosphoproteome in HCC. When we compared the data set with a recently published work by Xu et al., in which the same mouse tissue samples were analyzed on proteome and 276 proteins were observed with significant abundance changes between HCC and WT mice (33). We only found 36 proteins that were regulated accordingly at phosphorylation and protein level (32 with upregulation at phosphorylation and protein level, 4 with down-regulation at phosphorylation and protein level). Because 1099 proteins were significantly regulated on phosphorylation level in HCC, which suggests that regulation on proteome in HCC was weak but prominent on phosphoproteome, regulation on phosphorylation may be partially attributed to the abundance difference between HCC and WT mice. In addition, recently published data from large-scale of phosphoproteome analysis of tissues and cancer cell lines both showed a weak correlation between protein abundance and phosphorylation levels (34 -39). Taken together with our results, these studies also highlighted the significance of phosphorylation regulation in cells. To identify relevant cellular components (CC) of the regulated proteins, we performed a gene ontology (GO) analysis and found that plasma membrane, cytoskeleton, and cytosol represented the dominant proportions of GOCC terms (Fig. 1E), whereas comparison of the proteins with regulated and unregulated phosphosites revealed an overrepresentation of proteins localizing to the cytoskeleton, cell junction and plasma membrane (Fig. 1F), which implies the significantly aberrant phosphorylation regulation in these components in HCC.
Phosphorylation of proteins within the plasma membrane is considered as an important approach to link extracellular and intracellular signaling. The cytoskeleton associated proteins are often phosphorylated by MAPKs and CDKs kinases, which are stimulus responsive and are effectors of signaling cascades (40), 16% of regulated phosphoproteins in our data localized to the cytoskeleton in HCC.
Phosphorylation regulation is crucial for the activation or de-activation of kinases. We detected 214 kinases, representing 3.9% of observed phosphoproteins, of which one-third demonstrated dysregulated phosphorylation in HCC. To determine whether the kinases with aberrant phosphorylation regulation were more likely to be connected to other kinases in modularity than by chance, we built the protein interactions using the previously described DAPPLE algorithm. It uses a refined database of quality-controlled protein-protein interactions, and a network-building algorithm with quality thresholds optimized by permutation tests to construct the protein interactions network (22)(23)(24). With which we constructed a network of 83 connections comprising the dysregulated kinases with refined quality-controlled protein-protein interactions (supplemental Fig. S4A), which suggests that HBx might induce the wide dysregulation of kinase interactions in HCC.
To assess the pathways involved with HBx in HCC, we performed pathway enrichment analysis (41) using 1,099 phosphorylation regulated proteins as seeds. We observed that proteins involved in VEGF-, Wnt-, and Hippo-pathways were significantly regulated by phos-phorylation in HCC (supplemental Fig. S4B). For example, 15 proteins involved in the VEGF pathway were significantly regulated at the phosphorylation level, whereas 20 regulated proteins were involved in the WNT signaling pathway. Moreover, this analysis also highlighted the cytoskeleton remodeling and cell adhesion pathways, which implies key roles for regulated phosphosites in tumor cell migration and invasion. Although the invasion and migration processes of tumors have still not been clearly clarified, cellular cytoskeleton remolding and adhesion have been demonstrated to be essential in the migration and invasion of cancer cells (42,43), and changes in the relative levels of cytoskeleton related proteins can specifically influence cell migration and invasion (44). We observed phosphorylation dysregulation on kinases in several pathways, which might indicate that kinases activities are regulated directly and signaling cascade dysregulation is induced (supplemental Fig.  S5, see supplemental Table S6). Interestingly, we detected 20 phosphatases, representing 1/3 of all observed phosphatases, that were dysregulated in HCC, suggesting that phosphorylation regulation may affect their activities. We also observed the altered phosphorylation of 81 transcription factors in HCC, implying that phosphorylation signaling may modify their activity or promote or block the enzyme to read genes.
However, only a few observed phosphosites were consistent with previously reported kinase-substrate relations, based on comparison with the PhosphositePlus database (45) of the 1451 regulated phosphosites on 1099 proteins, only 79 phosphorylation sites (derived from 70 substrates) were annotated phosphorylation by specific kinases with corresponding 170 kinase-substrate relations in the PhosphositePlus database (supplemental Fig. S6). The extracted kinase-substrate rela-tions may partially reflect the phosphorylation dysregulation in HCC. One example is the phosphorylation of JNK2 (MAPK9) on Y185, leading to its activation (46,47), then JNK2 phosphorylates stathmin (STMN1) at sites S38 and S25, promoting microtubule stabilization and affecting the activity of this protein in cellular microtubule cytoskeleton remodeling. Nevertheless, the majority of phosphosites cannot be traced back to corresponding relations with kinases, only ϳ 5% of phosphosites are annotated with corresponding relations in Phos-phositePlus. These detected kinase-substrate relations, provide only limited clues regarding kinase mediated signaling pathways in HCC.
Kinase-Substrate Relations in HCC-To further investigate kinase-substrate relations in HCC, we interrogated our large data set of regulated phosphosites to identify potential dysregulated kinases in hepatocellular carcinoma. Dysregulation of kinase activities could relay signals to corresponding substrates based on the recognition of specific residues. Therefore, we carried out two primary strategies to explore kinasesubstrate relations in HCC, using the Motif-X and NetworKIN algorithms (27,29). The combination of linear motif analysis with Motif-X to identify the overrepresented sequences, and NetworKIN algorithm to predict the responsible kinases, have been widely applied in various studies to identify the kinases likely responsible for the observed phenotypes (48 -53). Those algorithms have been thoroughly described elsewhere (27)(28)(29)(30). The identification of amino acid enrichment around phosphosites can be valuable for the identification of potential kinase activities. When we analyzed regulated sites in greater detail with Motif-X, (27) we identified 10 kinase motifs that were significantly enriched at upregulated phosphosites in HCC ( Fig. 2A and 2B). Among the extracted motifs, prolinedirected kinase substrate motif (XpSPX) is overrepresented in HCC, which is phosphorylated by mitogen-activated protein kinases (MAPKs), including ERK1 and ERK2 kinases. We also observed the enrichment of sequences with basophilic motifs that associate to the AGC family kinases, including the RRXpS, RSXpS, RXXpS, KXXpS, RXpS motifs, which are phosphorylated by PKCs. And the RRXpS, RSXpS, RXXpS and RXpS motifs link to PKA kinase, which also belongs to the AGC family kinases, whereas the acidophilic motifs include pSDXE, pSXXD, and pSXXE, those are phosphorylated by acidophilic kinases (casein kinase I and II, CK1, and CK2). All together suggest that portions of these sites were targeted by one or both kinases. To achieve greater detail in our predictions regarding the kinases that are responsible for the observed phosphosites alterations, we employed the NetworKIN algorithm, which predicts kinase-substrate relations based on conserved linear motifs and the cellular network context (29,30). The NetworKIN algorithm first identifies a consensus kinase motif in the phosphoprotein sequence to assign to one or more kinase families, then the contextual information in STRING was used to identify upstream kinases, causing a higher accuracy of this algorithm (29,30). We obtained the FIG. 2. Kinase activity prediction in HCC. A, Sequence motifs for phosphosites that were upregulated in HCC; motif scores were annotated (according to Motif-X). B, Bar plot of motif fold changes for upregulated phosphosites (compared with the background IPI mouse database). C, Kinases predicted to be regulated in HCC; highlighted are kinases that are more active in HCC (red; p Ͻ 0.05) or in wide-type mice (blue; p Ͻ 0.05). kinase prediction score distribution for regulated phosphosites and for all other phosphosites in HCC. Then, we used the Kolmogorov-Smirnov test to determine whether individual kinases activities were amplified or depressed in HCC by comparing the distribution of upregulated or downregulated kinase-substrate relations with unregulated kinase-substrate relations for each kinase (54). Ultimately, we inferred significant amplification activity for 24 kinases and the depression of 8 kinases in HCC (p Ͻ 0.05; see Table S5).
Further, when we mapped these kinases to the human kinome tree (55), we obtained a dramatic view of the kinome activity landscape in HCC. Members of specific kinase families are significantly dysregulated in HCC (Fig. 2C, adapted from (55). The activities of six members of the AGC family, ten members of the TK family, five members of the STE family, two members of the CMGC family and PDK2 were associated with elevated activity in HCC (See supplemental Table S5); as expected, these kinases demonstrating elevated activity included PKCs, ROCK2, and MAPKs confirming the Motif-X findings indicating that basophilic-motifs were overrepresented. Among those kinases demonstrating elevated activities, one member of the PKCs (PKC), as well as five members of the TK family (ITK, Lyn, BLK, BRK, and SRM) and PDK2 have not yet been functionally described in HCC, whereas PKC␦ and PKC␤ were previously reported downregulated in HCC (56).
Changes in kinase phosphorylation may either lead to the activation or inactivation of theses kinases. Therefore, we investigated these kinases in greater detail. We extracted the kinases observed in the high confident 11,008 phosphosites and applied two-sample t test on them (cut-off of p value 0.05). Interestingly, we found that six phosphosites on three members of PKC family and four members of the SFK family were also observed to experience amplified levels of phosphorylation, specifically pT641 and pS660 on PKC␤II (prkcb, iso2), pS504, pS643 and pS652 on PKC␦ (prkcd), pS675 on PKC (prkh), pS206 on Fgr, pS34 and pS460 on HCK, pY316 on Lyn, and pS17 on Src. All phosphosites detected in PKCs localized in their C-terminal kinase domains. All PKCs have a conserved kinase core (activation loop) and C-terminal extension that contains two conserved phosphorylation sites: the turn motif and hydrophobic motif (57,58). pT641 (in turn motif) and pS660 (in hydrophobic motif) of PKC␤II were significantly upregulated in HCC, which was essential for activity, and have been both reported to have various effects on activity (58,59). Mutation of T641 of PKC␤II results in an inactive protein, and studies on PKC␤II have shown that phosphorylation at S660 affects its affinity for Ca 2ϩ and for phosphatidylserine (60). All three sites observed on PKC␦ locate in the kinase domain, including Ser504, Ser643, and Ser652. Ser643 of PKC␦ is an important autophosphorylation site for controlling its enzymatic activity and biological function (61). Ser504 of PKC␦ locates in the activation loop of PKC␦ (DFG-APE), phosphorylation on which can provide a biochemical read-out for kinase activity (62)(63)(64). For PKC (PRKCH), Ser675 (in C-terminal hydrophobic motif) was observed upregulated in HCC, which is equivalent to Ser662 of PKC␦, Ser660 of PKC␤II, and Ser657 on PKC␣, Ser662 unphosphorylated PKC␦ was shown has less than one-tenth of the activity (58,65). Based on founds on Ser657 on PKC␣, (66) Ser662 of PKC␦, (58,65) and S660 of PKC␤II, (60) mutation of the hydrophobic site to an alanine residue in PKC reveals that phosphorylation of this site stabilizes the enzyme, (57,60,67) therefore, Ser675 on the hydrophobic motif of PKC may be an autophosphorylation site that affects its kinase activity (57). Moreover, PKC␦ and PKC are also predicted to be activated in HCC. We next investigated the altered phosphorylation levels in SFKs. The phosphorylated amino acid residue (Ser206) localizes in the SH2 domain of Fgr; phosphorylation of which may modulates the binding affinity and specificity for SFKs (68,69). Ser 460 of HCK localizes in the kinase domain of this kinase; Tyrosine residue (Tyr316) of Lyn also localizes in the kinase domain of this protein which is equivalent to Tyr338 of Src and phosphorylation of this site was reported autophosphorylation by Src,(70) the phosphorylation of the domain is also reported autophosphorylated by SFKs, (69), and was predicted to be dominated by SFKs themselves in our analysis. The phosphorylation of Ser17 of Src, which is a specific residue in the N terminus of Src, will be discussed in greater detail in the following section. Elevated activity of tyrosine kinases may explain the overrepresentation of phosphorylation on protein tyrosine residues, phosphorylation on tyrosine residues in PIK3R1, FCER1G, PTK2, TGFB1I1, and CTNNA1 was imputed by SFKs in our NetworKIN analysis. Based on the observed phosphorylation changes to kinases and the analysis of predicted kinase activities described above, we concluded that PKCs and SFKs might be two predominant kinase families that are activated by HBx in HCC.
Further, several members of the AGC family including GRKs, together with CaMKI␣, and BCKDK were found to exhibit decreased activities in HCC.
Construction of a Phosphorylation Network in HCC-Aberrant kinase activity triggered a range of signal cascades. To better understand the phosphorylation landscape in HCC, we focused on the upregulated phosphosites and constructed a network comprising the top 3% of predicted relations for each kinase, only kinases with experimentally observed peptides were used. Thus, this distinct network comprised 117 proteins corresponding to 156 relations and demonstrated alterations in response to phosphorylation (Fig. 3, see supplemental Data S1). The rebuilt network revealed several kinase-substrate relations that have been described in previous studies. For example, Tyr76 on FCER1G was reported to be phosphorylated by Lyn in PhosphositePlus and this relation was also confirmed with predicted relation in our study, as was Tyr576 on PTK2 by Src, and Tyr185 on MAPK9 by MAP2K7. These confirmations also highlighted the high prediction accuracy of the methods we used. Further, many previously undiscovered kinase-substrate relations were present in the network. The dominant kinases and their substrates also reflect phosphorylation within a cellular context and exhibited connectivity between nodes. In contrast to the notion of linear signaling transmission from cell membrane receptors, with the gradual relay of signals to effectors in the cytoplasm and to transcription factors in the nucleus, we found that individual kinases are spatial components of cells. Dysregulated substrates include adherent and anchoring junction proteins on the cell membrane, several kinases in the cytoplasm, actin and microtubule proteins in the cytoskeleton, and transcription factors in the nucleus. The rebuilt network exhibits a cellular view of phosphorylation dysregulation in HCC.
When we examined the related biological processes of kinase substrates, we found that some kinases are in involved in certain biological processes, whereas some other individual substrates are more likely to be involved in crosstalk between kinases. Most substrates for the PDK2 and ROCK2 kinases are partially shared with MAP2K7, PKAs, and the PKCs kinases; further, several substrates of them are involved in cytoskeleton remodeling. Intriguingly, PKCs and SFKs are involved in crosstalk with other kinases, which complicated the signal network. In addition, to further dissect the networks, we built the protein interactions using the previously thoroughly described DAPPLE algorithm (22)(23)(24)(25). That generated a highly confident interconnected proteinprotein interaction network for these proteins (Fig. S7, see supplemental Data S2). We observed concentrated interactions between the SFKs and other kinases, particularly between the kinase Src and other nodes (See also SP material regarding Cytoscape), and most predominant kinases in HCC interacted with Src by direct interaction, which indicates a specific role for SFKs in HCC. These kinases may play a central role in the development of HCC.
Kinase Activities in HCC and Inhibitor Combinations Targeting Predominant Kinases-The above analysis indicates that the activities of the kinase members of the PKC, SFK families and ROCK2 are amplified. We hypothesized that this amplification is induced by HBx in principle, and dysregulated kinase signaling supposedly contributes to the progression of HCC. To examine this hypothesis, we constructed two HBx stable transfection cell lines, LO2-HBx and HepG2-HBx (Fig. 4a). The amplified activities of PKC, PKC␦, ROCK2 and Src induced by HBx were validated in these two cell lines by Western blotting (Fig. 4B and 4C). Further, the activities of PKC, PKC␦, and ROCK2 were confirmed in clinical HBV positive HCC samples (Fig. 4D).
In our above analysis, we found that ROCK2 and PKCs modified a large proportion of substrates involving multiple characteristics of tumor cells, including the regulation of adhesion junction and cytoskeleton remodeling, whereas SFKs are more likely to be involved in kinase crosstalk or concentrated in direct interactions. These kinases also share certain specific substrates or regulate one another at specific sites, demonstrating a complex network rather than a linear cascade signaling pathway in vivo (Fig. 3, supplemental Fig. S7). Therefore, these activated kinases may act as core kinases that predominant in the progression of HCC, and appropriate control of their activities might delay or abrogate the progression of HCC. In addition, Src has been shown to be a signaling hub mediating various pathways in cancers (71,72), and is a key node in our reconstructed topological network (Fig. S7). Small molecule inhibitors targeting the activation of kinase signaling pathways have been included in the arsenal of agents for the treatment of different tumors. Considering that numerous proteins that interact with Src are related to the progression of HCC, targeting Src might inhibit related kinase cascade signaling pathways. To investigate this hypothesis, we tested the effects of SFKsI (an inhibitor of SFK kinases) and monitored the activities of Src, PKC, PKC␦, and ROCK2. Intriguingly, the activities of these four kinases were all decreased upon treatment with SFKsI (Fig. 4E), indicating potential downstream roles for these kinases in Src mediated signaling in HCC.
Our analysis and validation indicated that SFKs and PKCs are two predominant kinase families in HCC and suggested that ROCK2 is a core kinase hub in HCC. Therefore, we tested a PKCs inhibitor (sotrastaurin) and a ROCK2 inhibitor (GSK429286A) alone or in combination with inhibitors specifically targeting SFKs (saracatinib or PP2) to explore their inhibitory effects on cell growth (Fig. 4F). Sotrastaurin and GSK429286A were ineffective as single agents, but both were significantly synergistic to inhibit cell growth in combination with saracatinib or PP2. However, the combination of sotrastaurin and GSK429286A had no synergistic effects on the inhibition of cell growth. Further, we tested PKCs inhibitor (sotrastaurin) and a ROCK2 inhibitor (GSK429286A) alone or in combination with inhibitors specifically targeting SFKs (saracatinib) on control cells to explore their inhibitory effects on cell growth (supplemental Fig. S8). As shown in supplemental Fig. S8, only GSK429286A was slightly synergistic to inhibit LO2 cell growth in combination with saracatinib at 72 h, whereas only sotrastaurin was slightly systematic to inhibit HepG2 cell growth in combination with saracatinib at 72 h. Taken together, these data indicate that Src acts as a signaling hub in HCC and, together with PKCs or ROCK2, promotes the proliferation of tumors in HCC, whereas other signaling pathways mediated by Src also determine and maintain cell growth, even under circumstances involving the inhibition of PKCs and ROCK2 activities.
Src Ser17 is Involved in Cell Migration and Affects ROCK2 Tyr256 Phosphorylation-Our results indicated that Src is one of the predominant kinases in HCC. Src is one of the most extensively investigated proto-oncogenes, and decades of research have revealed its role in cell proliferation, differentiation, apoptosis, migration, and invasion (73,74). Recent work revealed that Src represents a universal signaling hub and is involved in multiple mechanisms in RTK-targeting therapies (72,75,76). However, Src mutations are seldom found in cancers, indicating that its high activity in tumors might be associated with dysregulation at the protein or PTMs levels. Because the first study identifying the Src structure in 1997 (77), the regulatory mechanisms of the C-terminal tail of Src have been extensively studied, and recent studies have shown the lipid binding capabilities of the SH3 domains as well as the phosphorylation of the unique SFKs domain, which affects SFKs activity (78 -80). Within the N-terminal sequence of Src, which contains a specific residue of Ser17 that can distinguish Src from other SFKs, which was shown to be upregulated in HBx transgenic cell lines and clinical HBV positive HCC tissues(T) compared with adjacent nontumor (NT) tissues ( Fig. 5A and  supplemental Fig. S9A). Kaplan-Meier analysis showed that patients with HBV-positive HCC and high phosphorylation levels of pSrc (Ser17) were associated with shorter disease-free survival rates (Fig. 5B). Moreover, high phosphorylation levels of pSrc (Ser17) were also significantly associated with vascular invasion (p Ͻ 0.01, see supplemental Table S7).
To investigate the roles of Src Ser17 phosphorylation in liver cancer cells, we carried out a transwell assay to examine migration ability. We found that the cell migration ability increased 2-fold in wild-type Src compared with the control group or Src-S17A ( Fig. 5C and 5D). There was no difference between control cells and Src-S17A, which suggested that the phosphorylation of Src Ser17 was likely pivotal for cell migration, potentially ascribed to the activation of Src on Ser17 and Src downstream signaling.
To assess whether Src Ser17 mediates Src kinase activity, we individually overexpressed either wild-type Src or mutant Src-S17A in HepG2 and L02 cells. The contribution of endogenous Src in both cells was negligible because of very low expression levels compared with transfected constructs (Fig. 5E and 5F). The Src protein kinase is composed of a C-terminal tail whose phosphorylation status determines its tyrosine kinase activity. Src can be deactivated by phosphorylation at Tyr530 and activated by autophosphorylation at Tyr419, and the N-terminal SH4 domain, which contains myristoylation site, is important for membrane localization (74,81). Western blotting data showed that there were no significant differences in the phosphorylation status of on Tyr419 or Tyr530 in wild-type Src compared with Src-S17A in both cells. This indicates that phosphorylation at Ser17 on Src appears to exert no direct effects on phosphorylation regulation of Src kinase active sites (Tyr419, Tyr530).
Given that Src Ser17 is involved in cell migration, we further examined ROCK2, a kinase that is involved in tumor invasion and activated by phosphorylation of the Tyr 256 residue (82,83) and, was demonstrated to be a predominant kinase in HCC in our analysis and validation ( Fig. 2C and Fig. 4B-4F). Surprisingly, there was obvious upregulation of Tyr256 phosphorylation on ROCK2 in wild-type Src compared with Src-S17A or the control cell lines. The results suggested that Ser17 might affect the phosphorylation of Tyr256 on ROCK2 via other mechanisms such as cellular relocalization or protein-protein interactions (80). DISCUSSION Here, we present the deepest in vivo system-wide quantitative phosphoproteomic study of HCC to date, in which we analyzed phosphoproteome dysregulation. These were 1493 phosphosites observed with significant regulation in HCC and compared with a recent work by Xu et al. in which the same HCC and WT tissues were interrogated on proteome, only 36 proteins were found with accordingly regulation at protein and phosphorylation levels, suggesting that phosphorylation regulation in HCC may be weakly attributed to the corresponding protein abundance. Differential phosphorylation can reflect changes in protein abundance, as well as changes in a site phosphorylation. To distinguish these factors, Huttlin et al. presented the first systematic and comprehensive investigation on the proteome and phosphoproteome in nine mouse tissues and found that many proteins are regulated by phosphorylation independently of their expression (34). They further showed little correlation between protein abundance and phosphorylation levels, either for the entire data set or for individual proteins (34). In addition, recently published data from large-scale of phosphoproteome analysis of tissues and cancer cell lines both showed weak correlation-between protein abundance and phosphorylation levels (34 -39). Taken together with previous studies, our work highlighted the significance of phosphorylation regulation in HCC.
Our approach identified the dysregulated phosphorylation of various proteins in HCC, including kinases, phosphatases, TFs and others. One-third of the detected kinases (75 in 214 detected kinases) and phosphatases (20 in 60 detected phosphatases), together with 81 TFs, were observed to be dysregulated in HCC, indicating that phosphorylation events on these proteins may be essential for tumorigenesis and determine the progression of HCC. Although many proteins have been implicated in processes in HCC, the functions of the majority remain undefined, which expands our knowledge of phosphorylation mediated regulation in HCC. Alterations to cortical cytoskeleton organization are also likely to play a role in RTK signaling (84,85). Our analysis also highlighted the dysregulation of cytoskeleton remodeling and cell adhesion pathways in HCC, and we observed a high proportion of phosphorylated proteins localized to the cytoskeleton, which suggests that the cytoskeleton might not only be involved in migration and invasion but, may also serve as a phosphorylation signaling transduction hub.
Although ϳ10,000 phosphosites can be detected within a single MS-based experiment, relatively unexplored kinasesubstrate relations have impeded our ability to undertake a systematic exploration of phosphorylation networks. Identify-

FIG. 5. Src Ser17 is involved in cell migration and affects the phosphorylation of ROCK2 Tyr 256 in HCC.
A, Upper: phosphorylation of Src S17 and elevated protein levels after stable transfection with an HBx-expressing vector in liver cell lines. Lower: relative fold change in band intensities normalized against GAPDH from triplicate experiments. B, Src S17 phosphorylation levels correlate with disease-free survival in patients with HCC. C, HepG2 cells stably transfected with wild-type or mutant Src (S17A) plasmids exhibit different migration abilities. D, Comparison of wild-type Src with control cells or the mutant group using Student's t test, error bars indicate meanϮS.D., *p Ͻ 0.05. E, Expression of Src but not mutant Src (S17A) results in higher kinase activity and directly phosphorylates ROCK2 at Y256. F, Normalized fold changes in Western blotting band intensities (E).
ing the responsible kinases for distinct phosphosites is crucial for understanding essential signaling for normal cells or tumors. The latest version of NetworKIN (version 3.0) can predict substrates for 222 kinases, that is less than half of the human kinome, thus upstream kinases of some substrates cannot be inferred by current version NetworKIN. When matching phosphosites across species, the less conserved phosphosites may be eliminated in further analysis and lead to missing of information. Besides, some predictions from NetworKIN are based on indirect probabilistic contextual information, a substrate may be linked to kinases with indirect probabilistic association. In this procedure, the STRING transfers evidence between target species, this may introduce bias when using NetworKIN to analyze phosphoproteome data derived from mouse. However, the systematic analysis of the phosphoproteome and strict statistical analysis would minimize this issue. Nevertheless, during the last decade, several algorithms (28 -30, 86) have greatly facilitated our understanding of kinase-substrate relations in various studies (35,49,50,52,53,87). With the aid of NetworKIN algorithm, investigation of upstream kinases with dysregulated phosphosites via our data analysis strategy revealed the reprogramming of kinase activities in HCC, and, together with the reconstruction of a phosphorylation signaling topological network, we obtained a broad view of kinase centric phosphorylation signaling in HCC, which provides an opportunity to identify the core factors that contribute to the genesis and development of liver cancer.
PKCs, ROCK2, and SFKs are extensively involved in signaling in malignant tumors. Studies investigating PKC proteins in the early 1980s provided the first solid proof that kinases are involved in carcinogenesis (88,89). However, subsequent studies have uncovered functionally antagonistic roles for PKC isozymes, and this effect appears to be isozymes-dependent and cell type-dependent under certain circumstances (90). For instance, PKC is generally considered to be a pro-survival factor whereas PKC␦ is considered a growth inhibitor, and the activation of PKC␦ leads to apoptosis in hepatocellular carcinoma (91)(92)(93), whereas recent evidence indicated that PKCs are tumor suppressors, challenging the ingrained concept of a role for PKCs in cancer (94). The PKC␤ and PKC were reported downregulated in HCC in previous study (56). Based on our analysis and validation, the PKC was found overactivated in HBV-positive HCC for the first time. ROCK2, another AGC kinase that is significantly overexpressed in various cancers and is involved in tumor invasion in HCC (82). Targeting ROCK2 led to the inhibition of cholangiocarcinoma tumorigenesis and metastasis (95). SFKs comprise the largest group of nonreceptor tyrosine kinases, decades of research have revealed key roles for SFKs in many aspects of tumor progression, including tumor proliferation, survival, adhesion, and metastasis (50,73). Importantly, our results suggested that HBx can indeed increase the activities of these kinases both in cell lines and clinical liver cancer tissues; that also demonstrated the precision of our analysis. Although we only focused on partial kinases exhibiting amplified activities in HCC, the roles of other dysregulated kinases are worth further exploration.
Targeting oncogenes in cancers with single agents has largely failed, which may be attributable to the incomplete inhibition of activated kinase pathways or the reactivation of key downstream effectors via the rewiring of signaling networks by tumor cells, particularly the RAF-MEK-ERK pathway (97,98). Inhibitors targeting RAF kinases can lead to "paradoxical activation" of CRAF, even resulting in tumorigenesis (99,100). Another explanation is that, "pathway addiction" is more likely to be present in cancers than "oncogene addiction," in which tumors switch their signaling to develop resistance (101,102). Sorafenib is a multikinase inhibitor that has been adopted to treat advanced-stage HCC and is able to extend overall survival with a median ranging from 8 months to 11 months (96). However, in the past decade, there has been a lack of more effective agents to treat HCC, (10,11). Further, it is also difficult to identify the core factors that contribute to HCC progression, potentially resulting in missed opportunities to select more suitable targets for treatment.
Here, with the interrogation of kinase-substrate relations, we were able to determine the predominant kinases that might be suitable for targeted cancer therapies. In accordance with this, the kinase Src is a downstream node common to multiple resistance pathways in various tumors, and, with our strategy, we showed that Src is a crucial node and both PKCs and ROCK2 are predominant kinases in our reconstructed phosphorylation network in HCC. The combinations of kinase inhibitors targeting Src and PKCs or ROCK2 both achieved significant and synergistic inhibition of cell growth, indicating that targeting "pathway addiction" with inhibitor combinations a potential approach for cancer therapies, under which circumstances multiple kinases are activated and resistance could arise via alterations in kinase signaling. This also raises the prospect of specifically targeting tumors at multiple levels, interfering and preventing tumors from undergoing "pathway addiction," delaying tumor-related processes or even defeating tumors.
Moreover, we analyzed and demonstrated that phosphorylation at Src Ser17 likely affects the phosphorylation of Tyr256 on ROCK2 and promotes cell migration. Src is responsible for the downstream activation of ROCK2 on Tyr256, and Networ-KIN predicted that Src Ser17 can be phosphorylated by ROCK2 and PKCs, suggesting a potential mechanism involving Src activation via a positive feedback loop and elevation of the Src/ROCK2 signaling cascade, affecting cell migration, invasion, and tumor metastasis. Activation of Src by Ser17 phosphorylation potentially be a novel mechanism in HCC progression, and targeting Src Ser17 phosphorylation with small molecules may be a strategy to intervene in the progression, which must be investigated further. This also dem-onstrates that our approach could uncover detailed molecular mechanisms of kinases.
These data are a valuable resource for future research examining HBx-induced HCC. In principle, this approach can be extended to more complex networks involving temporal and spatial signaling transduction, to monitor signaling shifts in tumor progression under specific circumstances by tracing the signaling network of predominant regulated kinases. Importantly, the strategy we present is generally applicable to any type of tumor. With the identification of predominant kinases in pretreatment biopsies of various tumors, and corresponding suitable reference samples (i.e. para tissues), kinase activity can be targeted and extended to personalized medicine and the adoption of appropriate combinations of agents to benefit individual patients, freeing many from unnecessary suffering under ineffective protocols.