Transcriptomic signatures of progressive and regressive liver fibrosis and portal hypertension

Abstract

Persistent liver injury triggers a fibrogenic program that causes pathologic remodeling of the hepatic microenvironment (i.e., liver fibrosis) and portal hypertension. The dynamics of gene regulation during liver disease progression and early regression remain understudied. Here, we generated hepatic transcriptome profiles in two well-established liver disease models at peak fibrosis and during spontaneous regression after the removal of the inducing agents. We linked the dynamics of key disease readouts, such as portal pressure, collagen area, and transaminase levels, to differentially expressed genes, enabling the identification of transcriptomic signatures of progressive vs. regressive liver fibrosis and portal hypertension. These candidate biomarkers (e.g., Tcf4, Mmp7, Trem2, Spp1, Scube1, Islr) were validated in RNA sequencing datasets of patients with cirrhosis and portal hypertension, and those cured from hepatitis C infection. Finally, deconvolution identified major cell types and suggested an association of macrophage and portal hepatocyte signatures with portal hypertension and fibrosis area.

Keywords: Disease; Fibrosis; Integrative aspects of cell biology; Model organism; Transcriptomics.

Graphical abstract

Figure 1

Experimental design and key liver disease readouts (A) Two different murine fibrosis models were used to study their hepatic transcriptome during peak fibrosis and fibrosis regression. (B) Experimental timelines are shown. Portal pressure was measured, blood was sampled, and liver tissue was harvested at the respective endpoints. (C) Comparison of key liver disease readouts in the two models. Results of the Wilcoxon rank-sum test for the differences in the respective parameters after one (R1) and two (R2) weeks of regression to the peak (CIR) time point are indicated (n.s. = non-significant; ∗p < 0.05; ∗∗p < 0.01; n = 6 per group of each model). Data are represented as mean ± SEM. (D) Illustration of representative Picrosirius Red staining for quantification of collagen proportionate area are shown. Red = collagen, green = fast green tissue counterstain.

Figure 2

Transcriptomic characterization of genes associated with peak fibrosis and regressive fibrosis (A) PCA projection demonstrates a clear separation of HC and mostly of CIR, while the R1 and R2 groups cluster together. (B) Genes with 10% highest loading score for principal components 1 (PC1) and 2 (PC2). Unsupervised gene clustering denotes similar expression patterns. (C) Volcano plots for most differentially expressed genes between animals at peak fibrosis induction (CIR) and healthy control animals (HC) are shown separately for the CCl4 and TAA models. (D) Volcano plots for most differentially expressed genes between animals with regressive cirrhosis (combined R1+R2) vs. peak fibrosis induction (CIR) are shown separately for both models. Thresholds: log2FC > 1.5, p_adjusted <0.01. PC = principal component. See also Figures S2–S4, and Data S1.

Figure 3

Functional annotation of genes involved in CCl4 and TAA fibrosis progression and regression (A–C) Venn diagram representing overlaps in up- and downregulated genes (number and proportion of genes are shown) associated with fibrosis induction in the two fibrosis models. The functional analysis displays molecular signaling pathways for both CCl4/TAA models that are (B) up- or (C) downregulated in peak fibrosis induction (CIR). (D) Venn diagram displaying overlaps in genes (number and proportion of genes are shown) linked to fibrosis regression in the two models. (E and F) (E) Functional upregulation and (F) downregulation of molecular signaling pathways for both CCl4/TAA models during fibrosis regression. A hypergeometric test with MSigDB annotation was applied for the genes overlapped in regression. (G) The scores in genes perturbed by both conditions demonstrate a strong negative correlation (Pearson’s r is shown). The scores were scaled for visualization. (H) Signaling pathways from genes with the 20% highest score for CIR are shown. (I) Pathways from genes with the 20% highest regression score (CIR vs. R1+R2) are shown. See also Figures S5–S8, and Data S2.

Figure 4

Dynamic signature patterns following liver disease course using pseudotemporal analysis (A) Cluster 1: pathway analysis of genes upregulated at peak fibrosis followed by subsequent downregulated (vs. CIR) and approximation to baseline (HC) levels in regression. (B) Cluster 2: pathways analysis of genes downregulated at peak fibrosis followed by subsequent upregulation (vs. CIR) and approximation to baseline (HC) levels in regression. NES = normalized enrichment score. See also Figures S9–S11, and Data S3.

Figure 5

Linking gene signatures to key surrogate parameters of hallmark biological readouts (A) Gene modules, identified via unsupervised weighted co-expression network analysis (WGCNA), correlating with key biological readouts in the merged dataset. (B) Supervised clustering of top 40 markers associated with targeted minimum loss-based estimation (TMLE). (C) Prioritization via overlaps between TMLE markers for each biological readout and respective WGCNA module. PP = portal pressure. CPA = collagen proportionate area. See also Figures S12–S14, and Data S4.

Figure 6

Network analysis prioritizes biomarkers and shows their functional links to liver disease (A) Multilayer network constructed from selected transcriptional factors (purple, bottom layer), co-expressed biomarkers (blue, middle layer), and their protein-protein interactions (orange, top layer). (B–E) Biomarkers were scored to identify network cores with the feature’s strongest associations. Node shapes indicate gene dynamics between cirrhotic and regression groups. CIR = cirrhosis, CPA = collagen proportionate area, PP = portal pressure, R = fibrosis regression. See also Figure S15 and Data S5.

Figure 7

Transcriptomic overlaps and validation of the prioritized biomarkers in human liver disease datasets The comparisons were performed between (A) patients with histological cirrhosis (F4) and non-diseased; (B) patients with advanced fibrosis/cirrhosis (F3-F4) and non-diseased; (C) patients with portal hypertension and non-diseased; (D) patients with hepatitis C virus cured with direct-acting antiviral therapy and those prior treatment start, accordingly. The right panel indicates detection of the prioritized hub genes. PH = portal hypertension. CPA = collagen proportionate area. PP = portal pressure. DAA = direct-acting antiviral therapy. TAA = thioacetamide. CCl4 = carbon tetrachloride. See also Table S1.

Figure 8

Deconvolution analysis shows cell-specific signatures in mouse RNA-seq datasets and their relationship to PP and CPA (A) Cell type deconvolution of our study dataset using published single-cell datasets. The heatmap illustrates deconvoluted cell types with the most robust signatures. (B–E) We found a linear correlation of the respective cell type score to PP and CPA. Cases shown in bold are those in which we consider cell score contributing to parameter variance (R2 > 0.6). The color indicates the respective cell type’s low (gray) or high (blue for CCl4, orange for TAA) score value. HSC = hepatic stellate cells. PP = portal pressure. CPA = collagen proportionate area. p < 0.05; ∗∗p < 0.01; ∗∗∗p < 0.001. See also Figures S16–S19.