Integrative Analysis of Microarray Data to Reveal Regulation Patterns in the Pathogenesis of Hepatocellular Carcinoma.

BACKGROUND/AIMS: The integration of multiple profiling data and the construction of a transcriptional regulatory network may provide additional insights into the molecular mechanisms of hepatocellular carcinoma (HCC). The present study was conducted to investigate the deregulation of genes and the transcriptional regulatory network in HCC.
METHODS: An integrated analysis of HCC gene expression datasets was performed in Gene Expression Omnibus. Functional annotation of the differentially expression genes (DEGs) was conducted. Furthermore, transcription factors (TFs) were identified, and a global transcriptional regulatory network was constructed.
RESULTS: An integrated analysis of eight eligible gene expression profiles of HCC led to 1,835 DEGs. Consistent with the fact that the cell cycle is closely related to various tumors, the functional annotation revealed that genes involved in the cell cycle were significantly enriched. A transcriptional regulatory network was constructed using the 62 TFs, which consisted of 872 TF-target interactions between 56 TFs and 672 DEGs in the context of HCC. The top 10 TFs covering the most downstream DEGs were ZNF354C, NFATC2, ARID3A, BRCA1, ZNF263, FOXD1, GATA3, FOXO3, FOXL1, and NR4A2. This network will appeal to future investigators focusing on the development of HCC.
CONCLUSIONS: The transcriptional regulatory network can provide additional information that is valuable in understanding the underlying molecular mechanism in hepatic tumorigenesis.

INTRODUCTION

Hepatocellular carcinoma (HCC) is the most common and deadly form of liver cancer accounting for ~80% of adult primary liver cancer, and it is one result of underlying liver disease.1 The prevalence and the incidence of HCC have progressively increased in the world recently.2 Because of the coexistence of cirrhosis and inflammation, early diagnosis of HCC is plagued by lack of clinical research and reliable indicators. Moreover, the lack of specific symptoms in the early stages of HCC also contributes to the poor prognosis of the disease.

Considering that, many researchers have pay attention to the HCC-specific biomarkers for early diagnosis of HCC currently, and a number of biomarkers have been identified. The most commonly used serological biomarker is α-fetoprotein (AFP) for detection of HCC in clinical practice, which is a specific glycoprotein produced primarily by the fetal liver.3 However, due to its low sensitivity and specificity, the clinical diagnostic accuracy of AFP is unsatisfactory. It was then reported that combination assay of high-sensitivity des-γ-carboxy prothrombin and AFP-L3 can improve the detection rate of HCC.4 In addition, the midkine can serve as a useful marker in the diagnosis of AFP-negative HCCs and at a very early stage.5 Furthermore, several diagnostic biomarkers were identified continuously, such as Dickkopf-1,6 Golgi protein 73,7 Glypican-3,8 γ-glutamyl transferase,9 α-l-fucosidase,10 transforming growth factor β-1,11 IGFs,12 squamous cell carcinoma antigen,13 osteopontin,14 heat shock proteins,15 and so on, among which most of them can be used together with AFP for diagnosis of HCC. Actually, combination of more than one biomarker may improve the accuracy of HCC diagnoses.16

Even so, none of the biomarkers have been considered as the reliable indicator in the early HCC diagnosis, which mainly because the pathogenesis of HCC remains undetermined. In this study, we extracted the gene expression profiles of HCC from Gene Expression Omnibus (GEO) database. By comparing the global gene expression profiles between HCC and the normal tissues, a set of differentially expressed genes (DEGs) were identified and the differentially expressed transcription factors (TFs) were further extracted. Additionally, the HCC-specific transcriptional regulatory network was constructed, which may provide better clues on the underlying regulatory mechanisms of pathogenesis of HCC and therapeutic applications.

MATERIALS AND METHODS

1. Eligible datasets of HCC

GEO database is a public functional genomics data repository (http://www.ncbi.nlm.nih.gov/geo/).17 By online search, the gene expression profiles of HCC were obtained from the GEO database. The following key words were used “carcinoma, hepatocellular” [MeSH Terms] OR hepatoma [All Fields] AND “Homo sapiens” [porgn] AND “Expression profiling by array” [Filter]. The eligible datasets were included and downloaded for integrated analysis, which were obtained from microarray experiments on the gene expression of HCC and normal tissues. In the eligible datasets, the etiologies for HCC patients were as follows: GSE17548 (10 HBV-related+3 HCV-related+2 cryptogenic HCC), GSE44074 (17 HBV-related+17 HCV-related), GSE45436 (93 cryptogenic HCC), GSE46408 (6 cryptogenic HCC), GSE50579 (8 HBV-related+9 HCV-related+14 alcohol-related+30 cryptogenic HCC), GSE57957 (39 cryptogenic HCC), GSE60502 (18 cryptogenic HCC), and GSE62232 (10 HBV-related+9 HCV-related+33 alcohol-related+29 cryptogenic HCC).

2. Detection of DEGs

All the raw data were preprocessed via background correction and normalization. The limma (linear models for microarray analysis) package in R, one of the most commonly bioconductor packages, was used to analyze the differential expression between HCC and the normal tissues by t-test. The p-value and false discovery rate (FDR) were also obtained. Genes with FDR <0.01 were defined as DEGs in this study.

3. Functional annotation of DEGs

To better understand the biological functions of DEGs in the pathogenesis of HCC, functional enrichment of DEGs were analyzed using the web-based tools in Database for Annotation, Visualization and Integrated Discovery (DAVID).18 The enriched gene ontology (GO) terms19 and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways20 were identified. Only the GO terms and KEGG pathways with p<0.05 were taken into account as significantly enriched among the DEGs.

4. Construction of transcriptional regulatory network

TFs are essential for the regulation of gene expression, which can provide better clues on the underlying regulatory mechanisms. TRANSFAC (transcription factor database) is a manually curated database of eukaryotic TFs, their genomic binding sites and DNA binding profiles. Based on the DEGs in HCC from integrated analysis, we searched TRANSFAC for DEGs coded TFs and their targeted genes, and used TRANSFAC position weight matrix for gene promoter scanning21 to identify DEGs which have the binding sites of the TFs in the promoter regions. Then the transcriptional regulatory network was visualized using Cytoscape.22

5. RNA preparation and qRT-PCR

The peripheral blood samples were collected from five patients with HCC and five healthy volunteers, and the samples were immediately stored into vacuum EDTA anticoagulant tubes. All samples were obtained with permission, and the project was approved by our medical ethics committee for the relating screening, inspection, and data collection of the patients. Total RNA was extracted from the blood samples using the TRIzol® Reagent (Invitrogen, Carlsbad, CA, USA) following the manufacturer’s protocols. Total RNA (1 μg) was reverse-transcribed in 25 μL reactions using SuperScript® III Reverse Transcriptase (Invitrogen), and the resulted cDNA was used as template for real-time polymerase chain reaction (PCR). Real-time PCR was carried out in ABI 7500 real-time PCR system with a Power SYBR® Green PCR Master Mix (Invitrogen). The results were analyzed using 2−ΔΔct method. ACTIN gene was used as the endogenous control. The sequences of primers used for real-time PCR were listed in Supplementary Table 1.

RESULTS

1. Comparison of global gene expression profiles of HCC and normal tissue

In present study, eight eligible gene expression profiles of HCC were included, and the detailed information of datasets was displayed in Table 1. Totally, 347 cases of HCCs (45 HBV-related+38 HCV-related+47 alcohol-related+217 cryptogenic HCC) and 121 controls of normal liver tissues were enrolled in the integrated analysis. A set of 1,835 genes were regarded as having DEG by applying the selection criteria of FDR <0.01, among which 1,145 were upregulated and 690 downregulated, respectively. All the DEGs were listed in Supplementary Table 2.

Table 1

Characteristics of the Individual Studies

GEO ID	Control	Case	Platform	Year	Author
GSE17548	0	15	GPL570 Affymetrix Human Genome U133 Plus 2.0	2009	Yildiz et al.35
GSE44074	0	34	GPL13536 Kanazawa Univ. Human Liver chip 10k	2013	Ueda et al.36
GSE45436	41	93	GPL570 Affymetrix Human Genome U133 Plus 2.0	2013	Wang et al.
GSE46408	6	6	GPL4133 Agilent-014850 Whole Human Genome Microarray 4×44K	2013	Jeng et al.
GSE50579	7	61	GPL14550 Agilent-028004 SurePrint G3 Human GE 8×60K Microarray	2013	Neumann et al.37
GSE57957	39	39	GPL10558 Illumina HumanHT-12 V4.0 expression beadchip	2014	Mah et al.38
GSE60502	18	18	GPL96 Affymetrix Human Genome U133A Array	2014	Wang et al.39
GSE62232	10	81	GPL570 Affymetrix Human Genome U133 Plus 2.0	2014	Zucman-Rossi et al.

GEO, Gene Expression Omnibus.

2. Functional enrichment of DEGs

By subjecting the DEGs to enrichment analysis on gene sets based on GO (biological process, cellular component, and molecular function) and KEGG pathways, we observed that DEGs were significantly enriched in various terms (Table 2). For biological process, the DEGs were mainly enriched in oxidation-reduction process, response to external stimulus, response to stimulus and defense response. For cellular component, extracellular region and integral component of membrane were involved. For molecular function, receptor activity and oxidoreductase activity were implicated. Based on the canonical signaling pathways documented in KEGG, pathways on p53 signaling pathway, pathways in cancer, cell cycle, DNA replication and homologous recombination were significantly enriched (Table 3).

Table 2

Top 15 Most Significantly Enriched Gene Ontology Terms of Differentially Expressed Genes

GO ID	GO term	No. of genes	FDR
Biological process
GO:0044281	Small molecule metabolic process	175	5.01E-18
GO:0019752	Carboxylic acid metabolic process	107	3.42E-17
GO:0043436	Oxoacid metabolic process	110	1.05E-15
GO:0006082	Organic acid metabolic process	111	9.64E-16
GO:0055114	Oxidation-reduction process	91	1.73E-14
GO:0032787	Monocarboxylic acid metabolic process	66	5.40E-13
GO:0009605	Response to external stimulus	78	9.74E-13
GO:0044282	Small molecule catabolic process	47	5.50E-12
GO:0016054	Organic acid catabolic process	43	4.91E-12
GO:0046395	Carboxylic acid catabolic process	43	4.42E-12
GO:0006629	Lipid metabolic process	96	3.63E-11
GO:0050896	Response to stimulus	190	7.48E-10
GO:0006631	Fatty acid metabolic process	44	1.35E-09
GO:0006952	Defense response	66	3.47E-09
GO:0033993	Response to lipid	55	6.00E-09
Molecular function
GO:0004872	Receptor activity	67	7.67E-16
GO:0016491	Oxidoreductase activity	72	1.11E-13
GO:0038023	Signaling receptor activity	54	2.31E-12
GO:0004888	Transmembrane signaling receptor activity	48	2.58E-11
GO:0060089	Molecular transducer activity	71	5.81E-11
GO:0004871	Signal transducer activity	57	1.27E-07
GO:0005102	Receptor binding	62	2.30E-05
GO:0005125	Cytokine activity	16	2.91E-05
GO:0004930	G-protein coupled receptor activity	22	1.09E-04
GO:0048037	Cofactor binding	39	2.21E-04
GO:0004252	Serine-type endopeptidase activity	11	2.23E-04
GO:0008236	Serine-type peptidase activity	12	6.36E-04
GO:0017171	Serine hydrolase activity	12	5.87E-04
GO:0004497	Monooxygenase activity	11	6.44E-04
GO:0043177	Organic acid binding	21	9.70E-04
Cellular component
GO:0005576	Extracellular region	70	1.36E-14
GO:0016021	Integral component of membrane	198	7.85E-15
GO:0031224	Intrinsic component of membrane	166	5.96E-15
GO:0044425	Membrane part	262	7.70E-15
GO:0044459	Plasma membrane part	102	3.74E-13
GO:0031226	Intrinsic component of plasma membrane	67	7.56E-13
GO:0044421	Extracellular region part	205	6.79E-13
GO:0005887	Integral component of plasma membrane	65	3.95E-12
GO:0005615	Extracellular space	72	7.80E-12
GO:0043230	Extracellular organelle	153	1.78E-09
GO:1903561	Extracellular vesicle	153	1.62E-09
GO:0065010	Extracellular membrane-bounded organelle	152	2.50E-09
GO:0070062	Extracellular exosome	152	2.31E-09
GO:0031988	Membrane-bounded vesicle	165	1.54E-08
GO:0005886	Plasma membrane	127	1.74E-08

Go, gene ontology; FDR, false discovery rate.

Table 3

Top 15 Most Significantly Enriched Kyoto Encyclopedia of Genes and Genomes Pathways of Differentially Expressed Genes

KEGG term	Count	FDR
Cell cycle	35	1.06E-14
Valine, leucine and isoleucine degradation	20	1.29E-12
Fatty acid metabolism	16	7.88E-09
Propanoate metabolism	14	1.22E-08
DNA replication	14	4.34E-08
Oocyte meiosis	24	6.58E-08
p53 signaling pathway	18	1.86E-07
Progesterone-mediated oocyte maturation	19	1.35E-06
Pathways in cancer	41	1.88E-06
Lysine degradation	14	2.30E-06
Mineral absorption	14	3.40E-06
Pyrimidine metabolism	19	3.63E-06
Homologous recombination	10	5.36E-06
Tryptophan metabolism	12	1.03E-05
Pyruvate metabolism	11	4.26E-05

KEGG, Kyoto Encyclopedia of Genes and Genomes; FDR, false discovery rate.

3. Transcriptional regulatory network

Transcriptional regulatory network is a directed graph describing regulatory effect of TFs on the expression of target genes. Based on the database of TRANSFAC, 62 differentially expressed TFs were identified, in which 23 were upregulated and 39 were downregulated. The constructed regulatory network consisted of 872 TF-target interactions between 56 TFs and 672 DEGs in the context of HCC (Fig. 1). The top 10 TFs covering the most downstream DEGs were identified as crucial TFs involved in the development of HCC and listed in Table 4, including ZNF354C, NFATC2, ARID3A, BRCA1, ZNF263, FOXD1, GATA3, FOXO3, FOXL1, and NR4A2.

Fig. 1

The established transcriptional regulatory network of hepatocellular carcinoma. Red and green nodes denote upregulated and downregulated transcription factors (TFs), respectively. Blue nodes denote differentially expressed genes predicted to interact with the corresponding TFs.

Table 4

Top 10 Transcription Factors Covering the Most Downstream Differentially Expressed Genes Involved in the Development of Hepatocel-lular Carcinoma

Transcription factor	logFC	Up/down	Count	Genes
ZNF354C	−2.12E+00	Down	81	TMEM206,TNFRSF25,STEAP4,FANCB,FATE1,NOX4,RECQL5,H2AFZ,ZNF687,FTSJ3,ZNF335,SYT9,PAK7,AIM1L,PROSC,CCNE2,SHPK,TCTEX1D1,DVL2,CCNB2,ESR1,ZNF337,SEMA5B,TRAF5,ZKSCAN5,ERRFI1,ZMIZ2,CSPG4,ANKZF1,MAD2L1,SMARCD1,NUDT1,CYP4V2,ATP11C,PRDM7,ALDH9A1,NMB,GPC2,HJURP,SLC7A6,DSTYK,KIAA1522,DDX49,GOLGA3,TP53I3,EGR2,TACC3,TRIP13,TTC27,FAM64A,ZNF76,SEMA3F,KDSR,SLC26A6,SPTAN1,SCNM1,DSN1,PIAS3,DUSP1,SCRIB,ABCF2,FAM122A,CCR1,NAMPT,ZNF692,LOC642852,KIF20A,RASSF1,HOXD9,NCAPG,CENPI,BCHE,CLDN8,RBM34,MT1H,COBLL1,INCENP,CENPT,MEA1,C1RL,TIGD6
NFATC2	−2.21E+00	Down	78	SNHG12,EPHA2,KIF18A,DTX1,BHLHE40,C21orf91,KDSR,WDR62,NXPH4,LGI1,CSTF2,C17orf53,STC2,DUSP28,SARDH,GABRD,PXDNL,HIST2H3C,ELOVL1,GLMN,PSKH1,TDRKH,KIAA0907,GIGYF1,GADD45B,ANGPTL6,SQLE,PPAPDC1A,ZNF581,SLC26A2,C5orf46,G6PD,MAZ,MCM6,TRAF2,LMNB2,CHRM2,TROAP,CLEC4GP1,ZIC5,NAAA,HGS,ATAD5,R3HDM1,STK40,GABRP,PANK1,PMFBP1,CLEC4M,ADRA1A,MTBP,STX11,ZSCAN16,ZNF572,CYP26A1,C5orf34,ARL16,PLGLB2,HHIPL2,NEU1,FAM111B,ERCC6L,GPATCH4,TK1,TSC1,RRM2,AP4M1,OLFML2B,TRIM45,CXCR2P1,RTEL1,C2orf44,DCAF4L1,HSF2BP,FBP1,ZNF335,EDC3,IDUA
ARID3A	1.60E+00	Up	74	ANTXR2,PLEKHG2,NXPH4,EPRS,MYD88,B4GALT7,RBM3,EPHX4,PPBP,MBNL2,CDK16,DSE,ZNF517,N4BP2L1,DIAPH3,PLK1,HCFC1,METTL3,AGPAT1,TTL,NUPL2,ATG2A,SMARCD1,ANGPTL1,NGFR,KCNN2,MCM8,RET,TPM3,LILRB5,TBC1D13,MMAA,FOXH1,ZKSCAN5,DNASE1,F8,BMPER,SOCS3,UXS1,IL33,CCDC137,SLC7A8,SOCS2,PGP,NSUN5P2,ZWILCH,BMP10,GOLPH3L,MAP2K4,GHR,PEX11G,POLQ,TSC22D2,CDKN2C,CDK5RAP1,LOXL2,GOLGA6L9,ZEB2,CYR61,PPAP2B,FEN1,C9orf117,UNC119B,C6,NDST3,PBLD,CUL7,COX4I2,SPNS1,CCDC142,SLITRK3,ABAT,TBX15,NFRKB
BRCA1	9.55E−01	Up	64	LOC646762,ANTXR2,POLR2J4,CCDC64,GNAO1,PI4KB,CBFA2T3,TIGD7,GPC2,MAD2L1,KRI1,RDH16,GMNN,CEP68,TBC1D16,CYP1A2,PCNXL3,CDK1,SMG5,KIAA0195,C8orf44,NTF3,PXDNL,NFE2L2,DNAJB11,LMNA,C16orf59,CKS2,KIF19,EFNA3,ILF2,SEPSECS,TMEM81,FTSJ3,ERCC6L,MED20,SMARCA4,GLOD4,GABARAPL1,MTUS2,GHR,CAND2,ACSL1,PRKAB1,ZCCHC2,KIF18B,GTF2IP1,FARP2,AGL,LGI1,VKORC1,THAP8,PPP2CB,MTHFD2L,ILF3,CHTF18,RPL32P3,FAM149A,CAD,MSH5,GRK6,LRAT,SYT10,RDM1
FOXD1	2.46E+00	Up	52	PRC1,NPY1R,NR4A1,NUP37,NEK2,HMGCL,SPSB4,TTC39A,CCNE1,RPL23AP7,PMP2,CYP39A1,VIPR1,MAPT,AKR1C3,DBH,GLI4,YKT6,MTHFD1L,SLCO1C1,SIPA1L3,ZNF692,EXOC7,LMNA,C21orf91,TCF19,STARD5,SLITRK6,PEMT,CEP250,EPB41L4B,DNMT3A,CCT6A,HPS5,ST8SIA6,CCT3,OSBPL3,TIGD1,RNF130,PHACTR3,MTL5,TSEN15,CENPF,MTUS2,KIAA0556,HJURP,ANKRD55,AK3,KIAA0100,PTP4A1,SAE1,DLGAP5
ZNF263	3.40E−01	Up	52	HOXD3,EPHX2,SBF1,NVL,DUSP6,SLC7A11,AARS2,R3HDM1,CLEC4G,NEIL3,UBXN8,LNX2,MYBL2,ARHGAP39,ITGA6,TJP2,KCNJ16,B4GALT3,NTN4,EHD3,PLK1,ADAM15,AXL,TSPYL1,NAT9,GJC1,INTS8,DSN1,MYO1B,CDK5,SLC7A2,CLPX,GABRP,TGFBR3,HMGB2,C22orf23,PPAP2B,FOXO1,HGF,NFE2L2,LONRF1,GINS2,ESR1,SLC28A3,NSUN6,UXS1,ANGPTL6,ACAD11,ZNF668,DEPDC1B,AK3,NUDT17
FOXO3	−3.95E−01	Down	51	CALN1,XPO5,CKAP2,CLEC4G,DNMT3B,TMEM132C,LOC100133669,CCHCR1,PDLIM5,C14orf93,ZNF74,CAMK2G,RECQL5,DSTYK,C1RL,PRMT3,THNSL1,CPLX2,MCC,DSE,ITGB1BP1,DBH,PHLDA1,HOXA13,ZC3H3,ISOC1,ATF3,FLJ12825,CENPP,ZNF174,APOL6,ZNF142,TMEM164,RFWD2,CCNE2,PCK1,CENPW,CLEC4GP1,TP73,DUSP14,NXF3,SAE1,TRAM2,NEIL3,OTUB2,CETP,PDE2A,TGFBR3,EMR1,WHSC1,SETDB1
GATA3	−1.75E+00	Down	51	SHC1,LPL,ANTXR2,HMGXB3,FCGR2B,TMEM206,SAE1,CENPP,CREBL2,LDLR,KIF18B,MRO,DDX11,PELI2,SHBG,RNF130,CD200,P4HA2,NEIL3,CENPF,DNAJC5,GOLGA6L9,ACADSB,SLC39A13,HDAC11,C1RL,PEMT,TRAF5,OSTF1,TTC13,ROBO1,SKAP1,IP6K1,HMMR,KIF15,MS4A6A,SLC7A8,GPATCH1,SFRP1,VPS72,ARFGAP1,STIL,ACOT2,CYP2C19,KIF14,TTK,RYBP,VAX2,AGPAT1,POLD1,LIN7C
FOXL1	1.37E+00	Up	51	CLCN7,PDLIM5,INS-IGF2,IL13RA2,DNAJC5,IPO9,DSCAM,SGPP1,PNP,TIGD7,TMEM81,AQP3,DPF2,ZNF330,GNA14,KIFC1,SNRPC,LOC646762,PLAC8,ASPM,SLC39A14,RACGAP1,CCBE1,STAB2,SMG5,LOC642846,IGFBP4,PELI2,G6PD,ITGA9,TPPP2,DGKI,ATG2A,CENPW,HEBP1,CENPM,PABPC1L,SDHD,CYP2R1,ALAS1,KMO,EPS8L3,KIAA1462,CCNB2,ZNF696,ADK,CXCL2,COL4A2,SLC39A1,CDC25A,DOCK3
NR4A2	−2.29E+00	Down	43	CENPL,FBP1,DDX49,FOXO1,UBE2T,ABCD3,ALDH9A1,CROCC,CCR1,FGF23,IL13RA2,ANKS3,SETD4,FAM53B,PRAM1,PDSS2,NFRKB,TMEM184B,POU6F2,TMEM164,FCN3,MKI67,STX6,SCAMP3,STAB2,HMGA1,MYOM2,TNFRSF25,TIPARP,LOC81691,POLR2J4,CDC7,CRYGS,ATF3,GABARAPL1,NECAB3,DNASE1L3,CNOT3,GMNN,TMEM53,PIGV,PCNXL3,POLR1E

4. Validation of TF-target by qRT-PCR

To validate the findings in the integrated analysis, the peripheral blood samples were used and a pair of TF-target was selected, including ZNF263 and NEIL3. The results of qRT-PCR showed that the expression pattern of selected genes in HCC were similar with that in the integrated analysis (Fig. 2). ZNF263 and NEIL3 were both upregulated in the blood of HCC patient compared with the healthy volunteers.

Fig. 2

qRT-PCR validation of a transcription factor-target pair in peripheral blood collected from five patients with hepatocellular carcinoma (HCC) and five healthy volunteers. Actin was used as an internal reference gene for normalization. CON, control (normal tissues). *p<0.05.

DISCUSSION

HCC is a complex disease that involves various molecule interactions. Chronic infection of HBV or HCV is a major risk factor in the development of the HCC, independently from alcohol abuse and metabolic disease. Moreover, the HBV and HCV infection can cause the disarrangement in cellular pathways through an indirect and/or direct mechanism in liver injury. This study was designed to determine their common mechanism by integrating sufficient number of the HCC samples with various etiologies.23

In the present study, integrated analysis of eight microarray data of HCC led to a set of 1,835 DEGs (1,145 upregulated and 690 downregulated) in HCC compared with normal tissues. Functional annotation showed that DEGs were closely related to common pathways for cancers, including p53 signaling pathway and pathways in cancer. Previous study suggested that viral hepatitis infection was associated with cellular inflammation, oxidative stress, and DNA damage, that may lead to subsequent hepatic injuries such as chronic hepatitis, fibrosis, cirrhosis, and finally HCC.23 Our results revealed that cell cycle, DNA replication and homologous recombination were significantly enriched, which may be due to the DNA damage. Moreover, DEGs were significantly enriched in various GO terms, such as oxidation-reduction process, response to external stimulus, response to stimulus, defense response, and oxidoreductase activity, which may be due to the cellular inflammation and oxidative stress in the process of viral hepatitis infection.

Network analysis allows structured grouping of genes, and network construction is an important stage in the pathogenesis studies. A study of core network in HCC revealed that miRNAs mainly regulate biological functions related to mitochondria and oxidative reduction, while TFs mainly regulate immune responses, extracellular activity and the cell cycle.24 A recent study suggested that 86 crosstalks involving 52 pathways were identified through the DEGs between adjacent nontumor and HCC samples.25 Totally, 62 differentially expressed TFs were identified in this study, and the transcriptional regulatory network was constructed, which consisted of 872 TF-target interactions between 56 TFs and 672 DEGs in the context of HCC. The HCC-specific transcriptional regulatory network may help to better understand its underlying molecular mechanism of pathogenesis and tumorigenesis of HCC.

Based on the constructed transcriptional regulatory network, a set of crucial TFs caught our attention, which covered the most downstream DEGs, including ZNF354C, NFATC2, ARID3A, BRCA1, ZNF263, FOXD1, GATA3, FOXO3, FOXL1, and NR4A2. Previous studies reported that the forkhead box O (FOXO) TFs are involved in various cancer development including HCC, suggesting that FOXO factors function as tumor suppressors in a variety of cancers.26,27 As an important member of FOXO family, threonine 32 (Thr32) of FOXO3 is critical for TGF-β-induced apoptosis via Bim in HCC.28 Our results showed that FOXO3 was down-regulated in HCC compared with the normal tissues, which provided additional evidence for FOXO3 playing a role in HCC. Up to now, little was known about the function of ZNF263. The results of a recent study revealed that ZNF263 may be closely related to the stress- and age-related diseases, and it can have both positive and negative effects on transcriptional regulation of its target genes.29 Herein, we found that ZNF263 was one of the significantly upregulated TFs, indicating that it may appeal to future investigators to study this TF in the development of the many complex diseases including HCC.

To extract more information about the crucial TFs involved in HCC, the functions of their targets were further explored. NEIL3 (nei endonuclease VIII-like 3) is a kind of DNA glycosylase, which can initiate base excision repair by hydrolysing the N-glycosidic bond and releasing the damaged base.30 It was reported that the genetic abnormalities of NEIL3 may be related to hepatocarcinogenesis.31 We found that NEIL3 was upregulated in HCC compared with normal tissues. CLEC4G is a member of C-type lectins, which are important in various immune functions, including inflammation and immunity to tumor and virally infected cells. CLEC4G was predominantly expressed in liver and was expressed at very low levels or even undetectable in liver cancer tissue.32 CLEC4G was a top-listing downregulated gene, which suggested that disruption of expression of these C-type lectin proteins may have a role in the pathogenesis of HCC.33 HCC is a highly vascularized tumor and it can be amenable to antiangiogenic treatment modalities. The vascular remodelling and endothelial transdifferentiation are major pathogenic events in HCC development, in which Stabilin-2 may play the important role.34

The present study provided useful information on the transcriptomic landscape and to the constructed transcriptional regulatory network will be hopeful to better understand its underlying molecular mechanism in hepatic tumorigenesis. These findings shed light on several important TFs targets that may potentially drive hepatocarcinogenesis, and further functional characterization are needed to verify our findings.