Abstract
Alpha-fetoprotein (AFP) is a classic biomarker for hepatocellular
carcinoma (HCC). AFP-positive HCC (AFP^+ HCC) has been intensively
investigated; however, the genomic, transcriptomic and
microenvironmental characteristics of AFP-negative HCC (AFP^− HCC)
remain to be deciphered. Here we show that tumors display mild
differences in genetic alterations between AFP^− HCC and AFP^+ HCC
patients, while AFP^− HCC exhibits hyperactive arachidonic acid
metabolism. Furthermore, the transcription activity of androgen
receptor (AR) is significantly increased in AFP^− HCC and plays a
positive regulatory role in arachidonic acid metabolism and its
metabolite 11,12-epoxyeicosatrienoic acid (11,12-EET). The
tumor-derived 11,12-EET exhibits high affinity for EGFR that promotes
the migration and tube formation of endothelial cells in vitro.
Combination of lenvatinib and bicalutamide (an AR antagonist) enhances
the therapeutic efficacy for AFP^− HCC. Overall, we uncover the
AR-mediated hyperactive arachidonic acid metabolism in AFP^− HCC, and
reveal AR-11,12-EET-EGFR axis-induced angiogenesis, providing a
promising strategy of combined AR antagonist with lenvatinib for AFP^−
HCC treatment.
Subject terms: Tumour angiogenesis, Cancer metabolism, Cellular
signalling networks, Tumour biomarkers
__________________________________________________________________
Alpha-fetoprotein is a known biomarker for hepatocellular carcinoma.
Here, the authors analyse differences between AFP positive and negative
patients to identify alterations in androgen receptor signalling.
Introduction
Liver cancer is the third leading cause of cancer-related deaths, with
an overall low 5-year survival rate of 20.8%^[66]1,[67]2.
Hepatocellular carcinoma (HCC) accounts for 80–90% of liver cancer
worldwide^[68]3. Alpha-fetoprotein (AFP) is the most widely used
biomarker for HCC in clinical practice, primarily aiding in HCC
screening, diagnosis, prognosis prediction, and treatment efficacy
monitoring^[69]4. However, our knowledge on genetic features, molecular
characteristics and tumor ecosystem profiling of AFP-related HCCs is
very limited.
Extensive evidence has established that AFP exerts dual oncogenic
effects by enhancing tumor cell proliferation and suppressing apoptosis
in HCC^[70]5. In addition, AFP also plays pivotal roles in regulating
immune responses, for instance, it can induce natural killer (NK) cell
apoptosis, promote Treg differentiation and inhibit dendritic cell (DC)
maturation, ultimately repressing CD8^+ T cells and resulting in an
immunosuppressive microenvironment^[71]6. Recently, comprehensive
multi-omics approaches integrating whole-exome sequencing (WES), bulk
RNA-sequencing (RNA-seq) and single-cell RNA-sequencing (scRNA-seq)
have been adopted to reveal the multi-dimensional characteristics of
tumor immune ecosystem in AFP-positive HCC (AFP^+ HCC, serum AFP
level ≥ 20 ng/mL at diagnosis). The findings showed that
tumor-associated macrophages (TAMs) suppressed CD8^+ T cell function
through the SPP1-CD44 axis in the tumor microenvironment (TME) of AFP^+
HCC, thus providing a promising target for achieving a more favorable
efficacy in AFP^+ HCC treatment^[72]7. Therefore, AFP^+ HCC has been
widely explored in previous studies^[73]6–[74]10.
However, approximately 30–40% of HCC showed negative AFP serum level
(AFP^− HCC, serum AFP level < 20 ng/mL at diagnosis)^[75]11,[76]12. Of
note, several recent global clinical trials have indicated that AFP^−
HCC patients showed moderate treatment response to the first-line drug
lenvatinib or its analogues^[77]13–[78]15. Nevertheless, due to the few
studies on AFP^− HCC, the underlying mechanism remains largely unknown
without available interventions. Therefore, we aimed to comprehensively
decipher the genomic, transcriptomic and microenvironmental
characteristics of AFP^− HCC that helped guide the precise treatment
strategy for this subtype of HCC.
In the present study, we performed an integrated omics analyses
containing WES on 314 samples, scRNA-seq on 25 samples, bulk RNA-seq on
434 samples, single-nucleus ATAC-sequencing (snATAC-seq) on 4 samples,
methylation profiling by Reduced Representation Bisulfite Sequencing
(RRBS) on 26 samples and lipidomics on 34 samples to systematically
compare the genetic, transcriptional and microenvironmental differences
between AFP^+ HCC and AFP^- HCC. While no significant genetic
differences were observed between the two subtypes, AFP^− HCC exhibited
pronounced activation of arachidonic acid metabolism mediated by the
androgen receptor (AR), which induced substantial secretion of the
downstream metabolite 11,12-epoxyeicosatrienoic acid (11,12-EET) that
directly activated epidermal growth factor receptor (EGFR) on
endothelial cells (ECs), thereby leading to aberrant angiogenesis.
Thus, we have revealed an angiogenesis regulatory pathway controlled by
AR-11,12-EET-EGFR axis in AFP^− HCC, which cannot be effectively
inhibited by lenvatinib. This finding provided the scientific rationale
for combination treatment of lenvatinib and inhibition of AR-mediated
angiogenesis for AFP^− HCC.
Results
Profiling of cell diversity in the AFP^− HCC and AFP^+ HCC
To comprehensively elucidate the distinct features between AFP^− HCC
and AFP^+ HCC, we conducted integrated omics analyses including
whole-exome sequencing on 314 samples (157 tumor and 157 non-tumor
adjacent tissues) from 157 patients, scRNA-seq on 25 samples (17 tumor
and 8 non-tumor adjacent tissues) from 17 patients, bulk RNA-seq on 434
samples (434 tumor tissues) from 434 patients, snATAC-seq on 4 samples
(4 tumor tissues) from 4 patients, methylation sequencing on 26 samples
(26 tumor tissues) from 26 patients and lipidomics on 34 samples (14
tumor tissues and 20 blood plasma specimens) from 34 patients
(Fig. [79]1a and Supplementary Fig. [80]1a). Clinical information of
these patients was provided in Supplementary Data [81]1, with no
significant differences in HBsAg, antiviral treatment, liver cirrhosis
status, tumor number; modestly significant bias in gender, age and
tumor size (p < 0.05); and significantly lower microvascular invasion
(MVI), Barcelona Clinic Liver Cancer (BCLC) stages and
Edmondson-Steiner (E-S) grades in AFP^− HCC versus AFP^+ HCC
(p < 0.001).
Fig. 1. Profiling of genetic alteration and cell diversity in HCC.
[82]Fig. 1
[83]Open in a new tab
a Schematic overview of the research. The numbers of patients and
multiple experimental validation methods are provided. The lipidomics
datasets (n = 34) comprise two subsets: tumor tissue lipidomics
(n = 14) and blood plasma lipidomics (n = 20). b Genetic profiles and
associated clinical features of all 157 HCC patients with WES data. The
top panel shows mutation numbers followed by clinicopathological
features. The middle panel exhibits significantly mutated genes
calculated through MutSigCV. The bottom panel shows the CNV profiles of
top chromosomal lesions with the most significant q values. Different
alteration types and clinical features are signified with different
color codes. Significance of these variables between AFP^− HCC and
AFP^+ HCC was calculated by chi-square test or Fisher’s exact test,
which was provided in Supplementary Data [84]1. c The UMAP
visualization of 47 cell types from scRNA-seq data. d Scatter plot
showing positive correlation between serum AFP level and transcriptomic
AFP expression in scRNA-seq (left panel) and bulk RNA-seq (right panel)
data of FAH-SYSU cohort. Transcriptomic AFP expression in scRNA-seq
data was calculated as the average expression across tumor cells per
sample. R indicates the correlation coefficient calculated by Spearman
correlation test. The number of dots indicates the number of patients.
e Representative images of IHC staining in FFPE tissues (left panel)
and protein expression of AFP in AFP^− HCC and AFP^+ HCC (right panel).
The number of dots indicates the number of patients. Scale bars, 50 μm.
Data were shown as mean ± SD. **p < 0.01 by two-sided t-test. f Scatter
plot showing positive correlation between proteinic AFP expression and
available transcriptomic AFP expression. R indicates the correlation
coefficient calculated by Spearman correlation test. The number of dots
indicates the number of patients. Source data are provided as a Source
Data file.
To decipher the genomic profiles of AFP^− HCC, we performed WES on 157
tumors, including 44 AFP^− HCC and 113 AFP^+ HCC, and their paired
non-tumor adjacent tissues (NATs). Among them, 18,441 non-synonymous
somatic mutations (range 0–636, average 117.46 per sample) and 2776
indels (range 0–420, average 17.68 per sample) across all samples were
characterized. Consistent with previous findings^[85]16–[86]18, nine
recurrently mutated genes, including TP53, CTNNB1, AXIN1, ARID1A, and
TSC2 were identified (Fig. [87]1b), which showed no significance
between AFP^– HCC and AFP^+ HCC (Supplementary Data [88]2). By
identifying copy number variation (CNV) events and calculating GISTIC
scores, we observed mild differences in amplification and deletion
events between AFP^− HCC and AFP^+ HCC, with significant events
occupying 4% and 3%, respectively (Supplementary Fig. [89]1b, c). The
frequently amplificated and deleted events were listed in Supplementary
Data [90]3. In addition, no significant differences of immune editing
score, clonal expansion of recurrently mutated genes and COSMIC
mutation signatures were observed between the two subtypes
(Supplementary Fig. [91]1d–g). Therefore, these results indicated that
AFP^– HCC and AFP^+ HCC patients displayed similar genomic alterations.
Next, we aimed to decipher the transcriptomic differences between the
two subtypes by performing scRNA-seq on 6 AFP^− HCC and 11 AFP^+ HCC
tumor samples, with 8 paired NATs. The detailed clinical information of
these patients was provided in Supplementary Data [92]4. The baseline
of patients with scRNA-seq data was matched between AFP status groups
(AFP^− HCC vs. AFP^+ HCC). Notably, we also analyzed the paired WES
data from samples with matched scRNA-seq profiles (excluding one case
due to tissue unavailability). No significant differences were observed
between the two HCC subtypes in terms of gene mutations, CNVs, immune
editing scores, clonal expansion of recurrently mutated genes, or
COSMIC mutation signatures (Supplementary Fig. [93]2a–e and
Supplementary Data [94]5). In scRNA-seq data, total 210,545 single
cells were collected and were classified into 47 cell types based on
gene expression panorama (Fig. [95]1c, Supplementary Fig. [96]3a–d, and
Supplementary Data [97]6). Subsequently, CNVs were inferred to validate
the definition of tumor cells through inferCNV (Supplementary
Fig. [98]3e). We then further determined that the transcriptome level
of AFP in tumors was positively correlated with its serum levels in
both FAH-SYSU and TCGA-LIHC cohorts (Fig. [99]1d and Supplementary
Fig. [100]3f), and exhibited high expression in AFP^+ HCC tumor
samples, with no detectable expression in tissues obtained from healthy
liver, NATs of AFP^+ HCC, AFP^– HCC tumor and NATs of AFP^– HCC
(Supplementary Fig. [101]3g). Similar results were observed in staining
of FFPE specimen (Fig. [102]1e, f), suggesting that transcriptomic
expression of AFP could represent the serum AFP level and define the
two HCC subtypes.
AFP^– HCC tumors presented hyperactive arachidonic acid metabolism
Next, via differential gene expression analysis based on scRNA-seq
data, significantly altered genes (550 upregulated and 876
downregulated genes) in tumor cells from AFP^− HCC were identified
compared to AFP^+ HCC. The metabolism-related genes, such as CYP2A6,
FGGY, and CYP2C8, were significantly upregulated in AFP^− HCC, while
proliferation-related genes such as IGF2, EPCAM, and ARID3A, were
significantly upregulated in AFP^+ HCC (Fig. [103]2a). Consistently,
Gene Set Enrichment Analysis (GSEA) suggested that the AFP^– HCC tumors
enriched metabolism and cell adhesion-related pathways, such as
xenobiotic metabolic process, constitutive androstane receptor pathway,
arachidonic acid metabolism and extracellular matrix (ECM) receptor
interaction, a critical molecular feature of early-stage HCC^[104]19,
whereas the AFP^+ HCC tumors upregulated proliferation and
carcinogenesis-associated pathways, such as Myc targets, regulation of
cell division and non-canonical WNT signaling pathways (Fig. [105]2b).
These results were also supported by previous reports^[106]8. To gain
additional insights into the molecular differences of the two HCC
subtypes, we defined the 434 HCC patients by Hoshida’s classification
(Hoshida’s S1, S2, and S3) based on bulk RNA-seq data using the Nearest
Template Prediction (NTP) algorithm^[107]20,[108]21. S1 reflected
abnormal activation of the WNT and TGF-β signaling pathways; S2 was
characterized by proliferation features and activation of MYC and AKT
signaling, and exhibited the highest AFP expression among the three
subtypes; S3 showed well differentiation of tumor cells, with the
lowest AFP expression among the three Hoshida’s subtypes. The results
showed that AFP^− HCC presented similar features to the S3 subtype,
whereas AFP^+ HCC exhibited Hoshida’s S1/2 subtypes characteristics
(Supplementary Fig. [109]4a). Moreover, the differentially expressed
genes (DEGs) and GSEA results from the bulk RNA-seq data also verified
the hyperactive metabolic features in AFP^− HCC and the proliferative
properties in AFP^+ HCC (Supplementary Fig. [110]4b, c).
Fig. 2. Transcriptomic differences of tumor cells between AFP^– HCC and AFP^+
HCC.
[111]Fig. 2
[112]Open in a new tab
a Differentially expressed genes (DEGs) calculated from scRNA-seq data
between AFP^− HCC and AFP^+ HCC. b Pathways enriched in tumor cells in
AFP^– HCC and AFP^+ HCC from scRNA-seq data by GSEA. The p-values,
calculated by permutation test, of all listed pathways were less than
0.05. c The volcano plot from scRNA-seq data manifested differentially
metabolic activity of tumor cells in AFP^– HCC when compared to AFP^+
HCC. The color and dot size correspond to the log[2]FC in metabolic
activity between AFP^− HCC and AFP^+ HCC. d Correlation analysis
between metabolic pathways and transcriptomic AFP expression in
scRNA-seq data by Spearman correlation test. e Metabolic pathway
diagram of arachidonic acid metabolism and gene signature scores of
three branches in arachidonic acid metabolism. Box plot depicts the
median and interquartile range, and the lower and upper hinges denote
the 25–75% interquartile range (IQR), with whiskers extending up to a
maximum of 1.5 times IQR. f Scatter plot showing negative correlation
between transcriptomic AFP expression and CYP gene signature in
scRNA-seq data by Spearman correlation test. g The heatmap of genes
involved in arachidonic acid metabolism in AFP^– HCC and AFP^+ HCC,
colored by average expression value. h Relative expression of
CYP2C8/CYP1B1 detected by RT-qPCR in tumor tissues from 15 AFP^− HCC
and 15 AFP^+ HCC patients. Data were shown as mean ± SD. i The volcano
plot shows metabolites detected by targeted lipidomics. The
abbreviations ARA and DTA stand for arachidonic acid and
docosatetraenoic acid, respectively. j The expression levels of AFP
were detected by RT-qPCR and western blot in HCC cell lines. k ELISA
analysis revealed elevated levels of 11,12-EET in the CM derived from
AFP^low cell lines compared to AFP^high cell lines. Data were shown as
mean ± SD. *p < 0.05, **p < 0.01, ***p < 0.001 by Wilcoxon rank sum
test in (e, g, i), two-sided t-test in (h) and one-way ANOVA with
Tukey’s multiple comparisons test (j, k). Source data are provided as a
Source Data file.
To further delineate the metabolic features in AFP^− HCC, we employed
the scMetabolism algorithm^[113]22 to tumor cells from scRNA-seq data.
This analysis revealed a significant upregulation of the arachidonic
acid metabolism pathway in AFP^– HCC (Fig. [114]2c and Supplementary
Data [115]7), which exhibited a strong negative correlation with AFP
expression (Fig. [116]2d). These findings were further confirmed by the
GSEA results from scRNA-seq and bulk RNA-seq data, respectively
(Fig. [117]2b and Supplementary Fig. [118]4c). Compared to all other
cell types, tumor cells displayed the highest scores for arachidonic
acid metabolism and the strongest difference between AFP^− HCC and
AFP^+ HCC (Supplementary Fig. [119]4d). Based on the well-established
main branches of arachidonic acid metabolism, comprising the
cyclooxygenase (COX), lipoxygenase (LOX), and cytochrome P450 (CYP)
pathways^[120]23, we conducted differential expression profiling and
correlation analysis. These analyses demonstrated that CYP pathway,
which generates bioactive epoxyeicosatrienoic acids (EETs) and
hydroxyeicosatetraenoic acid (HETEs) to regulate cell proliferation,
survival, invasion, metastasis, and angiogenesis^[121]23, was the
significantly elevated signature in AFP^− HCC (Fig. [122]2e). Notably,
the CYP pathway was the only branch exhibiting a strong negative
correlation with AFP expression (Fig. [123]2f and Supplementary
Fig. [124]4e, f). This observation indicated that tumor cells in AFP^−
HCC were more likely to reshape their metabolic patterns through the
CYP pathway.
Elevation of 11,12-EET in AFP^– HCC tumor cells versus AFP^+ HCC
In agreement with the elevated activity of CYP pathway in AFP^− HCC,
the CYP family genes, including CYP2C8 and CYP1B1, were upregulated in
AFP^− HCC patients (Fig. [125]2g), and were confirmed by RT-qPCR assays
in tumor samples (Fig. [126]2h). We then further sought to determine
the key metabolite participating in the arachidonic acid metabolism in
AFP^− HCC. Through targeted lipidomics, we identified five
significantly upregulated lipid-related metabolites in AFP^– HCC
compared to that in AFP^+ HCC (Fig. [127]2i). We decided to study
11,12-EET because: (1) the top two metabolites (ARA and DTA) served as
the common substrates for COXs, LOXs, and CYPs as documented in prior
studies^[128]23–[129]26, whereas only the CYP pathway showed
significant activation in AFP^− HCC; (2) 11,12-EET was the specific and
direct metabolic product of the CYP pathway; (3) the upregulated CYP2C8
and CYP1B1 in AFP^− HCC were the key metabolic enzymes responsible for
11,12-EET generation^[130]23,[131]27; (4) elevation of 11,12-EET
inversely correlated with serum AFP levels (Supplementary
Fig. [132]4g); (5) moreover, 11,12-EET was reported to exhibit diverse
biological functions^[133]23 that might play an important role in tumor
biology.
To validate the expression relation between AFP and CYPs/11,12-EET, AFP
expression level in five HCC cell lines was further measured, where
PLC5, Hep3B and Huh7 were defined as AFP^high HCC cell lines; MHCC97H
and HCCLM3 belonged to AFP^low HCC cell lines (Fig. [134]2j). In line
with DEGs and lipidomics analysis, AFP^low cell lines displayed
elevated CYP2C8 and CYP1B1 expression (Supplementary Fig. [135]4h, i)
and higher 11,12-EET levels (Fig. [136]2k) in comparison to the
AFP^high cell lines, solidifying 11,12-EET as the pivotal metabolite in
AFP^– HCC.
Upregulated AR transcriptional activity in AFP^– HCC tumor cells
Transcription factors (TF) are key molecules that control gene
expression, and their activities determine the cell functions and
response to intricacies of environmental disturbances^[137]28. Thus, to
elucidate the potential regulatory mechanisms for the upregulation of
arachidonic acid metabolism in AFP^− HCC, we assessed the TF activity
across tumor cells in scRNA-seq data using the single-cell regulatory
network inference and clustering algorithm (SCENIC)^[138]29. The SCENIC
results showed that the TF activities of MSC, KLF9, EBF1, IRF6, and AR
were significantly upregulated in AFP^− HCC tumor cells (Fig. [139]3a
and Supplementary Data [140]8). Notably, AR transcriptional activity
showed the most substantial negative correlation with AFP expression
and positive association with CYP-related gene expression (Fig. [141]3b
and Supplementary Data [142]9). In agreement with this finding, AR was
further identified as the second most significant TF in AFP^– HCC based
on the bulk RNA-seq data (Supplementary Fig. [143]5a, b).
Fig. 3. Upregulation of AR transcription activity in AFP^− HCC tumor cells.
[144]Fig. 3
[145]Open in a new tab
a SCENIC analysis in scRNA-seq data shows TF activities in AFP^− HCC
and AFP^+ HCC tumor cells. b Spearman correlation analysis between
activity of TFs listed in (a), and AFP expression (colored in blue) as
well as CYP genes signature (colored in red) in tumor cells. c
Detection of AR expression using IHC in AFP^– HCC (n = 10) and AFP^+
HCC (n = 10). Scale bars, 50 μm. d Correlation between AFP and AR
expression by Spearman correlation test (n = 20). e Relative expression
of AR detected by RT-qPCR (AFP^− HCC vs. AFP^+ HCC = 15:15). f AR
expression in AFP^– HCC compared to AFP^+ HCC, as determined by bulk
RNA-seq data from FAH-SYSU (AFP^– HCC vs. AFP^+ HCC = 161:273) and
TCGA-LIHC cohort (AFP^– HCC vs. AFP^+ HCC = 142:127). g Concentrations
of testosterone and 5α-DHT detected by lipidomics using blood plasma
from 10 AFP^−/AFP^+ HCC. h Relative expression of AR detected by
RT-qPCR and western blot in HCC cells. i Relative concentration of
11,12-EET detected by ELISA assays in CM from HCC cells. j Relative
expression of CYP2C8 and CYP1B1 with or without 5α-DHT stimulation. k
Validation of the AR knockdown (shAR) in MHCC97H and HCCLM3 using
western blot. l 11,12-EET concentrations detected by ELISA assays in CM
from MHCC97H and HCCLM3 cells (shNC) compared to their shAR
counterparts. m Detection of 11,12-EET levels in AR-overexpressing
Hepa1-6 cells (oeAr) versus control cells (oeNC). Validation of AR
binding with CYP2C8 and CYP1B1 promoters by ChIP-PCR (n) and ChIP-qPCR
(o) in MHCC97H and HCCLM3. p Normalized snATAC-seq tracks of AFP, AR,
CYP2C8 and CYP1B1. Three biological replicates were employed (h–j, l,
m, o). Data were shown as mean ± SD (c, g–j, l, m, o). *p < 0.05,
**p < 0.01, ***p < 0.001 and ns stands for no significance by Wilcoxon
rank sum test (f, g), two-sided t-test in (c, e, i, j, m) and one-way
ANOVA with Tukey’s multiple comparisons test (h, l, o). Source data are
provided as a Source Data file.
To further examine the TF activity of AR in AFP^− HCC, we analyzed and
intersected the predicted target genes of AR from scRNA-seq and bulk
RNA-seq data, as shown in Supplementary Data [146]10 and [147]11,
respectively. We observed that these potential AR-target genes
exhibited lower methylation levels in AFP^− HCC, as revealed by RRBS
and TCGA-LIHC 450 K methylation sequencing, compared to those in AFP^+
HCC (Supplementary Fig. [148]5c, d). The reduced methylation of these
target genes indicated the function of AR acting as a TF in AFP^− HCC.
Moreover, AR nuclear expression was observed significantly higher in
AFP^− HCC tumors than that in AFP^+ HCC tumors, and was significantly
negatively correlated with AFP expression (Fig. [149]3c, d). Taken
together, these findings suggested that AR might be a core TF of AFP^−
HCC tumor cells.
As a male hormone receptor and major oncogenic driver in prostate
cancer^[150]30, AR was also highly expressed in other cancer types,
including breast cancer and HCC (Supplementary Fig. [151]6a). However,
the role of AR in HCC is controversial. Some studies have reported the
negative prognostic effect of AR in HCC patients^[152]31, while others
have demonstrated the anti-tumor functions of AR through inhibiting p38
and NF-κB pathways in HCC cells^[153]32. Therefore, we sought to
determine the regulatory functions of AR in AFP^− HCC. We observed a
tightly positive correlation of AR with CYP pathway genes in HCC
(Supplementary Fig. [154]6b). Interestingly, in addition to the
upregulation of TF activity, the transcriptomic levels of AR were also
upregulated in AFP^− HCC tumors (Fig. [155]3e, f). Generally, HCC is a
male-dominant cancer, and AR is highly expressed in male HCC
patients^[156]30,[157]33. Given the higher expression and activities of
AR in AFP^− HCC, we then further quantitatively analyzed plasma
concentrations of testosterone and 5α-dihydrotestosterone (5α-DHT), the
primary endogenous ligands of AR in HCC patients^[158]34. Notably,
these androgen levels exhibited significant elevation in AFP^− HCC
compared to AFP^+ HCC (Fig. [159]3g), further suggesting that AR might
be significantly activated in AFP^− HCC patients.
AR activates the CYP pathway to promote 11,12-EET production in AFP^− HCC
tumor cells
To explore whether AR participated in the activation of the CYP
pathway, we identified PLC5, Hep3B, Huh7 and Hepa1-6 as AFP^highAR^low
cell lines, while MHCC97H, HCCLM3 and RIL-175 were labeled as
AFP^lowAR^high cell lines by RT-qPCR and western blot assays
(Fig. [160]3h and Supplementary Fig. [161]6c, d). Next, AFP^lowAR^high,
but not the AFP^highAR^low cell lines, boosted the production of
11,12-EET and the expression of CYP2C8 and CYP1B1 in the presence of
5α-DHT stimulation (Fig. [162]3i, j). In consistent with this finding,
depletion of AR markedly repressed 11,12-EET production in
AFP^lowAR^high cells (Fig. [163]3k, l); whereas overexpression of AR
significantly promoted the production of 11,12-EET in AFP^highAR^low
cells (Fig. [164]3m and Supplementary Fig. [165]6e). These data
together indicated the critical role of AR in regulation of the CYP
pathway activity in AFP^− HCC.
Next, we aimed to elucidate whether AR could directly regulate the
expression of CYP2C8 and CYP1B1. By using JASPAR database^[166]35, we
predicted AR binding sites on the CYP2C8 and CYP1B1 promoters and
generated the primers. ChIP assays were then employed, and the data
confirmed that AR could bind with CYP2C8 and CYP1B1 promoters
(Fig. [167]3n). Besides, ChIP-qPCR results showed that AR depletion
could significantly attenuate its TF activity on CYP2C8 and CYP1B1
(Fig. [168]3o). Furthermore, AFP^− HCC tumors exhibited higher
chromatin accessibility signals of AR, CYP2C8, and CYP1B1, whereas
displayed lower accessibility signals of AFP in comparison to AFP^+ HCC
(Fig. [169]3p and Supplementary Data [170]12). Taken together, these
observations indicated that AR positively regulated the CYP pathway
activity in AFP^− HCC and promoted the production of 11,12-EET by
enhancing the transcription of CYP2C8 and CYP1B1.
Enrichment of tip-like endothelial cells in AFP^− HCC microenvironment
Accumulating evidence has suggested that metabolic alterations in tumor
cells modulated the composition and function of surrounding cells,
thereby reshaping the TME^[171]36,[172]37. To study the tumor
microenvironment of AFP^– HCC, we first compared the cellular
composition between AFP^− HCC and AFP^+ HCC by calculating the
percentage of each cell subtype relative to the total cell population
in each sample using the scRNA-seq data. The analysis revealed a
significant enrichment of the endothelial cell subtype EC_PLPP1 in
AFP^– HCC (Fig. [173]4a), which was further supported by the results
analyzed by the Milo algorithm^[174]38, a statistical framework that
assesses differential abundance by assigning cells to partially
overlapping neighborhoods on a k-nearest neighbor graph (Supplementary
Fig. [175]7a, b). Furthermore, the results of cellular composition were
also validated by xCell deconvolution based on bulk RNA-seq data
(Fig. [176]4b) and IHC staining of CD31, a specific biomarker of
endothelial cell^[177]39, in HCC tumor samples (Fig. [178]4c). Notably,
the CD31 level exhibited a negative correlation with AFP expression as
analyzed from the IHC data (Supplementary Fig. [179]7c). Altogether,
these results implied that the angiogenesis ability might be a distinct
signature of AFP^− HCC compared to AFP^+ HCC.
Fig. 4. Enrichment of endothelial cells in AFP^– HCC.
[180]Fig. 4
[181]Open in a new tab
a Cell proportion analysis of each cell cluster between AFP^− HCC and
AFP^+ HCC patients, as shown by boxplots and pie charts. P values are
displayed only for cell types with significant differences between the
two groups (AFP^− HCC vs. AFP^+ HCC = 6:11). b Endothelial cell
proportion deconvoluted by xCell using bulk RNA-seq from FAH-SYSU
cohort (AFP^− HCC vs. AFP^+ HCC = 161:273) and TCGA-LIHC cohort (AFP^−
HCC vs. AFP^+ HCC = 142:127), which showed higher infiltration of
endothelial cells in AFP^− HCC. c Representative images of IHC staining
in FFPE tissues (left panel) and microvascular density quantification
(right panel) in AFP^− HCC and AFP^+ HCC. The number of dots indicates
the number of patients (AFP^− HCC vs. AFP^+ HCC = 10:10). Data were
shown as mean ± SD. Scale bars, 50 μm. d The UMAP visualization of four
subtypes of endothelial cells (upper panel) and corresponding tissue
distribution (lower panel). e Marker gene expression for endothelial
cells. f Pathway enrichment analysis for endothelial cells using GO-BP
gene sets. g Functional scoring of endothelial cell phenotypes based on
tip-like and stalk-like gene signatures. *p < 0.05, **p < 0.01,
***p < 0.001, and ns stands for not significant by Wilcoxon rank sum
test in (a, b), two-sided t-test in (c) and Wilcoxon rank sum test with
FDR correction (g). Source data are provided as a Source Data file.
We thus further explored the distribution and function heterogeneity of
endothelial cells between AFP^– HCC and AFP^+ HCC using the scRNA-seq
data (Fig. [182]4d and Supplementary Fig. [183]7d). The EC_FTL with
high expression of FTH1 and FTL, was associated with metabolic pathways
such as cell respiration, oxidative phosphorylation, and iron ion
transport^[184]40 (Fig. [185]4e, f). EC_CLEC4G, with high expression of
canonical liver sinusoidal endothelial cell (LSEC) markers like CLEC4G,
CD14 and CLEC4M, showed preferential enrichment in NATs of AFP^+ HCC
and presented active pathways related to wound healing, regulation of
immune effector process, and receptor-mediated endocytosis
(Fig. [186]4e, f and Supplementary Fig. [187]7d). EC_PLPP1 with
expression of IGFBP3, CXCR4 and PLPP1 displayed preferential enrichment
in AFP^– HCC tumor tissues, and upregulated pathways related to
ameboidal-type cell migration, leukocyte migration, and endothelial
development (Fig. [188]4e, f and Supplementary Fig. [189]7d), likely to
be the tumor-associated endothelial cell. EC_ACKR1 with high expression
of ACKR1 and VWF was featured by the enriched pathways related to
regulation of cell-cell adhesion, positive regulation of T cell
activation, and T-helper cell differentiation (Fig. [190]4e, f), likely
to be the venous-related endothelial cell^[191]41.
Endothelial cells commonly exhibited two functional phenotypes in tumor
angiogenesis. The tip-like endothelial cells were characterized by cell
budding, migration, and formation of filopodia, while the stalk-like
endothelial cell presented proliferation, adhesion, and tight
connection abilities to ensure the stability of new
sprouts^[192]42,[193]43. Here, we observed that EC_PLPP1, enriched in
AFP^− HCC, exhibited a high tip-like signature, and significantly
upregulated tip-like endothelial cell-related pathways (Fig. [194]4g
and Supplementary Fig. [195]7e), indicating the enhancing angiogenesis
ability in AFP^− HCC.
Upregulated AR transcriptional activity promotes angiogenesis in AFP^− HCC
Next, we aimed to examine whether AR orchestrated the enrichment of
tip-like endothelial cells in AFP^− HCC. To this end, we performed
correlation analysis and found that the expression of nuclear AR was
significantly positively correlated with CD31 expression
(Fig. [196]5a). Furthermore, spatial transcriptomics analysis using a
public dataset revealed that AR^+ tumor spots were in closer proximity
to PECAM1^+ spots (namely endothelial cell spot) compared to AR^– tumor
spots (Supplementary Fig. [197]8a, b). These results implied that AR
might play a crucial role in regulating angiogenesis within the AFP^−
HCC TME.
Fig. 5. Activation of EGFR with AR-induced 11,12-EET from AFP^− HCC tumor
cells promoted aberrant angiogenesis.
[198]Fig. 5
[199]Open in a new tab
a Scatter plot showing positive correlation between proteinic nuclear
AR expression and proteinic CD31 expression. R indicates the
correlation coefficient calculated by Spearman correlation test. The
number of dots indicates the number of patients. b Sensorgrams for
detecting the binding affinity of 11,12-EET and EGFR using SPR assay. c
Representative images (left panel) and quantitative analysis (right
panel) demonstrating the pro-angiogenic effects of 11,12-EET on HUVEC
migration and tube formation compared to non-11,12-EET treatment group
(NC). A total of three biological replicates were employed. Data were
shown as mean ± SD. d Protein expression of VEGFR2, phosphorylated
VEGFR2 (p-VEGFR2), EGFR2 and phosphorylated EGFR2 (p-EGFR2) in HUVECs
when co-cultured with 11,12-EET compared to non-11,12-EET treatment
group (NC). Representative images (left panel) and corresponding
statistical results (right panel) of migration (e) and tube formation
(f) of HUVECs when co-cultured with the CM from shNC and shAR cell
lines, including MHCC97H and HCCLM3, respectively. A total of three
biological replicates were employed. Data were shown as mean ± SD. g
Representative images (left panel) and quantitative analysis (right
panel) of migration and tube formation of HUVECs when co-cultured with
CM from MHCC97H. The MHCC97H cells were pre-treated with 10 nM 5α-DHT
for 48 h, followed by separate treatments with either 6 μM lenvatinib,
10 μM gefitinib, 20 μM gefitinib or vehicle control for 24 h. A total
of three biological replicates were employed. Data were shown as
mean ± SD. Scale bars, 100 μm. *p < 0.05, **p < 0.01, ***p < 0.001 and
ns stands for no significance by two-sided t-test in (c) and one-way
ANOVA with Tukey’s multiple comparisons test (e–g). Source data are
provided as a Source Data file.
Due to the strong association between AR transcriptional activation and
11,12-EET production in AFP^− HCC, we interrogated whether AR-mediated
11,12-EET could promote sprouting angiogenesis by interacting with
tip-like endothelial cells. To this end, we first predicted the binding
affinity of 11,12-EET to the endothelial-specific pro-angiogenic
targets by Autodock^[200]44,[201]45, and found that 11,12-EET exhibited
the strongest binding affinity to EGFR among the candidate receptors
(Supplementary Fig. [202]8c). Moreover, the data from surface plasmon
resonance (SPR) assay confirmed that 11,12-EET could directly bind to
EGFR (Fig. [203]5b). As a result, EGFR expression and the downstream
signaling pathway were also significantly upregulated in EC_PLPP1 of
AFP^− HCC (Supplementary Fig. [204]8d, e). We then performed in vitro
assays to examine the interaction between AR-mediated 11,12-EET
production and endothelial cells. The 11,12-EET treatment was shown to
promote vascular migration and tube formation of human umbilical vein
endothelial cells (HUVECs) (Fig. [205]5c). Of note, the phosphorylation
level of EGFR, but not VEGFR2, was upregulated in the presence of
11,12-EET stimulation (Fig. [206]5d). Bulk RNA-seq data derived from
HUVECs treated with 11,12-EET also confirmed the activation of
EGFR-related signaling pathway, while VEGFR-related pathway remained
unaffected (Supplementary Fig. [207]8f, g). In addition, we cultured
HUVECs with the conditioned medium (CM) from AFP^lowAR^high cell lines
(MHCC97H and HCCLM3). The results demonstrated that vascular cell
migration and tube formation were robustly inhibited in the
AR-knockdown group (Fig. [208]5e, f). Notably, we further revealed that
gefitinib, an EGFR inhibitor, showed a significantly inhibitory effect
on HUVEC migration and tube formation than lenvatinib in a
dose-dependent manner (Fig. [209]5g), demonstrating that lenvatinib
showed mild inhibitory effect on tumor angiogenesis in the presence of
11,12-EET stimulation in AFP^− HCC. Overall, AR-11,12-EET-EGFR axis
played a critical role in promoting angiogenesis in AFP^− HCC.
AR modulates tumorigenesis and vascularization through 11,12-EET-EGFR axis in
vivo
To demonstrate the eliciting effects of AR on tumorigenesis and
vascularization in vivo, we established the subcutaneous HCC models
using two human cell lines (MHCC97H, defined as AFP^lowAR^high; Huh7,
defined as AFP^highAR^low) in immunodeficient NCG mice, and two murine
cell lines (RIL-175, defined as Afp^lowAr^high; Hepa1-6, defined as
Afp^highAr^low) in immunocompetent C57BL/6 mice, respectively.
Lentiviral transfection-mediated gene expression was performed to
generate AR-knockdown counterparts for AFP^lowAR^high cell lines and
AR-overexpression counterparts for AFP^highAR^low cell lines
(Fig. [210]3k and Supplementary Figs. [211]6e, [212]9a).
In NCG mice, AR knockdown significantly suppressed tumor burden,
manifested as reduced invasive growth, decreased tumor volume and
weight (Fig. [213]6a–e), accompanied with lower expression of the
proliferative marker Ki67 (Supplementary Fig. [214]9b). In contrast, AR
overexpression significantly fueled the tumor growth (Fig. [215]6f–i
and Supplementary Fig. [216]9c). We further examined the expression
profile of AR-11,12-EET-EGFR axis. The results showed that knockdown of
AR significantly inhibited the expression levels of CYP2C8 and CYP1B1
(Supplementary Fig. [217]9d) as well as the production level of
11,12-EET (Fig. [218]6j). The expression of downstream molecule EGFR on
endothelial cells (Fig. [219]6k) and the microvascular density
(Supplementary Fig. [220]9e, f) were also reduced upon AR knockdown,
while overexpression of AR generated significant contrast results
(Fig. [221]6j, [222]l and Supplementary Fig. [223]9g–i).
Fig. 6. Modulation of AR in tumorigenesis and vascularization in vivo.
[224]Fig. 6
[225]Open in a new tab
a Experimental design using NCG mice (n = 6 per group). Representative
images of tumors (b), tumor growth curve (c), tumor volume measurements
(d) and tumor weight measurements (e) in MHCC97H subcutaneous HCC mouse
models (shNC vs. shAR = 6:6). Representative images of tumors (f),
tumor growth curve (g), tumor volume measurements (h) and tumor weight
measurements (i) in Huh7 subcutaneous HCC mouse models (oeNC vs. oeAR =
6:6). j Relative concentrations of 11,12-EET measured by ELISA assays
in MHCC97H and Huh7 subcutaneous HCC tumors, and their corresponding
shAR and oeAR counterparts (n = 6 per group), respectively. Flow
cytometry analysis of CD45^-CD31^+EGFR^+ endothelial cell (EC)
fractions in MHCC97H (k) and Huh7 (l) subcutaneous HCC tumors (n = 6
per group). m Experimental design using C57BL/6 mice (n = 6 per group).
Representative images of tumors (n), tumor growth curve (o), tumor
volume measurements (p) and tumor weight measurements (q) in RIL-175
subcutaneous HCC mouse models (shNC vs. shAr = 6:6). Representative
images of tumors (r), tumor growth curve (s), tumor volume measurements
(t) and tumor weight measurements (u) from Hepa1-6 subcutaneous HCC
mouse models (oeNC vs. oeAr = 6:6). v Relative concentration of
11,12-EET detected by ELISA kit in RIL-175 and Hepa1-6 subcutaneous HCC
tumors, and their corresponding shAr and oeAr counterpart (n = 6 per
group), respectively. Flow cytometry analysis of EGFR^+ ECs fractions
in RIL-175 (w) and Hepa1-6 (x) subcutaneous HCC tumors (n = 6 per
group). Data were shown as mean ± SD. *p < 0.05, **p < 0.01,
***p < 0.001 by two-way ANOVA in (c, g, o, s) and two-sided t-test in
(d, e, h–l, p–q, t–x). Source data are provided as a Source Data file.
In C57BL/6 mice, similar to the results observed in NCG mice models
using human HCC cell lines, knockdown of Ar significantly repressed the
HCC growth and tumor angiogenesis while overexpression of Ar promoted
HCC progression and vascularization within the same species
(Fig. [226]6m–x and Supplementary Fig. [227]10a–h) Altogether, these
data revealed the important functions of AR in regulation of
tumorigenesis and vascularization in AFP^− HCC.
To further consolidate the regulatory link between AR and
11,12-EET-EGFR axis, we established an orthotopic HCC mouse model using
MHCC97H cells, comprising four experimental groups: NC (control), shAR
(AR knockdown), oeCYP2C8 (CYP2C8 overexpression), and shAR plus
oeCYP2C8 (AR knockdown with CYP2C8 rescue) (Supplementary
Fig. [228]11a, b). The data demonstrated that AR knockdown
significantly suppressed tumor growth (Supplementary Fig. [229]11c, d)
accompanied with reduced Ki67 (Supplementary Fig. [230]11e),
CYP2C8/CYP1B1 expression (Supplementary Fig. [231]11f) and 11,12-EET
production level (Supplementary Fig. [232]11g) as well as decreased
endothelial EGFR expression level (Supplementary Fig. [233]11h) and
microvascular density (Supplementary Fig. [234]11i, j), while CYP2C8
overexpression partially rescued the shAR-eliciting anti-tumor and
anti-angiogenesis effects (Supplementary Fig. [235]11). The in vivo
experimental data provided bona fide evidence supporting the regulatory
role of AR in control of 11,12-EET-EGFR axis that promoted angiogenesis
in AFP^− HCC.
Combining AR inhibitor bicalutamide with lenvatinib suppresses tumor
progression and angiogenesis in AFP^− HCC in vivo
To investigate the therapeutic potential of targeting the
AR-11,12-EET-EGFR for AFP^− HCC, we established an orthotopic HCC mouse
model by injecting MHCC97H cells into the liver lobes of NCG mice
(Fig. [236]7a). The control and treatment groups were set as follows:
(1) shNC groups including: (1) vehicle control (Ctrl); (2) gefitinib
monotherapy (GEF); (3) bicalutamide monotherapy (BIC); (4) lenvatinib
monotherapy (LEN); (5) bicalutamide plus gefitinib combination
(BIC + GEF); (6) bicalutamide plus lenvatinib combination (BIC + LEN).
(2) shAR groups including: (1) vehicle control (Ctrl); 2) gefitinib
monotherapy (GEF); (3) bicalutamide monotherapy (BIC); (4) lenvatinib
monotherapy (LEN); 5) bicalutamide plus gefitinib combination
(BIC + GEF); (6) bicalutamide plus lenvatinib combination (BIC + LEN).
Fig. 7. Combined inhibition of AR-mediated anti-angiogenesis with lenvatinib
repressed tumor progression in AFP^− HCC in vivo.
[237]Fig. 7
[238]Open in a new tab
a Schematic of the experimental design using NCG mice (n = 5 per
group). Representative gross images of MHCC97H orthotopic tumors (b)
and statistical analysis of tumor volume (c). The “+” indicates
treatment with gefitinib (GEF), bicalutamide (BIC) or lenvatinib (LEN),
while the “–” indicates treatment with vehicle control (n = 5 per
group). Data were shown as mean ± SD. d, e Representative images of IHC
staining (left panel) and microvascular density quantification (right
panel) in MHCC97H orthotopic tumors across different treatment groups
(n = 5 per group). Data were shown as mean ± SD. f A graphical summary.
In AFP^− HCC patients, TF activity of AR was elevated and promoted the
expression levels of CYP2C8 and CYP1B1, which enhanced arachidonic acid
metabolism and the production of the downstream metabolite 11,12-EET.
The 11,12-EET metabolite activated EGFR on endothelial cells, thereby
promoting angiogenesis. Targeting this pathway, a combination of AR
antagonists with lenvatinib might represent a promising therapeutic
strategy for AFP^− HCC. Scale bars, 50 μm. *p < 0.05, **p < 0.01,
***p < 0.001, and ns stands for no significance by one-way ANOVA with
Tukey’s multiple comparisons test (c, e). Source data are provided as a
Source Data file.
In the shNC group, GEF monotherapy showed negligible effects on tumor
burden compared to controls, whereas both LEN and BIC monotherapies
elicited moderate tumor suppression (Fig. [239]7b, c and Supplementary
Fig. [240]12a, b) accompanied by reduced angiogenesis (Fig. [241]7d,
e). Notably, the combination of BIC + GEF moderately reduced the tumor
growth, while the combination of BIC + LEN demonstrated the most
remarkable suppression of tumor burden (Fig. [242]7b, c and
Supplementary Fig. [243]12a, b).
In the shAR group, the vehicle controls (AR knockdown solo) exhibited
significantly lower baseline tumor burden and microvascular density
compared to shNC controls (Fig. [244]7b–e), further validating the role
of AR in regulation of tumor progression and vascularization. Among the
monotherapy groups, only LEN treatment significantly repressed tumor
burden compared to the shAR alone. Of note, both BIC + GEF and
BIC + LEN significantly damped down the tumor burden, with the
BIC + LEN showing more pronounced suppressive effects in tumorigenesis
and vascularization (Fig. [245]7b–e).
Furthermore, to solidify the synergistic effects of AR antagonist and
lenvatinib, we also established a Hepa1-6 orthotopic HCC mouse model
(oeNC vs. oeAR). The results showed that lenvatinib monotherapy
moderately inhibited tumor growth but this inhibitory effect was
significantly attenuated upon Ar overexpression (Supplementary
Fig. [246]12c, d). Besides, combination of lenvatinib with bicalutamide
significantly suppressed tumor burden, reduced Ki67 level and decreased
microvascular density in the oeAr group (Supplementary
Fig. [247]12c–f), further supporting the therapeutic potential of this
combined treatment approach.
Taken together, we uncovered the hyperactive metabolic activity of
arachidonic acid metabolism mediated by AR in AFP^− HCC, and shed
insights into the mechanism of AR-induced aberrant angiogenesis,
providing a promising therapeutic strategy of combining AR antagonists
(e.g., bicalutamide) and lenvatinib for AFP^− HCC (Fig. [248]7f).
Discussion
In the present study, we identified arachidonic acid metabolism as a
specifically activated metabolic pathway in AFP^− HCC, with AR being
its key regulatory factor by integrated omics analysis and multiple
experiments. Mechanistically, AR could promote the transcription of
CYP2C8 and CYP1B1 by binding to their promoter regions, thereby
enhancing the arachidonic acid metabolism that produced 11,12-EET,
which activated EGFR on endothelial cells and eventually fueled tumor
angiogenesis. Building on the potent roles of AR in metabolic
reprogramming^[249]46, including its regulation on lipid pathways of
prostate cancer^[250]47, and cholesterol metabolism^[251]48 as well as
glycolysis^[252]49 in HCC, here, we unveil the impact of AR-mediated
arachidonic acid metabolism on the reprogramming of AFP^− HCC TME,
providing AR as the target for AFP^− HCC treatment.
To our knowledge, the role of AR in HCC remains highly controversial.
AR is originally regarded as an oncogenic driver of HCC, maintaining
the stemness of HCC cells and promoting cell cycle
process^[253]30,[254]50,[255]51. However, parts of clinical trials have
not shown significant beneficial effects of anti-androgenic therapy in
HCC patients^[256]52, and some studies showed that hepatic AR
suppresses HCC metastasis^[257]32. These findings implied that hepatic
AR might have dual but opposing roles, resulting in a dilemma in
personalized treatment of HCC, and further fueling our interest in
potential AR-driven TME reprogramming of AFP^− HCC tumors.
Furthermore, our findings demonstrated that the AR-CYP2C8/CYP1B1 axis
promoted the generation of the downstream metabolite 11,12-EET, which
could activate the EGFR on endothelial cells, thus inducing a tip-like
phenotype in endothelial cells and ultimately leading to aberrant
angiogenesis. Currently, anti-angiogenic therapy targeting VEGFRs such
as lenvatinib is strongly recommended as the first-line systemic
treatment for the advanced HCC and is widely applied in the
clinic^[258]53. However, several recent global clinical trials have
indicated that AFP^− HCC patients showed moderate treatment response to
the first-line drug lenvatinib or its analogues^[259]13–[260]15.
Coincidentally, another study reported that HCC could be categorized
into three subtypes (C1, C2, and C3) based on metabolism-related
characteristic genes^[261]54. Among them, the C1 subtype exhibited
higher metabolic activity and lower AFP expression, conferring high
activity of angiogenesis-related characteristics^[262]54, which
exhibited resistant to cabozantinib, a second-line anti-angiogenic drug
targeting VEGFR2^[263]54. In this study, we employed various animal
models and demonstrated that AR^high HCC mice possessed higher tumor
burden and vascular density than those in AR^low HCC mice, while
combination therapy with bicalutamide and lenvatinib significantly
inhibited tumor burden and vascular density in AFP^− HCC.
Interestingly, both in vitro and in vivo findings suggested that AFP^−
HCC might promote angiogenesis through AR-induced metabolic
reprogramming and subsequent generation of 11,12-EET that binding to
EGFR on endothelial cells, which might lead to the insensitivity of
AFP^− HCC to canonical anti-angiogenic drugs.
In conclusion, our findings reveal the critical role of AR in
regulation of arachidonic acid metabolism in AFP^− HCC and provide
insights into the mechanism of AR-CYPs-11,12-EET-EGFR axis-induced
tumor angiogenesis. These findings not only disclose the link between
tumor metabolism and abnormal vascular growth but also highlight the
promising strategy of combining AR antagonists to enhance the efficacy
of lenvatinib therapy in AFP^− HCC.
Methods
Study approval
The study was approved by the Institutional Research Ethics Committee
of First Affiliated Hospital of Sun Yat-sen University (FAH-SYSU) under
Ethical Approval ID [2023]380. Informed consent, including to publish
clinical information potentially identifying individuals, was obtained
from all patients.
All animal experiments were performed in accordance with protocols
approved by the Institutional Animal Care and Use Committee, Sun
Yat-sen University (IACUC, SYSU). Approval numbers are
SYSU-IACUC-2024-000262 and SYSU-IACUC-2024-000264.
Patients and samples
A total of 448 treatment-naïve patients who underwent curative liver
resection at FAH-SYSU were included, with tumor samples and NAT
samples. This study has received approval from the Ethics Committee of
FAH-SYSU and adhered to relevant ethical standards. Clinical
information of patients was summarized in Supplementary Data [264]1 and
Supplementary Data [265]4. Two experienced pathologists confirmed the
diagnoses through histological examinations. The sequencing types
utilized in this study and their corresponding sample information were
provided in Supplementary Fig. [266]1A.
Exome and transcriptome capture, library preparation and sequencing
For exome capture and library preparation, total DNA was extracted from
snap-frozen tissues using the QIAGEN DNeasy Blood & Tissue Kit (Qiagen,
Hilden, Germany). Qualified DNAs from tissue samples were fragmented
into 200-300 bp fragments using Covaris technology with adapters
ligated to both ends. Next, the extracted DNA was amplified by
ligation-mediated PCR (LM-PCR), then purified and hybridized to
NimbleGen or Agilent human exome array for enrichment. Non-hybridized
fragments were washed away. For transcriptome capture and library
preparation, total RNA was extracted using the Trizol (Invitrogen,
USA). Poly(A) mRNA was isolated from total RNA using oligo-dT-coupled
beads. The mRNA was fragmented into short fragments by adding a
fragmentation buffer. Using these short fragments as templates, the
first-strand cDNA was synthesized with random hexamer primers. The
second-strand cDNA was then synthesized using buffer, dNTPs, RNase H,
and DNA polymerase I. The short fragments were purified using the
QIAQuick PCR extraction Kit (Qiagen, Hilden, Germany), dissolved in EB
buffer for end reparation and adding poly(A). Next, the short fragments
were ligated with sequencing adapters. Finally, exome and transcriptome
libraries were constructed and loaded onto the Illumina NovaSeq 6000
sequencing platform, generating paired-end reads of 150 base pairs.
Mutation calling and annotation
FASTQ files were filtered for adapter sequences and low-quality reads
using TrimGalore (v.0.6.1,
[267]https://github.com/FelixKrueger/TrimGalore). FastQC was then used
to evaluate the quality of each sequencing read. Qualified reads were
mapped to the human reference genome GRCh38 in paired-end mode using
BWA-mem2^[268]55 (v.2.2.1). Duplicate reads were marked using
Sambamba^[269]56 (v.0.8.2). Following GATK best practices,
Mutect2^[270]57 was then employed to detect confident somatic
mutations. ANNOVAR^[271]58 was used to annotate the final somatic
mutations. Finally, MutSigCV^[272]59 (v.1.3.5) was employed to identify
significantly mutated genes with all SNVs and INDELs used as input.
Genes with q-values < 0.05 were defined as significantly mutated. Tumor
mutation burden (TMB), defined as the number of somatic mutations per
megabase (Mb) in the coding region, was determined and quantified using
our WES data. To avoid overestimating TMB, mutations with population
frequencies higher than 0.05 in the Genome Aggregation Database and
driver mutations curated in COSMIC database were excluded.
Detection of CNVs
CNVkit^[273]60 (v.0.9.10) was employed to call CNVs, with tumor copy
number estimates derived from BAM files using batch processing and
segmentation methods. GISTIC^[274]61 (v.2.0.23) was applied to process
the results files from the previous step, identifying focal
amplifications and deletions in cytoband level with amplification
threshold set to 0.3, deletion threshold set to −0.3. The G-score was
calculated based on the amplitude of the CNVs and its frequency of
occurrence in the sample, reflecting the frequency of CNVs in the
segments. To compare CNV differences between AFP^− HCC and AFP^+ HCC
patients, we first identified frequently amplified and deleted events
across all patients (q value < 0.1). Next, we conducted a comparative
analysis of frequently amplified and deleted events across cytobands in
AFP^− HCC and AFP^+ HCC samples. For each cytoband, we calculated the
amplification or deletion frequency of CNVs in both groups and used
chi-square test or Fisher’s exact test with FDR adjustment to evaluate
the differences. Statistical significance was defined as FDR of less
than 0.05.
Microsatellite instability estimation
MSIsensor^[275]62 was used to detect frameshift mutations in
microsatellite regions and distinguish them as somatic mutations. In
brief, a chi-square test or Fisher’s exact test was employed to compare
the length distribution of microsatellite sequences repeated at each
locus between matched tumor and NAT samples. If the MSIsensor score of
a sample was greater than 3.5, it was considered microsatellite
instability (MSI); otherwise, it was classified as microsatellite
stability (MSS).
Cancer cell fraction analysis
PyClone^[276]63 (v.0.13.1) was used to assign somatic mutations to
clonal clusters and estimate their cellular prevalence. The cancer cell
fraction (CCF) of each somatic mutation was estimated based on the SNV
frequency and copy number profile. For each sample, the corresponding
tumor purity was used to adjust the estimated CCF.
Immune editing
Polysolver^[277]64 (v.1.0) was employed to determine the HLA typing for
each patient using the corresponding NAT sample. The non-synonymous
mutations that occurred in each sample were obtained from the processed
VCF file. Corresponding amino acid sequences were acquired from UniProt
database, with the amino acid sequences sliding around the mutation
site (9–11 residues) being selected as input. These mutated sequences,
along with the HLA typing information for each sample, were
collectively input into NetMHCpan^[278]65 (v.4.0). Peptide segments
with predicted binding affinities less than 500 nM were considered as
immunogenic and used for subsequent analysis. The same process was
applied to TCGA-LIHC cohort to construct a background condition.
Finally, the immune editing score was calculated according to the
literature published previously^[279]66,[280]67.
Mutation signature identification
SigProfilerExtractor^[281]68 (v.1.1.21) was used to de novo decompose
the mutation spectrum in the processed VCF files for all samples. The
decomposed signatures were compared to known signatures from the
Catalogue of Somatic Mutations in Cancer (COSMIC) database (v.3) using
cosine similarity (>0.8) to determine their corresponding types.
Wilcoxon rank sum test was applied to assess activity differences in
contributions to each signature between AFP^− HCC and AFP^+ HCC.
RNA sequencing data analysis
By filtering low-quality reads and removing sequencing adapters using
TrimGalore (v.0.6.1, [282]https://github.com/FelixKrueger/TrimGalore),
the FASTQ files of all samples were mapped to the human reference
genome GRCh38 using STAR^[283]69 (v.2.7.9) with default parameter. The
gene-level read counts and expression values in fragment per kilobase
per million mapped fragments (FPKM) units were obtained using the
RSEM^[284]70 (v.1.2.28) pipeline. Based on the quantified gene read
count, DESeq2^[285]71 was employed to assess the differentially
expressed genes between AFP^− HCC and AFP^+ HCC, with a significance
threshold set as adjusted p-value < 0.05 and |log[2]-transformed fold
changes (log[2]FC)| > 1. GSEA^[286]72 between AFP^− HCC and AFP^+ HCC
was conducted using pathway databases (Hallmark, KEGG, and REACTOME)
from the Molecular Signatures Database (MSigDB,
[287]https://www.gsea-msigdb.org/gsea/msigdb/index.jsp). In short, all
detected genes were ranked based on the log[2]FC between AFP^− HCC and
AFP^+ HCC for GSEA using the pre-ranked mode with default parameter
settings. The xCell^[288]73 tool was employed to estimate the
proportion of endothelial cells in treatment-naïve patients from
FAH-SYSU cohort and TCGA-LIHC cohort data.
Classification of HCC
To delineate the transcriptional profiles of AFP^− HCC and AFP^+ HCC,
we employed the NTP algorithm^[289]21 to predict the molecular subtypes
of HCC based on the previously published Hoshida’s
classification^[290]20,[291]21. NTP was utilized with default
parameters. Differences in Hoshida’s classification were performed by
chi-square test, based on methods described in previous
studies^[292]74–[293]76. Namely, we conducted a chi-square test and
calculated the ratio of observed to expected sample numbers for each
subtype (Eqs. [294]1 and [295]2).
[MATH:
X2= ∑fo−fe2fe
:MATH]
1
[MATH:
Ro/e=fo
f
e :MATH]
2
In the equation above,
[MATH: X2
:MATH]
represents the chi-square value. f[0] indicates the observed sample
counts, while
[MATH: fe
:MATH]
denotes the expected sample counts, which can be derived from
chi-square test.
[MATH:
Ro/e :MATH]
signifies the ratio of observed to expected sample counts for each
subtype in different groups. If
[MATH:
Ro/e :MATH]
> 1, we assumed that a subtype showed enrichment preference in a
specific group.
Methylation profiling by RRBS
Methylation profiling was performed using RRBS. Briefly, a total of
50 ng genomic DNA was digested with the MspI restriction enzyme for a
sufficient duration (1–3 h) to generate DNA fragments, and then ligated
with a methylated adapter containing complementary sticky end. The
ligated products were bisulfite-converted with sodium bisulfite using
the MethylCode^TM Bisulfite Conversion Kit (ThermoFisher, MECOV50),
PCR-amplified and purified using magnetic beads according to
manufacturer’s instruction. Purified library was quantified using
Qubit, and then diluted to a specified concentration suitable for
high-throughput sequencing. DNA library was conducted for sequencing
and methylation analysis.
Methylation analysis
FASTQ files were filtered for adapter sequences and low-quality reads
using TrimGalore (v.0.6.1,
[296]https://github.com/FelixKrueger/TrimGalore). FastQC was then used
to evaluate the quality of each sequencing reads. Subsequently, the
qualified reads were aligned to GRCh38 reference assembly and the
methylation state was determined using Bismark^[297]77 (v.0.24.1). In
order to investigate whether there were differences in the methylation
levels of specific genes between AFP^- HCC and AFP^+ HCC, we calculated
the methylation levels within the promoter regions, defined as the
0-1500 bp upstream of the transcription start site (TSS) of the target
genes.
Generation of scRNA-seq data
Fresh tumor samples were collected within 30 min after resection and
preserved in MACS Tissue Storage Solution (130-100-008, Miltenyi
Biotec). Subsequently, the samples were digested into a single-cell
suspension following these steps: (I) Physical slicing of the tumor
tissue using a blade, followed by digestion in RPMI 1640 containing 5%
FBS using Tumor Dissociation Kit (130-095-929, Miltenyi Biotec) and
DNaseI (DN25-100 mg, Sigma Aldrich) for 20 min at 37 °C; (II) The
digested mixture was filtered through 70 μm and 30 μm MACS strainers
(130-098-462; 130-098-458, Miltenyi Biotec). Cell suspension was
centrifuged at 400 x g for 6 min, treated with 1× Red Blood Cell Lysis
Buffer (00430054, Invitrogen) for 5 min at 4 °C. Cells were washed
twice with PBS, and AOPI staining was performed to assess the
concentration and viability of the single-cell suspension; (III) The
concentration of the cell suspension was adjusted to ~700–1300
cells/μL, followed by reverse transcription and cDNA amplification
according to the standard protocol (Single Cell 5’ Reagent Kits User
Guide) to prepare the library construction, and sequenced on the
Illumina NovaSeq 6000.
scRNA-seq raw data processing
The raw sequencing data were processed using CellRanger (v.6.1.2, 10x
Genomics). The human genome GRCh38 was used to generate expression
matrices. To eliminate contamination from ambient RNA present in
droplets, we applied SoupX^[298]78 (v.1.5.2) to each sample. Cells with
detectable expression of at least three genes were retained, while
cells with fewer than 200 detected genes were excluded. Additionally,
cells with gene counts exceeding 6000 and cells with a mitochondrial
gene percentage greater than 15% were removed to ensure high-quality
cells.
Normalization of the data was performed based on the raw counts using
the Seurat^[299]79 (v.4.1.0) R pipeline. The “vst” method in the Seurat
FindVariableFeatures function identified the top 2000 highly variable
genes for principal component analysis. A shared nearest neighbor graph
was constructed using the top 50 principal components, followed by cell
clustering using the Louvain algorithm. The resolution parameter was
set to 0.7. Harmony^[300]80 (v.1.0.3) pipeline was applied to remove
batch effects across all samples, encompassing both malignant and
non-malignant cells, using default parameters. Visualization of the
results was achieved using Uniform Manifold Approximation and
Projection (UMAP). Typical doublets were selected based on
co-expression of lineage-specific marker genes (e.g., T/B cells,
NK/myeloid cells) and subsequently removed.
Single cell CNV calling
To confirm the identification of tumor cells, we employed inferCNV
(v.1.14.2, [301]https://github.com/broadinstitute/inferCNV/wiki) to
calculate CNV score for each cell. The gene expression matrix of tumor
cells was extracted from the Seurat object, and CNVs were computed
using cancer-associated fibroblasts and endothelial cells as
references. The parameters used were “denoise” = TRUE, “hidden Markov