Shira Perez,Anat Lavi-Itzkovitz,Moriah Gidoni,Tom Domovitz,Roba Dabour,Ishant Khurana,Ateret Davidovich,Ana Tobar,Alejandro Livoff,Evgeny Solomonov,Yaakov Maman,Assam El-Osta,Yishan Tsai,Ming-Lung Yu,Salomon M. Stemmer,Izhak Haviv,Gur Yaari,Meital Gal-Tanamy
High-Resolution Genomic Profiling of Liver Cancer Links Etiology With Mutation and Epigenetic Signatures
Hepatocellular Carcinoma (HCC),Hepatitis C Virus (HCV),Cancer Epigenetics,Cancer Genomics,Regional Variation in Mutation Rates
Background & Aims Hepatocellular carcinoma (HCC) is a model of a diverse spectrum of cancers because it is induced by well-known etiologies, mainly hepatitis C virus (HCV) and hepatitis B virus. Here, we aimed to identify HCV-specific mutational signatures and explored the link between the HCV-related regional variation in mutations rates and HCV-induced alterations in genome-wide chromatin organization. Methods To identify an HCV-specific mutational signature in HCC, we performed high-resolution targeted sequencing to detect passenger mutations on 64 HCC samples from 3 etiology groups: hepatitis B virus, HCV, or other. To explore the link between the genomic signature and genome-wide chromatin organization we performed chromatin immunoprecipitation sequencing for the transcriptionally permissive H3K4Me3, H3K9Ac, and suppressive H3K9Me3 modifications after HCV infection. Results Regional variation in mutation rate analysis showed significant etiology-dependent regional mutation rates in 12 genes: LRP2, KRT84, TMEM132B, DOCK2, DMD, INADL, JAK2, DNAH6, MTMR9, ATM, SLX4, and ARSD. We found an enrichment of C->T transversion mutations in the HCV-associated HCC cases. Furthermore, these cases showed regional variation in mutation rates associated with genomic intervals in which HCV infection dictated epigenetic alterations. This signature may be related to the HCV-induced decreased expression of genes encoding key enzymes in the base excision repair pathway. Conclusions We identified novel distinct HCV etiology-dependent mutation signatures in HCC associated with HCV-induced alterations in histone modification. This study presents a link between cancer-causing mutagenesis and the increased predisposition to liver cancer in chronic HCV-infected individuals, and unveils novel etiology-specific mechanisms leading to HCC and cancer in general.
High-Resolution Genomic Profiling of Liver Cancer Links Etiology With Mutation and Epigenetic Signatures

Hepatocellular carcinoma (HCC) is a model of a diverse spectrum of cancers because it is induced by well-known etiologies, mainly hepatitis C virus (HCV) and hepatitis B virus. Here, we aimed to identify HCV-specific mutational signatures and explored the link between the HCV-related regional variation in mutations rates and HCV-induced alterations in genome-wide chromatin organization.

To identify an HCV-specific mutational signature in HCC, we performed high-resolution targeted sequencing to detect passenger mutations on 64 HCC samples from 3 etiology groups: hepatitis B virus, HCV, or other. To explore the link between the genomic signature and genome-wide chromatin organization we performed chromatin immunoprecipitation sequencing for the transcriptionally permissive H3K4Me3, H3K9Ac, and suppressive H3K9Me3 modifications after HCV infection.

Regional variation in mutation rate analysis showed significant etiology-dependent regional mutation rates in 12 genes: LRP2, KRT84, TMEM132B, DOCK2, DMD, INADL, JAK2, DNAH6, MTMR9, ATM, SLX4, and ARSD. We found an enrichment of C->T transversion mutations in the HCV-associated HCC cases. Furthermore, these cases showed regional variation in mutation rates associated with genomic intervals in which HCV infection dictated epigenetic alterations. This signature may be related to the HCV-induced decreased expression of genes encoding key enzymes in the base excision repair pathway.

We identified novel distinct HCV etiology-dependent mutation signatures in HCC associated with HCV-induced alterations in histone modification. This study presents a link between cancer-causing mutagenesis and the increased predisposition to liver cancer in chronic HCV-infected individuals, and unveils novel etiology-specific mechanisms leading to HCC and cancer in general.

SummaryWe show a novel distinct genomic signature for hepatitis C virus (HCV) etiology in hepatocellular carcinoma and identified a correlation between mutation loads and HCV-induced changes in epigenetic marks associated with chromatin decompaction. Thus, this study illustrates a direct link between HCV etiology and cancer. Cancer evolution is a 2-step dynamic process. First, variants (typically DNA sequence changes) occur randomly, and then the frequency at which the variant is retained in the cell population is determined stochastically or by selection forces. Accordingly, sequence changes are categorized into driver mutations, which commonly are attributed to selection, and passenger mutations, which are not thought to be affected by selection. Passenger mutations appear at low frequency and comprise most of the mutations found in the cancer genome., Because passenger mutations, unlike driver mutations, are not subjected to processes of selection,1, 2, 3 the genomic regional variation in mutation rates of passenger mutations was suggested to bear an imprint of the mechanisms that generated them.,, Specifically, regional variation in mutations rates of passenger somatic mutations in cancer are not distributed uniformly across the genome and their clustering mimics that of driver mutations.5, 6, 7 Recent genome-wide analysis of diverse genetic and epigenetic data sets have indicated that regional variations in mutation rates are linked to replication timing, transcriptional activity, proximity to the nuclear periphery, recombination rate, DNA accessibility, and epigenetic chromatin organization.7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17 Epigenetic chromatin organization, a key factor in shaping regional variation in mutations rates, is controlled by the cumulative effect of histone subunit post-translational modifications and base modifications, such as CpG methylation., These epigenetic factors determine the level of chromatin compaction and, therefore, the accessibility of DNA to both mutagenic agents and the counteracting function of DNA repair mechanisms. The complex interplay between these 2 processes results in a distinct mutational pattern. Some studies have shown that the highest somatic mutation loads in cancer genomes were observed in inaccessible compact chromatin regions, whereas the lowest were present in accessible chromatin domains, owing to more efficient DNA repair.,,,,20, 21, 22, 23 Furthermore, the epigenetic chromatin organization of the cell-of-origin of a tumor was observed to be the best correlated with its regional variation in mutations rates. Previous studies have supported a correlation between compact chromatin and decreased mutations, either owing to less accessibility to the mutagens or in cases of impaired DNA damage repair.,,, Therefore, a specific mutation signature in cancer may point to etiology-specific mechanisms of mutational processes and DNA repair activity and may provide a deeper insight into the disease and avenues of early detection and prevention. Hepatocellular carcinoma (HCC) is the fifth-most common cancer worldwide and the second-highest cause of cancer death in men, and its incidence has been increasing over the past few decades.24, 25, 26 HCC serves as a model for a diverse spectrum of cancers because it is induced by a number of well-recognized etiologic agents, mainly hepatitis C virus (HCV) and hepatitis B virus (HBV) infections, and specific chemical agents, such as alcohol and aflatoxin., Therefore, this disease provides a unique platform that enables us to test whether different etiologic agents dictate distinct mutational landscapes across the cancer genome. Indeed, recent studies have sequenced liver cancer genomes to identify an etiology-dependent genomic signature in HCC., However, to date, no HCV-specific mutational signature was found.,29, 30, 31, 32 HCV is the most common risk factor for HCC in the Western world. It is a major public health problem with more than 70 million people infected worldwide and at risk for developing life-threatening liver disease. HCV causes hepatitis, cirrhosis, liver failure, and HCC. Although the oncogenic effect of chronic HCV infections has been widely studied, the sequence of events that follow from chronic viral infection to the key drivers of this malignancy, somatic mutations, is still not understood. Recently, we and others have shown that HCV infection induces genome-wide epigenetic alterations that reprogram host gene expression and also remain persistent after cure of infection.35, 36, 37 Thus, here, we aimed to investigate the relationship between HCV as an environmental exposure factor, cancer risk, epigenetic markers, and regional variations in mutation rates of cancer somatic mutations in HCC, and in cancer in general. We show that HCC cases have a regional variation in mutation rates that is distinct between patients with no viral history, those with HCV infection, or those with HBV infection. We further show that HCV-related HCC regional variations in mutation rates were superimposed on HCV-infection driven hepatocyte epigenetic chromatin organization.

We first set out to evaluate whether an HCV-specific mutation signature could be detected in HCC. Most of the mutations in cancer are low-frequency passenger mutations, rather than high-frequency driver mutations. Because passenger mutations, unlike driver mutations, are not subjected to processes of selection, they were suggested to bear an imprint of the mechanisms that generated them. We hypothesized that low-frequency mutations, which most probably are passenger mutations, may be accumulated in hepatocytes during progression to HCC and may provide insight into the effect of HCV etiology on cancer evolution. Therefore, we sought to identify low-frequency mutations by using a targeted sequencing approach that enables much higher sequencing depth than whole-genome or whole-exome sequencing approaches used in other studies that aimed to identify high frequency driver mutations in HCC and did not identify HCV-specific signatures.,29, 30, 31, 32, To design a targeted sequencing platform for HCC, we selected a set of genes that most frequently are mutated in HCC based on previous reports29, 30, 31, 32,, and those found in Catalogue of Somatic Mutations in Cancer (COSMIC) (http://cancer.sanger.ac.uk/cancergenome/projects/cosmic) and the International Cancer Genome Consortium (http://www.icgc.org) databases. In all, we chose 224 genes whose mutation rate in HCC samples was found to be ≥4%. We designed an Agilent SureSelect Target Enrichment panel, which includes the coding exons of 216 genes and 150-bp intervals from an additional 8 genes containing HCC-related mutations (total size, 1.044 million bp). The HCC-specific genomic mutations were studied in 64 human HCC biopsy samples from 3 etiology groups: HBV-infected (13 biopsy specimens), HCV-infected (28 biopsy specimens), or nonviral (23 biopsy specimens) (Table 1). The SureSelect platform that we designed was applied to sequence all HCC samples. Genomic DNA, extracted from sections of formalin-fixed, paraffin-embedded (FFPE), microdissected to include tumor cell islands, was used to construct indexed next-generation sequencing (NGS) libraries, captured by the custom-designed interval (described earlier), and was subjected to NGS on an Illumina HiSeq2500 at an average depth of ×2000. Aligned reads were compared with the reference genome, and results of variant calling were used to generate a list of mutations (variant call format) for each sample. Data analysis showed high sequencing depth (mean, 2000 reads per base) compared with previous reports,, and close to 100% coverage of each gene in the SureSelect design shown in a box-and-whiskers plot of read-depth distribution across all samples. The mutations that were found for each sample were filtered for quality, depth, and strand bias. To avoid germline variations from entering our HCC somatic mutations list, we exploited the predicted frequency of a germline allele within the sequencing of a normal diploid cell (typically 33%–66% of the reads, for a heterozygote). We thus divided the mutations into low-frequency mutations (found in <25% of variants in the sample), midfrequency mutations (found in 25%–75% of variants in the sample), and high-frequency mutations (found in >75% of variants in the sample). We considered low-frequency mutations, from which we removed all variants using the known single-nucleotide polymorphisms (SNPs) available in the SNP database (https://www.ncbi.nlm.nih.gov/projects/SNP), to provide a list of passenger somatic mutations for all HCC samples. To validate that the list of low-frequency mutations minus dbSNP represent the somatic mutations in the sample, we first sequenced peripheral white blood cells from 3 patients, which we considered a normal reference, and compared variant calling. Figure 1A shows an average overlap of >90% between our results using the reference genome as the normal sample of the same patient vs the results obtained when we removed dbSNP from the low-frequency mutations. Correspondingly, the low-frequency mutations minus dbSNP of the 3 normal DNA samples obtained from blood gave rise to an empty list; thus, there are no variants in the normal tissue, compared with a median of 40 mutations per sample in the cancer tissue of the same patient, validating that the low-frequency mutations are indeed somatic mutations (Figure 1B). Furthermore, the overall number of mutations that we found in the 64 HCC cases analyzed with that method was, on average, 10 times higher than the number of somatic mutations reported in the COSMIC database within our interval (Figure 1C). Figure 1D presents an example of map distribution of identified mutations for the gene CLSTN2 in all samples. Figure 1 shows the abundance of novel mutations that we discovered (in red), of which only a few were identified previously in the COSMIC database (in blue). These 3 observations suggest that looking at low-frequency mutations, enabled by the high sequencing depth of our targeted sequencing approach, ensures that the mutations taken into the analysis are low-frequency somatic mutations and not SNPs.

Next, using Euclidean distance based on the number of mutations per gene in a sample, we analyzed the ratio between the distances within samples of HCV etiology and the distances between samples from other etiologies, in comparison with a randomized ratio. Figure 2A shows a P value heat map of the intergroup/intragroup distances calculated by a permutation test. When all mutations or only medium- and high-frequency mutations found in each gene were included in the analysis we observed nonsignificant distinct patterns among the etiology groups (P = .322 and .698, respectively). In contrast, when considering only low-frequency mutations we found that the P value for the inner distance between HCV samples was significant (P = .031), compared with nonviral and HBV samples. This result supports the concept that low-frequency mutations, which most likely are passenger mutations, may indicate a genomic signature of etiology in cancer. If we combined HCV with either 2 groups together, the difference of distances between groups was no longer significant, neither between HBV compared with HCV and nonviral etiologies, nor when nonviral etiology was compared with HCV and HBV etiologies combined. For the following analyses we included only low-frequency mutations. We next investigated whether specific genes are enriched in mutations in HCV etiology compared with other etiologies. We calculated the expected number of mutations in each gene according to its relative size, assuming that mutations are spread uniformly among genes. Then, we calculated the probability of observing the actual number of mutations for that gene for each etiology. The heatmap in Figure 2B shows the calculated P value for each gene, with a specific pattern of mutated genes for each etiology. The Venn diagram in Figure 2C shows the intersect between the significantly enriched mutated genes for each etiology in the Poisson test (false discovery rate, FDR, <0.1). Our results define a novel HCV-specific genomic signature consisting of 12 genes: LRP2, KRT84, TMEM132B, DOCK2, DMD, INADL, JAK2, DNAH6, MTMR9, ATM, SLX4, and ARSD. Interestingly, the number of genes in the HCV etiology signature we identified is significantly higher compared with the other HCC etiologies tested (Figure 2C). We evaluated whether our mutation signature was identified previously in publicly available HCC data. Data from The Cancer Genome Atlas Liver Hepatocellular Carcinoma (TCGA–LIHC) project (https://portal.gdc.cancer.gov/projects/TCGA-LIHC), containing clinical information and sequences from 375 HCC samples, integrating the sequencing data with clinical phenotype and viral information (https://xenabrowser.net/datapages), was retrieved. We used maftools (https://www.bioconductor.org/packages/release/bioc/html/maftools.html) for visualizations and pairwise and groupwise Fisher exact tests to find differentially enriched genes for HBV (140 samples), HBV–HCV (77 samples), HCV (54 samples), and nonviral (104 samples) etiology groups. The 12-gene signature we identified was not mutated significantly in HCV compared with the other groups, although all had a higher FDR in TCGA–LIHC data set. Because we selected only low-frequency mutations that were identified as a result of our substantial higher sequencing depth (mean, 2000 reads per base) compared with previous data,,, our data do not include high-frequency driver mutations as previous published data do. Therefore, it is expected that the passenger mutation signature we found was not identified in TCGA–LIHC data set. Moreover, the data in TCGA-LIHC are sampled from different populations compared with our cohort. These results show that the 12 gene HCV-related mutation signature is unique to this study. Integration of the data with survival information showed lower survival of HCC patients with mutated DMD, INDAL, DOCK2, JAK2, and LRP2 genes compared with wild type, in which only DOCK2 was significant with a hazard ratio of 2.92 and a P value of .000753 (Figure 2D). Interestingly, the mutation frequency in the DOCK2 gene we observed in the liver samples, which is 39 of 64 (number of samples with mutations in the gene/number of all samples) and overall 61% (19 of 28 in HCV, 7 of 13 in HBV, and 13 of 23 in nonviral), was almost 10 times higher than observed in TCGA–LIHC data (25 of 369; 6.8%). This observation is consistent with the finding that with our high-resolution approach we identified a much higher number of mutations than in previous data. Beyond the categoric, etiology-related mutation signature, the rate of mutations showed enrichment of specific gene sets. A ranked gene list was generated for each etiologic group in the following way: first, for each sample, a ranked gene list was created according to the number of mutations. Next, for each gene, the average rank of that gene in all samples was calculated. This ranked list was used for Gene Set Enrichment Analysis (GSEA) (http://software.broadinstitute.org/gsea/index.jsp). GSEA analysis found enrichment of different gene sets for different etiologies that are known to be involved in cancer biology pathways as indicated in Figure 3A and B. Two gene sets were identified as significantly unique to the mutation pattern in HCV etiology. The gene set MARTINEZ_RB1_AND_TP53_TARGETS_DN is enriched positively (normalized enrichment score, NES = 2.05) in HCV, and the gene set SCHAEFFER_PROSTATE_DEVELOPMENT_48HR_UP is enriched negatively (NES = -0.78).

Identifying base substitution biases in cancer may indicate the mutagenic process underlying the malignancy.,, Previous reports that studied high-frequency driver mutation in HCC identified a specific base substitution pattern for HBV etiology, and other nonviral etiologies; however, no specific signature was found for the HCV etiology. We examined the 192 possible trinucleotide mutations in each etiology, including the flanking bases 5′ and 3′ to each substitution. When all samples were analyzed, we identified a novel signature of predominant C->T substitutions in HCV compared with HBV etiologies (t test at FDR, <0.1) (Figure 4A and B). Because most of the samples from HCV etiology are from patients who also had cirrhosis (87%), we only analyzed samples with a cirrhosis background from all etiologies to eliminate a potential effect of cirrhosis on the mutation pattern. After selecting for cirrhosis patients, the novel signature of predominant C->T substitutions in HCV compared with HBV or nonviral etiologies was maintained (t test at FDR, <0.02) (Figure 4C and D). These results show that this signature is HCV-specific. Not only is the type of base substitution specific for a distinct mutagenic process, but also the context of the mutation within specific flanking bases. By analyzing only cirrhotic samples, we identified higher T_A or C_A context in C->T substitutions in HCV compared with HBV or nonviral etiologies, however, this result was not statistically significant owing to an insufficient number of mutations (Figure 4E, black for HCV, red for HBV, and green for nonviral). For HBV etiology, we identified a predominant C->A substitution compared with HCV etiology, as identified in a previous study, but this also was not statistically significant (Figure 4A and C). It is important to note that a potential artifactual source for the observed HCV-associated C>T mutation rate are the C>T transversions, known to result from exposure to formalin during sample collection. This potential artifact does not agree with the observed results. First, the HBV and nonviral cases (all cases without HCV infection history) serve as controls for the HCV cases because those were collected under the same standard operating procedure, and yet their C>T variant allele frequency was significantly lower. Furthermore, as expected, the observed C>T variant allele frequency in the non-HCV cases was within the range of the FFPE-derived C>T variant allele frequency observed across more than a thousand FFPE cases sequenced in our genome center.

The hypothesis we aimed to test is that etiology-related changes in epigenetic chromatin organization after virus infection dictate specific mutation profiles. Thus, we explored the interplay between the genomic signature we identified in HCV etiology to the epigenetic signature induced by HCV infection. First, we performed chromatin immunoprecipitation followed by next-generation sequencing (ChIP-seq) of 9 HCV-infected normal liver biopsy samples and 3 noninfected control normal liver samples. To compare these samples in terms of histone modifications, we used lysine 9 (H3K9Me3) as a marker for closed chromatin. We have selected this specific marker after a previous report showed that it revealed a liver-specific epigenetic profile and was highly correlated with regional variation in mutation rates. We identified 2303 differential peaks in HCV-infected compared with noninfected samples (nominal P < .05). Of these, 1166 peaks were up-regulated and 1137 were down-regulated in HCV-infected samples compared with controls (Figure 5A). Even higher numbers of peaks annotated to gene promoters (±transcription start site [TSS] 5 kB) were identified, as follows: a total of 3239 differentially modified genes (DMGs), of which 1378 DMGs were up-regulated and 1861 were down-regulated in HCV-infected samples compared with controls (Figure 5C). These results point to a genome-wide effect of HCV infection on the H3K9Me3 profile in the liver. A heat map of the top DMGs show clustering of control samples vs HCV-infected samples informative for HCV-specific epigenetic alterations (Figure 5B and D). Next, we evaluated whether the 12 highly mutated genes in HCV etiology correlated with peaks of the epigenetic marker H3K9Me3 of closed chromatin. We integrated our analysis of mutated genes in HCV and our differential peaks of H3K9Me3 ChIP-seq data. Then, we tested the degree of difference between the mean log2 fold change (log2FC) in our set of highly mutated genes in HCV, compared with the mean of randomized log2FC taken from the SureSelect design genes in the H3K9Me3 ChIP-seq experiment (Figure 6A). We found that after HCV infection, the mean log2FC of H3K9Me3, which is an epigenetic marker for closed chromatin, is significantly lower in the HCV-etiology–highly mutated genes than the mean log2FC of this epigenetic marker in random genes from all 224 genes in the targeted sequencing design (P = .021) (Figure 6B). These results show that, on average, the level of H3K9Me3 is lower in genomic regions containing the 12 HCV highly mutated genes. Therefore, the mutations in HCV-related HCC are more abundant in chromatin structure that is less compacted in HCV-infected livers vs noninfected livers, compared with mutations in the rest of the design genes. To validate these results, we designed an in vitro experiment in cell culture, in which HCV-infected and noninfected Huh7.5 cells (6 pairs: 3 control and 3 HCV) were treated with the mutagen ethyl methane sulfonate (EMS). To identify mutated genes, we then sequenced the treated cells using the targeted sequencing approach with the SureSelect Target Enrichment panel as described earlier for the liver biopsy specimens. We identified mutations that appeared only after EMS treatment. We generated a gene list that is ranked according to the number of mutations identified after EMS treatment in HCV vs control and P value. The top 5 up genes with the lowest P value and more mutated in HCV after EMS vs control were identified and included the following: ARID4B, PREX2, PAM, RB1, and MIR383. We also defined the top 5 down genes with the lowest P value and less mutated in HCV after EMS vs control, which included the following: PDE4DIP, NIPBL, RIMS1, TLK1, and EPHA6. Interestingly, the gene DOCK2, which was within the 12 HCV highly mutated gene signatures in the liver biopsy and showed lower survival in patients with a mutated gene (Figure 2D), also is included in the top 10 up genes in the cell culture experiment. We then aimed to identify potential correlation between changes in epigenetic chromatin organization after virus infection and the specific mutation profiles. Therefore, we used ChIP-seq data that we recently reported to allocate epigenetic changes after HCV infection in Huh7.5 cells that included histone H3 acetylated on lysine 9 (H3K9Ac) and trimethylated lysine 4 (H3K4Me3) as markers for opened chromatin, and H3K9Me3 as a marker for closed chromatin. Differential peaks for either H3K4Me3, H3K9Ac, or H3K9Me3 were filtered for the 224 sequenced genes in our design. We integrated our analysis of differentially mutated genes (up and down) in HCV-infected cells and our ChIP-seq data experiments on HCV-infected cells. We tested the degree of difference between the mean log2FC in our set of differentially mutated genes, compared with the mean randomized log2FC taken from the SureSelect design genes in each ChIP-seq experiment, as described earlier for liver biopsy specimens. We found that after HCV infection, the mean log2FC of H3K4Me3 and H3K9Ac, which are epigenetic markers for open chromatin, was higher in the highly mutated genes in HCV than the mean log2FC of these epigenetic markers in 5 randomly sampled genes of the 224 sequenced genes in the targeted sequencing design (Figure 6C and D). In contrast, the mean log2FC of H3K9Me3, which is an epigenetic marker for closed chromatin, was lower in the highly mutated genes in HCV infection than the mean log2FC of this epigenetic marker in 5 randomly sampled genes of the 224 sequenced genes in the targeted sequencing design (Figure 6C and D), as was found in the analysis of the integration of epigenetics and mutation signatures in the liver biopsy samples. The opposite results were obtained for the analysis of lowly mutated genes in HCV infection: the mean log2FC of H3K4Me3 and H3K9Ac was found to be lower in the lowly mutated genes in HCV and the mean log2FC of H3K9Me3 was higher in the lowly mutated genes in HCV infection than the mean log2FC of this epigenetic marker in all the genes in the design (Figure 6E and F). The combined probability for all 6 tests (up and down genes, each evaluated for 3 epigenetic markers) showed statistical significance (P = .016) between lower chromatin compaction induced by HCV (lower level of H3K9Me3 and higher levels of H3K4Me3 and H3K9Ac) and higher mutation load, and vice versa. Collectively, these results show that mutations in HCV infections are more abundant in chromatin structure that has been opened after HCV infection. To validate these results, we obtained data from ATAC sequencing of HCV-infected and noninfected Huh7.5 cells that assessed chromatin accessibility. Chromatin openness at gene bodies of highly and lowly mutated genes and for all of the other 224 genes was calculated. A significantly lower number of peaks was identified for the lowly mutated genes in HCV-infected cells compared with noninfected, while no difference was observed between the 2 groups in all of the other genes and the highly mutated genes (Figure 6G). These results validate that a lower mutation load is correlated with high chromatin compaction induced by HCV.

The finding that mutation loads in HCV-related HCC with chromatin regions that are decompacted after HCV infection are correlated is consistent with previous studies showing that high mutation loads in various cancers were observed in open chromatin regions in cases in which the action of DNA repair mechanisms are impaired.,,, Therefore, we hypothesized that this HCV etiology–specific pattern may be formed by the known down-regulation of DNA repair mechanisms after HCV infection. Because we observed higher C->T mutations in HCV-related HCC, we analyzed the base excision repair pathway, which is the primary mechanism that removes base lesions and usually repairs deamination of cytosine residues to uracil, which results in C->T mutations, and evaluated whether it is altered after HCV infection. We have subjected the outlier RNA-seq differentially expressed genes and H3K9Ac ChIP-seq DMG data in HCV-infected cells to gene ontology analysis using the MetaCore database from Thomson Reuters (version 6.11, build 41105, GeneGo, Thomson Reuters, New York, NY) using the pathway analysis tool. This tool enables us to perform gene network enrichment analysis from integration of RNA-seq and H3K9Ac ChIP-seq data. By using the RNA-seq and H3K9Ac ChIP-seq data, we generated ranked lists of DMGs for ChIP-seq signal and differentially expressed genes (for RNA-seq signal) (ranked from the most significant decrease to the most significant increase). The 2 gene lists were uploaded to the MetaCore database, which provided an integrated list of pathways changed both in epigenetic markers and transcription. Analysis of the base excision repair pathway genes presented in Figure 7 shows that both RNA levels and H3K9Ac levels generally are down-regulated by HCV infection.

This study provides evidence for a direct link between events induced after HCV infection and HCC, which usually occurs after many years of chronic HCV infection, and offers insight into etiology-induced cancer evolution. Although it is well known that HCV is one of the most common risk factor for HCC, one of the fundamental issues is whether HCV is directly carcinogenic, or whether the increased risk of cancer is related only to indirect mechanisms associated with the inflammation and cirrhosis it causes. Evidence suggests that both indirect host-mediated and direct HCV-mediated mechanisms contribute to HCV-induced HCC, but the direct link still is inconclusive (reviewed by Lemon and McGivern, Rusyn and Lemon, and Shlomai). It has been suggested that because both HCV and HBV induce chronic inflammation that contributes to HCC, significant differences in genomic profiles between these 2 etiologies would support an association between cancer and direct viral carcinogenesis. However, to date, a specific genomic signature has been reported in HBV-related HCC, but not in HCV-related HCC, except for an association with driver mutations in TERT and CTTNB1 genes.,29, 30, 31, 32, Here, we show a novel distinct genomic signature for HCV etiology in HCC based on passenger mutations. This genomic signature consists of a panel of 12 genes that are enriched significantly in mutations in HCV-related HCC compared with other etiologies. This passenger mutation signature was identified with high sequencing depth (mean, 2000 reads per base) compared with previous reports that considered driver mutations only. Because we analyzed only low-frequency mutations we were able to accurately identify somatic mutations without the need to analyze matched tumor surrounding tissue to subtract SNPs as in previous reports.29, 30, 31, 32,,49, 50, 51 Therefore, by using this approach we may have avoided excluding somatic mutations found in cancer tissue that also may have occurred in the surrounding tissue that appears normal as a result of the field effect reported in HCC,, which further increased the sensitivity and the overall mutation numbers we identified. Our findings show the significance of high-resolution sequencing and identify low-frequency mutations for detecting etiology-specific mutation signatures. Interestingly, the DOCK2 gene was found to be among the top mutated genes in the HCV-related HCC gene signature as well as in the cell culture experiment, and mutations in this gene correlated with lower survival in HCC patients. DOCK2 is known to be involved in cytoskeleton remodeling and to play a role in the proliferation and migration of lymphocytes. Although DOCK2 was identified as a prognostic factor in various cancers, its role in HCC only now is being acknowledged. Two recent studies showed a significant correlation between HCC and high mutation frequency in the DOCK2 gene. Analyzing HCC cohorts obtained from TCGA and the International Cancer Genome Consortium databases containing 695 patients in total found that DOCK2 was among the 50 most frequently mutated genes identified. Importantly, the mutation status in 4 genes, including DOCK2, was found to be a prognostic factor for the overall survival of HCC patients. In another recent study, mutations in the DOCK2 gene also correlated with a low survival rate and a poor prognosis in HCC. Importantly, we found a higher number of mutations, and, specifically, a novel signature of a significantly higher number of C->T mutations in HCV etiology compared with HBV and nonviral etiologies, that was not identified in prior studies.,, The C->T mutation signature was found more often when C was surrounded by T and A in HCV etiology compared with HBV and nonviral etiologies. This mutation signature matches the known TC(A|T)→(T|G) apolipoprotein B messenger RNA editing enzyme, catalytic polypeptide signature (APOBEC) signature identified in human cancer genomes., APOBECs are a family of cytosine deaminases that have a role in antiviral innate immune response (reviewed by Krokan et al), and have been shown to restrict several viruses by a hypermutation mechanism, including human papilloma virus and HBV. Interestingly, the APOBEC mutation signature was found in cervical cancer, which is known to be caused by human papilloma virus. Increased expression of APOBEC3G was shown in chronic HCV infections, suggesting a role for APOBEC in the HCV-specific mutation signature in HCC. The C->T transition also was identified in internal tumors as a signature mutation for reactive oxygen species (ROS), which is produced in response to oxidative stress. Although chronic inflammation contributes to the generation of ROS, it is unlikely that this is the mechanism that is accountable for the higher C->T mutations in HCV etiology because chronic inflammation is common in both chronic HCV and HBV infections. Supporting the direct effect of HCV on this genomic signature is substantial evidence for the induction of oxidative stress and ROS directly by HCV proteins. Both HCV core and NS5A proteins were reported to generate ROS in vitro and in vivo and to produce high hepatic levels of 8-oxo-20-deoxyguanosine, which generate mutations. Moreover, significantly higher hepatic levels of 8-oxo-20-deoxyguanosine were found in HCV patients than in HBV patients, correlating with the higher C->T mutations we observed in HCV-related HCC compared with HBV-related HCC. This implies an HCV-specific mechanism of mutagenesis in HCC. In the cell line experiment, on the other hand, we used the exogenous mutagen EMS to accelerate the mutation process, which may account for the pattern differences in the highest mutated genes we identified in the cell lines compared with HCC biopsy specimens. Previously, Polak et al showed that the epigenetic state of the cell of origin determines regional variations in mutation rates in cancer in tissue differentiation levels. Here, we expand this paradigm to environmental etiologic agents, which also redirect the mutation rate by reorganizing the chromatin. This was enabled by the high-resolution approach we used, looking at low-frequency mutations that are maintained in the cancer genome, reflecting the events that occurred after virus infection and throughout cancer initiation and progression. We further observed a correspondence between the level of chromatin compaction formed after HCV infection and the HCV-specific (passenger) mutation profile in HCC. We show that mutation loads are found in HCV-induced decompacted chromatin regions and vice versa, both in liver tissues and in vitro experiments. These findings counteract some studies that have provided evidence for higher mutation loads found in cancer genomes correlating with more compacted chromatin regions owing to lower accessibility of DNA repair machineries, but are in line with cases of DNA damage repair deficiencies.,,, Specifically, signatures associated with base excision repair deficiency as well as APOBEC signatures correlated with open chromatin. Therefore, we suggest that the mechanism that induces this specific signature may be related to the activation of the mutagenic process, possibly APOBEC and/or ROS, and the down-regulation of the base excision repair pathway by HCV, which usually repairs deamination of cytosine residues to uracil caused by APOBEC and ROS., Our striking observation that epigenetic alterations remain as epigenetic signatures after curing HCV with direct-acting antivirals (DAAs) suggests that because an epigenetic memory mechanism maintains these events, they may be reflected in the cancer genome, even if the virus is no longer present in the cells. These findings suggest that the cancer may arise from the HCV-infected cells because the presence of the virus in the cells is required to induce the epigenetic alterations. Collectively, these observations point to a link between the HCV-induced epigenetic signature and the HCV-specific mutation signature in HCC and opens a new direction for further study of HCV pathogenesis. It is important to note that the observations reported here were limited to the cohort size and Caucasian ethnicity used in this study, and should be validated in a larger sample size and in other ethnicities. In summary, we show a novel genomic signature for HCV in HCC that links epigenetic chromatin organization and regional variations in mutation rates. These findings enhance our understanding of etiologic effects on cancer evolution, and may open new avenues for HCC treatment and prevention. The potential plasticity of the epigenetic chromatin organization together with our observations highlight its role in shaping the mutation landscape in HCC on the background of chronic HCV infections, and supports the utilization of drugs that alter the epigenetic structure to revert the oncogenic effects of chronic hepatitis infections and to prevent HCC. Moreover, by looking at mutations that appear in a minority of patient cancer cells, we observed a larger mutation load per sample than in previous reports.,29, 30, 31, 32 Because the best biomarker of response to immunologic checkpoint inhibitors is mutation load, a higher mutational load in HCV-related HCC could support the development and use of immune-related therapy for subgroups of HCC. This possibility may explain the recent encouraging clinical data for immunologic checkpoint inhibitors in HCV-related HCC cases,70, 71, 72 but would require further clinical investigation.

For mutation signature analysis, FFPE blocks from liver biopsy specimens of HCC patients were obtained from the pathology archive of Rabin Medical Center in Israel. Samples included patients representing 3 etiologic groups: 13 samples from patients with HBV infection, 28 samples from patients with HCV infection, and 23 samples from patients with no history of hepatitis virus infection. All patients were of Caucasian ethnicity. Clinical data are described in Table 1. DNA purification was performed using Maxwell 16 FFPE Plus LEV DNA Purification Kit (Promega) from microdissected malignant regions, marked by a pathologist (early stage lesions). For ChIP-seq analysis, resected liver tissues were obtained from the Rabin Medical Center in Israel and the Kaohsiung Medical University in Taiwan. In total, 12 human liver biopsy specimens were collected from 9 HCV-infected livers and 3 control noninfected normal liver tissues. All samples were cryopreserved before being used for ChIP assays.

Huh-7.5 cells (a generous gift from Charles Rice, Rockefeller University) and Huh7/FT3-7 cells (a generous gift from Stanley M. Lemon, University of North Carolina at Chapel Hill) were cultured as described previously. To generate hepatocyte-like cells, Huh-7.5 cells were maintained in 2% human serum for 2–4 weeks before use, as previously described.,

HJ3–5 is a H77C/JFH1 chimeric virus containing core-NS2 proteins from genotype 1a and NS3-NS5B from genotype 2a (a generous gift from Stanley M. Lemon, University of North Carolina at Chapel Hill). Virus stocks were produced in Huh7/FT3-7 cells and viral titers were determined by focus forming unit (FFU) assay in Huh-7.5 cells, as described previously.

The optimal noncytotoxic concentration of EMS was evaluated by cell viability and proliferation XTT assay and was determined to be 0.64 mg/mL. Huh7.5 cells were infected by the HJ3–5 chimeric virus with a multiplicity of infection of 0.1 and passaged for 2 weeks until approximately 100% of the cells were HCV positive (evaluated by immunostaining for viral proteins). HCV-infected and control noninfected cells were treated with EMS for 4 days and DNA was purified using the QIAamp DNA minikit (Qiagen, Valencia, CA).

For detection of mutations in liver biopsy specimens and EMS-treated cell lines, DNA libraries were prepared using the NEBNext Ultra DNA Library Prep Kit for Illumina (NEB). The SureSelect DNA Design (Agilent Genomics) contains 19,400 probes against 828 target IDs comprising 5421 regions for 224 genes of interest, with a total region size of 1.044 million bp. The SureSelect Kit biotinylated 120 mer RNA library baits were hybridized with the multiplexed DNA libraries and the target-enriched genomic DNA library was captured by streptavidin beads. The captured DNA was polymerase chain reaction–amplified and sequenced. NGS was performed using Illumina HiSeq 2500 platforms with paired-end reads of 250 bp according to the manufacturer's instructions.

Paired-end reads were aligned to the human genome 19 hg19/GRCh37 reference genome using Burrows-Wheeler Aligner (BWA). Mutations in tumors were predicted by FreeBayes (https://arxiv.org/abs/1207.3907v2) and called variants were filtered further for read quality (Q, >20), minimum read depth (>200), and mutation frequency (<25%). An additional step of retrieving variants at a position where a variant already was identified, by comparing high-frequency variants with the identified variant, resulted in a list of additional low-frequency variants. These variants together with the filtered low-frequency variant data (<25%) were referred as the overall identified low-frequency somatic mutations data. Data were normalized by the length of each gene and by the number of mutations per sample.

To compare somatic mutation calling with and without a normal sample, 3 matched tumor–normal paired samples were used. For each of these samples, the normal sample (from blood) was compared with the reference genome, and the variants present at a frequency >25% were used to build a new reference genome on top of hg19 (after subtraction of the known SNPs available in the dbSNP database (https://www.ncbi.nlm.nih.gov/projects/SNP). This normal-reference then was used as a new reference for tumor sample variant calling. The low-frequency variants from this analysis (<25%) were compared with variants identified from directly comparing tumor samples with the hg19 reference genome (<%25).

The observed number of mutations was tested against the expected number in each gene under the assumption of a Poisson distribution and constant mutation rate for single nucleotide variants across genomic regions. Correction for multiple testing was performed by FDR. Highly mutated genes in HCV/HBV were defined as those with an FDR <0.1 for the Poisson test.

Sequencing and clinical data were obtained from TCGA–LIHC project containing 375 samples that included HBV (140 samples), HBV–HCV (77 samples), HCV (54 samples), and nonviral (104 samples) etiology groups (maf download https://portal.gdc.cancer.gov/projects/TCGA-LIHC), and clinical phenotype with viral information (https://xenabrowser.net/datapages). We performed maftools for visualizations. Using the clinicalEnrichment function command line, we performed pairwise and groupwise Fisher exact tests to find differentially enriched genes for etiology groups. Survival analysis was performed by the mafSurvive function command line.

A ranked gene list was created for each etiologic group as follows. First, for each sample, a ranked gene list was created according to the number of mutations. Second, for each gene, the average rank of that gene in all samples was calculated. This ranked list was used for GSEA (http://software.broadinstitute.org/gsea/index.jsp), with the following parameters: 10 permutations, c2 gene set collection (>4000 gene sets), and a minimum of 5 genes per gene set.

The number of mutations in the matrix of 192 possible mutations occurring in a trinucleotide context (C>A, C>G, C>T, T>A, T>C, T>G; somatic mutations presented in a pyrimidine context, and included the bases immediately 5' and 3' to each substitution) was determined for each sample. Frequencies obtained by dividing each mutation by the total number of mutations for that sample were used for a t test between etiologies.

Raw reads were trimmed using trimmomatic to remove adapter's low-quality bases. Next, paired-end reads were aligned to the hg19/human_g1k_v37 using BWA. Mutations were predicted by FreeBayes (https://arxiv.org/abs/1207.3907). Called variants were filtered further by variant call format filter with the following parameters: QUAL > 1 and QUAL / AO > 10 and SAF > 0 and SAR > 0 and RPR > 1 and RPL > 1. Insertions and deletions (INDELs) were removed and further filtering was based on coverage, SNPs with coverage below 200 were discarded, which resulted in 98 TrueSeq on-target SNPs being evaluated further. For each pair of pre– and post–EMS-treated cell lines (6 pairs: 3 control and 3 HCV), the number of mutations for each gene was set as the number of SNPs that appeared only after EMS treatment. This number was normalized according to the total number of SNPs after treatment in each pair of cell lines.

By applying makeComparison and aggregateGeneSet functions from the Qusage method on normalized counts, the t-distributions for each gene were defined, which were compared between control and HCV groups. Genes were sorted according to P value, and the top 10 genes were chosen: the top 5 up genes with the lowest P value and more mutated in HCV after EMS vs control, and the top 5 down genes with the lowest P value and less mutated in HCV after EMS vs control.

ChIP assays for resected liver tissues were performed using the iDeal ChIP seq kit for Histones (Diagenode) according to the manufacturer’s instructions. For each sample, 20–40 mg liver tissue was chopped, homogenized, and subjected to ChIP using H3K9Me3 (Active Motif) or a rabbit IgG antibody as a negative control (Diagenode). Libraries for NGS were prepared by using the NEBNext Ultra DNA Library Prep Kit for Illumina, and sequenced using the Illumina HiSeq 2500 platform with single-end reads of 50 bp according to the manufacturer's instructions. Inputs for each biopsy specimen also were sequenced. ChIP-reads (75 bp) were trimmed for adapters and quality. The sequences were aligned to the GRCh37 (hg19) reference genome using BWA-MEM (0.7.17, default parameters). Duplicate alignments were removed (Samtools). Peaks on gene promoters (±TSS 5 kB) were generated and counted using featureCounts (v1.6.0), with a minimum alignment quality threshold of 20. De novo peaks were generated for each sample using MACS2, with an input sample used to define the background (MACS2 - 2.1.1). Peaks were merged, and filtered for regions with a fold-enrichment greater than 4 that were used for differential analysis. Aligned reads within the consensus peak regions were counted using featureCounts (v1.6.0), with a minimum alignment quality threshold of 20 and annotated to hg19 Ensembl Transcripts release 75. All differential analysis of HCV-infected vs noninfected was performed using edgeR.

For integration of the mutation signature and H3K9Me3 ChIP-seq data from liver biopsy specimens, the differential peaks that were identified for H3K9Me3 ChIP-seq data were filtered for the 224 sequenced genes in our SureSelect design. Additional filters included peaks that were >10,000 bp away from these genes and -0.2 < log2FC < 0.2. Log2FC of the enriched genes was computed from the mean of the median values of these genes in ChIP and sampling the mean log2FC of the size of the enriched HCV genes 10,000 times. The P value refers to the expected direction of values for each ChIP-seq. For closed chromatin, the test was as follows: P (X < x). For integration of the mutation signature and H3K4Me3, H3K9Ac, or H3K9Me3 ChIP-seq data from the cell line experiment, a similar procedure was performed as was applied on ChIP-seq data from the liver biopsy specimen experiment, except the log2FC cut-off value was set to -0.1 < log2FC < 0.1. For the open chromatin markers the test was as follows: P (X > x); and for the closed chromatin marker the test was as follows: P (X < x). For integration of the mutation signature and ATAC-seq data from the cell line experiment, first, ATAC-seq data of HCV-infected and noninfected Huh7.5 cells were obtained from Hlady et al. Reads were mapped to the human genome assembly hg19 using Bowtie (version 1.1.2). Then, chromatin openness at gene bodies was calculated by quantifying the number of reads per kilobase per million mapped reads within a gene. Gene body coordinates were defined using RefSeq gene annotations taken from the University of California Santa Cruz Genome Browser database. Data visualization was performed using Rstudio (version 2021.09.2) (http://www.rstudio.com). A 2-way paired t test was used to assess the significance of differences in the ATAC-seq signal over genes between HCV-infected and noninfected cells.

The pathway analysis was generated from integrated RNA-seq and ChIP-seq data. Using the RNA-seq and H3K9Ac ChIP-seq data, we generated ranked lists of DMGs (for the ChIP-seq signal in the vicinity of a specific gene; for peaks within the gene body, both exons or introns and within TSS, up to 5-kb upstream of the TSS) and differentially expressed genes (for RNA-seq signal) (ranked from the most significant decrease to the most significant increase). The 2 gene lists were used for gene ontology analysis using the MetaCore database from the Thomson Reuters (version 6.11, build 41105, GeneGo; Thomson Reuters, New York, NY) pathway analysis tool. This pathway analysis tool was used to perform gene network enrichment analysis from integration of RNA-seq and H3K9Ac ChIP-seq data.

FreeBayes (https://arxiv.org/abs/1207.3907v2) was used to predict mutations in tumors. Mutations were normalized by the length of each gene and by the number of mutations per sample. Poisson distribution was used to calculate highly mutated genes in biopsy samples, with FDR <0.1 considered statistically significant. To identify highly and lowly mutated genes in cell lines after EMS treatment we applied the makeComparison and aggregateGeneSet functions from the Qusage method. Pairwise and groupwise Fisher exact tests to find differentially enriched genes for etiology groups from TCGA–LIHC project data were used. For calculation of frequencies of base substitutions between etiologies the t test was used, with FDR <0.02 considered statistically significant. For the ChIP-seq experiment, peaks on gene promoters (±TSS 5 kb) were generated and counted using featureCounts (v1.6.0). De novo peaks were generated for each sample using MACS2, with an input sample used to define the background (MACS2 - 2.1.1). Aligned reads within the consensus peak regions were counted using featureCounts (v1.6.0). All differential analysis of HCV-infected vs noninfected was performed using edgeR. Euclidean distance was used to analyze the ratio between the distances within samples of HCV etiology and the distances between samples from other etiologies, in comparison with a randomized ratio. Significance was tested by a permutation test.

FFPE blocks from liver biopsy specimens (for mutation signature analysis) of HCC patients were collected using protocols approved by the Rabin Medical Center Institutional Review Boards. Resected liver tissues (for ChIP-seq analysis) were obtained from the Rabin Medical Center in Israel and the Kaohsiung Medical University in Taiwan using protocols approved by the Institutional Review Boards. All human tissues were collected in accordance with the ethical standards of the Declaration of Helsinki. Written informed consent was received before participation. All authors had access to the study data and reviewed and approved the final manuscript.