|Home | About | Journals | Submit | Contact Us | Français|
Hepatocellular carcinoma (HCC) is the third-leading cause of cancer-related deaths worldwide. It is often diagnosed at an advanced stage, and hence typically has a poor prognosis. To identify distinct molecular mechanisms for early HCC we developed a rat model of liver regeneration post-hepatectomy, as well as liver cells undergoing malignant transformation and compared them to normal liver using a microarray approach. Subsequently, we performed cross-species comparative analysis coupled with copy number alterations (CNA) of independent early human HCC microarray studies to facilitate the identification of critical regulatory modules conserved across species.
We identified 35 signature genes conserved across species, and shared among different types of early human HCCs. Over 70% of signature genes were cancer-related, and more than 50% of the conserved genes were mapped to human genomic CNA regions. Functional annotation revealed genes already implicated in HCC, as well as novel genes which were not previously reported in liver tumors. A subset of differentially expressed genes was validated using quantitative RT-PCR. Concordance was also confirmed for a significant number of genes and pathways in five independent validation microarray datasets. Our results indicated alterations in a number of cancer related pathways, including p53, p38 MAPK, ERK/MAPK, PI3K/AKT, and TGF-β signaling pathways, and potential critical regulatory role of MYC, ERBB2, HNF4A, and SMAD3 for early HCC transformation.
The integrative analysis of transcriptional deregulation, genomic CNA and comparative cross species analysis brings new insights into the molecular profile of early hepatoma formation. This approach may lead to robust biomarkers for the detection of early human HCC.
Hepatocellular carcinoma (HCC) is the fifth most common cancer type, and is the third leading cause of cancer mortality worldwide [1,2]. Recent reports show that HCC is becoming more wide-spread and has dramatically increased in North America Western Europe and Japan [2-4]. Additionally, there is an increasing incidence of the disease among younger age groups that warrants further investigation [5,6].
Recently considerable attention has been placed on global gene expression studies as well as genomic aberrations in order to understand the pathogenesis of HCC, and to look for possible early markers of detection [7-13]. Although notable successes have been achieved, there still exist significant challenges due to the heterogeneous nature of HCC (and other cancers) as well as the complexity of the molecular pathogenesis of this disease. Depending on etiological and accompanying pathological conditions, such as viral infection, cirrhosis, inflammation, fibrosis and others, the HCC signature genes identified thus far vary considerably. Additionally, the study of tumor formation in the liver is difficult due to the continuous transcriptome changes that occur during regeneration after hepatectomy , as well as age related gene expression changes [15-17]. Similarly, cancer progresses through a series of histopathological stages during which genetic alterations accumulate, and a natural consequence of this are the dynamic changes in gene expression patterns that occur during hepatocellular carcinogenesis. Developing animal models of HCC provide an experimental ground for dissecting the genetic and biological complexities of human cancer and contribute to our ability to identify and characterize pathogenic modifications relevant to early stages of cancer development and progression [18,19]. The previous studies have used cross-species comparative genomics approach successfully to understand the molecular pathogenesis of various cancers [20-23]. Hence, combining cross-species comparative and/or functional genomics approaches with independent datasets from human and animal models of HCC along with genomic DNA copy number alterations enhances the ability to identify robust predictive markers for HCC [23-26].
Here we present a comparative and integrative functional genomics approach to find an early marker for HCC. We developed a rat model and analyzed the transcriptomes of early HCC versus regenerated liver and normal liver in both young and old age animals using a microarray of more than 27,000 annotated genes from Celera and public repositories. We then performed cross-species comparative genomics analysis to identify genes that are conserved in rat and human early HCCs by re-analyzing independent datasets for human early HCC microarray expression profiling data [8,27] and also comparing with the Stanford HCC microarray data . Finally, we performed an integrative analysis of DNA genomic copy number alterations (CNAs) and gene expression profiles (schematically outlined in Figure Figure1).1). Our findings include some genes already reported to be associated with human HCC, thus validating our approach. We also report many other novel genes which were not reported previously in liver cancers. Furthermore, we validated the high expression of eight potential biomarker genes from the blood of patients with early HCC using realtime RT-PCR. Our comparative and integrative genomics approach involving the integration of multiple high dimensional independent datasets may lead to robust biomarkers for the detection of early HCC.
We performed genome-wide gene expression profiling of 24 samples for early HCC, regenerated liver and normal liver of both young and old rats using Applied Biosystems Rat Genome Survey microarray which includes more than 27,000 annotated genes from Celera and public repositories. To find genes that were differentially expressed across three different "treatment" types (i.e. early HCC, regenerated, and normal), and two age groups (old and young), we performed two-factor ANOVA to look for variations due to treatment, age and their interactions. The ANOVA identified 432 and 4063 genes that were significantly modulated by treatment type and age with p < 0.01, respectively. In addition, we found 322 genes that showed a significant interaction of age and treatment effect (data not shown). The unsupervised two-dimensional hierarchical clustering as well as principal component analysis (PCA) using genes which varied significantly with the treatment effect clustered samples according to their treatment type for both old and young (Figure (Figure2A2A and and2B),2B), hence supporting the conclusion that gene expression profiles robustly reflected the histological classification.
The ANOVA identified 432 genes that showed significant expression differences due to treatment in both age groups (p-value < 0.01) were subjected to a template matching algorithm (TMA)  to identify HCC specific genes conserved across both age groups. We identified 96 up-regulated and 38 down-regulated genes specific to HCC (template-match p-values < 0.01) (Figure (Figure2C2C and Additional file 1).
The gene ontology and functional network analyses of HCC specific genes were performed using the Ingenuity knowledge base. The biological functions assigned to the dataset are ranked according to the significance of that biological function to the dataset. The enriched functional categories and diseases include carcinogenesis, cell cycle, immune response, cell morphology, cellular development, and growth and proliferation (Figure (Figure2D).2D). The PANTHER also revealed that signal transduction (p-value = 7.17 × 10-4), proteolysis (p-value = 1.59 × 10-3), cell motility (p-value = 3.91 × 10-3), immunity and defense (p-value = 1.07 × 10-2), and cell proliferation and differentiation (p-value = 4.9 × 10-2) were among the most enriched biological processes in the HCC specific genes. The most significantly altered pathways include p53, p38 MAPK, insulin/IGF pathway-protein kinase B signaling cascade, apoptosis, interleukin and integrin signaling pathways. The gene interaction network also corroborated with the altered pathways (Figure (Figure2E2E).
We found a significant interaction of age and treatment effect in our two-factor ANOVA analysis. Indeed, the hepatic transcriptome changes with ageing as in other cancer types, and age is a potential confounding factor embedded in gene expression profiles [15,16]. Therefore, we stratified samples as young and old cohorts, and identified HCC and regeneration specific genes using one-way ANOVA in each age group separately. The ANOVA identified 925 and 408 significantly dysregulated genes (up or down) due to the different treatment types in old and young animal groups (p-value < 0.02), respectively. The hierarchical clustering in both dimensions (samples and genes), as well as the PCA clearly separated samples based on the treatment type (Additional file 2), The gene expression clustering distance between the HCC group and the other two groups (regenerated and normal) was the greatest in both age groups (Additional file 2).
Early HCC signature genes in each age group were obtained by overlapping gene lists. Each circle in the Venn diagram represents the differential expression between two treatment types. We identified 80 genes and 100 genes specific to hepatoma in young and old, respectively (Figure (Figure3A3A and Additional file 3). As seen from the heatmap of HCC specific genes, these sets of genes were exclusively up/down regulated in the HCC group only (Figure (Figure3B3B and Additional file 3). The expression of Pbsn, Lum, Adam8, Ctse, Calb3, Fbn1, Agtpbp1, Prom1, Ela1, Tnfsf13, and Ap2b1 were significantly altered in both young and old rats with early HCC. The genes A2m, Cdh13, Mas1, Slack, Cidea, and Dcn were significantly dysregulated exclusively in HCC in the young; whereas Cxcl5, Lox, Slc25a2, Rmt1, and Nid2, were specific to HCC in old rats. We also identified regeneration specific genes in old and young rats in a similar approach (data not shown).
Early HCC signature genes in young animals were mainly associated with cancer, cell cycle, immune response, cellular function and maintenance, development, cell adhesion-mediated signaling, proteolysis, and signal transduction, whereas genes in old animals were highly associated with cancer, cellular movement, extracellular transport and import, cell adhesion, tissue development, cell morphology, cell-to-cell signaling and interaction. On the other hand, the regeneration specific genes were mainly associated with mRNA transcription and regulation, lipid metabolism, protein modification, protein phosphorylation, cell morphology, cellular development, small molecule biochemistry, and cellular growth and proliferation (Table (Table1).1). To gain more insights into HCC pathogenesis for young and old age early HCCs, we carried out gene interaction networks of HCC specific genes in young and old (Figure 3C, D and and3E3E and Additional file 3, respectively). The interaction networks highlight the important role of p53, p38 MAPK, ERK/MAPK, PI3K/AKT signaling, NF-κB and TGF-β pathways in early rat HCC.
To identify how many of our early rat HCC signature genes were conserved in early human HCCs, we re-analyzed two independently performed microarray datasets for early human HCC from Wurmbach et. al.  and Mas et. al. . The Wurmbach et. al.'s dataset composed of 19 early HCC patients and 10 normal controls, and the Mas et. al. dataset was composed of 16 cirrhotic livers with early HCC, 38 HCV associated HCC, and 19 normal livers. In addition, we also compared our signature genes with the OncodB.HCC database for 57 HCC patients from the Stanford HCC microarray data. The comparison of our rat early HCC signature genes (human orthologous) with the re-analyzed early human HCCs datasets revealed that 154 unique genes were conserved with early human HCCs (p < 0.001), and 35 of those were shared by all datasets analyzed (Table (Table2,2, Figure Figure1).1). We found that unsupervised clustering using the conserved signature gene list across species was sufficient to separate individuals in both Wurmbach's and Mas' human HCC samples as either early HCC patients or normal controls (data not shown). The gene interaction network analysis of the 154 signature genes indicated the importance of NF-κB, RAS and JNK activation in early hepatoma formation (Additional file 4). The network analysis of 35 cross-species conserved early HCC signature genes reveals the important roles of ERK/MAPK, PI3K/AKT, and TGF-β pathways (Figure (Figure4).4). It also indicated a potential critical regulatory role of MYC, ERBB2, HNF4A, and SMAD3 for malignant transformation to early HCC (Figure (Figure4B4B).
Various studies have reported chromosomal instability at chromosomal regions associated with many cancers, including human HCC copy number (CN) status [7,9,26,29,30]. These genomic modifications, which in part are reflected in changes in DNA copy number (CN), may alter the transcriptional control mechanism, and hence impact gene expression levels [31,32]. Hence, we compiled genes located in CNA regions reported in three independent genome copy number studies of human HCC [7,9,10]. The integration of 154 early HCC signature genes with the copy number data resulted in 75 genes that mapped to human CNA chromosomal loci genes including COL1A1, CCNA2, NFATC2, F2, DCK, MMP2, GJA1, VIM, LGALS3BP, and SP100. The interaction network of those genes further corroborated with the activation of the NF-κB, p38 MAPK, AP1, and JNK pathways (Figure (Figure5).5). We found that almost 50% of the 35 cross-species conserved signature genes that are common to different types of early HCCs analyzed were mapped to genomic locations within the CNA regions (Table (Table22).
As a validation of our results, we analyzed four independently performed microarray datasets for early human and rat HCCs [19,33-35] using the analysis procedure defined in the "Methods" section on the new datasets. The first validation dataset was from Chiang et.al.  using Affymetrix short oligo arrays. The dataset was composed of 91 HCV-related HCC tumor samples, of which 65 were very early or early stage disease, which we used in our re-analysis and comparison. The re-analyzed validation dataset showed a significant number of genes (p < 10-5) in common with our analysis results. In fact, more than 50% of our cross-species conserved genes were also differentially expressed in early HCC compared to normal controls in the validation dataset (Table (Table2).2). The significance of overlaps was calculated using hypergeometric distributional assumption  and p-values were adjusted using Bonferroni correction for multiple comparisons . In addition, unsupervised clustering was performed using our 35- gene signature to cluster the samples from Chiang et. al. We found that using our signature gene list was sufficient to separate individuals in Chiang et. al.'s study as either early HCC patients or normal controls (Additional file 5).
Moreover, we found a significant number of genes in common with the second dataset from Boyault et. al.  which consisted of 57 human HCCs and five samples of pooled non-tumorus tissues, and a third gene expression dataset from Liao et. al.  consisting of human HCC from various stages (we used expression data for only the early stage of the disease) (Table (Table2).2). Furthermore, we obtained consistent results with the DEN-induced HCC in rats from [19,38]. The IPA functional and network analysis of all validation datasets revealed a significant number of overrepresented functional categories and pathways in common with our results. Of note, cell death, cancer, cellular development, cellular growth and proliferation, organismal development, transport, and cell cycle came up as significantly enriched categories in both the validation datasets and our analyses. The interaction networks analyses of significantly dysregulated genes in validation datasets highlight the important roles of MYC, ERK/MAPK, AKT, NF-κB and TGF-β signaling pathways. The similar findings between our results and the independent validation sets argue against random chance accounting for the observed enrichment of these functional categories and pathways.
To confirm the microarray results by an independent method, we validated expression levels of six randomly selected differentially regulated genes (Pbsn, Cdh13, Lum, Nid2, Dcn, Slc22a5) in early rat HCC by realtime quantitative RT-PCR. A highly significant correlation existed between the microarray and realtime RT-PCR results (r = 0.97, p value < 0.001) (Figure (Figure6),6), thus demonstrating the reliability of our gene expression measurements. The selected genes and their interaction networks with other genes are shown (Figures 3C, D and and3E,3E, ,4A,4A, and Additional file 3).
To further validate the differential expression of potential biomarker genes using realtime RT-PCR from the blood of early HCC patients and healthy control subjects (ten subjects in each group), we selected eight genes (GJA1, VIM, IGFBP3, COL1A1, SP100, MMP2, LGALS3BP, and DPP4) among early HCC gene signature with CNA (denoted with asterisk in Table Table2)2) and were differentially regulated in at least one of the independent datasets. We confirmed a statistically significant increase in the expression of these biomarker genes in early HCC patients relative to healthy control subjects (p-value < 0.05) (Table (Table2,2, Figure Figure7);7); hence demonstrating the robustness of the cross-species integrated genomics procedure.
The present study sought to identify evolutionarily conserved inter-species biomarkers for early HCC differentiated from liver regeneration using integrative and cross-species comparative genomic approaches. The main contributions of this study are as follows: First, we developed a rat model of liver regeneration post-hepatectomy (return to quiescence), as well as liver cells undergoing malignant transformation and compared them to normal liver using a comprehensive microarray of 27,000 publicly available and Celera annotated rat genes. We included the liver regeneration in our model, as regeneration is a critical component in the surgical treatment of HCC, and frequently associated with HCC occurrence . Though liver cells can regenerate, they do not typically transform and lead to HCC . Therefore, an early HCC marker needs to be unique for the tumorigenic process and not overlap with the transcriptome changes that occur in regenerating or normal liver tissue. As ageing is also known to be a confounding factor embedded in gene expression profile data [15-17], we included age as a factor in our multi-factor statistical analysis, and identified age-specific differences in early HCC. Secondly, we performed cross-species comparative analysis to identify genes that are conserved in early rat HCCs and in multiple independently performed early human HCCs which would facilitate the identification of critical regulatory modules conserved across species in the expression profiles. Finally, we integrated genomic CNA data associated with the human HCCs with the transcriptomic profile, and performed validation analyses both in silico and with quantitative realtime RT-PCR (as schematically outlined in Figure Figure1).1). As CNAs have clear impact on expression levels in a variety of tumors [26,29,30], this dual strategy is very effective for interpreting the DNA and RNA level anomalies in cancer, in order to identify genes involved with tumor initiation and progression [24,26].
The validation analyses demonstrated great concordance of our results with other data sets using various microarray platforms. The ABI 1700 system has a unique approach in identifying dysregulated genes since it targets genes from both Celera and Public databases and utilizes chemiluminescently enhanced detection that is likely to determine relatively rare mRNAs. Also, our confirmatory quantitative realtime RT-PCR experiments displayed a strong correlation with the microarray results, adding to the validity of the present observations. This is in agreement with some recent studies showing a linear relationship for real-time and conventional reverse transcription and therefore validates the robustness of mRNA quantification using either microarrays or quantitative RT-PCR. Hence this allowed us to identify potential biomarkers for human early HCC and to gain further insight into the mechanism of early hepatoma formation.
We performed a two-step algorithm to identify early rat HCC signature genes: In the first step, a two-way ANOVA was performed including treatment and age as well as their interactions into our statistical model and we identified genes uniquely expressed in early HCC in both young and old rats (Figure (Figure2).2). In the second step, because the interaction between age and treatment was significant, we stratified our samples as young and old cohorts, and HCC specific genes were identified using two one-way ANOVA in each age group separately . Finally, the HCC specific genes in both young and old identified in the two steps were combined before performing the cross-species comparison (as detailed in "Material and Methods" section, and schematically shown in Figure Figure1).1). When comparing the early HCC group with the regeneration or normal cohorts to identify the differentially expressed genes we used a set of criteria: S/N ratio > 3 in > 50% of the samples, a p-value < 0.02 and absolute fold change > 1.8. These observations are consistent with the study of Guo et. al., in which gene lists ranked by fold change and filtered with non-stringent statistically significant tests were more reproducible across platforms than those generated through other analytical procedures [43,44]. In addition, as Ghosh et. al. discuss on combining data from multiple gene expression studies, if two studies independently discover that the same gene/protein to be differentially expressed, then the chance of error is significantly reduced .
The comparison of our signature genes with three different independently performed early human HCC microarray data sets revealed a significant number of early rat HCC genes (human orthologous) conserved across early human HCCs (p < 0.001). Indeed, many of those genes were related to cancer. For example, LUM, CCNA2, IGFBP3, HPX, COL1A1, SRPX, VIM, TGFBR1, DCN, MMP2, CD14, DCK, BIRC3, GJA1, LOX, SP100, PROM1 and CREB1 were known to regulate tumorigenesis, neoplasia, apoptosis, growth, differentiation and proliferation. Some of the most significantly activated canonical pathways included hepatic fibrosis/hepatic stellate cell activation (CD14, COL1A1, FGFR2, IGFBP3, MMP2, and TGFBR1), and p38 MAPK signaling pathways (CREB1, PLA2G2A, TGFBR1). The network analysis of early HCC signature genes indicated the activation of ERK/MAPK, PI3K/AKT, and TGF-β signaling pathway, and a potential critical regulatory role of MYC, ERbB-2, HNF4A, and SMAD3 for early HCC (Figures (Figures2,2, ,33 and and4).4). MAPKs are implicated in diverse cellular processes such as cell survival, differentiation, adhesion, and proliferation . The gene network analysis of differentially expressed genes further confirmed the altered pathways. Moreover, it also indicated the importance of NF-κB, RAS and JNK activation in early hepatoma formation (Figure (Figure55 and Additional file 4).The role of MYC in various types of carcinogenesis has been extensively investigated . Most recently, JNK1 activation, and increased expression of ErbB-2 were found to be associated with HCC . Thus, our current findings are consistent with previously performed independent cancer studies, including those for HCC. However, the novelty of our approach is that using comparative and integrative genomics, we provide evidence for the potential central role of these genes in the earliest phase of liver malignant transformation.
Our comparative genomics analysis resulted in a 35-gene cross-species conserved signature for all types of early HCCs. Over 70% of the conserved genes were associated with cancer according to the IPA knowledgebase, including LGALS3BP[50,51], VIM[52,53], DCN[54,55], IGFBP3, FGFR2, GJA1, SP100, DPP4, PROM1, BIRC3, MMP2, and COL1A1[64,65] (Table (Table2).2). Furthermore, using literature mining tools, such as MILANO , we found that almost 90% of our signature genes were reported to be cancer related.
There are areas of genomic instability reported in many cancers, including HCC, and some regions commonly exhibit either deletion or increased gene dosage, leading to changes in DNA copy number (CN) [9,26,29,30]. Integrating the gene expression with the CN data reveals the chromosomal regions with concordantly altered genomic and transcriptional status in tumors [24,32,67]. Hence, focusing on differentially expressed genes with concomitant altered DNA copy number may identify novel early HCC markers of malignant transformation and progression. The presence of altered DNA CN and LOH may contribute to cancer formation [30,31,68]. Therefore, the pattern of genomic modifications in a tumor represents a structural fingerprint that may include the transcriptional control mechanisms and locally impact gene expression levels [31,32]. We identified that more than 50% of our cross-species conserved early HCC signature genes were found to be copy number dependent (Table (Table22).
We found significant expression of LGALS3BP (Lectin, galactoside binding soluble 3 binding protein) and COL1A1 located on Chromosome 17q. The LGALS3BP is a 90-kD protein, designated serum protein 90 K that was found at elevated concentrations in the serum of patients with various types of breast, lung, colorectal, ovarian, and endometrial cancer [50,51]. It is a secreted glycoprotein that binds galectins, beta1-integrins, collagens, and fibronectin, and has some relevance in cell-cell and cell-extracellular matrix adhesion . Another gene which could be a potential biomarker for early HCC is dipeptidyl peptidase IV (DPP4). DPP4 is a serine protease, which plays an important role in immune regulation, signal transduction, and apoptosis. It has been shown that DPP4 may have a critical function in tumor progression in several human malignancies [60,70]. Matrix metalloproteinases (MMP) also are involved with early carcinogenic events, tumor growth, tumor invasion and metastasis [63,71,72]. Matrix metalloproteinases (MMPs) are zinc-dependent endopeptidases that cleave and degrade a wide spectrum of extracellular matrix components, and are involved with extracellular matrix remodeling during the process of tumor invasion and metastasis . Alterations in MMP expression and their endogenous inhibitor (TIMP) may contribute to HCC metastasis [71-73].
It is worth mentioning that, the gene "Probasin" (Pbsn), was significantly up-regulated (fold change > 15 in both old and young early rat HCCs). The high expression of Pbsn in our rat model was also confirmed with realtime RT-PCR. Pbsn is a member of the lipocalin family and has not yet been associated with HCC in rats and has no known human ortholog. However, it has been shown that Pbsn is highly expressed in prostate and implicated in both benign prostatic hyperplasia and prostate cancer [74-76] and taste bud tumorigenesis in rats . Also since the promoter of this gene exhibits strong androgen receptor-specific and tissue-specific regulation, Pbsn is proposed to be a potential candidate for targeted therapies for advanced prostate cancer .
We have also found significant expression of lumican (LUM) and decorin (DCN) in both early rat and human HCCs. LUM and DCN are members of a small leucine-rich proteoglycan (SLRP) family. Lumican has been shown to participate in the maintenance of tissue homeostasis and modulation of cellular functions including cell proliferation, migration, adhesion, and differentiation . Decorin has been reported to have a number of functions including suppressing cancer cell growth and metastasis andacting with extracellular matrix molecules to influence cell adhesion and fibril stability . DCN acts as a natural inhibitor of TGF and is considered to be a specific antagonist of EGFR . In addition, the altered expression of lumican and decorin has been associated with various human cancers including breast, pancreatic, lung, ovarian, melanoma, colorectal, osteosarcoma and ductal adenocarcinoma [54,78-82].
Genes whose protein products are released into the extracellular space would be ideal tumor markers for clinical applications, as it would be possible to detect these proteins in patients' biological fluids rather than through the use of invasive biopsies. Moreover, previous studies have found that cells derived from peripheral blood could be used to assess disease-associated gene signatures [83-89]. In our study, we confirmed the high expression of eight selected candidate biomarker genes (GJA1, VIM, IGFBP3, COL1A1, SP100, MMP2, LGALS3BP, and DPP4) by using realtime RT-PCR from the blood of early HCC patients. These genes and other potential biomarker genes identified through our integrated-comparative genomics approach (listed in Table Table2)2) and their encoded proteins will be further studied in a large cohort of patients to determine if they have a role in early HCC pathology and if they could be novel early HCC biomarkers detectable in biological fluids.
In summary, to our knowledge, this is the first study to examine HCC differentiated from regeneration in both old and young rats, and coupled with a cross-species comparative and integrative genomics approach to identify genes that could be potential biomarkers for early human HCC. The results of our study include the depiction of refined and delineated biological pathways differentially modulated in HCC that is built around TP53, p38 MAPK, ERK/MAPK, PI3K/Akt, NF-κB, TGF-β, MYC, and ERbB-2, including their target genes that were not previously implicated with early HCC. Our cross-species comparative and integrative genomics approach which involved integration of multiple high dimensional independent datasets has led to potentially robust biomarkers for the detection of early HCC. The signature genes that we identified could be considered as "evolutionarily conserved cross-species biomarkers for early HCC with genomic copy number alterations". Further studies are needed to identify if any of the potential biomarkers identified in this study can be readily and reproducibly detected in blood, urine or other bodily fluids. This could then form the basis of a useful diagnostic test for the detection of early HCC.
Male Sprague-Dawley rats were maintained at the King Fahad National Centre for Children's Cancer and Research Animal Facility. This facility is managed in accordance with AALAS regulations. Ten young adult (5 months) and ten old adult (21 months) animals were subjected to partial hepatectomy. Actual survival rates allowed for four animals in each group to be analyzed. The re-growth of one lobe of the liver was completed within one month, by which time the liver cells again became quiescent, which was confirmed by histological analysis. In parallel, separate animals were treated with diethylnitrosoamine (200 mg/kg), which was injected intraperitoneally to induce the formation of HCC. Once the early HCC formation became apparent within 2-4 weeks, the rats were sacrificed. In the partial hepatectomized animals, one unaffected lobe and the regenerated lobe of the liver were removed independently. In the carcinogenic treated animals the tumors were carefully dissected to avoid removing normal tissue. All tissues were snap-frozen and stored at -80°C until required for RNA isolation. Small pieces of tissue were removed for formalin fixation to be used for histological examination.
Twenty blood samples were collected for this study (10 early HCC and 10 from healthy controls). Histopathological classification of HCC and clinical staging of early HCC were performed according to International Working Party  as previously described . Patients diagnosed with the early HCC and healthy controls were recruited under an institutional review board-approved project (RAC# 2060040); all subjects provided written, informed consent before entry in the study. A total of 4 ml (in two separate PaxGene tubes) of whole blood samples were collected for each individual according to manufacturer's guidelines (QIAGEN Inc., Valencia, CA, USA). The total RNA isolation was performed using PreAnalytiX - PAXgene Blood RNA System (QIAGEN Inc.) by strictly following the manual and protocols provided with the kit-system.
Total RNA was isolated according to standard protocols. Quality Control of RNA was done using Bioanalyzer 2100 RNA 6000 NanoAssay and RNA above RIN = 8 was included to the study (Agilent Technologies, Santa Clara, USA). Rat Genome Survey Microarray (Applied Biosystems, Foster City, CA, USA) was utilized for microarray studies. cDNA synthesis, cRNA and labeling, chemiluminescence detection, image acquisition and analysis were performed following the manufacturer's protocols, guidelines and recommendations.
Images were auto-gridded and the chemiluminescent signals were quantified, then background subtracted using the Applied Biosystems 1700 Chemiluminescent Microarray Analyzer software v 1.1. For transcriptome analysis, detection thresholds were used following the manufacturer's recommendations. Detection threshold was set as S/N > 3 and quality flag < 5000. The microarray data were analyzed from 24 samples (2 samples were excluded for quality reasons). The open source R software package http://www.r-project.org and tools from the BioConductor project were used for normalization and determination of differentially expressed genes . Two-factor Analysis of Variance (ANOVA) was performed to include both "treatment" (HCC, regenerated and normal, which will be referred to as treatment in the remainder of the manuscript), as well as age (old and young) factor together with feature selection algorithm (also known as template matching(TMA))  to look for treatment as well as age specific variation. Significantly modulated genes specific to HCC were defined as those with ANOVA (treatment) p- value < 0.01, and TMA p-value < 0.01. Additionally, samples were stratified as young and old cohorts, and HCC specific genes identified using one-way ANOVA in each age group separately. When comparing HCC group with regenerated and normal controls to identify the differentially expressed genes specific to the HCC, we used a combination of three criteria. We considered genes that are "present" in at least half of the samples in either group. HCC specific genes were defined as those with in absolute fold change > 1.8 and p-value < 0.02. These observations are consistent with the study of Guo et. al., in which gene lists ranked by fold change and filtered with non-stringent statistically significant tests were more reproducible across platforms than those generated through other analytical procedures . A validation datasets were generated from three independent human HCC studies by Chiang et. al.  (GSE9843), which was composed of 91 HCV-related HCC tumor samples, of which 65 were very early or early stage disease (we used in our re-analysis only very early and early HCC datasets) and Boyault et. al.  (E-TABM-36) which consisted of 57 human HCCs and five samples of pooled non-tumorus tissues, and from Liao et. al.  (GSE 6222) consisting of various stages of HCC (we used expression data for only the early stage of the disease). Furthermore, we compared our results with the results from two independent studies [19,38] with the DEN-induced HCC in rats. The raw data was analyzed by using dChip and open source R/Bioconductor packages. The dChip outlier detection algorithm was used to identify outlier arrays, and probes "present" in at least 50% of the samples in either group were filtered. The data was normalized by the GC Robust Multi-array Average (GC-RMA) algorithm [93,94]. Unpaired t-tests were performed to determine significant differences in gene expression levels between patients and normal controls. The Hierarchical clustering using Pearson's correlation with average linkage clustering was performed by MeV 4.0 .
Information about genes participating in known biological process and pathways were derived by using DAVID Bioinformatics Resources, Expression Analysis Systematic Explorer (EASE), and PANTHER (Protein ANalysis THrough Evolutionary Relationships) Classification Systems . For each molecular function, biological process or pathway term, PANTHER calculates the number of genes identified in that category in both a list of differentially regulated genes and a reference list containing all the probe sets present on the chip and compares these results using the binomial test to determine if there are more genes than expected in the differentially regulated list . Over-representation is defined by p < 0.05. Statistical analyses were performed with the MATLAB software packages (Mathworks, Natick, MA, USA), R and Bioconductor and PARTEK Genomics Suite (Partek Inc, St. Lois, MO, USA).
Functional pathway, gene ontology and network analyses were executed using Ingenuity Pathways Analysis (IPA) 6.3 (Ingenuity Systems, Mountain View, CA). The differentially expressed signature gene lists for hepatoma and regeneration in different age groups were mapped to its corresponding gene object in the Ingenuity pathway knowledge base. These so-called focus genes were then used as a starting point for generating biological networks. A score was assigned to each network in the dataset to estimate the relevance of the network to the uploaded gene list. This score reflects the negative logarithm of the P that indicates the likelihood of the focus genes in a network being found together due to random chance. Using a 99% confidence level, scores of ≥2 were considered significant. Significances for biological functions or pathways in the signature genes for such functions or pathways compared with the ABI Rat Genome Survey Microarray as a reference set. A right-tailed Fisher's exact test was used to calculate a p-value determining the probability that the biological function (or pathway) assigned to that data set is explained by chance alone.
Human early HCC datasets from two independent studies by Mas et. al.  using Affymetrix HG-U133A 2.0 array, and Wumbach et. al.  using Affymetrix HG-U133 Plus 2.0 were re-analyzed. The raw data were analyzed using R/Bioconductor packages and Partek Genomics Suite (Partek Inc.). The data were normalized by the GC Robust Multi-array Average (GC-RMA) algorithm. Unpaired t-tests were performed to determine significant differences in gene expression levels between patients and normal controls. The cross mapping of Applied Biosystems Rat Genome Survey microarray probes were mapped to human orthologs through "AB1700 rat annotation spreadsheet" designed by Applied Biosystems on the basis of sequence identity. The transcripts present on both platforms (AB1700 and Affymetrix) were identified using Resourcerer . Genes within copy number altered regions based on three independent genome CNA studies of human HCC [7,9,10] were determined using NCBI MapViewer http://www.ncbi.nlm.nih.gov/mapview, and integrated those with the gene expression profiling data (Figure (Figure11).
Confirmatory realtime RT-PCR experiments were performed using the ABI 7500 Sequence Detection System (Applied Biosystems). 50 ng total RNA procured from the same microarray study samples were transcribed into cDNA using a Sensicript Kit (QIAGEN Inc., Valencia, CA, USA) under the following conditions: 25°C for 10 min, 42°C for 2 hrs, and 70°C for 15 min in a total volume of 20 μl. Six differentially expressed rat genes (Pbsn, Cdh13, Lum, Nid2, Dcn, Slc22a5) and eight human genes (GJA1, VIM, IGFBP3, COL1A1, SP100, MMP2, LGALS3BP, and DPP4) were selected and primers designed using Primer3 software. For the human samples, blood total RNA was utilized. After primer optimization, realtime PCR experiments were performed with 6 μl cDNA using Quantitech SyBr Green Kit (QIAGEN), employing GAPDH as the endogenous control gene. All reactions were conducted in triplicates and the data was analyzed using the delta delta CT method [101,102].
The authors declare that they have no competing interests.
DC and NK conceived the research problem and designed the methodology. MAC developed the rat model. MG provided suggestions for the rat model. NK and AA carried out microarray experiments. DC collected and processed the data, performed statistical and bioinformatics analyses, and drafted the manuscript together with NK and BHP. JQ and MMS provided suggestions for method design and commented on the manuscript. MAC, PTO, BHP, and AAQ provided general interpretation of results. All authors read and approved the final manuscript.
Selected HCC specific genes, conserved across both age groups (old and young), and significantly modulated with respect to regenerated and normal liver.
Comparison of expression profiles of HCC and regeneration within the same age group. (A, D) Heatmap of significantly dysregulated genes due to different treatment types in young and old, respectively. (B, E) Hierarchical clustering of samples separated based on treatment type in young and old, respectively. The gene expression clustering distance between the HCC group and other two groups (regenerated and normal) was the greatest in both age groups (C, F) Principle component analysis (PCA) which contained almost 76% of the variance in the data matrix clearly separated samples based on the treatment type in young and old, respectively.
Heatmap and gene interaction networks of HCC specific genes in the old age group. (A) Venn diagram characterizing differential gene expression between and specific to different treatment types (the HCC, the regenerated and the normal). The number of HCC specific genes, 100, is circled in black. (B) Heatmap of HCC specific genes exclusively dysregulated (up/down regulated) in the HCC group only. (C-E) Functional network analysis of HCC specific genes. Top three scoring gene interaction networks (with highest relevance scores) are shown. Nodes represent genes, with their shape representing the functional class of the gene product, and edges indicate biological relationship between the nodes (see legend in Figure Figure2).2). (F) Top network functions associated with three networks shown. An IPA score of three indicates that there is 1/1000 (score = -log (p-value)) chance that the focus genes are assigned to a network randomly.
The gene interaction network analysis early HCC signature genes that are conserved in rat early HCC and in either of multiple human early HCCs (A, B) The top two scoring gene interaction networks of 154 cross-species conserved signature genes indicated the importance of NF-κB, RAS and JNK activation in early hepatoma formation. Nodes represent genes, with their shape representing the functional class of the gene product, and edges indicate biological relationship between the nodes (see legend in Figure Figure55).
The unsupervised Principle Components Analysis (PCA) was performed using our 35-gene signature to cluster samples from independent validation dataset of Chiang et. al. Our signature gene list was sufficient to separate individuals in Chiang et al.'s study as either early HCC patients or normal controls.
We thank Maqbool Ahmad for providing help in RNA isolation, and Bedri Karakas for critically reading the manuscript. We would like to thank to King Faisal Specialist Hospital and Research Center for the financial support.