Wild-type C57BL/6J male mice were allowed to take food and water ad libitum. The colony rooms were maintained at a constant temperature and humidity with a 12:12 light/dark cycle. All animal protocols were approved by the Animal Care and Use Committee of the Institute of Basic Medical Sciences, Chinese Academy of Medical Sciences and Peking Union Medical College.
2. Sample preparation
Fourteen two-month-old mice were used in the young group and ten eighteen-month-old mice were used in the old group. Mice were euthanized by massive bloodletting from the orbital vessel after anesthesia with tribromoethanol. Next, the whole liver was detached from each mouse, immediately dissected, and stored separately. Samples acquired for phosphoproteome and transcriptome analyses were immediately frozen in liquid nitrogen and transferred to -80 ℃ until use. Samples for H & E and Masson analyses were quickly placed in 4% paraformaldehyde (PFA). Samples for oil red O analysis were appropriately embedded into optimal cutting temperature compound (OCT) and stored at -80 ℃ until use.
3. Morphological analysis
3.1 H & E and Masson staining
Liver samples were fixed in 4% PFA overnight. Fixed tissues were dehydrated by 75% ethanol for 4 h, 85% ethanol for 2 h, 90% ethanol for 2 h, 95% ethanol for 1 h, absolute ethanol for 30 min twice, ethanol-dimethylbenzene for 5 min, and dimethylbenzene twice for 10 min. Dehydrated tissues were embedded in paraffin and cut in 4 mm-thick sections. The paraffin embedded sections were successively placed in dimethylbenzene twice for 20 min each, absolute ethanol for 10 min twice, 95% ethanol for 5 min, 90% ethanol for 5 min, 80% ethanol for 5 min, 70% ethanol for 5 min, and washed by distilled water.
For H & E staining, hydrated sections were placed in hematoxylin solution for 3-8 min, 1% hydrochloric acid/ethanol differentiation solution for seconds, 0.6% ammonia for seconds, and eosin solution for 1-3 min.
For Masson staining, hydrated sections were processed according to the manufacturer’s protocol of the Masson staining kit (Wuhan Goodbio Technology Co., Ltd, G1006).
Stained sections were subsequently transferred into 95% ethanol for 5 min twice, absolute ethanol for 5 min twice, and dimethylbenzene for 5 min twice. Next, the sections were dried and sealed with neural gum. Pictures were taken with a Nikon Eclipse CI imaging system.
3.2 Oil Red O staining
Liver samples embedded in OCT were moved to a freezing microtome and cut into 8 mm-thick sections at -20 ℃. The sections were dried at room temperature for 10 min, fixed with 4% paraformaldehyde (PFA) for 15 min, and washed with phosphate-buffered saline (PBS) for 5 min three times. Sections were transferred into oil red O solution (G1016, Goodbio Technology Co., Ltd) for 10 min, followed by 75% ethanol for 2 s, and then washed with water for 1 min. Next, the sections were transferred to hematoxylin solution for 1 min, 1% hydrochloric acid/ethanol differentiation solution for 3 s, and 0.6% ammonia for 3 s, after which they were washed with water. Excess water was removed, and glycerin gelatin was used to seal the sections. Pictures were taken using a Nikon Eclipse CI imaging system (Japan).
4. Phosphoproteome and analysis
4.1 Protein Extraction and Digestion
Mouse liver tissue samples were ground in liquid nitrogen and sonicated with lysis buffer (9 M Urea, 10 mM Tris–HCl (pH 8.0), 30 mM NaCl, 50 mM IAA, 5 mM Na4P2O7, 100 mM Na2HPO4 (pH 8.0), 1 mM NaF, 1 mM Na3VO4, 1 mM sodium glycerophophate, 1% phosphatase inhibitor cocktail 2 (Sigma, St. Louis, USA), 1% phosphatase inhibitor cocktail 3 (Sigma, St. Louis, USA), 1 tablet of EDTA-free protease inhibitor cocktail (Roche, Basel, Switzerland) for every 10 mL of lysis buffer). Then the total lysate was centrifuged at 17,000 × g for 8 min at 4 ℃ to remove debris. The evaluation of protein concentration and subsequent in-gel digestion was performed as previously described . In brief, 1 mg proteins were used as starting material for each sample, the proteins were reduced with 5 mM DTT for 30 min at 45 °C and were alkylated with 20 mM iodoacetamide for 30 min at room temperature in the dark. Next, proteins were resolved on a 10% SDS-PAGE gel and running for a 0.8 cm length and then stained with Coomassie Blue G-250. The entire gel lane was sliced into 1 mm3 pieces and destained followed by in-gel digestion with 10 ng/μL Trypsin (Promega, Madison, WI, USA) at 37 °C incubation overnight.
4.2 Phosphopeptides Enrichment
Phosphopeptides were enriched by Ti4+-IMAC method  as previously reported with minor modifications. Briefly, peptide mixtures were resuspended in 500 μL loading buffer (5% ACN, 50 mM NH3HCO3) and acidified with 500 μL binding buffer (80% ACN, 6% TFA) followed by the addition of 30 mg Ti4+-IMAC beads. Then the mixtures were mixed wildly on Vortex Genie2 for 30 min and centrifuged at 17,000 × g for 6 min to remove the supernatant. The beads were washed with 1.8 mL wash buffer 1 (80% ACN, 6% TFA) once and 1.8 mL wash buffer 2 (80% ACN, 0.1 TFA) twice and all the wash buffers in each step were removed by centrifugation at 17,000 × g for 6 min. Afterwards the beads were resuspended in 1 mL elution buffer (10% NH3·H2O) and the mixture was vortexed for 15 min and sonicated in ice water for 15 min followed by centrifugation at 17,000 × g for 6 min. The eluted supernatant was collected. The beads were vortexed for 5 min in another 500 μL elution buffer and centrifuged to collect the supernatant. The eluted supernatants were combined and centrifuged at 21,000 × g for 8 min to further remove the beads. Then the supernatant was vacuum dried and the eluted phosphopeptides were stored at -80 ℃.
4.3 Stage-Tip Separation
The homemade C18 Stage-Tip  was used to separate the phosphopeptides in 3 fractions. The Stage-Tip was firstly activated with 40 μL methanol, then washed with 40 μL wash buffer (80% ACN in 10% NH3·H2O) twice and 40 μL loading buffer (10% NH3·H2O) twice. The phosphopeptides were resuspended in 40 μL loading buffer and loaded onto the Stage-Tip. Then the phosphopeptides were eluted in a sequential gradient of ACN (0% ACN, 2% ACN, 5% ACN, 8% ACN, 10% ACN, 20% ACN, 40% ACN, 50% ACN, 80% ACN) in 10% NH3·H2O buffer (20 μL/fraction). All the 9 fractions were combined to 3 fractions (0% ACN, 8% ACN, 40% ACN for fraction 1; 2% ACN, 10% ACN, 50% ACN for fraction 2; 5% ACN, 20% ACN, 80% ACN for fraction 3). All of these 3 fractions were lyophilized immediately.
4.4 LC-MS/MS Analysis and Database Search
The fractionated enriched peptides were eluted on a Thermo Fisher EASY-nLC 1200 liquid chromatography and analyzed by Orbitrap Fusion Lumos Tribrid MS (Thermo Fisher Scientific). Peptides were separated on a 16 cm column with a 75 μm inner diameter packed with Dr. Maisch GmbH reversed-phase material Reprosil-Pur 120 C18-AQ, 1.9 μm resin (SinoAmerican Proteomics, LLC). Separation was performed using 75 min runs at a flow rate of 300 nL/min through nonlinear gradient. The elution gradient was as follows: 3-10% B for 3 min, 10-24% B for 55 min, 24-32% B for 10 min, 32-90% B for 4 min, and 90% B for 3 min (Phase A: 0.1% FA and 2% ACN in ddH2O; Phase B: 0.1% FA in 80% ACN). The initial MS spectrum (MS1) was analyzed over a range of m/z 350-1500 with a resolution of 60,000 at m/z 400. The automatic gain control (AGC) was set as 4 × 105. The subsequent MS spectrum (MS2) was analyzed using data-dependent mode searching for the top 40 intense ions fragmented in the linear ion trap via high-energy collision dissociation (HCD) with the NCE set at 30. The MS2 resolution was 15,000 at m/z 400. The cycle time was 4s.
MS/MS raw files were processed in MaxQuant (version 18.104.22.168) against the uniprot Mus musculus database (downloaded at 2018.04.10, containing 16,972 sequences) for peptide identification, label-free quantification, and phosphosite localization. The parameters set for database searching were as follows: cysteine carbamidomethyl was specified as a fixed modification. Oxidation of methionine, phospho (STY) and protein N-acetylation and were set as variable modifications. The tolerances of precursor and fragment ions were set at 20 ppm. For digestion, trypsin was set as protease with utmost two missed cleavage permitted. False discovery rate was set as 1% at protein, peptide and modification level. Matched between runs and label free quantification were selected.
4.5 Data management
Proteins/peptides in the reverse decoy database and potential contaminant database were excluded. The localization probability of phosphorylation ranged from 0 to 1, and peptides with localization probabilities of > 0.75 were grouped into Class I and were selected for further analysis. Normalization was performed with Perseus (v.22.214.171.124)  by dividing the intensity by the median of each group. Using a two-sample t-test, a site with p-value < 0.05 was considered as differentially expressed. Upregulated and downregulated sites were determined as fold-change (FC, mean values of old mice/mean values of young mice) > 1 or < 1, respectively.
4.6 Functional enrichment analysis
Kyoto Encyclopedia of Genes and Genomes (KEGG)  enrichment analysis was performed using KOBAS . We chose hypergeometric test/Fisher’s exact test as the statistical method and QVALUE as the FDR correction method. KEGG terms with a corrected p-value of < 0.05 were considered as significantly enriched
4.7 Motif analysis
Phosphorylated peptides were submitted to MoMo modification motifs (http://meme-suite.org/tools/momo)  to get the sequence characteristics using the motif-x algorithm.
4.8 Kinase prediction and kinase-substrate interaction analysis
All the identified phosphosites were uploaded to NetworKIN  for kinase prediction following the algorithm’s instructions. Based on the Gene Set Enrichment Analysis (GSEA) principle , kinase activity were calculated and matched to the human kinome tree. Kinome tree modified courtesy of Cell Signalling Technology Inc. (www.cellsignal.com) and annotated using Kinome Render .
For strict kinase-substrate analysis, network among proteins with differentially expressed phosphosites was constructed according to the manually curated PhosphositePlus database . We used the cytoscape software to present the kinase-substrate interactions.
5. Transcriptome and analysis
RNA isolation, library preparation, and sequencing were performed by Novogene Bioinformatics Technology Co., Ltd (Tianjin, China). Briefly, a total of 3 mg RNA per sample was used as input material for RNA sample preparation. Firstly, ribosomal RNA was removed with the Epicentre Ribo-zeroTM rRNA Removal Kit (RZH1046, Epicentre, USA), and rRNA free residue was cleaned up by ethanol precipitation. Subsequently, sequencing libraries were generated using the rRNA-depleted RNA by NEBNext® UltraTM Directional RNA Library Prep Kit for Illumina® (NEBE7770, NEB, USA) following manufacturer’s recommendations. Fragmentation was carried out using divalent cations under elevated temperature in NEBNext First Strand Synthesis Reaction Buffer (5×). First strand cDNA was synthesized using random hexamer primer and M-MuLV Reverse Transcriptase (RNase H). Second strand cDNA synthesis was subsequently performed using DNA Polymerase I and RNase H. In the reaction buffer, dTTP was replaced by dUTP. Remaining overhangs were converted into blunt ends via exonuclease/polymerase activities. After adenylation of 3’ ends of DNA fragments, NEBNext Adapter with hairpin loop structure was ligated to prepare for hybridization. To select cDNA fragments of preferentially 150~200 bp in length, library fragments were purified with AMPure XP system (Beckman Coulter, Beverly, USA). Next, 3 mL USER Enzyme (NEB, USA) was used with size-selected, adaptor-ligated cDNA at 37 ℃ for 15 min followed by 5 min at 95 ℃ before PCR. PCR was performed with Phusion High-Fidelity DNA polymerase, Universal PCR primers, and Index (X) Primer. At last, products were purified (AMPure XP system) and library quality was assessed on an Agilent Bioanalyzer 2100 system. Clustering of the index-coded samples was performed on a cBot Cluster Generation System using TruSeq PE Cluster Kit v3-cBot-HS (Illumina) according to the manufacturer’s instructions. After cluster generation, the libraries were sequenced on an Illumina HiSeq 4000 instrument, and 150 base pair and paired-end reads were generated. For quality control, raw data in fastq format were first processed using Novogene Perl scripts. Clean data were obtained by removing reads containing adapters, reads containing poly-N, and low-quality reads from the raw data. In addition, the Q20, Q30, and GC contents of the clean data were calculated. All down-stream analyses were based on the clean data with high quality. RNA sequencing data were deposited in the Sequence Read Archive under the BioProject ID PRJNA609589.
Reference genome and gene model annotation files were downloaded from the Ensembl website (genome: ftp://ftp.ensembl.org/pub/release-97/fasta/mus_musculus/dna/Mus_musculus.GRCm38.dna.primary_assembly.fa.gz; gtf: ftp://ftp.ensembl.org/pub/release-97/gtf/mus_musculus/Mus_musculus.GRCm38.97.gtf.gz). HISAT2 (v2.0.5) was used to build the reference genome index and align paired-end clean reads to the reference genome. Then, StringTie (v1.3.3) was used to assemble the mapped reads of each sample and calculate FPKMs of coding genes. FPKM means fragments per kilo-base of exon per million fragments mapped, calculated based on the length of the fragments and reads count mapped to this fragment. Transcripts with FPKM values > 1 in over 50% of the samples in either group were considered validated. The edgeR R package (v3.243) provided statistical routines for determining differential expression in digital transcript or gene expression data using a model based on a negative binomial distribution. Transcripts with adjusted p-values < 0.05 were considered to be differentially expressed. Up- and down-regulated transcripts were determined based on the log2 fold-change (FC, generated by edgeR, mean values of old mice/mean values of young mice) > 0 or < 0, respectively.
6. Transcription factor (TF)-target interaction analysis
TRRUST collects 6,552 TF-target interactions for 828 mouse TFs . According to the trust_rawdata.mouse.tsv file downloaded on 2021/01/08, we selected TFs from phosphoproteome-identified proteins with differentially expressed phosphosites and targets from transcriptome-identified differentially expressed transcripts. TF-target interaction network was constructed by the Cytoscape software.