Single-cell sequencing reveals cell heterogeneity and aberrantly activated pathways associated with microvascular invasion in hepatocellular carcinoma

Abstract

Introduction:

Hepatocellular carcinoma (HCC) is the most common primary liver cancer, with microvascular invasion (MVI) identified as a major predictor of early recurrence. However, the intratumor cellular heterogeneity of MVI, the identification of pertinent biomarkers, and the role of intercellular signalling interactions in MVI progression are unclear. This study aims to explore these aspects using single-cell transcriptomic analysis.

Methods:

The present study utilized single-cell transcriptomic data from public databases to conduct an in-depth transcriptome analysis of tumour tissues and adjacent nontumor tissues from five patients with hepatocellular carcinoma, with a particular focus on samples from three patients exhibiting microvascular invasion. Bioinformatics tools were employed to analyze gene expression patterns and signalling pathways.

Results:

The findings indicated that MVI-positive malignant cells activate multiple signalling pathways to facilitate invasion and metastasis. Specific malignant cell subtypes strongly associated with MVI were identified, exhibiting distinctive gene expression patterns related to proliferation, invasion, and metabolic reprogramming of tumour cells. Further analysis revealed that the laminin and VEGF signalling pathways are crucial for remodelling the tumour microenvironment and angiogenesis associated with MVI. The MARCKSL1 gene was predominantly expressed in MVI-positive malignant cells and may contribute to MVI progression by interacting with the PTN signalling network. Additionally, MARCKSL1 is linked to tumour resistance to multiple anticancer drugs.

Discussion:

This study sheds light on the molecular characteristics and functional heterogeneity of MVI-associated malignant cell subpopulations. The single-cell transcriptome and bioinformatics analyses provided insights into the mechanisms driving MVI, potentially aiding the development of targeted diagnostic and therapeutic strategies. Future research should further validate the role of MARCKSL1 in MVI progression and explore its potential clinical applications.

1 Introduction

Hepatocellular carcinoma (HCC) is the third leading cause of cancer-related death worldwide, with the number of deaths and diagnoses expected to increase by more than 55% globally by 2040 (Marrero et al., 2018; Rumgay et al., 2022). Currently, the primary treatments are liver resection and liver transplantation. However, the recurrence rates remain high after treatment, with 5-year recurrence rates of 70% and 35% after liver resection and liver transplantation, respectively (Xu et al., 2019). In recent years, significant attention has been given to the role of microvascular invasion (MVI) in HCC. MVI, defined as the invasion of tumour cells into the spaces between vascular endothelial cells, including portal veins, hepatic arteries, and lymphatic vessels, is an independent risk factor for early postoperative recurrence and poor prognosis in HCC patients (Gouw et al., 2011). Notably, for patients with a solitary small HCC lesion less than 5 cm in diameter, the presence of MVI significantly reduces recurrence-free survival (RFS) and overall survival (OS) rates (Sheng et al., 2020; Hong et al., 2021; Xiong et al., 2022; Sumie et al., 2008). Thus, there is an urgent need for more specific molecular biomarkers with prognostic and therapeutic significance.

The rapid advancement of single-cell RNA sequencing (scRNA-seq) technology in recent years has revolutionized the understanding of cellular heterogeneity in various pathological tissues (Ramachandran et al., 2019; Kuppe et al., 2021). ScRNA-seq has led to significant discoveries in liver cancer research. Studies have shown that tumour-associated macrophages (TAMs) in liver cancer are closely linked to poor patient prognosis, and they have identified critical genes in the inflammatory response of TAMs, such as SLC40A1 and GPNMB (Ma et al., 2019; Zhang et al., 2019). Additionally, scRNA-seq has been used to map various immune cell subpopulations within liver cancer tissues, including T cells and dendritic cells. Each subpopulation plays a unique role in the liver cancer microenvironment. For example, LAMP3-positive dendritic cells mediate immune suppression, while TREM2-positive TAMs inhibit the infiltration of CD8⁺ T cells into tumour tissue (Zhang et al., 2019; Zheng et al., 2017; Tan et al., 2023).

Despite these findings, a comprehensive understanding of the expression profiles of malignant cells in hepatocellular carcinoma, particularly during the progression of MVI, is lacking, and their specific roles in tumours are unclear. The present study investigated the expression profiles of malignant cells in hepatocellular carcinoma, systematically classified these cells, and detailed the cellular heterogeneity associated with MVI, as well as the molecular biological features of specific malignant subpopulations. A machine learning approach was used to construct a prognostic model based on signature genes of malignant cells, which not only enhanced the prognostic utility of the signature gene but also identified a previously unreported molecule, namely, MARCKSL1. Further studies indicated that MARCKSL1 may promote the development of MVI through interactions with the PTN signalling network. The present findings suggested that MARCKSL1, a potential therapeutic target for hepatocellular carcinoma and MVI progression, may be crucial for improving therapeutic strategies and clinical outcomes, particularly for patients with MVI.

2 Materials and methods

2.1 Data collection

TCGA-LIHC data were sourced from UCSC Xena (https://xenabrowser.net/datapages/), and scRNA-seq data were obtained from the GEO database (accession number GSE242889). The dataset included samples from five patients diagnosed with HCC, three of whom had MVI, and two without MVI. For each patient, both HCC tissue and adjacent non-tumor tissue were collected. After stringent quality control, we analyzed a total of 42,070 cells, which included 4,070 high-quality sequenced cells. The dataset underwent single-cell sequencing using the Illumina NovaSeq 6,000 platform.

2.2 Data processing

For preprocessing, quality control, normalization, and dimensionality reduction clustering of the single-cell data, we employed Seurat v4.3.0. Potential doublets were identified and removed using DoubletFinder v2.0.3. The quality control criteria included the expression of at least 400 genes per cell and a mitochondrial gene threshold of 20%. Subsequent normalization, identification of highly variable genes, and dimensionality reduction clustering were performed using Seurat’s default parameters and standard workflow. Harmony v1.2.3 was utilized to integrate data from different samples. Cell cluster naming was conducted by aggregating marker genes from the literature and through manual annotations. The FindAllMarkers function was employed to identify differentially expressed genes between cell subgroups using the Wilcoxon test, with selection criteria of an adjusted p-value <0.05 and |log2FC| > 0.5.

2.3 Cell type identification

We performed differential expression analysis of all genes within cell clusters using Seurat’s FindAllMarkers function to identify marker genes for each cluster. Criteria for identifying marker genes included an adjusted P value <0.05, expression percentage >0.25, and |log2[fold change (FC)]| > 0.25. Subsequently, different cell clusters were identified and annotated using the singleR packages, based on the compositional patterns of the marker genes. These annotations were manually validated and corrected with reference to the CellMarker database.

2.4 Cell-cell communication

The ‘CellChat’ R package v1.1.3 was utilized to infer cellular communication within the tumor microenvironment based on receptor-ligand interactions (Gupta et al., 2021). This involved counting the number of links and collecting communication probabilities to calculate communication networks. The interaction time and the total strength of interactions between two arbitrary cell populations were visualized. Scatter plots were drawn to illustrate the primary sender (source of signals) and receiver cells (targets) within a two-dimensional space. This visualization helps identify the most significant contributors to the outgoing or incoming signals among a group of immune cells. We employed a pattern recognition approach, the global communication model, to discern how various immune cell types and signaling pathways operate in concert.

2.5 Pseudotemporal trajectory analysis

Cellular pseudotemporal trajectories were constructed using the Monocle two algorithm, an R package developed by Qiu et al. for single-cell trajectories (Huang et al., 2019). This algorithm employs machine learning techniques to reduce the high-dimensional expression spectrum into a low-dimensional space, organizing it into trajectories with branching points. Dynamic expression heatmaps were constructed using the plot_pseudotime_heatmap function. Integrated Machine Learning-Based Approach to Derive Feature Signatures: We integrated up to 10 machine learning algorithms, including Random Survival Forest (RSF), Elastic Net (Enet), Lasso, Ridge, Stepwise Cox, CoxBoost, Partial Least Squares Regression for Cox (plsRcox), Supervised Principal Components (SuperPC), Generalized Augmented Regression Model (GBM), and Survival Support Vector Machine (survival-SVM). Based on these approaches, a consensus model was generated. In total, 101 algorithm combinations were executed to match the predictive model based on the leave-one-out cross-validation (LOOCV) framework. The TCGA-LIHC dataset was divided into training and test datasets.

2.6 Immune infiltration assessment

We utilized the CIBERSORT algorithm to quantitatively assess the level of immune cell infiltration in patients with pancreatic adenocarcinoma (PAAD), exploring differences in cell abundance between high-risk and low-risk patient groups (Malumbres and Barbacid, 2009). Additionally, we calculated and analyzed the Pearson correlation between immune cell abundance and risk scores. To further delineate potential differences in immune function, we applied enrichment scores obtained from single-sample gene set enrichment analysis (ssGSEA) (He et al., 2021). Subsequently, we used the Wilcoxon test to compare immune function between the high-risk and low-risk groups.

2.7 Statistical analysis

All statistical analyses and data visualization were conducted using R software (version 4.1.3). Pearson’s correlation coefficient was used to assess the relationship between two continuous variables. For quantitative data, values were compared between subgroups using two-tailed, unpaired Student’s t-test or one-way analysis of variance (ANOVA) with Tukey’s multiple comparison test. A P-value <0.05 was considered statistically significant.

3 Results

3.1 Single-cell sequencing and cytotyping of cancerous and paracancerous tissues

The present study analysed 10 samples from the GSE242889 dataset, which included five patients with HCC. Among these patients, three patients exhibited MVI, while two patients did not exhibit MVI. HCC tissue and adjacent nontumor tissue were collected from each patient for single-cell RNA sequencing. After performing stringent quality control and removing duplicate cells, 42,070 cells were retained for analysis. These cells were categorized into 30 distinct clusters, and the origin of each cell was visualized (Figure 1A; Supplementary Figure S1A, B). The cell clusters were annotated using classical cell marker genes (Li et al., 2024), and they were classified them into 9 cell types, illustrating the distribution of cells in both cancerous and paracancerous samples (Figures 1B, C; Supplementary Figure S1C, D, average number of detected genes per defined cell type(Supplymentary Data.xlsx)). The following cell types were identified: myeloid cells with high LYZ expression; malignant cells with high TF expression; T cells with high CD3D and CD3G expression; B cells with high CD79A expression; NK cells with high NKG7 and KLRD1 expression; endothelial cells with high PECAM1, CLDN5, and FLT1 expression; mesenchymal cells with high ACTA2 expression; and haematopoietic progenitor cells (HPCs) with high EPCAM, KIT, MS4A2, and GATA2 expression (Figure 1D). The heatmap in Figure 1E displays the top 50 genes highly expressed in each cell type, highlighting the specificity of their compartmentalization. In this dataset, myeloid cells constituted the largest proportion of cells, followed by T cells (Supplementary Figures S1A–B). Figures 1F, G show the proportions of the cell types in different samples, in which the proportions of myeloid cells, malignant cells, and T cells were greater. Moreover, the malignant cells were almost exclusively derived from tumour tissues, demonstrating significant heterogeneity between samples.

FIGURE 1

3.2 Malignant tumour cells with MVI exhibit increased invasive and metastatic potential

The differences in the expression of various cell types were analysed between cancerous tissues and their adjacent noncancerous counterparts using the Wilcoxon rank sum test. Genes associated with lipid metabolism, such as APOA2, APOC3, APOC1, and APOE, were significantly upregulated in multiple cell types of hepatocellular carcinoma, indicating enhanced lipid metabolic activities (Figure 2A). Additionally, using the Augur algorithm in conjunction with the random forest model (Squair et al., 2021), the cell types that exhibited significant transcriptomic perturbations across different biological states were identified. Among the nine identified cell types, malignant cells showed the most significant transcriptomic changes (area under the curve (AUC) = 0.955), followed by endothelial cells (AUC = 0.866) and T cells (AUC = 0.851) (Figure 2B). In particular, the changes in endothelial cells were important. The DEGs in both states were analysed, and protein‒protein interaction (PPI) network mapping was performed using the STRING database (Figure 2C). Subsequent Gene Ontology (GO) and Kyoto Encyclopaedia of Genes and Genomes (KEGG) analyses of these genes revealed alterations in cellular connectivity and pathways related to oxidative phosphorylation (Figure 2D; Supplementary Figure S2).

FIGURE 2

The cell type distribution was compared between tumour tissues with MVI (1T_C21, 2T_C24, and 3T_C25) and those without MVI (4T_C29 and 5T_C36). There was a greater abundance of T cells, B cells, and NK cells in tumours with MVI than in those without MVI, suggesting a more active immune microenvironment in these samples (Figure 2E; Supplementary Figure S3A). Myeloid cells were the most numerous, and there was a significant preference for the distribution of different malignant cell subtypes based on MVI status (Supplementary Figures S3B–D). Furthermore, the Augur algorithm indicated that malignant cells experienced the greatest transcriptome perturbation in both MVI and non-MVI tumours (AUC = 0.914, Supplementary Figure S2F), underscoring a significantly active immune microenvironment (Supplementary Figures S2E, F). Gene set enrichment analysis (GSEA) on the differentially expressed genes in malignant cells from MVI and non-MVI tumours was performed. The anchoring junction, focal adhesion, and cell-substrate junction pathways were significantly activated in MVI-positive malignant cells, confirming their increased invasiveness and metastatic potential (Figures 2G–I). These findings highlighted the critical role of MVI in the progression of hepatocellular carcinoma and identified potential targets for future therapeutic strategies.

3.3 Gene expression heterogeneity of different malignant cell subtypes in microvascular invasion

A detailed reclassification of malignant cells in tumour tissues was performed, which identified nine distinct cell subtypes (Figure 3A). Analysis of the distribution of these subtypes in samples with and without MVI elucidated that the MCs_1 subtype was predominantly found in samples without MVI, whereas the MCs_2 to MCs_6 and MCs_9 subtypes were present in samples with MVI (Figures 3B, C). To further explore gene expression patterns, the gene expression scores were calculated for individual cells, which reflected their relative gene expression levels. These scores were used to employ an unsupervised clustering approach to visualize the unique gene expression patterns across different subtypes (Figure 3D).

FIGURE 3

The MET and TNFSF10 receptor tyrosine kinases, both highly expressed in subtype 2, were closely associated with the proliferation, survival, invasion, and metastasis of hepatocellular carcinoma cells. Aberrant activation of MET has been shown to correlate with progression and poor prognosis in hepatocellular carcinoma (Comoglio et al., 2018; Gupta et al., 2021; Huang et al., 2019). Moreover, high expression of CDK1 in subtype 5 is closely associated with the development of many cancers, including hepatocellular carcinoma (Malumbres and Barbacid, 2009). CYP3A4, which is highly expressed in subtypes 4 and 6, is the most abundant cytochrome P450 enzyme in the liver. Targeting CYP3A4 in oxidative metabolism has been reported to be an important strategy for enhancing the sensitivity of hepatocellular carcinoma cells to chemotherapeutic drugs (Figure 3E) (He et al., 2021; Özkan et al., 2021). Together, these findings revealed that each malignant cell subtype exhibits distinct functional characteristics. Subtypes with high expression of oncogenes are significantly correlated with MVI-positive samples, highlighting their potential role in the aggressive behaviour of these tumours.

3.4 Activation of angiogenic and mesenchymal transition pathways in malignant cell subtypes with MVI

To investigate the differences between malignant cell subtypes with MVI in hepatocellular carcinoma, copy number variation (CNV) analyses were conducted for subgroups of malignant cells with or without MVI. The inferCNV tool was used to perform a detailed assessment of intratumoral CNV heterogeneity in the subgroups. As anticipated, the CNVs exhibited significant heterogeneity among the different malignant cell subgroups (Figure 4A). Chromosomes four and nine showed significant copy number amplifications in some subgroups, whereas chromosomes 6 and 13 had significant deletions. The CNVs of each cell were visualized by a t-SNE plot (Figure 4B), and the quantitative CNV scores of each malignant cell subtype were visualized using a boxplot (Figure 4C). The CNV levels of the MCs_4 subtype were greater than those of the other subtypes, suggesting a more malignant phenotype. Similarly, the MCs_2 and MCs_6 subtypes also displayed elevated CNV levels.

FIGURE 4

Further analysis of the cancer signature gene sets (hallmark gene sets) from MSigDB was performed to assess the expression of each malignant cell subtype across different pathways (Figure 4D). The MCs_1 subtype significantly upregulated the angiogenesis and epithelial–mesenchymal transition (EMT) pathways in tumours with MVI, which indicated that a distinct gene expression pattern promoted tumour angiogenesis. In contrast, the MCs_4 and MCs_6 subtypes upregulated multiple pathways.

A differential expression analysis for MCs_4 and MCs_6 with or without MVI was conducted using limma. The Myc target pathway was significantly activated in both subtypes with MVI, suggesting rapid proliferation and metabolic reprogramming of tumour cells, which may contribute to their invasive and metastatic capacities in MVI. Moreover, the upregulation of the Wnt/Beta-Catenin and E2F Targets pathways further confirmed the high growth and invasive potential of these tumour cells (Figures 4E, F). Differences in the remaining malignant cell subpopulations are detailed in Supplementary Figure S4.

3.5 Single-cell trajectory analysis of MVI malignant cell subpopulations in hepatocellular carcinoma patients

To analyse the changes in different malignant cell subtypes with and without MVI, the Monocle software package was used to perform a time-series analysis of the nine malignant cell subtypes, which identified five distinct cell states. The state without MVI was identified as the starting point of the trajectory, while the state with MVI marked the endpoints of two differentiation trajectories (Figures 5A–C). The MCs_4 and MCs_6 cell subtypes predominantly focused on one of the trajectory endpoints, whereas the MCs_2 and MCs_3 cell subtypes focused on the other, suggesting that these subtypes may represent two distinct fates (Figure 5D). Further analysis explored changes in gene differentiation trajectories before and after node 1, and the top 30 and top 200 genes were identified (Figure 4E; Supplementary Figure S5A). By tracking the dynamic changes in the expression of genes along the trajectories (Figure 4F), two distinct expression patterns among the top 30 genes were identified, as determined by two statistical approaches for identifying DEGs. This finding was partially validated by BEAM analyses at node 2 (Supplementary Figures S5B, C).

FIGURE 5

Intersection analysis of these top 30 genes identified MET, TIMP1, and VTN as key genes for progression along the trajectory, with consistent findings using both statistical methods. These genes were highly expressed in the MCs_4, MCs_5, and MCs_6 cell subtypes, as well as in states four and 5 (Figure 4G; Supplementary Figure S5D). Research has shown that MET activates signalling pathways, such as the PI3K/AKT and RAS/MAPK pathways, promoting the proliferation, survival, and migration of hepatocellular carcinoma cells, as well as influencing vascular endothelial cells to promote neointimal formation (Bussolino et al., 1992; Bladt et al., 1995; Huh et al., 2004). TIMP1 regulates the remodelling of the extracellular matrix by inhibiting metalloproteinase activity, which helps tumour cells adapt to changes in the matrix, thereby enhancing cell survival by inhibiting apoptosis (Guccini et al., 2021; Justo and Jasiulionis, 2021), potentially affecting blood vessel formation during MVI. Vitronectin (VTN) promotes tumour cell migration and invasion by enhancing adhesion to the stroma and interactions with cell surface integrins, possibly supporting tumour nutrient and oxygen supply by influencing vascular endothelial function and promoting neointimal formation (Wei et al., 2018; Zhu et al., 2015). In conclusion, the present study highlighted the critical roles of MET, TIMP1, and VTN in the progression of MVI in hepatocellular carcinoma, underscoring their previously unreported contributions to this process (Figures 5E–F).

3.6 Immune cell dynamics and MVI-dependent alterations in HCC microenvironment

To comprehensively understand the alterations in the tumor immune microenvironment, we conducted a systematic analysis of lymphoid and myeloid cells in our dataset. Using unsupervised clustering methods, we initially classified lymphoid cells into 10 distinct subpopulations and annotated their cell types based on subgroup-specific highly expressed genes (Figures 6A, B). Cell abundance analysis revealed significant enrichment of proliferating lymphocytes (ProL) and IgL-expressing plasma cells (IgL-Plasma) in the tumor group, while the proportion of CD4⁺ T cells (CD4T) was notably decreased. Interestingly, in tumors positive for microvascular invasion (MVI), secretory B cells (secB) were significantly reduced compared to MVI-negative tumors, whereas CD4T showed an opposite trend, indicating a close relationship between MVI status and immune cell composition (Figure 6C). Transcriptomic analysis demonstrated that CD8⁺ T cells (CD8T) exhibited the most significant transcriptional changes in tumor versus normal tissue comparisons, while CD4T showed the most prominent alterations in MVI-related comparisons (Figure 6D). Detailed analysis revealed that genes related to oxidative phosphorylation and energy transfer were significantly upregulated in tumor CD8T, while cytoplasmic ribosome-related pathway genes were downregulated, suggesting metabolic reprogramming and functional state changes (Figure 6E). In MVI-positive samples, upregulated genes in CD4T were primarily enriched in cytoplasmic translation and focal adhesion pathways, indicating a potentially highly activated state (Figure 6F).

FIGURE 6

In the myeloid cell analysis, we identified seven distinct subpopulations and performed detailed analysis of their functional marker expression profiles (Figures 6G, H). The tumor group showed decreased proportion of monocytes (Mono) but increased proportions of tissue-resident macrophages (TRM) and activated macrophages (actMac) (Figure 6I). Notably, in MVI-positive tumors, the proportion of type 1 conventional dendritic cells (cDC1) was significantly lower than in MVI-negative tumors (Supplementary Figure S6I). Transcriptomic analysis indicated that TRM underwent the most significant transcriptional remodeling. Specifically, in tumor TRM, genes related to lipid metabolism and transport (such as APOA1) and nutrient metabolism and transport (such as RBP4, FABP5) were significantly upregulated, suggesting these cells may have acquired special metabolic functions and participate in nutrient redistribution within the tumor microenvironment. Meanwhile, the downregulation of antigen presentation molecules (HLA-A/B) and immunoglobulin (IGKC) suggested potentially suppressed immune functions, facilitating tumor immune escape (Figures 6J, K). Assessment of cell subpopulation functional states using MSigDB gene sets revealed that myeloid cell subpopulations in normal tissue generally exhibited low pathway activity levels, while differences in MVI status were reflected in the activation levels of multiple functional pathways, further supporting the association between microvascular invasion and immune cell functional changes (Figure 6L). These findings reveal the dynamic characteristics of immune cells in the hepatocellular carcinoma microenvironment and highlight the crucial role of microvascular invasion in reshaping the tumor immune microenvironment, providing new insights for understanding HCC progression mechanisms and developing immunotherapy strategies.

3.7 Upregulated laminin and VEGF signalling pathways in the MVI malignant subgroups

The three malignant cell subpopulations with high CNV, namely, MCs_2, MCs_4, and MCs_6, are potentially associated with the progression of MVI in hepatocellular carcinoma. To further explore the communication between these malignant cells and other cell types, the myeloid cells that exhibited the highest percentage of interaction were subdivided (Supplementary Figures S5E, F) and subjected to CellChat analysis (Supplementary Figure S6A). Significant intercellular communication was identified within all three malignant cell subpopulations (Supplementary Figures S6B–E). Notably, the MCs_2 subtype demonstrated the most extensive receptor‒ligand interactions among the three subtypes, engaging both the ApoA and complement signalling pathways, as well as acting as both a receiver and a sender of signals. The VTN, CypA, and VEGF signalling pathways were significantly activated in these subpopulations, underscoring the enhanced invasiveness and angiogenic capacity of these tumour cells (Figure 7A). Among the receptor‒ligand interactions, the laminin and VEGF signalling pathways were the most significantly activated, suggesting remodelling of the tumour microenvironment and the extracellular matrix.

FIGURE 7

By binding to cell surface integrins, laminin activates multiple downstream pathways that significantly influence cell attachment, migration, and survival. This mechanism is crucial for maintaining tissue architecture and promoting tumour invasiveness and metastasis (Rick et al., 2019). Single-cell data analysis revealed that two malignant cell subpopulations, namely, MCs_4 and MCs_6, were the primary emitters of laminin signalling, with the MCs_4 subtype delivering stronger signals, particularly to the endothelial cell population (Figure 7B). The MCs_4 subtype was the dominant signaller of endothelial cells, and all three malignant cell subpopulations were identified as key influencers of this pathway (Figures 7D, F). The expression of laminin signalling-related molecules is depicted in Figure 7H. Moreover, the VEGF signalling pathway, which is known to support tumour growth and metastasis by promoting angiogenesis through vascular endothelial cell proliferation, migration, and neointima formation (Carmeliet, 2005), was also highlighted. The present study demonstrated that the endothelial cells were the primary receivers of VEGF signalling, with the MCs_4 and MCs_6 subtypes acting as the main emitters (Figures 7C, E). All three malignant cell subpopulations were identified as potential effectors of VEGF signalling (Figure 7G). Analysis of the expression of molecules in the VEGF pathway (Figure 7I) revealed that VEGFA was highly expressed in the MCs_4 and MCs_6 subtypes, whereas VEGFB was predominantly expressed in the MCs_2 subtype.

The patterns of afferent and efferent signals across all malignant cell subpopulations were analysed using cophenetic- and silhouette-based methods (Supplementary Figures S6F, G). In all malignant cell subgroups, the MCs_5 subtype emerged as the most significant signalling efferent, followed by the MCs_2 and MCs_4 subtypes. Conversely, the MCs_2 subtype was the most significant signalling afferent, followed by the MCs_3 and MCs_5 subtypes. The MCs_2, MCs_4, and MCs_6 subtypes displayed similar afferent and efferent patterns, suggesting a potential association (Supplementary Figures S7B–D).

3.8 MVI malignant cell subpopulation-related traits predict patient survival

Data from TCGA-LIHC were utilized to fit matched prognostic data, and the data were screened for genes associated with patient prognosis using one-way Cox regression, which resulted in 673 genes (p < 0.0001). The genes specifically expressed in the MCs_2, MCs_4, and MCs_6 subtypes were intersected with those obtained by one-way Cox regression, which resulted in 452 genes. These models were fitted into 101 predictive models within the LOOCV framework. TCGA-HCC dataset was divided into training and validation sets at a 1:1 ratio, and the C-index was calculated for each model in both sets. The model combination of StepCox[both] + plsRcox yielded the highest average C-index (0.751), utilizing 174 features for survival prediction. However, the StepCox[backward] + Enet[alpha = 0.9] strategy, using only 27 features, achieved an average C-index of 0.722 (Figure 8A). Subsequently, an MVI malignant cell-related feature (MRS) model was constructed using the latter strategy.

FIGURE 8

In both the training and validation sets of TCGA-LIHC cohort, the low-risk group exhibited significantly longer OS (Figures 8B, D; OS, p < 0.0001). In the training set, the MRS model determined AUC values of 0.8 (71.72–87.48), 0.81 (73.53–87.82), and 0.85 (76.85–92.72) for 1-, 3-, and 5-year OS, respectively. Due to sample size limitations, the 1-year AUC was not calculated in the validation set, but the AUC values for 3-, 5-, and 7-year OS were 0.78 (63.19–84.47), 0.74 (56.93–93.97), and 0.75 (88.22-NA), respectively (Figures 8C, E). These results demonstrated that the MRS model, which is based on MVI malignant cell subtypes, is an effective predictive tool with satisfactory specificity and sensitivity. When combining the training and validation sets, the low-risk group consistently had significantly longer OS than the high-risk group across all samples (Figure 8F), with AUC values of 0.78 (70.27–86.29), 0.78 (71.64–83.36), and 0.79 (72.82–85.62) for at 1, 3, and 5 years, respectively (Figure 8G). Additionally, the CIBERSORT algorithm was used to analyse the proportion of infiltrating immune cells and the correlation between high- and low-risk samples in TCGA-HCC dataset (Figures 8H, I). There was a significantly greater proportion of infiltrated M0 macrophages in the high-risk group (P < 0.01) and a significantly lower proportion of mast cells in the resting state (P < 0.001).

3.9 MARCKSL1 promotes MVI progression through PTN-Related networks

The 27 signature genes in the MRS model were analysed, which identified MARCKSL1, a gene involved in various cellular processes, including cell migration, cytoarchitectural adjustments, and signalling. The role of MARCKSL1 has not been extensively demonstrated in hepatocellular carcinoma (Zhao et al., 2023; Chen et al., 2021; Yadav et al., 2024). Based on our preliminary findings, we hypothesized that MARCKSL1 plays a significant role in MVI and tumour progression in hepatocellular carcinoma. MARCKSL1 was primarily expressed in the MCs_2, MCs_4, and MCs_5 subtypes, indicating the potential association with MVI and cancer malignancy (Figure 9A). MARCKSL1-positive (MARCKSL1(+)) expression was almost exclusively found in tumours with MVI, while MARCKSL1-negative (MARCKSL1(−)) expression lacked this specificity (Figure 9B). Differential gene analysis between MARCKSL1(+) and MARCKSL1(−) malignant cells was performed, and the identified genes, along with differential ploidy, were used to conduct GSEA. This analysis revealed a significant increase in the number of anchoring junction entries in MARCKSL1(+) cells (NES = 2.64, p value = 0.001) (Figure 9C).

FIGURE 9

To further explore the interactions between malignant MARCKSL1(+) and MARCKSL1(−) cells and immune cells, CellChat tools were used to analyse all receptor‒ligand pairs, which indicated that MARCKSL1(+) cells delivered significant PTN signals to other cells (Figures 9D–F). PTN-related signals are known to promote angiogenesis and extracellular matrix remodelling through MMPs (Papadimitriou et al., 2022; Perez-Pinera et al., 2008; Mentlein, 2007). The present study demonstrated that MARCKSL1(+) cells were involved in other signalling pathways, including LAMC1 and PRSS3 signalling pathways, distinguishing these cells from MARCKSL1(−) cells (Figure 9G). Thus, these findings suggested that MARCKSL1 influences the progression of MVI and cancer malignancy through PTN-related signalling, which is supported by immunohistochemistry results from the Human Protein Atlas showing greater expression of MARCKSL1 in cancerous tissues than in adjacent noncancerous tissues (Figure 9H).

To further explore the therapeutic value of MARCKSL1 in cancer treatment, data from drug sensitivity analyses were collected using the NCI-60 cell line panel and RNA sequencing data from cell lines in the National Cancer Institute (NCI) database. In total, 75 clinical trials and 188 FDA-approved drugs were reviewed to investigate the relationship between MARCKSL1 expression levels and the IC50s of these drugs by calculating Pearson’s correlation coefficients (Figures 10A, B). The abundance of MARCKSL1 was correlated with tumour cell resistance to several drugs, including hypothemycin, vemurafenib, AP-26113, acrichine, erlotinib, and dabrafenib.

FIGURE 10

4 Discussion

The present study demonstrated significant differences in the transcriptomes of malignant cells associated with the presence and absence of MVI, particularly in terms of cell adhesion capacity and enhanced tumour invasion and metastasis. While previous research focused on the multicellular ecosystem and the roles of immune and stromal cells in MVI formation (Li et al., 2024), our study investigated the heterogeneity of MVI-associated malignant cells, revealing key molecular mechanisms. Specifically, we identified the upregulation of cell junction-related pathways (e.g., cell junction, NES = 2.32, P value = 0.001; structural molecule activity, NES = 2.49, P value = 0.001) and key signaling pathways such as VEGF and laminin signaling, which may promote tumour cell invasion along vessel walls (Leong et al., 2022; Shenoy and Lu, 2016). These findings suggest that malignant cells at the MVI stage may alter their cytoskeleton and adhesion patches to more efficiently cross the extracellular matrix and vessel walls. Additionally, we identified a subpopulation of malignant cells (MCs_1) in non-MVI hepatocellular carcinoma that exhibited activation of detoxification-related pathways (e.g., detoxification of inorganic compounds and copper ions), potentially enhancing resistance to chemotherapeutic agents (Lee et al., 2012). Overall, our study complements the multicellular ecosystem perspective by providing insights into the intrinsic mechanisms of malignant cells, offering a more comprehensive understanding of MVI complexity.

According to the inferCNV results, the MCs_4 and MCs_6 subtypes were the two most malignant cell subtypes. Comparison of the pathway changes in the MCs_4 and MCs_6 subtypes in two different biological states, namely, presence of MVI and absence of MVI, indicated that the Myc target pathway was significantly activated in both subtypes in the absence of MVI. Activation of Myc has been reported to promote the metabolic activity and proliferation of endothelial cells by enhancing their glycolysis and mitochondrial function, thereby effectively promoting vascular growth and expansion (Wilhelm et al., 2016). The present findings highlighted the impact of Myc-related networks on MVI. Moreover, significant enhancements in the G2M checkpoint and E2F target pathways were also found in two subgroups of malignant cells in MVI-present tumours compared to those in MVI-negative tumours. The G2/M checkpoint is involved in cell cycle control, whereas E2F is a key element in the regulation of the cell cycle and is associated with DNA replication and cell division. The upregulation of these pathways may be the molecular basis for the rapid proliferation and invasiveness of MVI tumour cells. Additionally, the present findings indicated downregulation of the apoptosis, coagulation, and hypoxia pathways in MVI-present tumours, suggesting that malignant MVI tumour cells may promote their survival and proliferation by inhibiting natural cell death mechanisms, adjusting interactions with the microenvironment, and adapting to hypoxic conditions. These findings offer new insights into how malignant MVI cells adapt and promote tumour growth, invasion, and metastasis by regulating key biological pathways.

In the present study, malignant cells were classified based on MARCKSL1 expression into the MARCKSL1(+) and MARCKSL1(−) cells. In addition to affecting the PTN signalling pathway, MARCKSL1(+) cells exhibited enrichment of GDF15-TGFBR2 in immune cells, a feature that distinguishes them from MARCKSL1(−) cells. TGF-β is known to play a critical role in early embryonic development and adult homeostasis in vivo (Xu et al., 2018). In cancer, overexpression of TGF-β is closely associated with metabolic disorders, dysfunction, epithelial–mesenchymal transition, immune deficits, and cancer progression (Su et al., 2020; Liu et al., 2021). Xu and colleagues reported that the TGF-β-associated pathway, through the regulation of FOXC1, affects tumour EMT, thereby promoting MVI (Xu et al., 2012). These findings suggested that the TGF-β signalling pathway is one of the mechanisms through which MARCKSL1 promotes MVI and the malignant progression of hepatocellular carcinoma, a connection that has not yet been reported and warrants further investigation. Additionally, the interaction of LAMC1 (laminin C1) with integrin subunits (e.g., ITGA9, ITGB1, ITGA6, and ITGA1) involves extracellular matrix interactions and cell adhesion. These interactions are also specific to MARCKSL1(+), which functions as a receptor‒ligand pair and has been reported to be critical for tumour cell migration, particularly in interactions with vascular endothelial cells. These findings suggested that this pathway is also involved in the regulation of MVI by MARCKSL1.

The present study utilized single-cell analysis techniques to investigate the functional switch of malignant cells in the progression of MVI in hepatocellular carcinoma and its impact on the immune microenvironment. For the first time, the present study reported and highlighted the potential impact of the laminin- and VEGF-related pathways on MVI. MARCKSL1, a molecule not yet widely reported, was identified to be associated with PTN-related signalling pathways. High expression of this gene was closely associated with poor patient prognosis, suggesting that it may become an important prognostic marker and therapeutic target.

5 Limitations of the study

The present study had several limitations. First, the single-cell and prognostic data were obtained from public databases and involved a limited number of patient samples; therefore, the present findings need to be validated in larger clinical samples. Second, the present study did not include in vitro or in vivo experimental validation. Future studies will confirm the present results with additional experiments and further investigate the role of MARCKSL1 in the progression of MVI in hepatocellular carcinoma.

References

• 1

  BladtF.RiethmacherD.IsenmannS.AguzziA.BirchmeierC. (1995). Essential role for the c-met receptor in the migration of myogenic precursor cells into the limb bud. Nature376 (6543), 768–771. 10.1038/376768a0

  • CrossRef
  • Google Scholar

• 2

  BussolinoF.Di RenzoM. F.ZicheM.BocchiettoE.OliveroM.NaldiniL.et al (1992). Hepatocyte growth factor is a potent angiogenic factor which stimulates endothelial cell motility and growth. J. Cell Biol.119 (3), 629–641. 10.1083/jcb.119.3.629

  • CrossRef
  • Google Scholar

• 3

  CarmelietP. (2005). VEGF as a key mediator of angiogenesis in cancer. Oncology69 (Suppl. 3), 4–10. 10.1159/000088478

  • CrossRef
  • Google Scholar

• 4

  ChenZ.ZhangW.SelmiC.RidgwayW. M.LeungP. S. C.ZhangF.et al (2021). The myristoylated alanine-rich C-kinase substrates (MARCKS): a membrane-anchored mediator of the cell function. Autoimmun. Rev.20 (11), 102942. 10.1016/j.autrev.2021.102942

  • CrossRef
  • Google Scholar

• 5

  ComoglioP. M.TrusolinoL.BoccaccioC. (2018). Known and novel roles of the MET oncogene in cancer: a coherent approach to targeted therapy. Nat. Rev. Cancer18 (6), 341–358. 10.1038/s41568-018-0002-y

  • CrossRef
  • Google Scholar

• 6

  GouwA. S.BalabaudC.KusanoH.TodoS.IchidaT.KojiroM. (2011). Markers for microvascular invasion in hepatocellular carcinoma: where do we stand?Liver Transpl.17 (Suppl. 2), S72–S80. 10.1002/lt.22368

  • CrossRef
  • Google Scholar

• 7

  GucciniI.RevandkarA.D'AmbrosioM.ColucciM.PasquiniE.MosoleS.et al (2021). Senescence reprogramming by TIMP1 deficiency promotes prostate cancer metastasis. Cancer Cell39 (1), 68–82.e9. 10.1016/j.ccell.2020.10.012

  • CrossRef
  • Google Scholar

• 8

  GuptaG.Al-MalkiW. H.KazmiI.ThangaveluL.GuptaP. K.JhaN. K.et al (2021). The role of HGF/MET in liver cancer. Future Med. Chem.13 (21), 1829–1832. 10.4155/fmc-2021-0128

  • CrossRef
  • Google Scholar

• 9

  HeX.SunH.JiangQ.ChaiY.LiX.WangZ.et al (2021). Hsa-miR-4277 decelerates the metabolism or clearance of sorafenib in HCC cells and enhances the sensitivity of HCC cells to sorafenib by targeting cyp3a4. Front. Oncol.11, 735447. Published 2021 Jul 26. 10.3389/fonc.2021.735447

  • CrossRef
  • Google Scholar

• 10

  HongS. B.ChoiS. H.KimS. Y.ShimJ. H.LeeS. S.ByunJ. H.et al (2021). MRI features for predicting microvascular invasion of hepatocellular carcinoma: a systematic review and meta-analysis. Liver Cancer10 (2), 94–106. 10.1159/000513704

  • CrossRef
  • Google Scholar

• 11

  HuangX.GanG.WangX.XuT.XieW. (2019). The HGF-MET axis coordinates liver cancer metabolism and autophagy for chemotherapeutic resistance. Autophagy15 (7), 1258–1279. 10.1080/15548627.2019.1580105

  • CrossRef
  • Google Scholar

• 12

  HuhC. G.FactorV. M.SánchezA.UchidaK.ConnerE. A.ThorgeirssonS. S. (2004). Hepatocyte growth factor/c-met signaling pathway is required for efficient liver regeneration and repair. Proc. Natl. Acad. Sci. U. S. A.101 (13), 4477–4482. 10.1073/pnas.0306068101

  • CrossRef
  • Google Scholar

• 13

  JustoB. L.JasiulionisM. G. (2021). Characteristics of TIMP1, CD63, and β1-integrin and the functional impact of their interaction in cancer. Int. J. Mol. Sci.22 (17), 9319. Published 2021 Aug 27. 10.3390/ijms22179319

  • CrossRef
  • Google Scholar

• 14

  KuppeC.IbrahimM. M.KranzJ.ZhangX.ZieglerS.Perales-PatónJ.et al (2021). Decoding myofibroblast origins in human kidney fibrosis. Nature589 (7841), 281–286. 10.1038/s41586-020-2941-1

  • CrossRef
  • Google Scholar

• 15

  LeeC.RaffaghelloL.LongoV. D. (2012). Starvation, detoxification, and multidrug resistance in cancer therapy. Drug Resist Updat15 (1-2), 114–122. 10.1016/j.drup.2012.01.004

  • CrossRef
  • Google Scholar

• 16

  LeongS. P.NaxerovaK.KellerL.PantelK.WitteM. (2022). Molecular mechanisms of cancer metastasis via the lymphatic versus the blood vessels. Clin. Exp. Metastasis39 (1), 159–179. 10.1007/s10585-021-10120-z

  • CrossRef
  • Google Scholar

• 17

  LiK.ZhangR.WenF.ZhaoY.MengF.LiQ.et al (2024). Single-cell dissection of the multicellular ecosystem and molecular features underlying microvascular invasion in HCC. Hepatology79 (6), 1293–1309. 10.1097/HEP.0000000000000673

  • CrossRef
  • Google Scholar

• 18

  LiuS.RenJ.Ten DijkeP. (2021). Targeting TGFβ signal transduction for cancer therapy. Signal Transduct. Target Ther.6 (1), 8. Published 2021 Jan 8. 10.1038/s41392-020-00436-9

  • CrossRef
  • Google Scholar

• 19

  MaL.HernandezM. O.ZhaoY.MehtaM.TranB.KellyM.et al (2019). Tumor cell biodiversity drives microenvironmental reprogramming in liver cancer. Cancer Cell36 (4), 418–430.e6. 10.1016/j.ccell.2019.08.007

  • CrossRef
  • Google Scholar

• 20

  MalumbresM.BarbacidM. (2009). Cell cycle, CDKs and cancer: a changing paradigm. Nat. Rev. Cancer9 (3), 153–166. 10.1038/nrc2602

  • CrossRef
  • Google Scholar

• 21

  MarreroJ. A.KulikL. M.SirlinC. B.ZhuA. X.FinnR. S.AbecassisM. M.et al (2018). Diagnosis, staging, and management of hepatocellular carcinoma: 2018 practice guidance by the American association for the study of liver diseases. Hepatology68 (2), 723–750. 10.1002/hep.29913

  • CrossRef
  • Google Scholar

• 22

  MentleinR. (2007). Targeting pleiotropin to treat osteoarthritis. Expert Opin. Ther. Targets11 (7), 861–867. 10.1517/14728222.11.7.861

  • CrossRef
  • Google Scholar

• 23

  ÖzkanA.StolleyD. L.CressmanE. N. K.McMillinM.DeMorrowS.YankeelovT. E.et al (2021). Tumor microenvironment alters chemoresistance of hepatocellular carcinoma through CYP3A4 metabolic activity. Front. Oncol.11, 662135. Published 2021 Jun 28. 10.3389/fonc.2021.662135

  • CrossRef
  • Google Scholar

• 24

  PapadimitriouE.MourkogianniE.NtenekouD.ChristopoulouM.KoutsioumpaM.LamprouM. (2022). On the role of pleiotrophin and its receptors in development and angiogenesis. Int. J. Dev. Biol.66 (1-2-3), 115–124. 10.1387/ijdb.210122ep

  • CrossRef
  • Google Scholar

• 25

  Perez-PineraP.BerensonJ. R.DeuelT. F. (2008). Pleiotrophin, a multifunctional angiogenic factor: mechanisms and pathways in normal and pathological angiogenesis. Curr. Opin. Hematol.15 (3), 210–214. 10.1097/MOH.0b013e3282fdc69e

  • CrossRef
  • Google Scholar

• 26

  RamachandranP.DobieR.Wilson-KanamoriJ. R.DoraE. F.HendersonB. E. P.LuuN. T.et al (2019). Resolving the fibrotic niche of human liver cirrhosis at single-cell level. Nature575 (7783), 512–518. 10.1038/s41586-019-1631-3

  • CrossRef
  • Google Scholar

• 27

  RickJ. W.ChandraA.DalleO. C.NguyenA. T.YagnikG.AghiM. K. (2019). Fibronectin in malignancy: cancer-specific alterations,protumoral effects,and therapeutic implications. Semin. Oncol.46 (3), 284–290. 10.1053/j.seminoncol.2019.08.002

  • CrossRef
  • Google Scholar

• 28

  RumgayH.ArnoldM.FerlayJ.LesiO.CabasagC. J.VignatJ.et al (2022). Global burden of primary liver cancer in 2020 and predictions to 2040. J. Hepatol.77 (6), 1598–1606. 10.1016/j.jhep.2022.08.021

  • CrossRef
  • Google Scholar

• 29

  ShengX.JiY.RenG. P.LuC. L.YunJ. P.ChenL. H.et al (2020). A standardized pathological proposal for evaluating microvascular invasion of hepatocellular carcinoma: a multicenter study by LCPGC. Hepatol. Int.14 (6), 1034–1047. 10.1007/s12072-020-10111-4

  • CrossRef
  • Google Scholar

• 30

  ShenoyA. K.LuJ. (2016). Cancer cells remodel themselves and vasculature to overcome the endothelial barrier. Cancer Lett.380 (2), 534–544. 10.1016/j.canlet.2014.10.031

  • CrossRef
  • Google Scholar

• 31

  SquairJ. W.SkinniderM. A.GautierM.FosterL. J.CourtineG. (2021). Prioritization of cell types responsive to biological perturbations in single-cell data with Augur. Nat. Protoc.16 (8), 3836–3873. 10.1038/s41596-021-00561-x

  • CrossRef
  • Google Scholar

• 32

  SuJ.MorganiS. M.DavidC. J.WangQ.ErE. E.HuangY. H.et al (2020). TGF-β orchestrates fibrogenic and developmental EMTs via the RAS effector RREB1. Nature577 (7791), 566–571. 10.1038/s41586-019-1897-5

  • CrossRef
  • Google Scholar

• 33

  SumieS.KuromatsuR.OkudaK.AndoE.TakataA.FukushimaN.et al (2008). Microvascular invasion in patients with hepatocellular carcinoma and its predictable clinicopathological factors. Ann. Surg. Oncol.15 (5), 1375–1382. 10.1245/s10434-008-9846-9

  • CrossRef
  • Google Scholar

• 34

  TanJ.FanW.LiuT.ZhuB.LiuY.WangS.et al (2023). TREM2+ macrophages suppress CD8+ T-cell infiltration after transarterial chemoembolisation in hepatocellular carcinoma. J. Hepatol.79 (1), 126–140. 10.1016/j.jhep.2023.02.032

  • CrossRef
  • Google Scholar

• 35

  WeiF.WuY.TangL.HeY.ShiL.XiongF.et al (2018). BPIFB1 (LPLUNC1) inhibits migration and invasion of nasopharyngeal carcinoma by interacting with VTN and VIM. Br. J. Cancer118 (2), 233–247. 10.1038/bjc.2017.385

  • CrossRef
  • Google Scholar

• 36

  WilhelmK.HappelK.EelenG.SchoorsS.OellerichM. F.LimR.et al (2016). FOXO1 couples metabolic activity and growth state in the vascular endothelium. Nature529 (7585), 216–220. 10.1038/nature16498

  • CrossRef
  • Google Scholar

• 37

  XiongY.CaoP.LeiX.TangW.DingC.QiS.et al (2022). Accurate prediction of microvascular invasion occurrence and effective prognostic estimation for patients with hepatocellular carcinoma after radical surgical treatment. World J. Surg. Oncol.20 (1), 328. Published 2022 Sep 30. 10.1186/s12957-022-02792-y

  • CrossRef
  • Google Scholar

• 38

  XuX.ZhangH. L.LiuQ. P.SunS. W.ZhangJ.ZhuF. P.et al (2019). Radiomic analysis of contrast-enhanced CT predicts microvascular invasion and outcome in hepatocellular carcinoma. J. Hepatol.70 (6), 1133–1144. 10.1016/j.jhep.2019.02.023

  • CrossRef
  • Google Scholar

• 39

  XuX.ZhengL.YuanQ.ZhenG.CraneJ. L.ZhouX.et al (2018). Transforming growth factor-β in stem cells and tissue homeostasis. Bone Res.6, 2. 10.1038/s41413-017-0005-4

  • CrossRef
  • Google Scholar

• 40

  XuZ. Y.DingS. M.ZhouL.XieH. Y.ChenK. J.ZhangW.et al (2012). FOXC1 contributes to microvascular invasion in primary hepatocellular carcinoma via regulating epithelial-mesenchymal transition. Int. J. Biol. Sci.8, 1130–1141. 10.7150/ijbs.4769

  • CrossRef
  • Google Scholar

• 41

  YadavV.JenaM. K.ParasharG.ParasharN. C.JoshiH.RamniwasS.et al (2024). Emerging role of microRNAs as regulators of protein kinase C substrate MARCKS and MARCKSL1 in cancer. Exp. Cell Res.434 (2), 113891. 10.1016/j.yexcr.2023.113891

  • CrossRef
  • Google Scholar

• 42

  ZhangQ.HeY.LuoN.PatelS. J.HanY.GaoR.et al (2019). Landscape and dynamics of single immune cells in hepatocellular carcinoma. Cell179 (4), 829–845.e20. 10.1016/j.cell.2019.10.003

  • CrossRef
  • Google Scholar

• 43

  ZhaoY.XieX.TianL.LiuF.SunY.LuH.et al (2023). MARCKSL1 interacted with F-actin to promote esophageal squamous cell carcinoma mobility by modulating the formation of invadopodia. Cancer Med.12 (3), 3299–3312. 10.1002/cam4.5079

  • CrossRef
  • Google Scholar

• 44

  ZhengC.ZhengL.YooJ. K.GuoH.ZhangY.GuoX.et al (2017). Landscape of infiltrating T cells in liver cancer revealed by single-cell sequencing. Cell169 (7), 1342–1356.e16. 10.1016/j.cell.2017.05.035

  • CrossRef
  • Google Scholar

• 45

  ZhuW.LiW.YangG.FuC.JiangG.HuQ. (2015). Vitronectin silencing inhibits hepatocellular carcinoma in vitro and in vivo. Future Oncol.11 (2), 251–258. 10.2217/fon.14.202

  • CrossRef
  • Google Scholar