Patients with HCC (n=4) were enrolled in the present study, and their tumor and adjacent non-tumor tissue samples were collected between May 2017 and July 2018 at the First Affiliated Hospital of Zhengzhou University (Zhengzhou, China). The distance between the adjacent non-tumor and tumor tissues was 1–5 cm. All patients were males and had hepatitis B virus infection. The mean age of the patients was 43 years old (SD, 2.8; age range, 41–46 years). The pathological grading of all 4 cases was moderate differentiated HCC according to the World Health Organization HCC grading classification (38). The diagnosis of HCC was based on the histological results of liver surgical resections. All patients provided written informed consent, and the study was approved by the Ethics Committee of Zhengzhou University (Zhengzhou, China).

To successfully extract proteins from liver tissues, the samples were first ground using liquid nitrogen and subsequently transferred to 5-ml centrifuge tubes. Subsequently, lysis buffer (8 M urea, 1% Protease Inhibitor Cocktail and 1% phosphatase inhibitor cocktail; MedChemExpress) was added at a ratio of 4:1 to the tissue. Subsequently, the lysate was sonicated for 3 sec with a 5 sec interval for three times using a high intensity Sonotrode set (Ningbo Scientz Biotechnology Co., Ltd.) at 4°C and then centrifuged at 12,000 × g for 10 min at 4°C. The supernatant was collected and the protein concentration was determined using a BCA protein assay kit (Beyotime Institute of Biotechnology).

Disulfide bonds of proteins were reduced using dithiothreitol (Sigma-Aldrich; Merck KGaA) at a final concentration of 5 mM for 30 min at 56°C. Proteins were alkylated using 11 mM iodoacetamide (Sigma-Aldrich; Merck KGaA) for 15 min at room temperature in the dark. Protein solutions were then diluted with 100 mM triethylammonium bicarbonate buffer (TEAB; Sigma-Aldrich; Merck KGaA) to an urea concentration <2 M. Trypsin was added at a 1:50 trypsin-to-protein mass ratio, and the sample was incubated overnight at 37°C. Subsequently, a second aliquot of trypsin was added at a 1:100 trypsin-to-protein mass ratio for 4 h at 37°C.

The digested peptides were desalted using a Strata X C18 SPE column (Phenomenex, Inc.) and were lyophilized in a vacuum centrifuge at 300 × g for 6 h at room temperature. Lyophilized peptides were reconstituted in 0.5 M TEAB buffer and were labeled using a TMT kit according to the manufacturer's protocol (cat. no. 90406; Thermo Fisher Scientific, Inc.). TMT reagents were dissolved in acetonitrile and then added to the sample. The labeling reaction was conducted for 2 h at room temperature. Subsequently, the labeled peptides were pooled and desalted, followed by lyophilization in a vacuum centrifuge, as aforementioned.

The phosphopeptides were fractionated via high-pH reverse-phase high-performance liquid chromatography (HPLC; Agilent 1260 HPLC; Agilent Technologies, Inc.) using a Thermo BetaSil C18 column (5-µm particles; 10-mm ID; 250-mm length; Thermo Fisher Scientific, Inc.). The column temperature was 35°C. The detection wavelength was 214 nm. The sample quantity was 1 ml. The mobile phase A was water and acetonitrile (98:2, v/v). The mobile phase B was water and acetonitrile (2:98, v/v). The elution gradient was as follows: 0–5 min (5% mobile phase B); 5–10 min (5-8% mobile phase B); 10–67 min (8-32% mobile phase B); 67–69 min (32-95% mobile phase B); 69–80 min (95% mobile phase B); 80–85 min (95-5% mobile phase B). The flow rate was 1 ml/min. The fractionated peptides were first suspended in loading buffer (50% acetonitrile and 6% trifluoroacetic acid) and then incubated with pre-treated immobilized metal ion affinity chromatography (IMAC) microspheres (J&K Scientific Ltd.). Phosphopeptide-bound IMAC microspheres were washed once with the loading buffer, followed by three washes with the washing buffer (30% acetonitrile and 0.1% trifluoroacetic acid). Subsequently, the phosphopeptides were eluted using 10% ammonium hydroxide, desalted with C18 ZipTips (EMD Millipore) and dried by vacuum centrifugation for LC-MS/MS analysis.

The enriched phosphopeptides were analyzed using an EASY-nLC 1000 UPLC system coupled in line with an Orbitrap Fusion™ mass spectrometer (Thermo Fisher Scientific, Inc.) via a nanospray ionization source in positive mode. Each sample was loaded onto a reversed-phase analytical column (75-µm inner diameter; 15-cm length) packed in-house in solvent A (0.1% formic acid) and resolved at a 400-nl/min flow rate using the following gradient: 6–23% solvent B (0.1% formic acid in 98% acetonitrile) for 26 min, 23–35% solvent B for 8 min and 35–80% solvent B for 3 min, then holding at 80% solvent B for the last 3 min. The electrospray voltage was set to 2.0 kV. The nitrogen gas temperature was 320°C and the nebuliser pressure was 15 psi. MS was operated in the data-dependent acquisition mode. The full scan MS1 spectra range was 350–1,800 m/z. The resolution in the Orbitrap was 60,000. The MS2 scan was measured at a resolution of 30,000 at an m/z of 100. The 20 most intense precursor ions were selected and fragmented by higher-energy collision-induced dissociation with a normalized collision energy at 35%. Automatic gain control was set at 5×10⁴ to determine the target value. The signal threshold was set to 5,000 ions/sec and the maximum ion injection time was set to 200 msec. Dynamic exclusion was set at 15 sec to avoid the repeated detection of the same fragment ion peaks. The mass spectrometry proteomics data were deposited in the ProteomeXchange Consortium via the Proteomics Identification Database (39) (dataset no. PXD013934). For the quantitative analysis, statistical significance was determined using unpaired Student's t-tests. P<0.05 and a differentially expressed ratio >1.5 or <0.67 was considered to indicate a statistically significant difference.

The LC-MS/MS data were analyzed using Maxquant v.1.5.2.8 (40). The data were searched against the SwissProt Human database (https://www.uniprot.org/taxonomy/?query=Homo+sapiens+%28Human%29+&sort=score) (20,130 entries) concatenated with its sequence-reversed database. Contaminants frequently observed were also added. Trypsin/P was the proteolytic enzyme, and two missed cleavages were allowed. Carbamidomethylation of cysteine was defined as a fixed modification, and oxidation of methionine, acetylation of N-terminal protein and phosphorylation of serine, threonine and tyrosine were specified as variable modifications. The minimum peptide length was set to 7 amino acid residues, and the maximum number of modifications per peptide was set to 5. Quantification was performed using TMT 10plex (Thermo Fisher Scientific, Inc.). The mass tolerance for precursor ions was set to 20 ppm in the first search and to 5 ppm in the main search, while for fragment ions, the mass tolerance was set to 0.02 Da. Regarding the level of proteins and peptide spectrum matches, the maximum false discovery rate was adjusted to 1%.

The UniProt-Gene Ontology annotation (GOA) database (http://www.ebi.ac.uk/GOA/) was used to detect the significantly enriched Gene Ontology (GO) terms. Since some identified phosphoproteins were not annotated on the UniProt-GOA database, the InterProScan software (version 5.47-82.0; http://interproscan-docs.readthedocs.io/en/latest/Introduction.html) based on the protein sequence alignment method was used. Subsequently, the identified phosphoproteins were automatically assigned to three categories of biological processes, cellular components and molecular functions. Wolfpsort software (version 0.1; http://wolfpsort.hgc.jp/) was used to predict subcellular localization. The InterPro database (https://www.ebi.ac.uk/interpro/) was used to analyze the functional domains of differentially expressed proteins. Pathway enrichment analysis was conducted using the Kyoto Encyclopedia of Genes and Genomes (KEGG) database (https://www.kegg.jp/). A two-tailed Fisher's exact test was used to filter the enrichment result to a corrected P<0.05. The protein-protein interaction (PPI) network was obtained using the Search Tool for the Retrieval of Interacting Genes/Proteins (STRING; http://www.string-db.org/) database. To determine the significance of the differentially expressed phosphoproteins against all the identified phosphoproteins, a two-tailed Fisher's exact test was used in each category, and a corrected P<0.05 was considered to indicate a statistically significant difference. Motifs of the phosphorylated sites for serine, threonine and tyrosine residues were analyzed using Motif-X (41). Peptide sequences were centered on each phosphorylation site and extended to 10 amino acids upstream and downstream of the site. The minimum occurrence was set to 20, and the significance was set to 0.000001. The new phosphorylation sites found in the present study were identified by comparing the dataset in the current study with the phosphoSitesPlus database (https://www.phosphosite.org/).

Tumor and adjacent non-tumor tissue samples collected from 6 different patients with HCC were used to further verify the bioinformatic results using western blotting. All patients were males and their mean age was 53 years (SD, 7.5; age range, 46–67 years). The tissues were also collected between May 2017 and July 2018 at the First Affiliated Hospital of Zhengzhou University. All patients provided written informed consent. Proteins were extracted from tissues on ice using RIPA buffer (Beyotime Institute of Biotechnology) supplemented with a protease inhibitor cocktail and phosphatase inhibitors for 30 min at 4°C. Protein concentration was determined using the BCA method. A total of 20 µg/lane protein was separated via 10% SDS-PAGE and transferred onto polyvinylidene difluoride membranes. The membranes were blocked at room temperature for 2 h with 5% BSA (Beyotime Institute of Biotechnology), followed by incubation overnight at 4°C with primary monoclonal antibodies against minichromosome maintenance complex component 2 (MCM2; 1:1,000; cat. no. 10513-1-AP; ProteinTech Group, Inc.), phosphorylated (p-)MCM2 (1:10,000; cat. no. ab109270; Abcam), vimentin (1:1,000; cat. no. 10366-1-AP; ProteinTech Group, Inc.), p-vimentin (1:1,000; cat. no. 13614; Cell Signaling Technology, Inc.) and GAPDH (1:10,000; cat. no. 60004-1-Ig; ProteinTech Group, Inc.), with the latter used as the loading control. After washing 3 times with TBS-Tween (0.1% Tween-20) for 5 min at 37°C, the membranes were incubated with HRP-conjugated Affinipure goat anti-mouse or anti-rabbit secondary antibodies for 2 h at room temperature (1:10,000; cat. nos. SA00001-1 or SA00001-2, respectively; ProteinTech Group, Inc.). Proteins were visualized using enhanced chemiluminescence detection reagents (Beyotime Institute of Biotechnology). The software used for densitometry was Quantity One 4.6.6 (Bio-Rad Laboratories, Inc.).

GraphPad Prism 5.0 (GraphPad Software, Inc.) was used for statistical analysis. The results of the experiments were expressed as the mean ± SD (n=6). Paired Student's t-test analysis was used to determine the statistically significant differences. P<0.05 was considered to indicate a statistically significant difference.

To identify the phosphoproteomic changes that result from HCC, four pairs of tumor and adjacent non-tumor tissues were resected from patients diagnosed with HCC. These tissues were frozen immediately after resection and stored in liquid nitrogen until lysis to maintain potentially unstable phosphorylation status. The proteins from HCC and normal tissues were extracted and then digested using trypsin. The tryptic peptides were labeled with TMT isobaric mass tags for multiple and precise quantification of phosphorylation levels, and then mixed together in equal quantities, followed by fractionation of the peptides using high-pH HPLC. The phosphopeptides were subsequently enriched by IMAC, and analyzed via LC-MS/MS. The distribution of mass errors was ~0, and all were <10 ppm (Fig. S1A). The majority of the identified peptides were between 8 and 20 amino acids in length, which is in accordance with the pattern of trypsin-digested peptides (Fig. S1B). These results indicated that the sample preparation method met the standards of LC-MS/MS for quantitative proteomics and that the MS data were accurate. The LC-MS/MS quantitative method then determined whether the peptides present in the samples were phosphorylated (Fig. 1), with MCM2, ADH4 and TP53BP1 used as examples.

The phosphoproteomic analysis of HCC and normal liver tissues resulted in the identification of 6,129 total unique phosphorylation sites from 2,687 proteins, which were predicted via Maxquant. A total of 4,780 phosphorylation sites from 2,209 proteins were quantified and were considered for further analysis (Table SI). Among these, single, double, triple and quadruple or higher phosphorylated peptides represented 55, 19, 11 and 15% of the total phosphopeptides, respectively (Fig. 2A). The observed phosphorylation sites were divided into four classes according to their phosphorylation site localization probability score calculated using Maxquant: 6,129 phosphorylation sites were classified as class I with a specific residue with high confidence (≥0.75); 1,008 were class II sites (0.5-0.75); 989 were class III sites (0.25-0.5); and 24 were class IV sites (<0.25). In total, 75% of all phosphorylation sites were classified as class I sites and were considered to be accurately quantified and localized, and were therefore further analyzed (Fig. 2B). The non-redundant phosphorylation sites in class I comprised 5,446 phosphorylated serine residues (89%), 655 phosphorylated threonine residues (11%) and 28 phosphorylated tyrosine residues (<1%) (Fig. 2B). Comparing the dataset in the current study with the phosphoSitesPlus database, 317 phosphoproteins and 138 phosphorylation sites were identified for the first time in the present study (Fig. 2C). For the quantitative analysis, statistical significance was determined using unpaired Student's t-tests. P<0.05 and a differentially expressed ratio >1.5 or <0.67 was considered to indicate a statistically significant difference. The difference between the HCC and normal tissue groups of each phosphosite was plotted with the corresponding P-value, as shown in Fig. 2D. A total of 533 differentially expressed phosphosites from 376 proteins were identified in HCC, specifically 74 phosphosites that were upregulated and 459 phosphosites that were downregulated compared with in normal tissues (Table SII).

Considering that tyrosine phosphorylation comprised <1% of all the localized phosphorylation sites observed, only phospho-serine motifs and phospho-threonine motifs were analyzed using Motif-X. The motif enrichment heat map of all the upregulated and downregulated phosphoproteins is shown in Fig. S2. In total, 53 phospho-serine motifs and 8 phospho-threonine motifs were identified (Table SIII). The majority of phospho-serine motifs were categorized into three classes according to the literature (42): Acidic, proline-directed and basic motifs. Proline-directed and acidic motifs accounted for 15 and 22 of the 53 identified phospho-serine motifs, respectively, while 12 of the 53 phospho-serine motifs belonged to the basic motifs class (Table SIII). The proline-directed motifs were predicted to be associated with MAPK, CDK, ERK1/2, PKA and PKC (Fig. 3A). The acidic motifs were predicted to be associated with CK2 (Fig. 3B), while the basic motifs were predicted to reflect activation of PKA and PKC (Fig. 3C).

All the identified proteins were subjected to GO annotation in order to assess differences in biological processes, molecular functions and cellular components to explain the biological effects of phosphoproteins. The phosphoproteins identified via LC-MS/MS were involved in various biological processes (Fig. 4A), including ‘cellular process’ (14%), ‘single-organism process’ (12%), ‘biological regulation’ (12%) and ‘metabolic process’ (9%). The identified phosphoproteome in the differentially expressed phosphoproteins was also classified according to molecular functions (Fig. 4B), revealing that most phosphoproteins were classified into the categories of ‘binding’ (56%), ‘catalytic activity’ (20%), ‘structural molecule activity’ (7%) and ‘molecular function regulator’ (7%). The cellular components of phosphoproteins were determined as ‘cell’ (24%), followed by ‘organelle’ (22%), membrane-enclosed lumen (13%), macromolecular complex (12%) and ‘membrane’ (11%) (Fig. 4C). Subcellular localizations of differentially expressed phosphoproteins were predicted and classified using wolfpsort software. The sub-cellular protein annotation for the differentially expressed phosphoproteins revealed that the majority of phosphoproteins were localized to the nucleus (53%), followed by cytoplasm (33%), mitochondria (9%) and plasma membrane (3%) (Fig. 4D).

Functional enrichment analysis revealed that the target phosphoproteins were enriched in numerous processes. GO enrichment analysis of molecular functions revealed that the phosphorylated proteins were significantly enriched in processes involved in binding activity and structural molecule activity. GO enrichment analysis based on cellular components illustrated that the phosphorylated proteins were enriched in various cellular components, such as ‘anchoring junction’, ‘extracellular organelle’ and ‘extracellular exosome’. GO enrichment analysis of biological processes demonstrated that the phosphorylated proteins were mainly enriched in biosynthetic and metabolic processes (Fig. 5A and B). KEGG enrichment analysis revealed that the upregulated targets were enriched in eight pathways (Fig. 5C), while the downregulated targets were enriched in 11 pathways (Fig. 5D). The upregulated phosphorylated proteins were most significantly enriched in the processes of ‘RNA transport’, ‘DNA replication’ and ‘tyrosine metabolism’, while the downregulated phosphorylated proteins were most significantly enriched in the processes of ‘arrhythmogenic right ventricular cardiomyopathy’, ‘ABC transporters’ and ‘glycolysis/gluconeogenesis’. These pathways were involved in genetic information processing, metabolism, environmental information processing, cellular processes and human diseases. Additionally, a domain enrichment analysis was performed. The identified proteins were annotated using the InterProScan software and the InterPro domain database. There were 9 domains that matched the upregulated proteins, such as ‘alcohol dehydrogenase, N-terminal’, ‘Gas2-related domain’ and ‘NAD(P)-binding domain’ (Fig. 5E). Additionally, there were 25 domains that matched the downregulated proteins, such as ‘ABC transporter type 1, transmembrane domain’, ‘calponin homology domain’ and ‘riboflavin synthase-like beta-barrel’ (Fig. 5F).

To further understand the cellular processes regulated via phosphorylation in HCC, a PPI network of the phosphoproteins was established (Fig. 6A). All the upregulated and downregulated phosphoproteins were uploaded to STRING, presenting a global view of the networks of phosphoproteins in HCC. It is well known that protein phosphorylation is closely associated with cell signal transduction. Therefore, STRING was used to construct a PPI network to illustrate the association between the identified phosphoproteins and signal transduction pathways in HCC (Fig. 6B). The PPI network was composed of 173 phosphoproteins participating in the hypoxia-inducible factor-1 (HIF-1), mTOR and MAPK signaling pathways, as well as in cell cycle and apoptosis, which may be involved in networks associated with cancer functions, directly or indirectly. There were 4 upregulated and 22 downregulated phosphoproteins, while LMNB2 included both upregulated and downregulated phosphorylation sites (Fig. 6B). The network indicated phosphoprotein dysregulation in HCC, thus providing information on the initiation and progression of HCC, and revealing potential targets for therapeutic intervention.

To verify the protein phosphorylation data obtained via LC-MS/MS, western blotting was performed. As shown in Fig. 7, MCM2 exhibited an increased phosphorylation level in HCC, while vimentin exhibited a decreased phosphorylation level in HCC compared with in normal tissues. Therefore, western blotting confirmed the reliability of the LC-MS/MS data.