A genomic enhancer signature associates with hepatocellular carcinoma prognosis

Abstract

Background & aims: Lifestyle and environmental-related exposures are important risk factors for hepatocellular carcinoma (HCC), suggesting that epigenetic dysregulation significantly underpins HCC. We profiled 30 surgically resected tumours and the matched adjacent normal tissues to understand the aberrant epigenetic events associated with HCC.

Methods: We identified tumour differential enhancers and the associated genes by analysing H3K27 acetylation (H3K27ac) chromatin immunoprecipitation sequencing (ChIP-seq) and Hi-C/HiChIP data from the resected tumour samples of 30 patients with early-stage HCC. This epigenome dataset was analysed with previously reported genome and transcriptome data of the overlapping group of patients from the same cohort. We performed patient-specific differential expression testing using multiregion sequencing data to identify genes that undergo both enhancer and gene expression changes. Based on the genes selected, we identified two patient groups and performed a recurrence-free survival analysis.

Results: We observed large-scale changes in the enhancer distribution between HCC tumours and the adjacent normal samples. Many of the gain-in-tumour enhancers showed corresponding upregulation of the associated genes and vice versa, but much of the enhancer and gene expression changes were patient-specific. A subset of the upregulated genes was activated in a subgroup of patients' tumours. Recurrence-free survival analysis revealed that the patients with a more robust upregulation of those genes showed a worse prognosis.

Conclusions: We report the genomic enhancer signature associated with differential prognosis in HCC. Findings that cohere with oncofoetal reprogramming in HCC were underpinned by genome-wide enhancer rewiring. Our results present the epigenetic changes in HCC that offer the rational selection of epigenetic-driven gene targets for therapeutic intervention or disease prognostication in HCC.

Impact and implications: Lifestyle and environmental-related exposures are the important risk factors of hepatocellular carcinoma (HCC), suggesting that tumour-associated epigenetic dysregulations may significantly underpin HCC. We profiled tumour tissues and their matched normal from 30 patients with early-stage HCC to study the dysregulated epigenetic changes associated with HCC. By also analysing the patients' RNA-seq and clinical data, we found the signature genes - with epigenetic and transcriptomic dysregulation - associated with worse prognosis. Our findings suggest that systemic approaches are needed to consider the surrounding cellular environmental and epigenetic changes in HCC tumours.

Keywords: Cancer prognosis; Epigenetics; Hepatocellular carcinoma; Multi-omics; Personalised medicine.

Graphical abstract

Fig. 1

Enhancer rewiring is partly attributable to genomic mutations. (A) Pie charts showing proportions of different genomic elements in all enhancers (left), gain-in-tumour enhancers (middle), and lost-in-tumour enhancers (right). (B) Types of enhancers and the total number detected. (C) Number of enhancer loci detected across the adj.normal (top) and tumour (bottom) tissues of different numbers of patients. (D) Stacked percentage bar chart of types of changes in the TAD boundaries between different pairs of samples, analysed using TADCompare. N vs. N, adj.normal to adj.normal; T vs. T, tumour to tumour; inter-patient, between samples of different patients; intra-patient, between samples of the same patient. (E) Heatmap of the type of mutations at mutation hotspots identified from MutEnricher. The left heatmap shows mutation hotspots at the gain-in-tumour enhancer regions, whereas the right heatmap shows mutation hotspots at the lost-in-tumour enhancer regions. Each row represents a mutation hotspot, corresponding to an enhancer locus. Barplots next to each heatmap represents the number of patients from which the mutation hotspot was detected. A mutation hotspot was considered truncal (red) if it was detected in all tumours of a patient, whereas mutations occurring in certain tumour sectors of a patient were considered branch mutations (blue). (F) Enriched somatic mutations at the most frequent mutation hotspot detected in (E). Numbers in parentheses represent the number of tumour sectors in which the mutation was detected from. The enhancer, ‘chr5:1295700’ (denoted in grey bar), overlaps with the TERT promoter (denoted in blue bar) and is a gain-in-tumour enhancer. adj.normal, adjacent normal; TAD, topologically associated domain; UTR, untranslated region.

Fig. 2

Differential enhancers in HCC tumours. (A) Schematic diagram of the overall project design. (B) Schematic diagram of how differential enhancers were identified (see Patients and methods). (C) Emission and transition probability matrix from ChromHMM. (D) H3K27ac signal profiles at the gained and lost enhancers. Red, normal; blue, tumour. (E) Heatmap of H3K27ac signal from adj.normal and tumour samples at differential enhancers, normalised to the respective control samples. Heatmaps are centred at the midpoint of each enhancer and extended upstream and downstream by 5 kb. Gained-in-tumour (top) and lost-in-tumour (bottom) enhancer regions are shown. The same heatmaps for all samples can be found in Fig. S1. adj.normal, adjacent normal; AFP, alpha-foetoprotein; ChIP-seq, chromatin immunoprecipitation sequencing; H3K27ac, H3K27 acetylation; HCC, hepatocellular carcinoma; RNA-seq, RNA sequencing; scRNA-seq, single-cell RNA sequencing.

Fig. 3

Heterogeneity in the H3K27ac signal at differential enhancers and the associated gene expression changes. (A) Heatmap of H3K27ac ChIP-seq signal at the gain-in-tumour enhancers. Enhancers shown on the rows were ordered using hierarchical clustering. H3K27ac ChIP-seq samples are shown on the columns, ordered by the eClusters identified by K-means clustering. (B) Heatmap of H3K27ac ChIP-seq signal at the lost-in-tumour enhancers. (C) PCA plot of the samples shown in (A) and (B), grouped by the eClusters. (D) Total H3K27ac signal at the differential enhancers (left, gain-in-tumour; right, lost-in-tumour), grouped by the eClusters. (E) Proportion of eClusters assigned to adj.normal and tumour tissues of each patient (top and middle), and the recurrence status for the corresponding patients with tumour tissues sequenced (bottom). All bar plots are grouped by the Edmondson grade of the patients’ tumour sample. (F) Boxplot of disease-free survival for patients shown in (E), grouped by the eClusters. (G) Gene set enrichment analysis results for the genes associated with the gain-in-tumour enhancers (top) and lost-in-tumour enhancers (middle and bottom). (H) Proportion of differentially expressed genes for genes associated with the gain-in-tumour enhancers (top) and with the lost-in-tumour enhancers (bottom). Results obtained from the conventional differential expression testing. (I) Heatmap of categorised expression changes matrix from patient-specific differential expression testing. Top shows the gain-in-tumour associated genes, whereas the bottom shows the lost-in-tumour associated genes. Colours of blue, light yellow, and red each represent the expression changes of upregulation, no change, and downregulation, respectively. adj.normal, adjacent normal; ChIP-seq, chromatin immunoprecipitation sequencing; eCluster, enhancer cluster; FDR, false discovery rate; GOBP, Gene Ontology Biological Process; H3K27ac, H3K27 acetylation; PCA, principal component analysis.

Fig. 4

Patient-dependent H3K27ac signal at the associated enhancers and expression changes of SOX4 (left) and GPC3 (right). (A) H3K27ac ChIP-seq signal at the enhancers associated with SOX4. Differential enhancer status of each enhancer is denoted in yellow. Chromatin loops were inferred from Hi-C/HiChIP libraries and are coloured based on the frequency of detection in Hi-C/HiChIP libraries. (B) Similar plot as A for GPC3. (C) Normalised H3K27ac coverage at the differential distal enhancer of SOX4, for the paired adj.normal (blue, bottom) and tumour (red, top) tissues of 6 patients. (D) Similar plot to (C) for the differential promoter enhancer of GPC3. (E) Per-patient differential expression testing results for SOX4 across 90 patients. Normalised gene counts for each adj.normal (blue) and tumour (pink) samples are shown on top, whereas the bottom plot shows the log₂ fold change values for those detected as differentially expressed (blue bars, upregulated in tumour; green bars, downregulated in tumour; grey, not significantly different). In both plots, each column represents a patient. (F) Similar plot to (E) for GPC3. adj.normal, adjacent normal; ChIP-seq, chromatin immunoprecipitation sequencing; H3K27ac, H3K27 acetylation.

Fig. 5

Gain-in-tumour associated genes show enrichment for foetal liver genes. (A) Gene set enrichment analysis result of gain-in-tumour enhancer-associated genes against cell type signature gene set C8 from MSigDB. (B) Heatmap of H3K27ac ChIP-seq signal for the gain-in-tumour enhancers associated with the foetal liver-expressed genes, ordered by the eCluster. (C) Hierarchical clustering of different liver tissues based on the amount of cells expressing foetal liver and gain-in-tumour enhancer-associated genes. Genes were further filtered based on the availability in the scRNA-seq data. (D) Selected pool of genes that showed significant correlation between its expression level and the serum AFP level. More genes are shown in Fig. S6. adj.normal, adjacent normal; AFP, alpha-foetoprotein; ChIP-seq, chromatin immunoprecipitation sequencing; eCluster, enhancer cluster; H3K27ac, H3K27 acetylation; HCC, hepatocellular carcinoma; scRNA-seq, single-cell RNA sequencing.

Fig. 6

The patient-specific upregulation pattern of gain-in-tumour enhancer-associated genes shows prognostic values. (A) Heatmap of patient-specific differential expression status for genes upregulated in at least 30% of the patients. The bottom heatmap shows the effect of co-clustering, and the columns are coloured by the assigned patient strata. (B) Disease-free survival plot and the risk table for the two patient clusters identified from co-clustering (Tarone–Ware test, p = 0.038). (C) Distribution of clinical parameters that were significantly different between the two patient clusters. The full breakdown of the clinical tables is shown in Table 1. (D) The cumulative event plots for patients with different Edmondson grades, grouped by the patient clusters. (E) Volcano plot of differential expression testing between PC0 and PC1 from RNA-seq data. (F) Heatmap of tumour-normal difference expression values of the epigenetic oncofoetal genes from TCGA-LIHC RNA-seq data, grouped by the patient clusters identified from co-clustering. (G) Progression-free survival plot and the risk table for the two patient clusters identified from TCGA-LIHC data (Tarone–Ware test, p = 0.035). PC, patient cluster; RNA-seq, RNA sequencing; TCGA-LIHC, The Cancer Genome Atlas Liver Hepatocellular Carcinoma.