Liver cirrhosis represents the end stage of chronic liver disease characterized by extensive fibrosis and destruction of normal hepatic architecture.¹ It is a major cause of morbidity and mortality worldwide, resulting from various etiologies including chronic viral hepatitis, alcohol abuse, metabolic disorders, and autoimmune diseases.^2-4 The pathogenesis of cirrhosis is complex and driven by multiple cellular and molecular factors that ultimately lead to inflammation, hepatocellular damage, and accumulation of extracellular matrix proteins.^5-8

A growing body of evidence indicates that immunological mechanisms play a central role in the initiation and progression of liver fibrosis in cirrhosis.^7-9 Among various immune cell populations, CD8+ T cells have emerged as major contributors to persistent hepatic inflammation and injury.^10-12 Under normal conditions, CD8+ T cells provide protection against infections and tumors by recognizing and eliminating aberrant cells.^13,14 However, in chronic liver disease, dysregulated CD8+ T cell responses can paradoxically amplify hepatocellular damage and fibrogenesis through cytotoxic effects and inflammatory cytokine production.¹⁵ Despite their significance, the precise molecular pathways governing CD8+ T cell-mediated effects in cirrhosis pathogenesis remain poorly defined.

Programmed cell death protein 1 (PD-1) is an immune checkpoint receptor that is upregulated on exhausted CD8+ T cells during chronic inflammation.^16-18 Recent studies have showed that PD-1 may critically control CD8+ T-cell functions in cirrhosis.^19,20 PD-1 signaling is thought to disable T cell responses in cirrhosis, leading to failure in immune surveillance.²¹ Hence, the complex interplay between CD8+ T cells and regulatory genes like PD-1 needs deeper investigation in cirrhosis immunopathology.

In this study, we aimed to delineate the heterogeneity of CD8+ T cells and dissect their functional contributions in liver cirrhosis through integrated transcriptomic and single-cell sequencing approaches. By profiling gene expression patterns in purified CD8+ T cells and related immune cell populations, we aimed to uncover previously unrecognized cell subsets, activation states, and regulatory genes that drive cirrhosis progression. These findings may provide greater mechanistic insights into cirrhosis pathogenesis while identifying potential diagnostic biomarkers and therapeutic targets.

1. Microarray data collection and differential expression analysis

The raw microarray dataset GSE14323 was downloaded from the GEO database (https://www.ncbi.nlm.nih.gov/geo/). The dataset included gene expression profiles of 153 patients with liver cirrhosis and 40 healthy controls. The R package limma²² was utilized to identify differentially expressed genes (DEGs) between cirrhotic and healthy groups. Genes with fold change >1.5 and adjusted p-value <0.05 were defined as DEGs.

2. Immune cell infiltration analysis

The abundances of immune cells in cirrhosis and healthy groups were estimated by the R package MCPcounter,²³ which calculates enrichment scores for immune cell populations based on gene expression profiles. Wilcoxon rank-sum test was used to determine statistical differences in immune cell proportions between the two groups. p<0.05 was considered significant.

3. Weighted gene co-expression network analysis

Weighted gene co-expression network analysis (WGCNA) was performed using the WGCNA R package²⁴ to identify modules of highly correlated genes across the whole transcriptome dataset GSE14323. A weighted adjacency matrix was constructed by calculating biweight midcorrelation between gene expression profiles. The matrix was subsequently transformed into a topological overlap matrix and used for average linkage hierarchical clustering to identify gene modules. Module eigengenes were computed and correlated with immune cell abundances to determine immune-related modules. Hub genes in key modules were identified for downstream analysis.

4. Functional enrichment analysis

Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses on module hub genes were carried out using the R package clusterProfiler (4.0).²⁵ Over-represented terms with adjusted p-value <0.05 were considered statistically significant.

5. Protein-protein interaction network analysis

The online database Search Tool for the Retrieval of Interacting Genes (https://cn.string-db.org/) was utilized to construct protein-protein interaction (PPI) networks for module hub genes. The networks were visualized using Cytoscape software²⁶ and densely connected modules were identified from the networks using the MCODE plugin²⁷ in Cytoscape.

6. Consensus clustering analysis

The ConsensusClusterPlus package²⁸ in R was applied to perform consensus clustering on the gene expression dataset GSE63898 downloaded from the GEO database. This allowed determining the optimal number of stable clusters based on the relative change in area under the cumulative distribution function curve.

7. Calculating Cirrhosis Disease Scores

Initially, three gene sets related to cirrhosis, namely C0023890, C0023893, and C1623038, were collected from the DisGeNET²⁹ database (https://www.disgenet.org/home/). Genes that were upregulated in cirrhosis in the GSE14323 dataset were then selected from each set to characterize the severity of cirrhosis (Supplementary Table 1). Subsequently, the Gene Set Variation Analysis (GSVA) enrichment scores for these three gene sets (C0023890, C0023893, and C1623038) were calculated in the GSE14323 and GSE63898 datasets using the GSVA package in R.³⁰

8. Single-cell RNA sequencing data pre-processing

The single-cell RNA sequencing dataset SCP2154 was obtained from the Single Cell Portal (https://singlecell.broadinstitute.org/single_cell). Only T cells and natural killer (NK) cells from five cirrhotic and five healthy samples (Supplementary Table 2) were retained for further analysis. CD8+ T cells were identified by the expression of CD3, CD8 and T-cell receptor genes. Cells with less than 500 detected genes were considered low-quality and filtered out.

9. Single-cell RNA-seq data dimensionality reduction and clustering

The R package Seurat (version 2.3.4) was utilized for dimensionality reduction, clustering, differential expression analysis, and trajectory analysis of the filtered single-cell RNA-seq data. The linear dimension reduction technique principal component analysis was applied followed by non-linear dimension reduction t-distributed stochastic neighbor embedding to visualize single-cell heterogeneity. Unsupervised graph-based clustering using Leiden algorithm³¹ was performed to identify cell subpopulations.

10. Differential expression analysis of single-cell RNA-seq data

Wilcoxon rank-sum test implemented in the Seurat package was used to identify marker genes for each single-cell cluster that showed differential expressions compared to other clusters. Genes with adjusted p-value <0.05 and log2 fold change >0.25 were considered significantly differentially expressed.

11. Trajectory analysis

Trajectory analysis was performed in Scanpy (v1.7.2) using the sc.tl.paga function based on the Leiden clustering results to infer developmental relationships between cell subtypes. The trajectories were visualized using the sc.pl.paga function.

12. Slingshot cell differentiation trajectory analysis

The R package slingshot³² was utilized for cell trajectory analysis. The expression matrix of all CD8+ T cells was used as input counts, with the eight annotated subtypes of CD8+ T cells serving as clustering information to construct a SingleCellExperiment object. T naive cells were selected as the starting point of the trajectory for cell trajectory inference.

13. Gene expression density analysis

The Scanpy package³³ was utilized to calculate local gene expression density for each cell in the reduced dimensionality space using tl.embedding_density function. The densities were visualized by the pl.embedding_density function.

14. Identification of cluster-specific gene signatures

To identify gene signatures representing each cell subtype, the tl.rank_genes_groups function in Scanpy was used to rank genes based on their expression levels and specificities within each Leiden cluster. The top 100 marker genes in each cluster were defined as the signature gene sets for downstream analysis.

15. Gene set enrichment analysis

Using the signature gene sets for each cell subtype, enrichment scores were calculated for the cirrhotic samples from GSE63898 using the R package GSVA.³⁰ The scores were correlated with cirrhotic status to evaluate functional associations.

16. Functional GO network enrichment analysis

The online tool Metascape (http://metascape.org) was used for GO and KEGG pathway enrichment analysis of signature genes. GO networks were constructed in Metascape and densely connected modules were extracted using MCODE plugin in Cytoscape.

1. Differential analysis and WGCNA analysis revealed key modules and genes related to immunity in liver cirrhosis

Microarray data related to liver cirrhosis from the GEO database (GSE14323) were analyzed. Initially, 301 DEGs were identified between patients with liver cirrhosis and healthy individuals, including 240 upregulated and 61 downregulated DEGs in liver cirrhosis (Fig. 1A, Supplementary Table 3). Further analysis showed that the most significantly upregulated DEGs in patients with liver cirrhosis included CXCL9, CXCL10, and those related to HLA class II histocompatibility antigen and other cytokines (Fig. 1B). Immuno infiltration analysis revealed a significant increase in the proportions of immune cells such as activated CD8+ T cells, activated CD4+ T cells, NK cells, natural killer T (NKT) cells, and myeloid-derived suppressor cells in the liver cirrhosis group (Fig. 1C). Additionally, based on WGCNA, a co-expression network analysis of all genes in GSE14323 was performed, resulting in the identification of 13 co-expression gene modules. To further investigate immune-related genes in liver cirrhosis, the correlation between the infiltration level of each immune cell and the co-expression modules was analyzed (Fig. 1D). The results showed that the lightgreen and darkorange modules exhibited the highest correlation with immune cell infiltration, with the lightgreen module showing the strongest correlation with activated CD8+ T-cell infiltration (r=0.84, p=9e-17). The lightgreen module included 843 co-expressed genes, while the darkorange module included 1429 co-expressed genes. Furthermore, combining the co-expressed genes with the previously analyzed DEGs between patients with liver cirrhosis and healthy individuals in GSE14323, 110 and 38 key module-related genes were identified in the lightgreen and darkorange modules, respectively (Fig. 1E and F).

Figure 1. Differential expression and weighted gene co-expression network analysis (WGCNA) identified key modules and genes related to immune infiltration in liver cirrhosis. (A) Volcano plot showing 301 differentially expressed genes (DEGs) between patients with cirrhosis and healthy controls. (B) Top 20 upregulated DEGs including HLA and chemokine genes. (C) Immune cells with significantly increased proportions in patients with cirrhosis compared with healthy controls. (D) Module-trait correlation analysis between module eigengenes and infiltration levels of immune cells. (E) Venn diagram showing overlap between the aforementioned DEGs and co-expressed genes in the lightgreen module. The blue circle represents DEGs, and the red circle represents the module co-expressed genes. (F) Venn diagram showing overlap between the aforementioned DEGs and co-expressed genes in the darkorange module. The blue circle represents DEGs, and the red circle represents the module co-expressed genes.

2. Functional analysis and network analysis revealed 19 key hub genes related to liver cirrhosis

To further analyze the functions of the identified key genes, we performed a GO-based network analysis, and a KEGG-based functional enrichment analysis on the 148 module-related genes from the lightgreen and darkorange modules. GO network analysis indicated that these genes were associated with functions such as MHC class II protein complex, cell activation, positive regulation of cell migration, and positive regulation of response to external stimulus (Fig. 2A). This suggests that these key genes may influence immune cell activation, migration and infiltration. KEGG functional analysis revealed that these genes were mainly associated with cell adhesion molecules, the PI3K-Akt signaling pathway, antigen processing and presentation, and Th17 cell differentiation, indicating their potential high expression in the cell membrane or extracellular matrix (Fig. 2B). Furthermore, a PPI network analysis was conducted on the 148 genes, which revealed extensive interactions and may indicate a high functional correlation among these genes. Additionally, using the MCODE algorithm, a subnetwork with the highest MCODE score was identified from the constructed PPI network (Fig. 2C). This subnetwork comprised 19 genes, with HLA-DRA, PTPRC (CD45), CD74, and GZMA being the most crucial genes in the module. Due to the high consistency and close functional relationship among these genes, we defined these 19 genes as key genes influencing immune cell infiltration in patients with liver cirrhosis. Except for the granule enzyme genes and cytokine genes, the majority of the 19 genes belonged to cell surface molecules, such as CD2, CD74, and PTPRC. This might suggest a possible association between these genes and immune cell differentiation and activation, as well as their direct involvement in liver cirrhosis.

Figure 2. Functional enrichment and network analysis of 148 module-related genes revealed 19 key hub genes related to liver cirrhosis. (A) Gene Ontology network analysis of the 148 module-related genes showing association with immune cell activation and migration. (B) Kyoto Encyclopedia of Genes and Genomes functional analysis of the 148 module-related genes indicating enrichment in cell adhesion, antigen presentation and cytokine pathways. (C) Protein-protein interaction network analysis and identification of 19 key hub genes.

3. Impact of the 19 identified key genes on CD8+ T-cell infiltration and cirrhosis severity

To further verify the influence of these 19 key genes on immune cells in liver cirrhosis, a consensus clustering analysis was conducted on 148 cirrhotic patients in the GSE63898 dataset based on the expression levels of these genes. These 19 genes divided the 148 patients into two main groups: group C1 comprising 70 patients, and group C2 with 78 patients (Fig. 3A). To evaluate the differences in the severity of cirrhosis between the two groups, we used GSVA scores based on three self-constructed cirrhosis-related gene sets to characterize the severity in cirrhotic patients. Initially, the efficacy of C0023890, C0023893 and C1623038, the three Cirrhosis Disease Scores, was validated in the GSE14323 dataset by their significant elevation in cirrhotic patients (Fig. 3B). Then, the Cirrhosis Disease Scores were applied to assess the severity in groups C1 and C2 in the GSE63898 dataset, revealing higher scores in group C1, suggesting a more severe cirrhosis in this group (Fig. 3C). Differential analysis between groups C1 and C2 revealed 440 significant DEGs, including genes related to HLA class II histocompatibility antigens, CXCL11, CXCL13, cytokine-related genes, and proteasome-related genes (Fig. 3D, Supplementary Table 4). Further analysis of immune cell infiltration levels showed significant differences between the two groups in activated CD8+ T cells, activated CD4+ T cells, and NK cells (Fig. 3E). Functional analysis of the DEGs indicated their association with T-cell activation, interferon signaling, immune response activation, antigen processing and presentation, and NF-kB signaling regulation (Fig. 3F). This may suggest significant functional differences related to these 19 key genes in cirrhotic patients. Network analysis and MCODE subnetwork mining of the DEGs identified nine key hub genes: CCL19, GABBR1, CXCL10, CCR7, CXCL11, CCR5, CCL5, CXCL13, and CXCL9 (Fig. 3G). Considering the functions of immune cells and the roles of the identified key genes in immune cell infiltration, activated CD8+ T cells may be a critical cell type affecting the progression of cirrhosis. Moreover, previous studies have demonstrated that these nine key genes are highly associated with the function of CD8+ T cells, such as CXCL13.

Figure 3. Key hub genes influencing the infiltration of CD8+ T cells in liver cirrhosis. (A) Consensus clustering based on hub gene expression segregates patients into two subgroups, group C1 and group C2. (B) Validation of Cirrhosis Disease Scores based on key cirrhosis pathogenesis gene sets in the GSE14323 dataset. (C) Validation of Cirrhosis Disease Scores according to the GSE63898 dataset. (D) Top 20 differentially expressed genes (DEGs) between groups C1 and C2. (E) Levels of immune infiltration of activated CD8+ T cells, activated CD4+ T cells, and natural killer (NK) cells between groups C1 and C2. (F) Functional enrichment analysis of DEGs between groups C1 and C2. (G) Protein-protein interaction network analysis revealed nine key chemokine hub genes.

4. Tex cells might represent a key T-cell subtype in the progression of liver cirrhosis

Based on prior analysis, activated CD8+ T cells may be a crucial immune cell type affecting the progression of liver cirrhosis. We therefore conducted further analysis on CD8+ T cells. The single-cell dataset SCP2154 of human liver cirrhosis was acquired from the Single Cell Portal database, which encompassed 328,783 cells. Key cell subtypes, namely T cells and NK cells, were selected from the dataset (Fig. 4A). The dataset was then screened to include cells from five cirrhotic patients and five healthy controls. CD3D+CD3E+CD8A+ T cells were finally chosen from these 10 samples, excluding cells with gene counts lower than 500. This resulted in a dataset of 4,878 CD8+ T cells, each with gene counts ranging from 500 to 2,500 (Supplementary Fig. 1A). All CD3D+CD3E+CD8A+ T cells were then subjected to Leiden clustering, categorizing the cells into nine clusters (Fig. 4B). We observed considerable heterogeneity of CD8+ T cells in terms of cell population and subgroup distribution among samples (Fig. 4C). Further analysis revealed that clusters 4 and 9 were predominantly detected in cirrhotic patients, while other subgroups mainly presented in healthy individuals (Fig. 4D). Analysis of the expression levels of eight known T-cell subtype marker genes in the nine clusters showed that cells in cluster 1 highly expressed CD69 and CXCR6, markers of resident memory CD8+ T cells (Trm); cells in cluster 2 expressed FOS, IL7R, JUN, KLRB1, etc., which were markers of memory CD8+ T cells (Tmm); cluster 3 cells expressed CCR7, SELL, TCF7, and LEF1, markers of naive CD8+ T cells (Tnaive); cluster 4 expressed H4C3, TOP1, TUBB, STMN1, NAPC11, markers of proliferating CD8+ T cells (Tpro); cluster 5 expressed CX3CR1, FGFBP2, FCGR3A, markers of effector CD8+ T cells (Tef); clusters 6 and 7 expressed GZMK, GZMA, IFNG, CCL4, CCL5, etc., markers of effector memory CD8+ T cells (Tem); cluster 8 expressed XCL1, IFNG, etc., markers of cytotoxic CD8+ T cells (Tct); and cluster 9 expressed PDCD1, LAG3, CXCL13, LAYN, HAVCR2, markers of exhausted CD8+ T cells (Tex) (Fig. 4E). Based on the gene expression patterns of the nine subgroups, they were annotated as eight T-cell subtypes (Fig. 4F and G), with significant differential expression of marker genes (Fig. 4F). CD8A was highly expressed across all cells and subgroups (Fig. 4H, Supplementary Fig. 1C), while CD4 and CD79A were lowly expressed through all subgroups (Supplementary Fig. 1C). The proportion of each CD8+ T-cell subtype in each sample was then calculated. We found a significant reduction in the proportion of Tem cells in cirrhotic patients, with their Tmm and Tct cells decreased, but Tex cells increased (Fig. 4I, Supplementary Table 5). This may suggest a significant shift in the immune microenvironment of cirrhotic patients, especially in the types of CD8+ T cells, with Tem and Tex cells possibly playing a key role.

Figure 4. Identification of CD8+ T-cell subtypes through single-cell RNA sequencing data. (A) T cells and NK cells were selected for further analysis. (B) Leiden clustering identified nine subtypes of CD8+ T cells. (C) Distribution of CD8+ T cells in five cirrhotic patients and five healthy individuals. (D) Comparison of CD8+ T-cell distribution between cirrhotic patients and healthy individuals. (E) Expression levels of eight known T-cell subtype marker genes across nine CD8+ T-cell subtypes, including resident memory CD8+ T cells (Trm), memory CD8+ T cells (Tmm), naive CD8+ T cells (Tnaive), proliferating CD8+ T cells (Tpro), effector CD8+ T cells (Tef), effector memory CD8+ T cells (Tem), cytotoxic CD8+ T cells (Tct), and exhausted CD8+ T cells (Tex). (F) Expression patterns of marker genes in the eight CD8+ T-cell subtypes. (G) Annotation of nine clusters into eight subtypes based on the expression levels of known T-cell subtype marker genes. (H) Expression distribution of CD8A across all CD8+ T cells. (I) Proportions of Tem, Tex, Tmm, and Tct cells in each cirrhotic and healthy individual.

A further analysis of the proportions of eight CD8+ T-cell subtypes in cirrhotic and healthy individuals revealed that Tex and Tpro cells were increased in cirrhotic patients, whereas all other CD8+ T-cell subtypes decreased (Fig. 5A and B). Analysis of DEGs across each CD8+ T-cell subtype identified 537, 260, 231, 605, 213, 401, 123, and 339 DEGs in Trm, Tmm, Tnaive, Tpro, Tef, Tem, Tct, and Tex, respectively (Supplementary Table 6). Signal intensity calculations based on DEGs in Tem and Tex cells were then performed (Fig. 5C). KEGG enrichment analysis of DEGs in Tem and Tex cells suggested their primary involvement in pathways like oxidative phosphorylation, cell adhesion molecules, T-cell receptor signaling pathway, chemokine signaling pathway, and NOD-like receptor signaling pathway (Fig. 5D). PPI analysis revealed key DEGs in Tem cells as CD8A, CD247, IFNG, CCL5, TNF, etc. (Fig. 5E), whereas in Tex cells, key DEGs included TOX, LAG3, CXCL13, TNFRSF9, TYROBP, HAVCR2, and PDCD1, indicating a pronounced exhaustion phenotype in Tex cells in cirrhosis (Fig. 5F). Further cellular trajectory analysis showed that Tex cells were at the center of CD8+ T-cell differentiation, with strong associations with Tpro, Tct, and Tnaive cells, which suggested that these cells may ultimately differentiate into Tex cells in cirrhosis (Fig. 5G). Thirteen genes including RGS1, CD2, COTL1, DUSP4, CXCL13, etc., were significantly overexpressed in Tex cells, suggesting their potential association with CD8+ T-cell exhaustion in cirrhosis (Fig. 5G). The high expression of these exhaustion-related genes may indicate that Tex represents an overly activated subtype of CD8+ T cells.

Figure 5. Analysis of gene expression characteristics in CD8+ T-cell subtypes. (A) Proportions of eight CD8+ T-cell subtypes in cirrhotic and healthy individuals. (B) Proportions of each CD8+ T-cell subtype in cirrhotic patients and healthy individuals. (C) Density distribution of significantly overexpressed characteristic genes in Tem and Tex cells. (D) Kyoto Encyclopedia of Genes and Genomes enrichment analysis of differentially expressed genes (DEGs) in Tem and Tex cells. (E) Protein-protein interaction (PPI) network of specific DEGs in Tem cells. (F) PPI network of specific DEGs in Tex cells. (G) Partition-based graph abstraction (PAGA) trajectory analysis based on eight CD8+ T-cell subtypes, showing the proportion of each cell type in diseases, patient distribution, and 13 genes highly expressed in Tex cells.

5. Role of CXCL13+ Tex cells in cirrhosis progression and immunological landscape

Previous analyses indicated high expression of the CXCL13 gene in Tex cells. Initially, we analyzed the potential differentiation trajectories of eight CD8+ T-cell subtypes, and identified lineage5 as the terminal differentiated state of CD8+ T cells in cirrhosis (Fig. 6A). CXCL13 showed specific overexpression at the differentiation endpoint of Tex cells (Fig. 6B), with RGS1 also exhibiting a gradual increase during Tex cell differentiation (Fig. 6B). A differential analysis between CXCL13+ and CXCL13- CD8+ T cells showed significant upregulation of genes such as CXCL13, DUSP4, PDCD1, LYST, RGS1, RGS2, and RGS10 in CXCL13+ CD8+ T cells (Fig. 6C, Supplementary Table 7). Gene feature analysis of the CXCL13+ CD8+ T-cell subcluster revealed a high overlap with CXCL13 expression (Fig. 6D), suggesting that CXCL13+ Tex cells might represent the most specialized subtype within Tex cells, potentially playing a crucial role in cirrhosis.

Figure 6. CXCL13+ Tex cells represented a key subtype associated with cirrhosis. (A) Cell differentiation trajectory inference based on eight CD8+ T-cell subtypes, where lineage5 indicated that Tex cells may be the terminal state of CD8+ T cells in cirrhosis. (B) Specific overexpression of CXCL13 in Tex cells. (C) Expression pattern of RGS1 in CD8+ T cells, which also showed high expression in Tex cells. (D) Differentially expressed genes between CXCL13+ and CXCL13- CD8+ T cells. (E) Gene Ontology network enrichment analysis of characteristic genes in CXCL13+RGS1+ CD8+ T cells. (F) Protein-protein interaction network analysis of characteristic genes in CXCL13+RGS1+ CD8+ T cells. (G) Expression distribution of significantly overexpressed genes in CXCL13+RGS1+ CD8+ T cells.

KEGG enrichment analysis of characteristic genes in CXCL13+ Tex cells revealed associations with pathways such as cell adhesion molecules, PD-L1 expression and PD-1 checkpoint, T cell receptor signaling, apoptosis, necrosis, hypoxia-inducible factor (HIF)-1 signaling, glycolysis, and pyruvate metabolism (Supplementary Fig. 2). GO network analysis of CXCL13+ Tex cell characteristic genes indicated their main involvement in functions like regulation of leukocyte activation, cytokine signaling in the immune system, leukocyte activation, and T cell modulation (Fig. 6E). Further network analysis identified key genes in CXCL13+ Tex cells (Fig. 6F), including PDCD1, TNFRSF9, LAG3, TOX, HLA-DRB5, HLA-DMA, HLA-DPB1, HLA-DRB1, and so on. The expression patterns of these genes not only showed high abundance of Tex cells, but also exhibited specific overexpression in CXCL13+ Tex cells (Fig. 6G). In summary, we hypothesized that CXCL13+ Tex may play a role in the immunosuppression associated with cirrhosis.

To further analyze the potential role of CXCL13+ Tex in liver cirrhosis, we conducted an in-depth analysis of cells within the Tex subgroup (Fig. 7A). The distribution maps of Tex cells revealed distinct differences and separation between cirrhotic and healthy individuals. A differential analysis between Tex cells from cirrhotic patients (cTex) and healthy individuals (hTex) showed significant upregulation of genes such as CXCL13, CCL3L1, TIGIT, CD27, FAM118A, TNFRSF9, LYST, DUSP4, and PDCD1, but downregulation of genes like EIF5A, LGALS1, S100A10, XCL1, GLIPR1, CD160, IFITM1, NCR3, KLRG1, and IL7R in cTex (Fig. 7B, Supplementary Table 8). Further analysis of the expression distribution of these genes demonstrated significant progressive changes between cirrhotic patients and healthy individuals, implying the importance of these genes (Fig. 7C and D). Moreover, we detected CXCL13 expression in cirrhotic patients but not the healthy group, suggesting that CXCL13 might be a key driver gene leading to the exhausted phenotype of cTex. Subsequent GO network enrichment analysis of significantly upregulated DEGs in cTex indicated their association with functions like regulation of leukocyte activation, positive regulation of immune response, cell activation, and regulation of T cell differentiation (Fig. 7E). Network analysis of the DEGs suggested that genes such as CD74, HLA-DRA, CCR7, CD40LG, and CXCR4 might also be crucial regulators inducing the exhausted phenotype of cTex (Fig. 7F).

Figure 7. Differences in Tex cells between cirrhotic patients and healthy individuals. (A) Distribution of Tex cells in cirrhotic patients and healthy individuals. (B) Differentially expressed genes (DEGs) in Tex cells between cirrhotic patients and healthy individuals. (C) Expression distribution of significantly overexpressed differential genes in cirrhotic patients. (D) Expression distribution of significantly overexpressed differential genes in healthy individuals. (E) Gene Ontology network enrichment analysis of DEGs in Tex cells between cirrhotic patients and healthy individuals. (F) Protein-protein interaction network analysis of DEGs in Tex cells between cirrhotic patients and healthy individuals.

6. Evaluating the impact of CD8+ T-cell subtypes on cirrhosis and liver cancer prognosis

To validate the functional role of characteristic genes in the eight CD8+ T-cell subtypes, GSVA was performed on the gene sets of Trm, Tmm, Tnaive, Tpro, Tef, Tem, Tct, and Tex cells in the GSE14323 and GSE63898 datasets. The results indicated a significant increase of Tex, Tef, Tpro, and Trm in the cirrhotic samples of the GSE14323 dataset, with Tex being the most prominent cell type (p=2.84e−14) (Fig. 8A). Further comparison of cell scores in the C1 and C2 groups of the GSE63898 dataset revealed that Tex, Tef, Tpro, and Trm were significantly elevated in the C1 group of GSE14323, with Tex again being the most prominent (p=0.003) (Fig. 8B). This suggested that Tex may facilitate disease progression in patients with cirrhosis. Subsequent correlation analysis between the CXCL13+ Tex cell signature score and the Cirrhosis Disease Score C0023890 showed a strong positive correlation (Fig. 8C). Validation of the CXCL13+ Tex cell signature score in the GSE14323 and GSE63898 datasets demonstrated its significant elevation in cirrhosis patients and the C1 group (Fig. 8D). Based on these analyses, we concluded that CXCL13+ Tex cells represent the most critical subtype in the progression of cirrhosis, with CXCL13+ being a key cellular factor in this subtype. Analyzing single-cell data from the liver cancer dataset GSE140228, CXCL13 expression was detected mainly in exhausted CD8+ T cells and Treg cells. We further validated the prognostic efficacy of the CXCL13+ Tex cell signature in the TCGA-LIHC (Cancer Genome Atlas Liver Hepatocellular Carcinoma) dataset through lasso regression, modeling 162 expressed genes in TCGA-LIHC. A prognostic model of 16 genes was derived (p=1.54e-10) (Fig. 8E and F), indicating the association between CXCL13+ Tex cell signature and the poorer prognosis in liver cancer. These results suggested that CXCL13+ Tex may be an important CD8+ T-cell subtype affecting the progression of liver cirrhosis and liver cancer, and might also be key in responding to immune checkpoint therapy.

Figure 8. Validation of CD8+ T-cell subtypes in Bulk-seq data and the impact of CXCL13+CD8+ T cells on prognosis. (A) Validation of Gene Set Variation Analysis (GSVA) scores based on gene features of Tef, Tex, Tpro, and Trm cells in the GSE14323 dataset. (B) Validation of GSVA scores based on gene features of Tef, Tex, Tpro, and Trm cells in the GSE63898 dataset. (C) Correlation analysis between the GSVA score of the CXCL13+CD8+ T-cell gene signature and the cirrhosis disease score. (D) Validation of the GSVA score of the CXCL13+CD8+ T-cell gene signature in the GSE14323 and GSE63898 datasets. (E, F) Validation of the prognostic value of CXCL13+CD8+ T-cell gene features in The Cancer Genome Atlas Liver Hepatocellular Carcinoma (TCGA-LIHC) dataset.

Through transcriptomic immune infiltration analysis, our study discovered a significant increase in the proportion of activated CD8+ T cells, activated CD4+ T cells, NK cells, NKT cells, and myeloid-derived suppressor cells in patients with liver cirrhosis. Each cell type, with unique immunological functions, might critically impact cirrhosis development. Activated CD8+ T cells, essential effector T cells that eradicate virus-infected hepatocytes via cytotoxicity, reflect persistent hepatic inflammation and cellular damage in cirrhosis.³⁴ Simultaneously, activated CD4+ T cells could trigger an overactive inflammatory response through cytokine secretion, thus advancing cirrhosis.³⁵ Furthermore, increases in NK and NKT cells may signify changes in the hepatic immune environment under hepatic damage and chronic inflammation.^36,37 Lastly, the escalation in myeloid-derived suppressor cells, an immunosuppressive cell population, might limit excessive immune responses, preventing further hepatic damage, but could also obstruct effective antiviral immune responses, potentially maintaining chronic viral infection.^38,39 In summary, the amplified proportions of these immune cells in cirrhosis could be linked to sustained inflammatory responses, changes in the hepatic immune environment, and immunosuppressive responses aiming to deter further hepatic injury. Their collective activities and interactions may be key in propelling cirrhosis development, offering crucial insights into the immune pathophysiology of cirrhosis and potential novel therapeutic avenues.

Further examination revealed CCL19, GABBR1, CXCL10, CCR7, CXCL11, CCR5, CCL5, CXCL13, and CXCL9 as potentially crucial genes influencing cirrhosis. Our data indicated that high expressions of these genes were associated with increased risks of cirrhosis. CCL19, CXCL10, CXCL11, CCL5, CXCL13, and CXCL9, which are chemokines involved in immune cell trafficking, likely exacerbate inflammation and fibrosis by recruiting immune cells to the liver, perpetuating inflammatory response and subsequent tissue damage.^40-45 Meanwhile, CCR5 and CCR7, receptors for the above chemokines, may influence the migration patterns of immune cells, contributing to inflammation and fibrosis.^46,47 Gamma-aminobutyric acid receptor B1 (GABBR1), a receptor for GABA, has known anti-inflammatory properties;⁴⁸ however, its role in cirrhosis remains unclear and requires further exploration. In summary, high expressions of these genes might increase cirrhosis risk by influencing immune cell trafficking, promoting liver inflammation, fibrosis, and potentially altering GABA-mediated signaling. These findings add nuance to the understanding of the genetic underpinnings of cirrhosis, potentially highlighting novel therapeutic targets.

CD8+ T cells are pivotal in controlling viral infections through the production of effector molecules such as granzyme B, interferon-gamma, and tumor necrosis factor-alpha.^49-51 These T cells eradicate antigen-expressing cells during acute infections, facilitating the formation of long-lived memory T cells.^49-51 However, in chronic liver infections like hepatitis B virus and hepatitis C virus, the development of memory T cells is hindered due to T cell exhaustion.^52,53 In our study, we observed a marked decrease in resident memory CD8+ T cells (Trm), memory CD8+ T cells (Tmm), and effector memory CD8+ T cells (Tem) in patients with cirrhosis, while exhausted CD8+ T cells (Tex) exhibited an opposite trend. Our results suggest a higher proportion of Tex cells in patients with non-alcoholic steatohepatitis-induced cirrhosis, which were potentially linked to prolonged hepatitis B virus infection, and may lead to CD8+ T cell exhaustion.

Previous research indicates that exhausted T cells are functionally and phenotypically heterogeneous. Our findings implicate RGS1+ CD8+ T cells (Tex) as crucial in cirrhosis, primarily impacting the regulation of T cell activation and cytokine signaling. RGS1, a regulator of G-protein signaling, influences T cell migration and homeostasis.⁵⁴ Elevated expression of RGS1 may lead to increased recruitment or retention of these specific CD8+ T cells in the liver, exacerbating local inflammation and injury.⁵⁵ In summary, the role of RGS1+ Tex in cirrhosis is attributed to their involvement in T cell activation regulation, cytokine signaling, and hepatic tissue damage, offering new insights into immune dysregulation in cirrhosis and identifying potential therapeutic targets. Moreover, our analysis highlights CXCL13+ Tex cells as a key subpopulation within Tex cells, significantly impacting pathways like PD-L1 expression and PD-1 checkpoint, T cell receptor signaling, apoptosis, HIF-1 signaling, and glycolysis. The modulation of PD-L1 expression and PD-1 checkpoint pathway may reflect an exhaustion state of cells, commonly seen in chronic liver diseases like cirrhosis.^56,57 CXCL13+ Tex cells could alter T cell activation via T cell receptor signaling, contributing further to immune dysregulation.⁵⁸ Their effects on HIF-1 signaling and glycolysis suggest metabolic reprogramming under chronic inflammation and hypoxia, conditions prevalent in cirrhosis.^59,60 Thus, the role of CXCL13+ Tex cells in cirrhosis is multifaceted, involving increased immune cell trafficking, potential immune exhaustion, T cell activation changes, apoptosis, and metabolic reprogramming. These findings deepen our understanding of cirrhosis immunopathology in cirrhosis and suggest potential therapeutic interventions.

Finally, we demonstrated elevated expression of PDCD1, HAVCR2, TIGIT, TNFRSF9, and DUSP4 in CXCL13+ Tex cells, indicative of CD8+ T cell exhaustion.^61,62 This exhaustion hampers the standard immune surveillance and clearance functions of CD8+ cells. Therapeutically reversing this state could potentially halt or reverse cirrhosis progression. PDCD1 encodes PD-1, an immune checkpoint that suppresses T cell inflammatory activity, thereby regulating immune responses and self-tolerance.^61-63 Similarly, HAVCR2 and TIGIT, expressed on T cells, may inhibit CD8+ T cell activation.^61,64 In cirrhosis, the upregulation of these genes, coupled with DUSP4, signals T cell exhaustion, contributing to an ineffective immune response in ongoing liver injury and exacerbating cirrhosis. Targeting these immune checkpoints, such as with checkpoint blockade therapy, might reactivate T cells, slowing or reversing cirrhosis progression.

In conclusion, liver cirrhosis involves complex immunological mechanisms driving persistent inflammation and progressive fibrosis. Through an integrated analysis of transcriptomic and single-cell profiling data, our study elucidates the heterogeneity and functional dynamics of CD8+ T cells in cirrhosis pathogenesis. We demonstrate a significant decline of memory CD8+ T cell subsets like Tem, paralleled by an expansion of dysfunctional, exhausted Tex cells in cirrhotic livers. Tex cells display upregulation of multiple inhibitory receptors like PD-1, TIGIT and LAG3, indicating a state of unresponsiveness. We identify CXCL13+ Tex cells as a distinct subset strongly associated with advanced fibrosis. This terminally exhausted population exhibits dysregulation of pathways related to T cell migration, activation, metabolic adaptation and immunoregulation. An exhausted gene signature comprising PDCD1, TIGIT and other regulators in CXCL13+ Tex cells correlates with cirrhosis severity. Overall, our integrated analysis defines specialized CD8+ T cell states underlying cirrhosis progression, with CXCL13+ Tex cells reflecting immune dysfunction. Key genes and pathways governing the accumulation and exhaustion of Tex cells present potential diagnostic biomarkers and immunotherapeutic targets in liver cirrhosis. Further research should explore approaches to counteract Tex-mediated immune suppression and regain CD8+ T cell effector functions. Reinvigorating exhausted CD8+ T cells through immune checkpoint blockade could help constrain inflammation and fibrosis associated with cirrhosis.

This work was supported by Research of Key Techniques and Application of Liver-Kidney Organ Transplantation (202302AA310018).

No potential conflict of interest relevant to this article was reported.

Study concept and design: R.J., T.Y. Data acquisition: J.L. Data analysis and interpretation: Z.K., D.L., Z.S. Drafting of the manuscript: Z.K., D.L., Z.S. Critical revision of the manuscript for important intellectual content: R.J., T.Y., D.L. Statistical analysis: Z.Y., L.W., G.Y., M.Y. Obtained funding: R.J. Administrative, technical, or material support; study supervision: T.Y. Approval of final manuscript: all authors.

Supplementary materials can be accessed at https://doi.org/10.5009/gnl230345.