Broad Transcriptomic Impact of Sorafenib and Its Relation to the Antitumoral Properties in Liver Cancer Cells

  Citation: Contreras, L.; Rodríguez-Gil, A.; Muntané, J.; de la Cruz, J. Broad Transcriptomic Impact of Sorafenib and Its Relation to the Antitumoral Properties in Liver Cancer Cells. Cancers 2022,14, 1204. https://doi.org/10.3390/ cancers14051204 Academic Editor: Hiroyuki Tsuchiya Received: 4 February 2022 Accepted: 21 February 2022 Published: 25 February 2022 Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations. Copyright: © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/). cancers Article Broad Transcriptomic Impact of Sorafenib and Its Relation to the Antitumoral Properties in Liver Cancer Cells Laura Contreras 1,2, Alfonso Rodríguez-Gil 1,3,4 , Jordi Muntané1,4,5,*,† and Jesús de la Cruz 1,2,*,† 1Instituto de Biomedicina de Sevilla, Hospital Universitario Virgen del Rocío/CSIC/Universidad de Sevilla, E-41013 Seville, Spain; [email protected] (L.C.); ar[email protected] (A.R.-G.) 2Departamento de Genética, Facultad de Biología, Universidad de Sevilla, E-41012 Seville, Spain 3Centro de Investigación Biomédica en Red de Cáncer (CIBERONC), E-28029 Madrid, Spain 4Departamento de Fisiología Médica y Biofísica, Universidad de Sevilla, E-41009 Sevilla, Spain 5Centro de Investigación Biomédica en Red de Enfermedades Hepáticas y Digestivas (CIBEREHD), E-28029 Madrid, Spain *Correspondence: [email protected] (J.M.); [email protected] (J.d.l.C.); Tel.: +34-955-923-122 (J.M.); +34-923-126 (J.d.l.C.) † These authors contributed equally to this work. Simple Summary: Hepatocellular carcinoma (HCC) is the fourth most frequent cause of cancerrelated mortality worldwide. While ablation, resection and orthotopic liver transplantation are indicated at an early stage of the disease, Sorafenib (Sfb) is the current most administrated first-line treatment for advanced HCC, even though its therapeutic benefit is limited due to the appearance of resistance. Deep knowledge on the molecular consequences of Sfb-treatment is essentially required for optimizing novel therapeutic strategies to improve the outcomes for patients with advanced HCC. In this study, we analyzed differential gene expression changes in two well characterized liver cancer cell lines upon a Sfb-treatment, demonstrating that both lines responded similarly to the treatment. Our results provide valuable information on the molecular action of Sfb on diverse cellular fundamental processes such as DNA repair, translation and proteostasis and identify rationalization issues that could provide a different therapeutic perspective to Sfb. Abstract: Hepatocellular carcinoma (HCC) is one of the most frequent and essentially incurable cancers in its advanced stages. The tyrosine kinase inhibitor Sorafenib (Sfb) remains the globally accepted treatment for advanced HCC. However, the extent of its therapeutic benefit is limited. Sfb exerts antitumor activity through its cytotoxic, anti-proliferative and pro-apoptotic roles in HCC cells. To better understand the molecular mechanisms underlying these effects, we used RNA sequencing to generate comprehensive transcriptome profiles of HepG2 and SNU423, hepatoblastoma- (HB) and HCC-derived cell lines, respectively, following a Sfb treatment at a pharmacological dose. This resulted in similar alterations of gene expression in both cell lines. Genes functionally related to membrane trafficking, stress-responsible and unfolded protein responses, circadian clock and activation of apoptosis were predominantly upregulated, while genes involved in cell growth and cycle, DNA replication and repair, ribosome biogenesis, translation initiation and proteostasis were downregulated. Our results suggest that Sfb causes primary effects on cellular stress that lead to upregulation of selective responses to compensate for its negative effect and restore homeostasis. No significant differences were found specifically affecting each cell line, indicating the robustness of the Sfb mechanism of action despite the heterogeneity of liver cancer. We discuss our results on terms of providing rationalization for possible strategies to improve Sfb clinical outcomes. Keywords: hepatocellular carcinoma cell line; gene ontology; RNA-Seq; RNA synthesis; Sorafenib Cancers 2022,14, 1204. https://doi.org/10.3390/cancers14051204 https://www.mdpi.com/journal/cancers Cancers 2022,14, 1204 2 of 18 1. Introduction Hepatocellular carcinoma (HCC) is the most common type of primary liver cancer in adults [ 1 ]. HCC has been described as the sixth most common neoplasia and the fourth most frequent cause of cancer-related mortality in men and women worldwide, being the main cause of death in people with cirrhosis (see [ 2 ] and references therein). The prognosis of HCC is strictly determined by the staging of the tumor and the hepatic function of the patients. The Barcelona-Clinic Liver Cancer Classification (BCLC) is currently a useful method for staging patients and recommending treatments depending on the number and size of the tumor nodules, liver function, the presence of vascular invasion, and extrahepatic metastasis [ 3 , 4 ]. While ablation, resection, and orthotopic liver transplantation are indicated at very early and early stages of the disease [ 5 , 6 ], Sfb is the current most administered drug for the treatment of patients advanced HCC [ 6 – 8 ]. Unfortunately, Sfb only provides very limited clinical benefit to patients, the five-year survival rate remaining extremely low [ 7 , 8 ]. Thus, there is a need to develop novel therapies that alone or in combination with Sfb could improve the outcomes for patients with advanced HCC [2,6,9]. Sfb is a multiple tyrosine kinase receptor inhibitor such as the vascular endothelial growth factor receptors 2 and 3 (VEGFR2 and VEGFR3), the platelet-derived growth factor receptor beta (PDGFR-ß), FLT3 and c-KIT, as well as the Raf kinases, which are integral components of the Ras/Raf/mitogen-activated protein (MAP)/extracellular signal-regulated kinase (ERK) kinase (MEK)/ERK signaling cascade [ 10 , 11 ]; thus, Sfb also downregulates the phosphoinositide-3 kinase (PI3K)/Akt/mammalian target of rapamycin (mTOR) that regulates fundamental cellular processes [ 12 ]. Sfb exerts a potent anti-proliferative and pro-apoptotic activity against HCC cells [ 13 , 14 ]. We and others have shown that this is likely due to the generation of a potent endoplasmic reticulum (ER) stress that leads to the sequential induction of autophagy and apoptosis processes in hepatoblastoma (HB) and HCC cell lines and tumor-derived xenograft mice models [ 13 , 15 ]. However, many of the molecular events by which Sfb exerts its antitumor activity remain unclear. Previous studies have examined the expression level of numerous transcripts in selected HB and HCC cell lines, which differ in their genetic characteristics, treated with Sfb using DNA microarray technology [ 10 ]. The aim of the present study was to further determine the impact of Sfb in HB and HCC cells using high-throughput RNA-Seq technology. The study was performed in two different liver cancer cell lines differing in differentiation stage and the expression of p53: HepG2 (well-differentiated HB cell line; wild-type p53) and SNU423 (moderately differentiated HCC cell line; in-frame p53 gene truncation of amino acids 126–132) [ 16 ]. Herein, we infer, from our RNA-Seq analysis, the specific pathways that seem to be activated and inhibited upon Sfb treatment and compare the obtained results with previously published ones obtained using microarray analysis of gene expression. Both, similarities and differences in the transcriptome with and without Sfb in the two cells lines were discussed. 2. Materials and Methods 2.1. Cell Lines, Culture Conditions, and Sorafenib Treatment HepG2 and SNU423 cell lines were obtained from the American Type Culture Collection (ATCC/LGC Standards, S.L.U., Barcelona, Spain). As described earlier, cell lines were selected according to their origin, cellular differentiation stage, and p53 genetic status: HepG2 (wildtype p53. HB-8065 ™ ) [ 17 ] and SNU423 (in-frame p53 gene deletion. CRL-2238) [ 18 ]. Both cell lines were negative for mycoplasma contamination. Cells were cultured in minimal essential medium (MEM) with Earle’s balanced salts with L-glutamine (ref E15–825, PAA Laboratories Inc., Toronto, ON, Canada) supplemented with 10% fetal bovine serum (FBS, F7524, Sigma Aldrich, Lot No. 022M3395, endotoxin < 0.2 EU/mL), 1% sodium pyruvate (ref. S11–003, PAA Laboratories Inc.), 1% non-essential amino acids (ref. M11–003, PAA Laboratories Inc., Toronto, Canada) and penicillin–streptomycin solution (100 U/mL–100 µ g/mL; ref. P11–010, PAA Laboratories Inc., Toronto, Canada); cells were grown in culture flasks at 37 ◦ C in a humidified incubator with 5% CO 2 until reaching a density of 100,000 cells/cm 2 . Sfb was Cancers 2022,14, 1204 3 of 18 added at the concentration of 10 µ M at 24 h after plating, and lysates were obtained after 12 h treatment. Sorafenib (Sfb, ref. FS10808; Carbosynth Ltd., Compton, UK) was dissolved in dimethyl sulfoxide (DMSO) as a stock solution (100 mM). 2.2. Cell Proliferation The measurement of bromodeoxyuridine (BrdU) incorporation was completed as exactly described in [ 13 ] using a commercial kit (ref. 11 647229001, Roche Diagnostics, Mannheim, Germany). Cells were seeded in 96-well plates at low density (15,000 cells/cm 2 ). Two hours before cell harvesting, 20 µ L of 10 µ M BrdU was added to the cultures. DNA was denaturalized with 200 µ L FixDenat solution included in the commercial kit for 30 min at room temperature. After removal, cells were incubated with 100 µ L of monoclonal anti- BrdU antibody HRP conjugated for 90 min at room temperature. Afterwards, cells were washed with PBS buffer (137 mM NaCl, 2.7 mM KCl, 10 mM Na 2 HPO 4 , 1.8 mM KH 2 PO 4 , pH 7.4) and incubated with 100 µ L revealing solution for 15 min at room temperature. Absorbance at 370 nm (A 370 ) was measured using an Infinite 200 PRO microplate reader (Tecan, Männedorf, Switzerland). Cell cycle progression was assessed by flow cytometric analysis. For this, cells were seeded 24 h before treatment to a final confluence of 70% in a 6-well plate. Afterwards, cells were treated with 10 µ M Sfb for 12 h and then harvested and fixed in 70% ethanol in PBS buffer overnight at 4 ◦ C. After this, cells were resuspended in PBS buffer and incubated with 0.5 mg/mL RNase A during 1 h at 37 ◦ C. Propidium Iodine (PI) was added to a final concentration of 0.05 mg/mL for 20 min at room temperature. Finally, cells were filtered to avoid aggregation and cell cycle progression was assessed using a FACSCantoTM Flow Cytometer and analyzed using the FACSDiva software (BD Biosciences, Allschwil, Switzerland). 2.3. RNA Extraction Total RNA was extracted from each sample using a RNeasy mini kit according to the manufacturer’s instructions (QIAGEN, Hilden, Germany). RNA from each sample was then stored at −80 ◦C before further analyses. 2.4. mRNA Library Preparation Libraries were prepared using RNA of cells treated either with the vehicle or with Sfb. Total RNA was prepared as described above. Then, concentration and quality of the RNA were assessed with Qubit (Qubit ™ DNA HS assay, Thermo Fisher Scientific, Waltham, MA, USA) and a 2100 Bioanalyzed Nano Chip (Agilent Technologies Genomics, Santa Clara, CA, USA), respectively. RNA Integrity Number (RIN) values were > 9 in all RNA samples. Polyadenylated RNA was isolated from the total RNA using NEBNext Oligo d(T) 25 beads (New England Biolabs, Ipswich, MA, USA) according to the manufacturer’s instructions. Samples were normalized to an equivalent concentration of 67.3 ng/ µ L and prepared for RNA Ilumina Sequencing. For HepG2 cell lines, three biological replicates were obtained, while for SNU423 cells only two biological replicates were collected. 2.5. RNA Sequencing and Data Analyses RNA sequencing was performed with the NextSeq500 Mid-Output and 2 × 75 pb length parameters (paired-end). RNA-Seq data were first filtered using the FASTQ Toolkit v1.0.0 program and then analyzed using the BaseSpace Onsite v3.22.91.158 from Illumina. Only genes that were upregulated or downregulated with a p-value < 0.5 and [log 2 (fold changes)] ≥± 0.5 were selected. Data presented in this study has been submitted to the Gene Expression Omnibus database under the accession number GSE186280. 2.6. Gene Set Enrichment Analysis and Over-Representation Analysis To identify functional categories significantly affected by the Sfb treatment, Gene Set Enrichment Analysis (GSEA) [ 19 , 20 ] and Over-Representation Analysis (ORA) were performed using the Reactome database (https://reactome.org; last accessed 3 February Cancers 2022,14, 1204 4 of 18 2022) [ 21 ] downloaded from the Molecular Signatures Database v. 7.1 (https://www.gseamsigdb.org; last accessed 3 February 2022). For GSEA, genes were ranked according to the p-value and the direction of change (up or downregulated) obtained from the RNA-Seq data. For the combined analysis of both cell lines, the product of the p-value of each cell line was used for the ranking. GSEA was performed using the GSEA software with the default parameters for pre-ranked lists. ORA was performed using a Fisher’s exact test [ 22 ], selecting the genes with a p-value lower than 0.001 in the RNA-Seq data, either in each cell line or in both cell lines. 2.7. Real-Time Quantitative PCR To confirm mRNA-Seq results, real-time quantitative PCR (RT-qPCR) was performed on a set of selected genes. Total RNA was obtained as previously described and RNA samples were treated with DNase I (Promega, Madison, WI, USA) to remove all traces of DNA. After this treatment, DNase I was inactivated by incubation of the samples at 65 ◦ C for 10 min. RNA was then reverse transcribed using SuperScript ™ III First Strand Synthesis for RT-PCR according to the manufacturer’s instruction (Invitrogen, Walthman, MA, USA) and random hexamer primers (Roche, Switzerland). Each reaction contained 250 ng of RNA in a total volume of 25 µ L. RT-qPCR was performed using SYBR ® Green Premix Ex Taq ™ 2X (Takara, Japan) and primer specifics of each transcript. The ribosomal RNA 28S (28S rRNA) was used as an internal normalization control. Primers pair used for the RT-qPCR were as follows: β -actin mRNA: 5 0 TCCCTGGAGAA GAGCTACGA3 0 (forward) and 5 0 AGGAAGGAAGGCTGGAAGAG3 0 (reverse); BAX mRNA: 5 0 TCCACCAAGAAGCTGAGCGAG3 0 (forward) and 5 0 GTCCAGCCCATGATGGTTCT3 0 (reverse); BIM mRNA: 5 0 CACCAGCACCATAGAAGAA3 0 (forward) and 5 0 ATAAGG AGCAGGCACAGA3 0 (reverse); BIRC3 mRNA 5 0 ATGCTTCTGTTGTGGCCTGA3 0 (forward) and 5 0 ACTCTGAACGAATCTGCAGCT3 0 (reverse); BOP1 mRNA: 5 0 CTGATTCACC AGCTGAGCC30(forward) and 50GACGCCACCAACAGGAAG30(reverse); CEBPβmRNA: 5 0 AAACTCTCTGCTTCTCCCTCTGC3 0 (forward) and 5 0 CTGACAGTTACACGTGGGTTGC3 0 (reverse); CPEB4 mRNA 5 0 CACCAACACCCTCCTCTTCC3 0 (forward) and 5 0 TTCAGGGG CGTTATTCCACC3 0 (reverse); DUSP1 mRNA: 5 0 CCTGAGTACTAGCGTCCCTG3 0 (forward) and 5 0 CAGGTACAGAAAGGGCAGGA3 0 (reverse); EIF4E2 mRNA: 5 0 TGAAAGAT GATGACAGTGGGGA3 0 (forward) and 5 0 CTGATTCTTGTCTCGTTCCGT3 0 (reverse); EPOP mRNA 5 0 AGTTTTCGGGGTGACAGTCC3 0 (forward) and 5 0 AGATGGAAGGAGGCAGG GAT3 0 (reverse); FEN1 5 0 TGGGGTCAAGAGGCTGAGTA3 0 (forward) and 5 0 GTGGATCCC TTGGGTTCTGG3 0 (reverse); GADD45B 5 0 TGGGAAGGTTTTGGGCTCTC3 0 (forward) and 5 0 TCCAGCGTCATGTTGCAATTATA3 0 (reverse); IDI1 5 0 AACCACCTCGACAAGCAACA3 0 (forward) and 5 0 TGTTCTCGTTCAGGTGACAA3 0 (reverse); PCNA 5 0 AAAGTCCAAAGTC AGATCTGGTC3 0 (forward) and 5 0 ACTGCATTTAGAGTCAAGACCCT3 0 (reverse); PHB mRNA: 5 0 TCAACATCACACTGCGCATC3 0 (forward) and 5 0 ATAGTCCTCTCCGATGCTG G3 0 (reverse); SMAD7 mRNA: 5 0 CCCCTCCTCTCCCTCATCAA3 0 (forward) and 5 0 GGCTGG CAGGAAGGGAATAA3 0 (reverse); TPI mRNA: 5 0 GGACTCGGAGTAATCGCCTG3 0 (forward) and 5 0 TGTTGGGGTGTTGCAGTCTT3 0 (reverse); VEGFA 5 0 CCATCCAATCGAGACC CTGG3 0 (forward) and 5 0 CTCCAGGCCCTCGTCATTG3 0 (reverse); 28S rRNA: 5 0 CAAAGCG GGTGGTAAACTCC30(forward) and 50TTCACGCCCTCTTGAACTCT30(reverse). 2.8. Western Blot Analysis Protein extracts were obtained by lysing cell pellets at 100 ◦ C for 10 min in 2 × Laemmli buffer (125 mM HCl-Tris, pH 6.8, 4% SDS, 0.02% bromophenol blue, 20% glycerol, 200 mM DTT). Cellular extracts were then sonicated in a Bioruptor (Diagenode, Liège, Belgium) for 1 min at high intensity. Protein extracts were subjected to 10% SDS-PAGE and transferred to nitrocellulose membranes (Amersham TM Protran ® 0.45 µ m, GE Healthcare Chicago, IL, USA). The membranes were blocked for 1 h with 5% skim milk in TTBS (15 mM HCl-Tris, pH 7.5, 200 mM NaCl, 0.2 M NaCl, 0.1% (v/v) Tween-20), followed by incubation with primary anti-NDUFS1 (Abcam, Cambridge, UK, ab169540, 1:10,000 dilution), anti-NDUFS2 Cancers 2022,14, 1204 5 of 18 (Abcam, ab103024, 1:2000 dilution), anti-NDUFV2 (Abcam, ab183715, 1:2000 dilution) and anti-GAPDH (Santa Cruz, sc-47724, 1:1000 dilution) at 4 ◦C overnight. After washing with TTBS buffer, the membranes were incubated with HRP-conjugated secondary antibody (Bio-Rad, Hercules, CA, USA) for 1 h at a 1:5000 dilution at room temperature. Proteins were detected using an enhanced chemiluminescence detection kit (Pierce TM Super-Signal West Pico, Thermo Fisher Scientific, Waltham, MA, USA) in a ChemiDoc TM Touch Imaging System (Bio-Rad) and the relative intensity value quantified with the Image Lab software provided with this system. 2.9. Statistical Analyses Statistical analyses were performed with the Prism 6.01 software (GraphPad, San Diego, CA, USA). Data were generated from several repeats (at least three) of different biological replicates (at least three). Mean ± S.D. were represented in the different graphs. To determine significance, Student’s tests for unpaired samples with confidence intervals of 95% were computed. Significance between conditions were indicated with the symbols *p< 0.05 , ** p< 0.01, *** p< 0.001, and **** p< 0.0001. Regression plots and determination of Pearson coefficients and p-value were performed using the R software (Institute for Statistics and Mathematics, Vienna, Austria). A Venn diagram was computed using EulerAPE software [23]. 3. Results and Discussion 3.1. Transcriptional Changes Caused by Sorafenib in Hepg2 Hepatoblastoma and SNU423 Hepatocellular Carcinoma Cell Lines Sfb is an inhibitor of several kinases involved in tumor cell proliferation and angiogenesis, including Raf, VEGFR and PDGFR [ 10 , 24 ]. Nowadays, Sfb is one of the most used molecular targeted drugs for the treatment of advanced inoperable HCC with significant but unfortunately modest anticancer results [ 7 , 8 ]. To better understand the biological consequences of the treatment, two different model liver cancer cell lines (HepG2 and SNU423 cell lines) were treated with Sfb 10 µ M for 12 h and their transcriptomes analyzed by RNA-Seq and compared to those of untreated cells grown in the same conditions. It has been previously shown that at this concentration, Sfb induces apoptosis and significantly reduces cell proliferation in cultures of the two liver cancer cell lines [16]. Principal component analyses demonstrated strong consistency between repeats of each sample and clearly separated data from treated vs. untreated cells in both cell lines (data not shown). Numerous changes were observed in both cell lines when treated with Sfb (see Table S1). When we took into consideration changes in gene expression with an established p-value lower than 0.001, our RNA-Seq analysis identified 2140 and 1347 differentially expressed genes (DEGs) upon Sfb treatment in HepG2 and SNU423 cells, respectively. From these, 1166 and 677 were upregulated while 974 and 670 were downregulated, respectively, Table S2 displays the list of DEGs (p-value < 0.001) and highlights those genes with log 2 (FC) either higher than 1.5 or lower than − 1.5 in Sfbtreated HepG2 (265 vs. 298 genes) and SNU423 (264 vs. 175 genes) cells, respectively. All these results are summarized in Figure 1A. To verify whether the RNA-Seq results were valid, we first checked a set of genes using RT-qPCR, chosen in a random way. In HepG2 cells, we checked ACTB ( β -actin), BIRC3, which encodes an inhibitor of apoptosis, CEBP β , which encodes for a leucinezipper transcriptional factor regulator of, among others, genes involved in the immune and inflammatory responses, CPEB4, a RNA binding protein that regulates activation of UPR, DUSP1, which encodes for the Dual Specificity Protein Phosphatase 1, GADD45B, the growth arrest and DNA damage-inducible protein, SMAD7, which encodes for a nuclear protein involved in the inhibition of the TGF-beta receptor and VEGFA, the vascular endothelial growth factor A as genes showing a positive log 2 (FC), and BAX, a BCL2 family member with a role as a mitochondrial apoptotic activator, the gene encoding the translation initiation factor eIF4E2, EPOP encoding for an scaffold protein, FEN1 encoding for a DNA Cancers 2022,14, 1204 6 of 18 nuclease, IDI1 encoding for a peroxisomal enzyme, PCNA, which encodes for a cofactor of the replicative DNA polymerase and the gene for the triosephosphate isomerase (TPI) as genes showing a negative log 2 (FC). In all cases, a significant and similar trend was obtained in the levels expression of the genes tested in Sfb-treated versus non-treated cells by RT-qPCR or RNA-Seq analysis (Figure S1); thus, our RT-qPCR analysis was consistent with the obtained RNA-Seq data. In SNU423 cells, we checked BIM, which encodes for a proapoptotic protein, BOP1 encoding for a factor involved in the assembly of 60S ribosomal subunits BIRC3,CPEB4,DUSP1,GADD45B, and SMAD7 as gene showing a positive log 2 (FC), and eIF4E2, EPOP,IDI1,PHB (Prohibitin), which encodes for a protein proposed to play an antiproliferative role, PCNA and TPI as genes showing a negative log 2 (FC). Similarly, gene expression levels obtained using RNA-Seq and RT-qPCR methods were in good agreement each other with only small variation in the magnitude of expression (Figure S2). Cancers 2022, 14, x 6 of 18 Figure 1. Differentially expressed genes (more than two-fold) in HepG2 and SNU423 cells upon Sfb treatment (10 µM, 12 h) detected using RNA-Seq. (A) Venn diagram depicting the unique and shared upregulated and downregulated genes upon treatment. The p-value cut-off was of 0.001. (B) A correlation analysis showing the relation between the RNA-Seq data of the HepG2 versus SNU423 cells upon the Sfb treatment. The regression lines for all genes, and those with p-values lower than 0.001 for HepG2, SNU423 or both cell lines are shown with the corresponding linear regression equations and Pearson coefficients. To verify whether the RNA-Seq results were valid, we first checked a set of genes using RT-qPCR, chosen in a random way. In HepG2 cells, we checked ACTB (β-actin), BIRC3, which encodes an inhibitor of apoptosis, CEBPβ, which encodes for a leucine-zip- per transcriptional factor regulator of, among others, genes involved in the immune and inflammatory responses, CPEB4, a RNA binding protein that regulates activation of UPR, DUSP1, which encodes for the Dual Specificity Protein Phosphatase 1, GADD45B, the growth arrest and DNA damage-inducible protein, SMAD7, which encodes for a nuclear protein involved in the inhibition of the TGF-beta receptor and VEGFA, the vascular endothelial growth factor A as genes showing a positive log2(FC), and BAX, a BCL2 family member with a role as a mitochondrial apoptotic activator, the gene encoding the translation initiation factor eIF4E2, EPOP encoding for an scaffold protein, FEN1 encoding for a DNA nuclease, IDI1 encoding for a peroxisomal enzyme, PCNA, which encodes for a cofactor of the replicative DNA polymerase and the gene for the triosephosphate isomerase (TPI) as genes showing a negative log2(FC). In all cases, a significant and similar trend was obtained in the levels expression of the genes tested in Sfb-treated versus non-treated cells by RT-qPCR or RNA-Seq analysis (Figure S1); thus, our RT-qPCR analysis was consistent with the obtained RNA-Seq data. In SNU423 cells, we checked BIM, which encodes for a proapoptotic protein, BOP1 encoding for a factor involved in the assembly of 60S ribosomal subunits BIRC3, CPEB4, DUSP1, GADD45B, and SMAD7 as gene showing a positive log2(FC), and eIF4E2, EPOP, IDI1, PHB (Prohibitin), which encodes for a protein proposed to play an antiproliferative role, PCNA and TPI as genes showing a negative log2(FC). Similarly, gene expression levels obtained using RNA-Seq and RT-qPCR methods were in good agreement each other with only small variation in the magnitude of expression (Figure S2). Interestingly, when the RNA-Seq data were globally analyzed, a strong overlap between the set of genes misregulated by Sfb in HepG2 and SNU423 cells was observed. Figure 1. Differentially expressed genes (more than two-fold) in HepG2 and SNU423 cells upon Sfb treatment (10 µ M, 12 h) detected using RNA-Seq. ( A ) Venn diagram depicting the unique and shared upregulated and downregulated genes upon treatment. The p-value cut-off was of 0.001. (B) A correlation analysis showing the relation between the RNA-Seq data of the HepG2 versus SNU423 cells upon the Sfb treatment. The regression lines for all genes, and those with p-values lower than 0.001 for HepG2, SNU423 or both cell lines are shown with the corresponding linear regression equations and Pearson coefficients. Interestingly, when the RNA-Seq data were globally analyzed, a strong overlap between the set of genes misregulated by Sfb in HepG2 and SNU423 cells was observed. Indeed, changes in gene expression altered by Sfb treatment were highly positively correlated in both cell lines (Figure 1B). The linear regression for the common DEG genes shows a very high correlation (y = 1.055 × − 0.0158), with a value of the Pearson Correlation coefficient as high as r= 0.9294. This result indicates that despite differences among the set of DEGs in each liver cancer cell line, the overall response to the Sfb treatment is very similar between both cell lines. Moreover, no special features regarding gene, ORF, 5 0 UTR or 3 0 UTR length, or GC content were found among the DEGs upon a Sfb treatment (data not shown). Focusing on HepG2 cells, which correspond to the cell line routinely studied in our laboratory, when genes are analyzed individually (see Figure 2A), those showing the high- Cancers 2022,14, 1204 7 of 18 est positive log 2 (FC) with the lowest p-value were mainly stress-response genes such as INHBE,ATF3,NUPR1,PPP1R13B,PPP1R15A (also known as GADD34), DDIT3,DUSP1, DUSP8, or ERN1, as well as genes such as PDGFA (Platelet Derived Growth Factor subunit A) and VEGFA (Vascular Endothelial Growth Factor A), which correspond to growth factors that activate the PDGFR and VEGFR receptors, respectively, in charge of transducing their extracellular signals into the cell and described as direct targets of Sfb [ 10 , 24 ]. Their upregulation by Sfb might be related to secondary cell response as a consequence of the inhibition of PDGFR and VEGFR receptors. Consistently, most abovementioned genes regulate important cellular processes leading to inhibition of cell growth and proliferation and/or induction of apoptosis. As examples, INHBE is a member of the TGF-beta (transforming growth factor-beta) superfamily of protein genes, which encodes a preproprotein that requires proteolytical processing to generate the inhibin beta E subunit. Inhibins downregulate different cellular processes, among them cell proliferation and apoptosis [ 25 , 26 ]. In agreement with our data, it has been reported that INHBE is upregulated under conditions of endoplasmic reticulum (ER) stress, in which it is involved in the Sfb signaling [ 13 , 27 ]. In addition, DUSPs dephosphorylate many key signaling molecules, including MAPKs, leading to the reduction in the duration, magnitude, or spatiotemporal profiles of the activities of MAPKs [ 28 ]. ATF3, a gene for a transcriptional factor that functions in general adaptive responses and whose expression is induced by various stimuli including ER stress [ 29 ], is also upregulated upon the Sfb treatment. NUPR1 is a transcription regulator that induces autophagy and apoptosis through upregulation of ER stress-related factors including DDIT3 [ 30 ]. DDIT3, also known as CHOP or GADD153, is a central transcriptional factor, also induced by ER stress, that triggers apoptosis through inhibition of BCL2 and upregulation of BIM, which regulate BAX-BAK-mediated mitochondrial outer membrane permeabilization [ 31 , 32 ]. We have also previously shown that Sfb induced a sustained and progressive increase in CHOP expression by Western blot analysis [ 13 ]. Finally, it is worth mentioning another example of an upregulated gene that is ERN1, which encodes the transmembrane kinase IRE1, which functions as a general sensor of unfolded proteins during the activation of the unfolded protein response (UPR) upon ER stress [ 33 ]. Additional DEGs showing upregulation corresponded to RAB42 and EMMOD2 related to GTP-binding or GTPase activation signaling transduction or FBXO9 and TRIM56 related to a ubiquitin ligase function. On the other hand, individual genes showing both lowest negative log 2 (FC) and p-value corresponded to chaperones, and co-chaperone proteins that alleviate the tendency of pre-existing proteins from aggregation or help the folding of nascent ones [ 34 ]. Thus, HSPA1A,HSPA8,DNAJA1 and HSPA1B are among the genes showing the most significant downregulation (Figure 2A). This result is, however, paradoxical since the ER stress caused by Sfb leads to the activation of UPR, thus, requiring further clarification. Similar results to these found for HepG2 cells were also uncovered for SNU423 cells when treated with Sfb (Figure 2B), indicating the robustness of the cellular responses to this drug. Next, we performed enrichment analyses of differentially expressed genes using the REACTOME database to identify the biological significance of the genes affected by the Sfb treatment common to both liver cancer cell lines of this study. Figure 3A summarizes the 20 over-represented REACTOME terms which are significantly associated with the upregulated genes. This analysis suggested that the upregulated DEGs were most significantly categorized into functional groups related to circadian regulation, group of pathways related to phosphorylated eIF2-alpha translation factor, FOXO-mediated transcription, metabolism of diverse lipids, events of membrane trafficking, its regulation by RAB GTPases and UPR, transport of amino acids and inorganic cations, biogenesis of mitochondria, apoptosis mediated by BH3-only proteins; strikingly, the activation of the signaling by receptor tyrosine kinases, which is the route that Sfb mainly inhibits, is among the categories that are also upregulated; this might be interpreted as a cellular response to the Sfb inhibition in order to regain cellular homeostasis. Figure 3B shows the 20 over-represented REACTOME terms which are downregulated by Sfb. These mainly include the routes of synthesis of cholesterol and steroid, DNA replication, metabolism of Cancers 2022,14, 1204 8 of 18 RNA, and cell cycle. Other processes such as translation and formation of the mitochondrial respiratory complex I were also downregulated; again, as a possible adaptation mechanism, Sfb upregulates the transcriptional activation of mitochondrial biogenesis (Figure 3A). Thus, we conclude that Sfb both induces the activation of the expression of genes belonging to pathways involved in signal transduction, that regulate cell metabolism, belong to the integrated stress response and trigger apoptosis and negatively affects the expression of genes mainly involved in cell growth, cell cycle control, and proliferation. Cancers 2022, 14, x 8 of 18 Figure 2. Differentially expressed genes in HepG2 (A) and SNU423 (B) cells upon Sorafenib treatment (10 µM, 12 h) detected using RNA-seq. The Volcano plot illustrates selected transcripts that show significant misregulation. The p-values are plotted in a −log10 format as y-values and fold change of the mRNAs are plotted in a log2 format as x-values. Representative transcripts that showed significant differential expression are indicated (blue). The rest of genes are shown in grey; density of genes is represented by a darker grey intensity. Next, we performed enrichment analyses of differentially expressed genes using the REACTOME database to identify the biological significance of the genes affected by the Sfb treatment common to both liver cancer cell lines of this study. Figure 3A summarizes the 20 over-represented REACTOME terms which are significantly associated with the upregulated genes. This analysis suggested that the upregulated DEGs were most significantly categorized into functional groups related to circadian regulation, group of pathways related to phosphorylated eIF2-alpha translation factor, FOXO-mediated transcription, metabolism of diverse lipids, events of membrane trafficking, its regulation by RAB GTPases and UPR, transport of amino acids and inorganic cations, biogenesis of mitochondria, apoptosis mediated by BH3-only proteins; strikingly, the activation of the signaling by receptor tyrosine kinases, which is the route that Sfb mainly inhibits, is among the categories that are also upregulated; this might be interpreted as a cellular response to the Sfb inhibition in order to regain cellular homeostasis. Figure 3B shows the 20 overrepresented REACTOME terms which are downregulated by Sfb. These mainly include the routes of synthesis of cholesterol and steroid, DNA replication, metabolism of RNA, and cell cycle. Other processes such as translation and formation of the mitochondrial respiratory complex I were also downregulated; again, as a possible adaptation mechanism, Sfb upregulates the transcriptional activation of mitochondrial biogenesis (Figure 3A). Thus, we conclude that Sfb both induces the activation of the expression of genes belonging to pathways involved in signal transduction, that regulate cell metabolism, belong to the integrated stress response and trigger apoptosis and negatively affects the expression of genes mainly involved in cell growth, cell cycle control, and proliferation. Figure 2. Differentially expressed genes in HepG2 ( A ) and SNU423 ( B ) cells upon Sorafenib treatment (10 µ M, 12 h) detected using RNA-seq. The Volcano plot illustrates selected transcripts that show significant misregulation. The p-values are plotted in a − log 10 format as y-values and fold change of the mRNAs are plotted in a log 2 format as x-values. Representative transcripts that showed significant differential expression are indicated (blue). The rest of genes are shown in grey; density of genes is represented by a darker grey intensity. 3.2. Sorafenib Inhibits Cell Growth Several studies have clearly shown that treatment with Sfb reduced cell viability and promoted cell death in HCC cell lines (e.g., [ 13 ] and references therein). Consistently, Figure 4A shows a substantial decrease in cell proliferation upon the treatment of Sfb in our growing conditions, which is compatible with the cytotoxic activity of this drug. Interestingly, our RNA-Seq analysis found, among the different functional groups of downregulated DEGs, those related to cell growth, ribosome biogenesis, translation, cotranslation protein targeting to membrane, cell cycle DNA replication and repair, cell cycle checkpoints, and sterol biosynthesis. All showed high negative values of normalized enrichment score (NES) in the GSEA analysis, thus reflecting a significant enrichment of these pathways in the list of predominantly inhibited genes upon Sfb treatment in HepG2 cells (Table S3 and Figure 3A). Cancers 2022,14, 1204 9 of 18 Cancers 2022, 14, x 9 of 18 Figure 3. Overrepresentation analysis (ORA) of commonly upregulated (A) and downregulated (B) genes in both HepG2 and SNU423 cells treated with Sorafenib (10 µM, 12 h) using the REACTOME database. The 20 categories with better p-values as well as a few selected ones (e.g., signaling by receptor tyrosine kinase, translation, complex I biogenesis) are shown. The categories of related biological processes are labeled in the same color. The number of genes with p-value < 0.001 in the different categories (n) is depicted. The p-values were −log10 transformed. 3.2. Sorafenib Inhibits Cell Growth Several studies have clearly shown that treatment with Sfb reduced cell viability and promoted cell death in HCC cell lines (e.g., [13] and references therein). Consistently, Figure 4A shows a substantial decrease in cell proliferation upon the treatment of Sfb in our growing conditions, which is compatible with the cytotoxic activity of this drug. Interestingly, our RNA-Seq analysis found, among the different functional groups of downregulated DEGs, those related to cell growth, ribosome biogenesis, translation, co-translation protein targeting to membrane, cell cycle DNA replication and repair, cell cycle checkpoints, and sterol biosynthesis. All showed high negative values of normalized enrichment score (NES) in the GSEA analysis, thus reflecting a significant enrichment of these pathways in the list of predominantly inhibited genes upon Sfb treatment in HepG2 cells (Table S3 and Figure 3A). Figure 4. Effects of Sorafenib on cell proliferation and cell cycle. (A) BrdU incorporation assay on HepG2 cells treated or not with 10 µM Sorafenib. Results are expressed as the mean ± S.D. of three Figure 3. Overrepresentation analysis (ORA) of commonly upregulated ( A ) and downregulated (B) genes in both HepG2 and SNU423 cells treated with Sorafenib (10 µ M, 12 h) using the REACTOME database. The 20 categories with better p-values as well as a few selected ones (e.g., signaling by receptor tyrosine kinase, translation, complex I biogenesis) are shown. The categories of related biological processes are labeled in the same color. The number of genes with p-value < 0.001 in the different categories (n) is depicted. The p-values were −log10 transformed. Cancers 2022, 14, x 9 of 18 Figure 3. Overrepresentation analysis (ORA) of commonly upregulated (A) and downregulated (B) genes in both HepG2 and SNU423 cells treated with Sorafenib (10 µM, 12 h) using the REACTOME database. The 20 categories with better p-values as well as a few selected ones (e.g., signaling by receptor tyrosine kinase, translation, complex I biogenesis) are shown. The categories of related biological processes are labeled in the same color. The number of genes with p-value < 0.001 in the different categories (n) is depicted. The p-values were −log10 transformed. 3.2. Sorafenib Inhibits Cell Growth Several studies have clearly shown that treatment with Sfb reduced cell viability and promoted cell death in HCC cell lines (e.g., [13] and references therein). Consistently, Figure 4A shows a substantial decrease in cell proliferation upon the treatment of Sfb in our growing conditions, which is compatible with the cytotoxic activity of this drug. Interestingly, our RNA-Seq analysis found, among the different functional groups of downregulated DEGs, those related to cell growth, ribosome biogenesis, translation, co-translation protein targeting to membrane, cell cycle DNA replication and repair, cell cycle checkpoints, and sterol biosynthesis. All showed high negative values of normalized enrichment score (NES) in the GSEA analysis, thus reflecting a significant enrichment of these pathways in the list of predominantly inhibited genes upon Sfb treatment in HepG2 cells (Table S3 and Figure 3A). Figure 4. Effects of Sorafenib on cell proliferation and cell cycle. (A) BrdU incorporation assay on HepG2 cells treated or not with 10 µM Sorafenib. Results are expressed as the mean ± S.D. of three Figure 4. Effects of Sorafenib on cell proliferation and cell cycle. ( A ) BrdU incorporation assay on HepG2 cells treated or not with 10 µ M Sorafenib. Results are expressed as the mean ± S.D. of three independent experiments. Statistical significance was analyzed using Student’s test (**** p< 0.0001). ( B ) Cell cycle profiles obtained by FACS of HepG2 cells treated or not with 10 µ M Sorafenib. For measuring DNA content, cells were stained with propidium iodide. G1, S and G2/M phases are indicated. The percentage of cells at the different phases are shown in the histogram. Results are expressed as the mean ± S.D. values of three independent experiments. Statistical significance was analyzed using Student’s test (* p< 0.05; ** p< 0.01). Regarding the cell cycle and DNA replication categories, we found that the levels of different E2F isoforms and of the kinase CDK4 were significantly and mildly downregulated, respectively. In agreement with a reduction in the function of different E2F isoforms, the most severe one of E2F2, Sfb causes a clear delay at the S and G2 phases of the cell cycle in HepG2 cells (Figure 4B). It should be emphasized that a selective CDK4/6 inhibitor, Pal- Cancers 2022,14, 1204 16 of 18 13. Rodríguez-Hernández, M.A.; González, R.; de la Rosa, A.J.; Gallego, P.; Ordoñez, R.; Navarro-Villarán, E.; Contreras, L.; Rodríguez-Arribas, M.; González-Gallego, J.; Álamo-Martínez, J.M.; et al. Molecular characterization of autophagic and apoptotic signaling induced by sorafenib in liver cancer cells. J. Cell. Physiol. 2018,234, 692–708. [CrossRef] [PubMed] 14. Fernando, J.; Sancho, P.; Fernández-Rodríguez, C.M.; Lledó, J.L.; Caja, L.; Campbell, J.S.; Fausto, N.; Fabregat, I. Sorafenib sensitizes hepatocellular carcinoma cells to physiological apoptotic stimuli. J. Cell. Physiol. 2012,227, 1319–1325. [CrossRef] 15. Gauthier, A.; Ho, M. Role of sorafenib in the treatment of advanced hepatocellular carcinoma: An update. Hepatol. Res. 2013 ,43, 147–154. [CrossRef] 16. Rodríguez-Hernández, M.A.; Chapresto-Garzón, R.; Cadenas, M.; Navarro-VillarÁn, E.; Negrete, M.; Gómez-Bravo, M.A.; Victor, V.M.; Padillo, F.J.; Muntané, J. Differential effectiveness of tyrosine kinase inhibitors in 2D/3D culture according to cell differentiation, p53 status and mitochondrial respiration in liver cancer cells. Cell Death Dis. 2020,11, 339. [CrossRef] 17. López-Terrada, D.; Gunaratne, P.H.; Adesina, A.M.; Pulliam, J.; Hoang, D.M.; Nguyen, Y.; Mistretta, T.A.; Margolin, J.; Finegold, M.J. Histologic subtypes of hepatoblastoma are characterized by differential canonical Wnt and Notch pathway activation in DLK+ precursors. Hum. Pathol. 2009,40, 783–794. [CrossRef] 18. Park, J.G.; Lee, J.H.; Kang, M.S.; Park, K.J.; Jeon, Y.M.; Lee, H.J.; Kwon, H.S.; Park, H.S.; Yeo, K.S.; Lee, K.U.; et al. Characterization of cell lines established from human hepatocellular carcinoma. Int. J. Cancer 1995,62, 276–282. [CrossRef] 19. Subramanian, A.; Tamayo, P.; Mootha, V.K.; Mukherjee, S.; Ebert, B.L.; Gillette, M.A.; Paulovich, A.; Pomeroy, S.L.; Golub, T.R.; Lander, E.S.; et al. Gene set enrichment analysis: A knowledge-based approach for interpreting genome-wide expression profiles. Proc. Natl. Acad. Sci. USA 2005,102, 15545–15550. [CrossRef] 20. Mootha, V.K.; Lindgren, C.M.; Eriksson, K.F.; Subramanian, A.; Sihag, S.; Lehar, J.; Puigserver, P.; Carlsson, E.; Ridderstrale, M.; Laurila, E.; et al. PGC-1alpha-responsive genes involved in oxidative phosphorylation are coordinately downregulated in human diabetes. Nat. Genet. 2003,34, 267–273. [CrossRef] 21. Jassal, B.; Matthews, L.; Viteri, G.; Gong, C.; Lorente, P.; Fabregat, A.; Sidiropoulos, K.; Cook, J.; Gillespie, M.; Haw, R.; et al. The reactome pathway knowledgebase. Nucleic Acids Res. 2020,48, D498–D503. [CrossRef] [PubMed] 22. Morine, M.J.; McMonagle, J.; Toomey, S.; Reynolds, C.M.; Moloney, A.P.; Gormley, I.C.; Gaora, P.O.; Roche, H.M. Bi-directional gene set enrichment and canonical correlation analysis identify key diet-sensitive pathways and biomarkers of metabolic syndrome. BMC Bioinform. 2010,11, 499. [CrossRef] [PubMed] 23. Micallef, L.; Rodgers, P. eulerAPE: Drawing area-proportional 3-Venn diagrams using ellipses. PLoS ONE 2014,9, e101717. 24. Molina-Ruiz, F.J.; González, R.; Rodríguez-Hernández, M.A.; Navarro-Villarán, E.; Padillo, F.J.; Muntané, J. Antitumoral activity of Sorafenib in hepatocellular carcinoma: Effects on cell survival and death pathways, cell metabolism reprogramming, and nitrosative and oxidative stress. Crit. Rev. Oncog. 2016,21, 413–432. [CrossRef] 25. Brüning, A.; Matsingou, C.; Brem, G.J.; Rahmeh, M.; Mylonas, I. Inhibin beta E is upregulated by drug-induced endoplasmic reticulum stress as a transcriptional target gene of ATF4. Toxicol. Appl. Pharmacol. 2012,264, 300–304. [CrossRef] 26. Phillips, D.J.; Woodruff, T.K. Inhibin: Actions and signalling. Growth Factors 2004,22, 13–18. [CrossRef] 27. Bernard, D.J.; Chapman, S.C.; Woodruff, T.K. Mechanisms of inhibin signal transduction. Recent Prog. Horm. Res. 2001 ,56, 417–450. [CrossRef] 28. Chen, H.F.; Chuang, H.C.; Tan, T.H. Regulation of Dual-Specificity Phosphatase (DUSP) ubiquitination and protein stability. Int. J. Mol. Sci. 2019,20, 2668. [CrossRef] 29. Lu, D.; Chen, J.; Hai, T. The regulation of ATF3 gene expression by mitogen-activated protein kinases. Biochem. J. 2007 ,401, 559–567. [CrossRef] 30. Cai, D.; Huang, E.; Luo, B.; Yang, Y.; Zhang, F.; Liu, C.; Lin, Z.; Xie, W.B.; Wang, H. Nupr1/Chop signal axis is involved in mitochondrion-related endothelial cell apoptosis induced by methamphetamine. Cell Death Dis. 2016,7, e2161. [CrossRef] 31. Wolfgang, C.D.; Chen, B.P.; Martindale, J.L.; Holbrook, N.J.; Hai, T. gadd153/Chop10, a potential target gene of the transcriptional repressor ATF3. Mol. Cell. Biol. 1997,17, 6700–6707. [CrossRef] 32. Verfaillie, T.; Salazar, M.; Velasco, G.; Agostinis, P. Linking ER stress to sutophagy: Potential implications for cancer therapy. Int. J. Cell Biol. 2010,2010, 930509. [CrossRef] [PubMed] 33. Chen, Y.; Brandizzi, F. IRE1: ER stress sensor and cell fate executor. Trends Cell Biol. 2013,23, 547–555. [CrossRef] [PubMed] 34. Kim, Y.E.; Hipp, M.S.; Bracher, A.; Hayer-Hartl, M.; Hartl, F.U. Molecular chaperone functions in protein folding and proteostasis. Annu. Rev. Biochem. 2013,82, 323–355. [CrossRef] [PubMed] 35. Bollard, J.; Miguela, V.; Ruiz de Galarreta, M.; Venkatesh, A.; Bian, C.B.; Roberto, M.P.; Tovar, V.; Sia, D.; Molina-Sánchez, P.; Nguyen, C.B.; et al. Palbociclib (PD-0332991), a selective CDK4/6 inhibitor, restricts tumour growth in preclinical models of hepatocellular carcinoma. Gut 2017,66, 1286–1296. [CrossRef] [PubMed] 36. Smits, V.A.J.; Cabrera, E.; Freire, R.; Gillespie, D.A. Claspin-checkpoint adaptor and DNA replication factor. FEBS J. 2019 ,286, 441–455. [CrossRef] [PubMed] 37. Maréchal, A.; Zou, L. DNA damage sensing by the ATM and ATR kinases. Cold Spring Harb. Perspect. Biol. 2013 ,5, a012716. [CrossRef] 38. Krokan, H.E.; Bjørås, M. Base excision repair. Cold Spring Harb. Perspect. Biol. 2013,5, a012583. [CrossRef] 39. Richly, H.; Schultheis, B.; Adamietz, I.A.; Kupsch, P.; Grubert, M.; Hilger, R.A.; Ludwig, M.; Brendel, E.; Christensen, O.; Strumberg, D. Combination of sorafenib and doxorubicin in patients with advanced hepatocellular carcinoma: Results from a phase I extension trial. Eur. J. Cancer 2009,45, 579–587. [CrossRef] Cancers 2022,14, 1204 17 of 18 40. Yang, X.D.; Kong, F.E.; Qi, L.; Lin, J.X.; Yan, Q.; Loong, J.H.C.; Xi, S.Y.; Zhao, Y.; Zhang, Y.; Yuan, Y.F.; et al. PARP inhibitor Olaparib overcomes Sorafenib resistance through reshaping the pluripotent transcriptome in hepatocellular carcinoma. Mol. Cancer 2021,20, 20. [CrossRef] 41. Yu, W.; Gu, K.; Yu, Z.; Yuan, D.; He, M.; Ma, N.; Lai, S.; Zhao, J.; Ren, Z.; Zhang, X.; et al. Sorafenib potentiates irradiation effect in hepatocellular carcinoma in vitro and in vivo. Cancer Lett. 2013,329, 109–117. [CrossRef] [PubMed] 42. Wang, C.; Wang, H.; Lieftink, C.; du Chatinier, A.; Gao, D.; Jin, G.; Jin, H.; Beijersbergen, R.L.; Qin, W.; Bernards, R. CDK12 inhibition mediates DNA damage and is synergistic with sorafenib treatment in hepatocellular carcinoma. Gut 2020 ,69, 727–736. [CrossRef] [PubMed] 43. In der Stroth, L.; Tharehalli, U.; Günes, C.; Lechel, A. Telomeres and Telomerase in the Development of Liver Cancer. Cancers 2020,12, 2048. [CrossRef] 44. Zhou, T.; Lv, X.; Guo, X.; Ruan, B.; Liu, D.; Ding, R.; Gao, Y.; Ding, J.; Dou, K.; Chen, Y. RACK1 modulates apoptosis induced by sorafenib in HCC cells by interfering with the IRE1/XBP1 axis. Oncol. Rep. 2015,33, 3006–3014. [CrossRef] [PubMed] 45. Shi, Y.H.; Ding, Z.B.; Zhou, J.; Hui, B.; Shi, G.M.; Ke, A.W.; Wang, X.Y.; Dai, Z.; Peng, Y.F.; Gu, C.Y.; et al. Targeting autophagy enhances sorafenib lethality for hepatocellular carcinoma via ER stress-related apoptosis. Autophagy 2011 ,7, 1159–1172. [CrossRef] [PubMed] 46. Fang, C.; Xie, H.; Zhao, J.; Wang, W.; Hou, H.; Zhang, B.; Zhou, D.; Geng, X. eIF4E-eIF4G complex inhibition synergistically enhances the effect of sorafenib in hepatocellular carcinoma. Anticancer. Drugs 2021,32, 822–828. [CrossRef] 47. Honma, Y.; Shimizu, S.; Takehara, T.; Harada, M. Sorafenib enhances proteasome inhibitor-induced cell death via inactivation of Akt and stress-activated protein kinases. J. Gastroenterol. 2014,49, 517–526. [CrossRef] 48. Jiang, C.; Xu, R.; Li, X.X.; Zhou, Y.F.; Xu, X.Y.; Yang, Y.; Wang, H.Y.; Zheng, X.F.S. Sorafenib and Carfilzomib synergistically inhibit the proliferation, survival, and metastasis of hepatocellular carcinoma. Mol. Cancer Ther. 2018,17, 2610–2621. [CrossRef] 49. Chen, K.F.; Yu, H.C.; Liu, T.H.; Lee, S.S.; Chen, P.J.; Cheng, A.L. Synergistic interactions between sorafenib and bortezomib in hepatocellular carcinoma involve PP2A-dependent Akt inactivation. J. Hepatol. 2010,52, 88–95. [CrossRef] 50. Rodríguez-Hernández, M.A.; de la Cruz-Ojeda, P.; López-Grueso, M.J.; Navarro-Villarán, E.; Requejo-Aguilar, R.; Castejón- Vega, B.; Negrete, M.; Gallego, P.; Vega-Ochoa, A.; Victor, V.M.; et al. Integrated molecular signaling involving mitochondrial dysfunction and alteration of cell metabolism induced by tyrosine kinase inhibitors in cancer. Redox Biol. 2020 ,36, 101510. [CrossRef] 51. Fiume, L.; Manerba, M.; Vettraino, M.; Di Stefano, G. Effect of sorafenib on the energy metabolism of hepatocellular carcinoma cells. Eur. J. Pharmacol. 2011,670, 39–43. [CrossRef] 52. Tesori, V.; Piscaglia, A.C.; Samengo, D.; Barba, M.; Bernardini, C.; Scatena, R.; Pontoglio, A.; Castellini, L.; Spelbrink, J.N.; Maulucci, G.; et al. The multikinase inhibitor Sorafenib enhances glycolysis and synergizes with glycolysis blockade for cancer cell killing. Sci. Rep. 2015,5, 9149. [CrossRef] [PubMed] 53. Zhang, C.; Liu, Z.; Bunker, E.; Ramírez, A.; Lee, S.; Peng, Y.; Tan, A.C.; Eckhardt, S.G.; Chapnick, D.A.; Liu, X. Sorafenib targets the mitochondrial electron transport chain complexes and ATP synthase to activate the PINK1-Parkin pathway and modulate cellular drug response. J. Biol. Chem. 2017,292, 15105–15120. [CrossRef] [PubMed] 54. Rodríguez-Hernández, M.A.; de la Cruz-Ojeda, P.; Gallego, P.; Navarro-Villarán, E.; Stankova, P.; del Campo, J.A.; Kucera, O.; Elkalaf, M.; Maseko, T.E.; Cervinkova, Z.; et al. Dose-dependent regulation of mitochondrial function and cell death pathway by sorafenib in liver cancer cells. Biochem. Pharmacol. 2020,176, 113902. [CrossRef] 55. Paech, F.; Mingard, C.; Grunig, D.; Abegg, V.F.; Bouitbir, J.; Krähenbühl, S. Mechanisms of mitochondrial toxicity of the kinase inhibitors ponatinib, regorafenib and sorafenib in human hepatic HepG2 cells. Toxicology 2018,395, 34–44. [CrossRef] 56. Li, Y.; Xia, J.; Shao, F.; Zhou, Y.; Yu, J.; Wu, H.; Du, J.; Ren, X. Sorafenib induces mitochondrial dysfunction and exhibits synergistic effect with cysteine depletion by promoting HCC cells ferroptosis. Biochem. Biophys. Res. Commun. 2021 ,534, 877–884. [CrossRef] 57. Fan, G.; Wei, X.; Xu, X. Is the era of sorafenib over? A review of the literature. Ther. Adv. Med. Oncol. 2020 ,12, 1758835920927602. [CrossRef] 58. Liu, G.; Kuang, S.; Cao, R.; Wang, J.; Peng, Q.; Sun, C. Sorafenib kills liver cancer cells by disrupting SCD1-mediated synthesis of monounsaturated fatty acids via the ATP-AMPK-mTOR-SREBP1 signaling pathway. FASEB J. 2019 ,33, 10089–10103. [CrossRef] [PubMed] 59. Yin, F.; Feng, F.; Wang, L.; Wang, X.; Li, Z.; Cao, Y. SREBP-1 inhibitor Betulin enhances the antitumor effect of Sorafenib on hepatocellular carcinoma via restricting cellular glycolytic activity. Cell Death Dis. 2019,10, 672. [CrossRef] 60. Eun, H.S.; Cho, S.Y.; Lee, B.S.; Seong, I.O.; Kim, K.H. Profiling cytochrome P450 family 4 gene expression in human hepatocellular carcinoma. Mol. Med. Rep. 2018,18, 4865–4876. [CrossRef] 61. Dong, X.C.; Chowdhury, K.; Huang, M.; Kim, H.G. Signal transduction and molecular regulation in fatty liver disease. Antioxid. Redox Signal. 2021,35, 689–717. [CrossRef] [PubMed] 62. Che, L.; Chi, W.; Qiao, Y.; Zhang, J.; Song, X.; Liu, Y.; Li, L.; Jia, J.; Pilo, M.G.; Wang, J.; et al. Cholesterol biosynthesis supports the growth of hepatocarcinoma lesions depleted of fatty acid synthase in mice and humans. Gut 2020,69, 177–186. [CrossRef] 63. Di Paolo, G.; De Camilli, P. Phosphoinositides in cell regulation and membrane dynamics. Nature 2006 ,443, 651–657. [CrossRef] [PubMed] 64. Reppert, S.M.; Weaver, D.R. Coordination of circadian timing in mammals. Nature 2002,418, 935–941. [CrossRef] 65. Bass, J.; Takahashi, J.S. Circadian integration of metabolism and energetics. Science 2010,330, 1349–1354. [CrossRef] [PubMed] Cancers 2022,14, 1204 18 of 18 66. Reinke, H.; Asher, G. Circadian clock control of liver metabolic functions. Gastroenterology 2016 ,150, 574–580. [CrossRef] [PubMed] 67. Jacobi, D.; Liu, S.; Burkewitz, K.; Kory, N.; Knudsen, N.H.; Alexander, R.K.; Unluturk, U.; Li, X.; Kong, X.; Hyde, A.L.; et al. Hepatic Bmal1 regulates rhythmic mitochondrial dynamics and promotes metabolic fitness. Cell Metab. 2015 ,22, 709–720. [CrossRef] 68. Sahar, S.; Sassone-Corsi, P. Metabolism and cancer: The circadian clock connection. Nat. Rev. Cancer 2009 ,9, 886–896. [CrossRef] 69. Fuchs, B.C.; Fujii, T.; Dorfman, J.D.; Goodwin, J.M.; Zhu, A.X.; Lanuti, M.; Tanabe, K.K. Epithelial-to-mesenchymal transition and integrin-linked kinase mediate sensitivity to epidermal growth factor receptor inhibition in human hepatoma cells. Cancer Res. 2008,68, 2391–2399. [CrossRef] 70. Hirschfield, H.; Bian, C.B.; Higashi, T.; Nakagawa, S.; Zeleke, T.Z.; Nair, V.D.; Fuchs, B.C.; Hoshida, Y. In Vitro modeling of hepatocellular carcinoma molecular subtypes for anti-cancer drug assessment. Exp. Mol. Med. 2018,50, e419. [CrossRef] 71. Clayton, N.S.; Ridley, A.J. Targeting Rho GTPase signaling networks in cancer. Front. Cell Dev. Biol. 2020,8, 222. [CrossRef] 72. Bruix, J.; Gores, G.J.; Mazzaferro, V. Hepatocellular carcinoma: Clinical frontiers and perspectives. Gut 2014 ,63, 844–855. [CrossRef] [PubMed] 73. Paech, F.; Bouitbir, J.; Krähenbaühl, S. Hepatocellular toxicity associated with tyrosine kinase inhibitors: Mitochondrial damage and inhibition of glycolysis. Front. Pharmacol. 2017,8, 367. [CrossRef] [PubMed] 74. Asghar, U.; Meyer, T. Are there opportunities for chemotherapy in the treatment of hepatocellular cancer? J. Hepatol. 2012 ,56, 686–695. [CrossRef] [PubMed] 75. Regan-Fendt, K.; Li, D.; Reyes, R.; Yu, L.; Wani, N.A.; Hu, P.; Jacob, S.T.; Ghoshal, K.; Payne, P.R.O.; Motiwala, T. Transcriptomicsbased drug repurposing approach identifies novel drugs against Sorafenib-resistant hepatocellular carcinoma. Cancers 2020 , 12, 2730. [CrossRef] [PubMed] 76. Ye, L.; Mayerle, J.; Ziesch, A.; Reiter, F.P.; Gerbes, A.L.; De Toni, E.N. The PI3K inhibitor copanlisib synergizes with sorafenib to induce cell death in hepatocellular carcinoma. Cell. Death Discov. 2019,5, 86. [CrossRef] 77. Tang, Y.; Chen, J.; Li, J.; Zheng, Y.; Zhong, X.; Huang, S.; Chen, B.; Peng, B.; Zou, X.; Chen, X. Pristimerin synergistically sensitizes conditionally reprogrammed patient derived-primary hepatocellular carcinoma cells to sorafenib through endoplasmic reticulum stress and ROS generation by modulating Akt/FoxO1/p27(kip1) signaling pathway. Phytomedicine 2021,86, 153563. [CrossRef] 78. ¸Sirin, N.; Elmas, L.; Secme, M.; Dodurga, Y. Investigation of possible effects of apigenin, sorafenib and combined applications on apoptosis and cell cycle in hepatocellular cancer cells. Gene 2020,737, 144428. [CrossRef]