Quantitative Proteomic Study Unmasks Fibrinogen Pathway in Polycystic Liver Disease

(1) Background: Polycystic liver disease (PLD) is a heterogeneous group of congenital disorders characterized by bile duct dilatation and cyst development derived from cholangiocytes. Nevertheless, the cystogenesis mechanism is currently unknown and the PLD treatment is limited to liver transplantation. Novel and efficient therapeutic approaches are th6us needed. In this context, the present work has a principal aim to find novel molecular pathways, as well as new therapeutic targets, involved in the hepatic cystogenesis process. (2) Methods: Quantitative proteomics based on SWATH–MS technology were performed comparing hepatic proteomes of Wild Type and mutant/polycystic livers in a polycystic kidney disease (PKD) murine model (Pkd1cond/cond;Tam-Cre−/+). (3) Results: We identified several proteins altered in abundance, with two-fold cut-off up-regulation or down-regulation and an adjusted p-value significantly related to hepatic cystogenesis. Then, we performed enrichment and a protein–protein analysis identifying a cluster focused on hepatic fibrinogens. Finally, we validated a selection of targets by RT-qPCR, Western blotting and immunohistochemistry, finding a high correlation with quantitative proteomics data and validating the fibrinogen complex. (4) Conclusions: This work identified a novel molecular pathway in cystic liver disease, highlighting the fibrinogen complex as a possible new therapeutic target for PLD.


Introduction
Polycystic liver disease (PLD) is a heterogeneous group of congenital disorders characterized by bile duct dilatation and cyst development derived from cholangiocytes (epithelial cells from bile duct). The numerous cysts spread and progress throughout the liver parenchyma compromising its function [1]. PLDs are inherited in dominant or recessive form and can occur in isolation (isolated polycystic liver disease or PLD) or, more commonly, as an extra-renal manifestation of autosomal dominant adult polycystic kidney disease (PKD) [2]. In autosomal dominant polycystic kidney disease (ADPKD),~85% of the patients develop liver cysts throughout hepatic parenchyma and its severity can range from a few cysts to a severe hepatic cystogenesis [3,4]. Meanwhile, in autosomal recessive polycystic kidney disease (ARPKD), liver cysts could develop, but the principal hepatic complication is a defective remodeling of the ductal plate with congenital hepatic fibrosis [5,6]. For both, PLD is the principal extrarenal manifestation and requires clinical management and treatment [2,5]. The principal clinical approach for PLD patients is to reduce their clinical symptoms. Management of bile duct complications is the main clinical feature to be controlled, and the final treatment strategy is usually liver transplantation [7], since there is no effective pharmacological treatment for PLD [2,5]. Multiple molecular pathways have been associated with the disease [8], such as increased fluid secretion [9,10], proliferation [11], fibrosis [12,13], ciliary defects [14] and others [7,8]; but the main mechanism undergoing cystogenesis remains to be elucidated. The identification of new potential therapeutic targets is urgently required. All these studies were based on molecular studies using different classical molecular tools: 3D cell culture [7,8], immunohistochemistry [7][8][9]11,12], ELISA [9,11], transmission electron microscopy (TEM) [12], immunocytochemistry [11], RT-qPCR (quantitative real-time PCR) [7,8,10,11] and Western blot [7,8,[10][11][12]. Approaching the study of disease from new perspectives, well-established preclinical animal models, and taking advantage of new technologies, can be a benefit in the discovery of new molecular mechanisms of disease to translate molecular discoveries into new therapies.
Liquid chromatography coupled to tandem mass spectrometry (LC-MS/MS) has been the most commonly used technology for the characterization of disease proteomes as well as the identification of possible disease biomarkers or drug screening [15,16]. An emerging strategy named SWATH-MS (sequential window acquisition of all theoretical fragment-ion spectra-mass Spectrometry) allows a reproducible, reliable, high-throughput quantification with up to thousands of proteins analyzed per biological sample, and with the possibility of simultaneously performing a large number of assays [16,17]. For this technology, the samples are digested with trypsin and analyzed by liquid chromatography coupled to a tandem mass chromatography works in the so-called Data-Independent Acquisition (DIA) mode [16] (see Material and Methods). Next-generation proteomics SWATH-MS technology allows quantifying the proteome and identifying differences in protein abundance in different samples and conditions.
In the present study, we performed a differential proteomic quantitative analysis based on SWATH-MS technology in polycystic liver disease from mutant animal models. We report a list of proteins with statistical significance of differential protein abundance, which were characterized by different enrichment analysis and data curation. Finally, we validated the SWATH data by identifying and highlighting the fibrinogen complex as a new mechanism of disease and possible new therapeutic target for PLD.

Murine Model
In this study, we used a murine Pkd1 conditional-knockout animal model, the C57/BL6 k Pkd1 cond/cond ;Tam-Cre [18,19]. Cre-mouse model was used as Wild Types (WT), and the Cre + as Mutant (KO). Genotypes were confirmed by standard PCR, following the conditions and primers described [18]. Deletion of Pkd1 was always at postnatal days 10 and 11 (p10 and p11) with a single intraperitoneal injection of tamoxifen (10 mg/40 g) (Sigma ® , St. Louis, MO, USA No. T5648), diluted in corn oil (Sigma-Aldrich ® , St. Louis, MO, USA No. C8267), both days in the mother of the litters to be induced. The animals were euthanized at postnatal day 30 (p30). The groups were established by different individuals and organized randomly (randomization) and keeping 1:1 ratio between males and females in all groups. Anatomical characteristic features were used to identify the sex of the animals. The mice were housed in a pathogen-free facility (SPF) in accordance with the established conditions of the University of Santiago.

Protein Extraction and Digestion
At the p30 sacrifice point, the dissected livers were stored at −80 • C. The whole liver was ground in liquid nitrogen using a liquid nitrogen-cooled mortar (DD Biolab ® , Barcelona, Spain No. 088763), with half of the tissue assigned to protein extraction and the other half to RNA extraction. Protein extracts were prepared with RIPA lyses buffer (10 Mm TRIS, 5 mM EDTA, 150 mM NaCl, 0.1% SDS, 1% TRITON X-100 and 1% sodium deoxycholate). 1% protease and phosphatase inhibitors (Sigma ® No. P8340 y No. P0044) were also added. This tissue lysate underwent a centrifugation process (for 30 min at 4 • C and 14,000 rpm). The supernatant was then recollected and quantified by the Bradford protein assay (Bio-Rad ® No. 5000001). Protein aliquots (100 µg) were concentrated in a single band in a 10% sodium dodecyl sulphate-polyacrylamide gel electrophoresis (SDS-PAGE), cut and submitted for manual digestion as previously described [20]. Finally, the extracted peptides were dissolved in 0.1% formic acid for further analysis. 4 µL of each sample (over 4 µg of protein) were analyzed by shotgun data-dependent acquisition (DDA) approach by micro-LC-MS/MS. The samples were separated by the Ekspert nLC425 micro-LC system (Eksigen ® , Dublin, CA, USA) using a YMC-TRIART C18 trap column with a 3 mm particle size and 120 Å pore size (YMC Technologies, Teknokroma, Barcelona, Spain) and a column Chrom XP C18 150 mm × 0.30 mm, 3 mm particle size and 120 Å pore size (Eksigen ® ) at a flow rate of 10 µL/min. The solvents used were solvent A (water, 0.1% formic acid) and solvent B (ACN, 0.1% formic acid). The gradient was from 5% to 95% B for 30 min, 5 min at 90% B and, finally, other 5 min at 5% B for column equilibration, for a total time of 40 min. The mass spectrometer coupled was a hybrid quadrupole-TOF mass spectrometer, 6600 (SCIEX ® , Framingham, MA, USA) operating with a data dependent acquisition system in positive ion mode. Using the mass spectrometer, a 250 ms survey scan was performed at 400 to 1250 m/z followed by MS/MS experiments at 100 to 1500 m/z (25 ms of acquisition time) for a total cycle time of 2.8 s. Fragmented precursors were added to the dynamic exclusion list for 15 s, any ion with charge +1 was excluded from the MS/MS analysis. MS raw file and database searches were combined and performed using ProteinPilot software v.5.0.1. (SCIEX ® , Framingham, MA, USA) using a mouse-specific UniProt Swiss-Prot database. It is necessary to specify iodoacetamide cysteine alkylation as fixed modification and trypsin digestion. The false discovery rate (FDR) was set to 1 for peptides and proteins with a confidence score greater than 99% [21].

Generation of the References Spectral Library
The spectral library was performed by a DDA method as described above, but instead each independent sample, a pool of each group was performed with 4 µL of each sample (over 4 µg of protein) were run in order to obtain more representative protein identification in the spectral library. All raw files were launched together into the Uniprot database to obtain the spectral library using the protein pilot conditions mentioned in the previous section.

Quantification by SWATH and Data Analysis
The SWATH-MS acquisition was performed using a DIA method. We used 4 groups with 4 biological replicates per group and 3 technical replicates per sample (see Table S1 for more details). The tissues were analyzed and, in all cases, compared in Wild Type (from Pkd1 cond/cond ;Tam-Cre − mouse) and Mutant (from Pkd1 cond/cond ;Tam-Cre + mouse) conditions. 4 µg of protein per sample and technical replicates were run using a SWATH method. For each set of samples, the width of the 100 variable windows was optimized according to the ionic density found in the library DDA runs using a Sciex SWATH variable window spreadsheet. Therefore, the SWATH method was based on a cycle of repetitions that consisted of the acquisition of 100TOF MS/MS scans (400 to 1500 m/z, high sensitivity Biomedicines 2022, 10, 290 4 of 21 mode, 50 ms acquisition time) of isolation windows of sequential overlapping precursor variables with widths (1 m/z overlap) covering the 400 to 1250 m/z mass range with a previous TOF MS scan (400 to 1500 m/z, 50 ms acquisition time) for each cycle. Total cycle time was 6.3 s.

Data Analysis
All proteins in the ion library that were identified by ProteinPilot with an FDR below 1% were quantified by the SWATH method. Once individual samples were acquired, the spectral alignment and targeted data extraction were performed by PeakView v.2.2. (SCIEX ® , Framingham, MA, USA) matching the reference spectral library as described below [22][23][24][25][26]. The retention times of the peptides that were selected for each protein, were realigned in each run according to the iRT peptides present in each sample and were eluted along the whole-time axis. The extracted ion chromatograms were then generated for each ion of the selected fragment; the peak areas for the peptides were obtained by adding the peak areas of the ions of the corresponding fragments. PeakView computed an FDR and a score for each assigned peptide according to the chromatographic and spectra components; only peptides with an FDR below 5% were used for protein quantitation. Similarly, to obtain the areas for protein quantification as a function of signal intensity, up to 10 peptides per protein and seven fragments per peptide were selected. All shared and modified peptides were excluded from processing. Five-minute windows and 30 ppm widths were used to extract the ion chromatograms.
Integrated area peaks (processed ".mrkvw" files from PeakView, 22 December 2017) were exported directly to MarkerView software (SCIEX ® , Framingham, MA, USA) for relative quantitative analysis. The export generates three files that contain quantitative information about individual ions, the summed intensity of different ions for a particular peptide and the summed intensity of different peptides for a particular protein [27][28][29][30]. MakerView uses processing algorithms that accurately find chromatographic and spectral peaks directly from raw SWATH data. MarkerView data alignment compensates for minor variations in both mass and retention time values, ensuring that identical compounds in different samples are accurately compared to each other. To control the possible uneven loss of samples in the different samples during the sample preparation process, we performed a TAS (total area sum) normalization [31][32][33]. Version R 3.6.2 (https://www.r-project.org/, 12 February 2020) and a set of packages has been used for the analysis [34]. Volcano plots were used to graphically visualize the results that were generated by plotting the log(2)-fold changes for all proteins identified against their −log(10) p-value. Unsupervised cluster analysis heatmap figures of proteins with significant differences according to SWATH-MS analysis were generated to detect different cluster of proteins in different studied groups.

Signaling Pathway Analysis
Different functional analysis was performed using differentially expressed proteins in order to evaluate the most relevant interaction networks, related pathways or cellular components and other associated proteins or pathways. UniProt protein accession numbers (https://www.uniprot.org/, 28 August 2021) were used. Proteins with differences in protein abundance were subjected to: String ® (https://string-db.org/, 25 October 2021) using Gene Ontology (GO) terms and FunRich ® (http://funrich.org/index.html, 21 October 2021) to study protein networks and biological association between proteins; and to Reactome ® (https://reactome.org/, 25 October 2021) to study the pathways involved in relation or interaction with significantly expressed proteins.

Histology and Measurement of Cystic Index
The largest lobe of the liver of the animals was always fixed in paraformaldehyde (4% paraformaldehyde buffered solution, pH = 7, Labbox No. FORM-D0P-10K) overnight at 4 • C. The livers were then dehydrated with xylene and ethanol and finally embedded in For quantifications, slides of longitudinal sections of the largest lobe of liver were always taken as representative images and, at least six different sections per animal, visualized under an Olympus BX51 microscope connected to an Olympus Camera DP70, were used for quantification of cystic index and number of cysts. The quantification of both parameters was calculated with Automatic Cyst Recognition software tool called CystAnalyser (https://citius.usc.es/transferencia/software/cystanalyser, 1 July 2020) as previously described [35].

Clinical Chemistry
Blood was collected before euthanasia, and serum was obtained using BD Vacutainer SST II Advance tubes. Serum alanine transaminase (ALT) and alkaline phosphatase (ALP) were measured using the ADVIA 2400 Chemistry System (Siemens Healthineers, Erlangen, Germany).

Real-Time Quantitative PCR
One microgram of total liver RNA was treated with DNase I (Invitrogen ® , Carlsbad, CA, USA, No. 18068015). Then, reverse transcription was performed with SuperScript ® III First-Strand Synthesis System kit (Invitrogen ® , No. 18080051), and gene expression analysis by quantitative real-time PCR using FastStart Universal SYBR GREEN Master (ROX) (Roche ® No. 04913914001) in a Mx3005P system (Agilent Technologies ® , Santa Clara, CA, USA). The primers used for this study were designed with Primer 3 Plus (http://tiny.cc/pf928y, 11 March 2020) and are based on the cDNA reference sequences published in Ensembl genome browser (http://www.ensembl.org/index.html, 11 March 2020). The full sequences of all primers used in the study are found in Table S2. The following program was used for all complete PCR reactions: a denaturation step with 30 s at 95 • C; 45 cycles of 30 s at 95 • C, 1 min at 60 • C and 30 s at 72 • C; followed by a cooling step. As we have previously carried out [36], each sample (n = 6 for each group) was run in triplicate in each experiment. The results were analyzed through the absolute quantification method. Consequently, a standard curve with known concentrations was used in the assay and, in addition, each sample was normalized to a constantly expression gene. In this study, the Gapdh gene (glyceraldehyde-3-phosphate dehydrogenase) was used as the housekeeping gene.

Western Blotting
Western blotting assays were performed on whole liver lysates. Total protein content was determined using Bradford protein assay (Bio-Rad ® , Hercules, CA, USA, No. 5000001) and the same amount of each protein (20-30 µg) was boiled at 95 • C for 5 min and separated in 10% SDS-PAGE under reducing conditions. The PageRuler TM Prestained Protein Ladder (Thermo Fisher ® , Waltham, MA, USA, No. 26616) was used as the protein molecular weight standard. Subsequently, the proteins were transferred to nitrocellulose membranes (Bio-Rad ® No. 1620115). To avoid nonspecific binding, the membranes were blocked for 1 h at room temperature using a 2.5% BSA solution (NZYTech ® , Lisboa, Portugal, No. MB046) for phosphorylated proteins, or SuperBlock Blocking Buffer (Thermo Fisher ® , Waltham, MA, USA, No. 37515). Subsequently, the membranes were incubated overnight at 4 • C with primary antibodies and next washed in 0.3% Tween-20 in PBS. Finally, the nitrocellulose membranes were incubated for 1 h at room temperature with the HRP-linked secondary antibodies. Signals were developed with the SuperSignal™ West Pico PLUS Chemiluminescent Substrate kit (Thermo Fisher ®

Statistics Analysis
Data are presented as means ± SEM for all cases (bar or point plots). First, normal distribution was confirmed using the Kolmogorov-Smirnov test. Two-tailed t-test was used to determine significance of differences between two groups. p < 0.05 was considered statistically significant. Furthermore, ** corresponds to p-values < 0.01 and *** to p-values < 0.001. All datasets were analyzed using GraphPad Prism software (version 9, https://www.graphpad.com/, 13 October 2021).
For proteomic studies, a Student's t-test (using MarkerView software, 22 December 2017) was performed to compare the groups in pairs in order to identify proteins, which were significantly differentially represented using adjusted p-value < 0.05 (multiple correction FDR method) and fold change >2 as cut-off. Unsupervised multivariate statistical analysis using principal component analysis (PCA) was performed to compare the data across the samples.

Kidney and Liver Does Not Follow the Same Cystic Mechanisms
Previous studies showed that early inactivation of the Pkd1 gene (before p12) triggers rapid development of polycystic disease (cystic window), whereas inactivation after p14 leads to late cyst formation after 4-5 months, with a mild phenotype (non-cystic window) [19]. Using the same orthologous model of ADPKD (Pkd1 cond/cond ;Tam-Cre) [18,19], we investigated whether this window of renal development corresponded in time with a similar window of hepatic development. To draw a direct comparison with our previous studies [37], mice were induced with tamoxifen to inactivate the Pkd1 gene at postnatal days 14 (p14) and 12 (p12) and sacrificed at the age of 30 days (p30). Interestingly, we observed a cystic phenotype at both time points with a different degree of disease (mild and severe), suggesting that the kidney and liver have independent cystic developmental windows (Figure 1a), and/or possible different timings and progression of disease. No differences were found between males and females at this slaughter age between both windows of inactivation (n = 6 was used for each group of mutants). The p14 mutant animals presented a milder phenotype than p12 mutant group with lower values of cystic index, number of cysts and hepatic function according to ALP serum value (Figure 1b-d). This is the first study showing different developmental windows/mechanisms for Pkd1 liver and kidney deficiency. Those molecular differences will have to be determined in future studies. rapid development of polycystic disease (cystic window), whereas inactivation after p14 leads to late cyst formation after 4-5 months, with a mild phenotype (non-cystic window) [19]. Using the same orthologous model of ADPKD (Pkd1 cond/cond ;Tam-Cre) [18,19], we investigated whether this window of renal development corresponded in time with a similar window of hepatic development. To draw a direct comparison with our previous studies [37], mice were induced with tamoxifen to inactivate the Pkd1 gene at postnatal days 14 (p14) and 12 (p12) and sacrificed at the age of 30 days (p30). Interestingly, we observed a cystic phenotype at both time points with a different degree of disease (mild and severe), suggesting that the kidney and liver have independent cystic developmental windows (Figure 1a), and/or possible different timings and progression of disease. No differences were found between males and females at this slaughter age between both windows of inactivation (n = 6 was used for each group of mutants). The p14 mutant animals presented a milder phenotype than p12 mutant group with lower values of cystic index, number of cysts and hepatic function according to ALP serum value (Figure 1b-d). This is the first study showing different developmental windows/mechanisms for Pkd1 liver and kidney deficiency. Those molecular differences will have to be determined in future studies. Bars represent means ± SEM in all cases. p < 0.05 by Student's t-test (two-tail) was considered as a significant result. ns represents not significance. * p < 0.05, ** p < 0.01, *** p < 0.001.

Shotgun and SWATH-MS Proteomic Analysis in PLD
As referred in Figure 2, we performed differential proteome analysis at both cystic stages (p12-severe cystic disease-and p14-mild cystic disease-) of Mutants (KO) in comparison to Wild Type (WT) animals, in order to identify and characterize relevant molecular mechanisms undergoing liver cystogenesis and disease progression.

Shotgun and SWATH-MS Proteomic Analysis in PLD
As referred in Figure 2, we performed differential proteome analysis at both cystic stages (p12-severe cystic disease-and p14-mild cystic disease-) of Mutants (KO) in comparison to Wild Type (WT) animals, in order to identify and characterize relevant molecular mechanisms undergoing liver cystogenesis and disease progression. Figure 2. Illustrated workflow scheme. We used an ADPKD murine model (Pkd1 cond/cond ;Tam-Cre mice) in which the genetic inactivation of Pkd1 was induced by tamoxifen administration at postnatal day 10 (p10) and p11 (rapid disease progression) or p15 and p16 (delayed disease progression) according to the developmental switch for renal cystogenesis, defining the two groups of the study named p12 and p14 respectively [19]. For both groups, we used an n equal or greater than six Wild Type (WT) and Mutant (KO) individuals in each condition, and all animals were sacrificed at p30. We studied the differential proteome of WT and KO livers in both severities of liver cystic disease (Figure 1), using both quantitative proteomic SWATH-MS and shotgun data-dependent acquisition DDA-MS analysis. Of the proteins with a significant change in abundance, we used bioinformatics Figure 2. Illustrated workflow scheme. We used an ADPKD murine model (Pkd1 cond/cond ;Tam-Cre mice) in which the genetic inactivation of Pkd1 was induced by tamoxifen administration at postnatal day 10 (p10) and p11 (rapid disease progression) or p15 and p16 (delayed disease progression) according to the developmental switch for renal cystogenesis, defining the two groups of the study named p12 and p14 respectively [19]. For both groups, we used an n equal or greater than six Wild Type (WT) and Mutant (KO) individuals in each condition, and all animals were sacrificed at p30. We studied the differential proteome of WT and KO livers in both severities of liver cystic disease (Figure 1), using both quantitative proteomic SWATH-MS and shotgun data-dependent acquisition DDA-MS analysis. Of the proteins with a significant change in abundance, we used bioinformatics tools to detect novel therapeutic targets and pathways involved in disease. Finally, we validated those targets using first in silico (DDA vs SWATH analysis) and finally in vivo strategies, as well as RT-qPCR, Western blotting and immunohistochemistry. p means postnatal day.

Protein Expression Pattern in Hepatic Cystogenesis
First, we performed proteomic mass spectrometry analysis by shotgun data-dependent acquisition DDA-MS or DDA. This analysis allows characterizing the proteome from complete and individual samples, providing the presence/absence of proteins with high Biomedicines 2022, 10, 290 9 of 21 sensitivity but without providing their quantification. We exclusively considered all the proteins expressed differentially in more than five independent samples from a total of six samples, for p14 and p12 stages of Wild Type and mutant individuals, respectively ( Figure S1). At p14, 1 WT-exclusive (or down-regulated for the disease) and 7 MUTexclusive (up-regulated) proteins were identified, while at p12, eight WT-exclusive and six MUT-exclusive proteins were observed ( Figure S1).
Next, we used a method that allowed us to quantify differentially expressed proteins. We performed a SWATH-MS quantitative proteomic analysis on the very same samples that we used for DDA, obtaining~1500 proteins per sample and more than 2000 proteins in each study group. First, we performed total area sum (TAS) normalization to evaluate that our samples followed a normal distribution pattern, as shown in Figures S2 and S3. Next, we studied the protein with significant differences between Wild Type and Mutant groups. Table 1 shows the proteins with a significant change in abundance with a two-fold increase (up-regulated) or decrease (down-regulated), and a significant adjusted p-value (according to the parametric Student's t-test) in WT vs MUT samples at p14 and p12. In total, six proteins showed significant differences in protein abundance between the WT-p14 and MUT-p14 groups (five up-regulated proteins and one down-regulated proteins). Between WT-p12 and MUT-p12, 26 proteins (20 proteins up-regulated and 6 proteins down-regulated) showed significant differences (Figure 3a and Table 1). Graphically, these variations can be observed through volcano plots, which were generated by plotting the log(2)-fold changes for all proteins identified against their −log(10) p-value (Figure 3b). Significant differences among protein abundance levels can be visualized in the heat map with individualized values according to the areas of the spectral library, and its clustering level in the cluster heat map analysis (Figures 3c and S4). These clusters detected by the SWATH-MS analysis help us to identify those qualitatively separated samples, in order to be considered for the quantitative analysis ( Figure S5). An unsupervised multivariate statistical analysis was performed using principal component analysis (PCA) to compare the data between samples ( Figure S6). The heatmap and PCA data demonstrated reproducibility among the sample triplicates as the three replicates closely clustered in both analyses. We can observe both in the PCA and in the heatmap cluster analysis ( Figures  S5 and S6) how some samples did not group correctly, overlapping between both clusters. Nonetheless, we can see the tendency to separate data in Wild Type and Mutant clusters. The lower number of proteins with significant differences in p14 vs p12 could be justified by the different degree of severity (Figure 2), suggesting that the number of proteins and pathways increases based on the severity and disease status of the cystic phenotype.

Clustering, Pathway Enrichment and Protein-Protein Interaction Analysis of Gene Expression in PLD
To establish the function related to the differentially expressed proteins identified by SWATH-MS analysis, we used several bioinformatic tools and forms of analyses. Gene Ontology (GO) terms and pathways analyses were performed using FunRich [39,40], String [41] and Reactome [42]. The proteins were sorted in FunRich gene enrichment analysis. Regarding the biological process, the most enriched group of proteins for p14 group (mild cystic) were related to fatty acid transport and prostaglandin biosynthesis, and for p12 (severe cystic) group with fibrinogen/fibrin activity, cell-cell/matrix adhesion and metabolic process (Figure 4a, upper panel). ANXA2 and the fibrinogen complex were the two most enriched cellular components in p14 and p12 groups, respectively. Furthermore, the extracellular space was the cellular component with the highest percentage of proteins in both groups (Figure 4a, lower panel). Similarly, the same analysis was established with string analysis, identifying different metabolic processes as the main biological process, and the extracellular space as the most enriched cell component; consistent with GO results ( Figure S7a). Finally, we repeated the analysis with the Reactome platform, finding that metabolism and immune system were the most enriched pathways ( Figure S7b).     Protein-protein or cluster analysis based on string analysis was also performed in order to investigate possible protein-protein interactions (interactions with experimental evidence). Numerous protein interactions were identified in different groups ( Figure S8), but one cluster at p14 and another at p12 were the most relevant based on the number of interactions between proteins and the level of expression (Figure 4b). The p14-cluster contains annexin A2 protein (ANXA2_P07356) associated with several cellular functions, such as fibrinolysis, cell motility (in epithelial cells), protein related to actin cytoskeleton and cell matrix interactions [43], interacting with two other proteins involved in the intracellular lipid transport (FABP1_P12710 and FABP5_Q05816). Interestingly, another relevant group of proteins was identified in the p12-cluster with a very strong protein-protein interaction between several fibrinogens (FGL1_Q71KU9, FIBB_Q8K0E8, FIBA_E9PV24 and FIBGG_Q8VCM7). Fibrinogens are involved in hepatocyte growth, and they mediate blood platelet spreading, interstitial collagen and fibrotic lesions [44]. Furthermore, additional proteins were found to interact with fibrinogens ( Figure S8), such as serum amyloid p-component (SAMP_P12246), hemopexin (HEMO_Q91X72) or hydroxyacid oxidase 2 (HAOX2_Q9NYQ2). Consistently with string results, ANXA2 and fibrinogen complexes were also the most enriched cellular components identified by FunRich analysis (Figure 4a). Interestingly, when comparing the data obtained by DDA and SWATH analysis, we confirmed that several proteins identified by SWATH were also identified by DDA analysis, including several related to fibrinogen complex proteins (SAMP/Apcs and S10A9/S100a9; Table S3).

Validation of SWATH-MS Analysis Unmasks the Fibrinogen Complex as a New Molecular Mechanism Related to PLD
Based on these data, we selected from the SWATH-MS analysis ( Figure S9) the altered proteins related to the fibrinogen complex at p12 and p14 stages for in vivo validation with RT-qPCR, Western blotting and immunohistochemistry, in order to establish their confirmation as possible targets of disease. Eight out of 10 of these selected proteins were validated at mRNA level by RT-qPCR ( Figure 5); one up-regulated at the 14-stage (ANXA2_P07356) and at p12-group (SAMP_P12246, FGL1_Q71KU9, ILK_O55222, S10A9_P31725, FIBB_Q8K0E8, HEMO_Q91X72, FIBA_E9PV24 and FIBG_Q8VCM7), and one down-regulated at the p12 group (HAOX2_Q9NYQ2). Interestingly, gene expression of the two proteins with more distance or fewer interactions to the fibrinogen group (ILK and S10A9) were the only ones that were not validated (Figure 5c).
Up-regulation of the fibrinogen complex (ANXA2, SAMP, FIBB, FIBA, FIBG and HAOX2 proteins) was further confirmed by Western blotting (WB) analysis ( Figure 6). Overall, gene expression and WB were highly consistent with SWATH-MS proteomic data. Finally, we also performed immunohistochemistry validation of the fibrinogen complex in mutant and polycystic liver samples (Figure 7). Fibrinogen staining was strongly upregulated in cyst-lining epithelial cells and bile duct dilatations (Figure 7a). These results unmask, for the first time, the fibrinogen complex as a possible pathway related to hepatic cystogenesis and as a possible future therapeutic target for PLD.
Overall, gene expression and WB were highly consistent with SWATH-MS proteomi data. Finally, we also performed immunohistochemistry validation of the fibrinogen com plex in mutant and polycystic liver samples (Figure 7). Fibrinogen staining was strongly up-regulated in cyst-lining epithelial cells and bile duct dilatations (Figure 7a). These re sults unmask, for the first time, the fibrinogen complex as a possible pathway related to hepatic cystogenesis and as a possible future therapeutic target for PLD. ; down-regulated p12 target hydroxyacid oxidase (HAOX2) (f). n = 6 was used for each protein group. GAPDH or β-ACTIN were used as housekeep ing genes. Bars represent means ± SEM. Student's t-test with two-tail was used and a value of p 0.05 was considered significant. * p < 0.05, ** p < 0.01, *** p < 0.001. (HAOX2) (f). n = 6 was used for each protein group. GAPDH or β-ACTIN were used as housekeeping genes. Bars represent means ± SEM. Student's t-test with two-tail was used and a value of p < 0.05 was considered significant. * p < 0.05, ** p < 0.01, *** p < 0.001. . n = 6 was used for each group. The immunohistochemistry analysis showed an up-regulation of fibrinogen through cystic epithelia. The samples corresponded to p12 group. Scale bars represents 100 μm. Bd means bile duct. Bars represent means ± SEM. Student's ttest with two-tail was used and a value of p < 0.05 was considered significant. *** p < 0.001.

Discussion
Polycystic liver disease (PLD) is commonly an extrarenal manifestation associated with autosomal dominant polycystic kidney disease (ADPKD). PLD is characterized by progressive dilation of the bile ducts and progression of multiple cysts, occupying at least half the volume of the liver parenchyma. Several of the current therapeutic strategies are based on surgical and pharmacological procedures to improve the symptoms of the disease. For this reason, the treatment of PLD has been ineffective so far, and the only curative option is liver transplantation [45]. This is mainly because the key intrinsic molecular mechanism of cystogenesis remains unknown. The main approaches are focused on the cAMP signaling pathway and somatostatin analogues, and several new targets (for example, histone deacetylase 6, Cdc25A phosphatase, PPAR-γ and matrix metalloproteases) have been evaluated in preclinical studies, but need to be tested clinically [8]. So, an emerging interest in understanding molecular pathways and developing new therapeutic strategies is increasing [2,8]. The purpose of present and novel SWATH-MS work was characterized by the proteomic changes in polycystic livers and provided new insights and therapeutic targets for PLD.
In previous studies, we have identified in an orthologous model of ADPKD (Pkd1 cond/cond ;Tam-Cre) a developmental window for kidney cystogenesis, suggesting that timing of secondary events may influence the severity of cystic kidney disease. Using the same strategy, we investigated whether this window of renal development corresponded in time with a similar window of hepatic development. Interestingly, we observed that cystic liver phenotype does not follow the same timing, suggesting that kidney and liver have possible different mechanisms undergoing cystogenesis. Further studies should be focused on identifying the common and different mechanisms of disease in both organs, and identifying the unknown liver developmental window for cystogenesis. For further analysis, we used a next generation quantitative proteomic approach SWATH-MS to perform differential proteome analysis of cystic liver disease in a mild and more severe PLD phenotype within the same kidney developmental window. Our data suggest that in advanced stages of disease progression, the number of pathways deregulated are increased, probably due to secondary effects of the disease. In the same way, these results reinforce the idea that cyst growth passes with two big stages: cyst initiation and cyst progression [46]. Interestingly, both proteome analyses from a mild (p14 group) and a severe (p12 group) cystic status identified protein complexes of fibrinolysis directly related to liver . n = 6 was used for each group. The immunohistochemistry analysis showed an up-regulation of fibrinogen through cystic epithelia. The samples corresponded to p12 group. Scale bars represents 100 µm. Bd means bile duct. Bars represent means ± SEM. Student's t-test with two-tail was used and a value of p < 0.05 was considered significant. *** p < 0.001.

Discussion
Polycystic liver disease (PLD) is commonly an extrarenal manifestation associated with autosomal dominant polycystic kidney disease (ADPKD). PLD is characterized by progressive dilation of the bile ducts and progression of multiple cysts, occupying at least half the volume of the liver parenchyma. Several of the current therapeutic strategies are based on surgical and pharmacological procedures to improve the symptoms of the disease. For this reason, the treatment of PLD has been ineffective so far, and the only curative option is liver transplantation [45]. This is mainly because the key intrinsic molecular mechanism of cystogenesis remains unknown. The main approaches are focused on the cAMP signaling pathway and somatostatin analogues, and several new targets (for example, histone deacetylase 6, Cdc25A phosphatase, PPAR-γ and matrix metalloproteases) have been evaluated in preclinical studies, but need to be tested clinically [8]. So, an emerging interest in understanding molecular pathways and developing new therapeutic strategies is increasing [2,8]. The purpose of present and novel SWATH-MS work was characterized by the proteomic changes in polycystic livers and provided new insights and therapeutic targets for PLD.
In previous studies, we have identified in an orthologous model of ADPKD (Pkd1 cond/cond ;Tam-Cre) a developmental window for kidney cystogenesis, suggesting that timing of secondary events may influence the severity of cystic kidney disease. Using the same strategy, we investigated whether this window of renal development corresponded in time with a similar window of hepatic development. Interestingly, we observed that cystic liver phenotype does not follow the same timing, suggesting that kidney and liver have possible different mechanisms undergoing cystogenesis. Further studies should be focused on identifying the common and different mechanisms of disease in both organs, and identifying the unknown liver developmental window for cystogenesis. For further analysis, we used a next generation quantitative proteomic approach SWATH-MS to perform differential proteome analysis of cystic liver disease in a mild and more severe PLD phenotype within the same kidney developmental window. Our data suggest that in advanced stages of disease progression, the number of pathways deregulated are increased, probably due to secondary effects of the disease. In the same way, these results reinforce the idea that cyst growth passes with two big stages: cyst initiation and cyst progression [46]. Interestingly, both proteome analyses from a mild (p14 group) and a severe (p12 group) cystic status identified protein complexes of fibrinolysis directly related to liver cystic disease. At the milder stages of the disease, we identified ANXA2 as one of the interesting proteins to validate as a potential therapeutic target for PLD. ANXA2 is a calcium-regulated membrane-binding protein produced by a wide range of cell types, including epithelial, dendritic, trophoblast, and tumor cells, as well as monocytes and macrophages. The main role is to maintain cell surface proteolytic activity, but also fulfills a range of intracellular functions, including exocytosis, endocytosis, membrane repair [47] and maintenance of adherent-like intercellular junctions [48]. Its numerous functions contribute to fibrinolysis, regulation of inflammation and immune system activation, and tissue injury and repair; as a result, ANXA2 dysfunction has been implicated in multiple human diseases, for example, extracellular matrix remodeling and hepatic fibrosis [49,50], which is a common feature in disease progression of PLD [8]. Interestingly, at the severe stage of liver disease, our proteome analysis discovered a number of different proteins related to the fibrinogen complex. During tissue and vascular injury, fibrinogen is converted enzymatically by thrombin to fibrin, and during fibrinolysis, fibrin is degraded by the main enzyme plasmin (the active form of plasminogen). Fibrinogens are associated with MAPK [51,52] and integrins [53], two pathways that have been shown to be effective as potential therapy for ADPKD in preclinical studies [54,55]. Even the animal model of one of the overexpressed fibrinogens (FIBA) has been shown to have cystic disease [56]. Moreover, the increased activity of fibrinogen-clusters in PLD offers a new perspective to treat the disease in comparison to the current experimental drugs [57]. Consequently, the decrease in this cluster could support a new mechanism of actions for drugs in which fibrinogen family proteins would be the therapeutic targets. A greater example would be the blockade of ERRγ (estrogen related receptor γ), an orphan nuclear receptor [58,59], which modulates fibrinogen levels in hypofibrinogenemia states caused by diet-induced obesity, diabetes mellitus type 2, liver injury and alcohol-induced oxidative stress [60,61]. Specifically, inverse agonists, such as GSK5182 (CAS: 877387-37-6) [60][61][62][63][64][65][66] and DN200434 [67,68], would offer such beneficial action through the decrease of fibrinogen levels.
In summary, we identified several novel targets and pathways involved in physiopathology of PLD by quantitative proteomic SWATH-MS technology. We emphasized the novel therapeutic opportunity presented by fibrinogens and fibrinolysis (and the other related identified proteins), mainly due the fact that all validations with significant differences carried out are related to the fibrinogen complex. Nevertheless, there is potential to confirm this in future preclinical studies as well as more functional studies to characterize the exact role of our target candidates in PLD. In conclusion, this work has created a new mechanism and opportunity for future research in PLD physiopathology, leading to possible new therapeutic approaches of the disease.

Summary
We reported a novel proteomic quantitative study based on SWATH-MS technology comparing proteomes of Wild Type and polycystic livers. Several novel pathways and targets were identified, expanding the knowledge about molecular mechanisms in the PLD field and highlighting the fibrinogen complex as a possible new therapeutic target.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/biomedicines10020290/s1, Figure S1: Qualitative DDA mass spectrometry analysis in the different groups of the study, Figure S2: TAS (Total Area Sums) normalization plots p14 group, Figure S3: TAS (Total Area Sums) normalization plots of p12 group, Figure  S4: Heatmap figure of unsupervised cluster analysis of proteins with significant differences according to SWATH-MS analysis, Figure S5: Unsupervised heatmap cluster analysis of all proteins detected by SWATH-MS, Figure S6: Unsupervised PCAs from SWATH-MS data, Figure S7: GO enrichment and pathway analysis of protein significantly regulated in hepatic cystogenesis according to String and Reactome, respectively, Figure S8: Full protein-protein interaction map according String, Figure  S9: Box-plots of quantitative proteomic SWATH-MS data of selected targets for validation, Table S1: ID Analyzed samples, distribution in each group and samples pooling for the acquisition of shotgun runs to build the spectral library, Table S2: Primers sequence list used for RT-qPCR analysis, Table  S3: List of proteins presented in 5 out of 6 independent samples of different groups after qualitative DDA analysis.

Informed Consent Statement: Not applicable.
Data Availability Statement: The data presented in this study are available in this article and this Supplementary Material. In addition, The mass spectrometry proteomics data have been deposited to the ProteomeXchange Consortium via the PRIDE [69] partner repository with the dataset identifier PXD031253.