Identification of Key Genes Affecting Results of Hyperthermia in Osteosarcoma Based on Integrative ChIP-Seq/TargetScan Analysis.

BACKGROUND The purpose of this study was to research the effects of hyperthermia on osteosarcoma (OS) by integrating the Chromatin Immunoprecipitation with the generation sequencing (ChIP-Seq) and TargetScan analysis of heat shock transcription factor 1 (HSF1). MATERIAL AND METHODS The HSF1 ChIP-seq dataset of GSE60984 was downloaded from the Gene Expressed Omnibus (GEO) database. The HSF1-binding sites were screened by MACS2 in OS cells after 10 and 20 min of hyperthermia, and they were annotated using the ChIPseeker package. The overlapped genes were selected out when HSF1-binding sites were located in the promoter region. The Database for Annotation, Visualization, and Integrated Discovery (DAVID) was used to perform Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis of the overlaps. The miRNA-gene pairs of the overlaps were screened out via TargetScan, and the miRNA-gene-regulated network was constructed by Cytoscape software. RESULTS 1880 and 1283 genes of promoter regions were obtained in the osteosarcoma cells after 10 and 20 min of hyperthermia, respectively, and 889 of them were overlapped. The overlapped genes were enriched in 122 GO terms and 3 KEGG pathways. There were 13 657 pairs involved in the miRNA-gene regulated network of the overlaps. CONCLUSIONS Some biomarkers were identified for OS treated with hyperthermia. Moreover, some GO terms (regulation of programmed cell death and regulation of cell death) and pathways (p53 signaling pathway, methane metabolism, and viral myocarditis) might be involved.

Background

Osteosarcoma (OS) is a common histological form of bone cancer, which is most prevalent in children and adolescents [1,2]. It was reported that OS was the eighth most common childhood cancer, accounting for about 20% of all primary bone cancers and 2.4% of all malignancies in pediatric patients [1]. Symptoms include pain at night, overt localized swelling, and bones that fracture easily. The causes may lie in OS stem cells and radiotherapy [3,4]. Various treatments, such as wide-excision surgery of tumors, radiotherapy, chemotherapy, and neoadjuvant chemotherapy, are usually adopted in OS [1]. Although significant improvements in the long-term outcome of OS patients have been made, it is still not satisfactory [5,6]. The best reported 10-year survival rate was 92%, and the protocol was an aggressive intra-arterial regimen that individualizes therapy based on arteriographic response [7]. The 3-year survival rates ranged from 50–75%, and 5-year survival rate ranged from 60–85% [8,9]. In addition, 30–40% children died of it, and 25–50% patients subsequently developed metastatic disease, which was the major cause of death [10]. Recently, some biomarkers for OS were selected, such as CD20, ECM, ITGA, and PMP22 [11,12].

Hyperthermia therapy is a treatment approach in which a particular area or the whole body is heated above normal temperatures to achieve therapeutic effects. Hyperthermia has been used in cancer treatment since 1989, and currently it is used either for ablation purposes as an alternative to surgery, or, less frequently, in combination with chemotherapy or radiation therapy [13]. It was reported that hyperthermia can activate various systemic anti-tumor immune responses, and leads to methamphetamine-induced toxicity [13,14]. The principle may be as follows: heat alters membrane characteristics, leading to modification in cell morphology, intracellular sodium and calcium levels, and membrane potential [15-17]. Some studies reported that hyperthermia can induce apoptosis and reduce migration of OS. In this study, we aimed to screen gene expression changes of OS induced by hyperthermia via integrating the Chromatin Immunoprecipitation with generation sequencing (ChIP-Seq) and TargetScan analysis, and further explored the related potential biomarkers and mechanism.

Material and Methods

The ChIP-seq dataset

The ChIP-seq dataset of GSE60984 [18] was downloaded from the Gene Expressed Omnibus (GEO) database (). In this dataset, there were 3 samples of OS cell: one was treated with both anti-HSF1 antibody and heat shocked for 10 min, another was treated with anti-HSF1 antibody and was heat shocked for 20 min, and the last was treated with neither antibody nor hyperthermia (the control sample). The sequencing platform was Illumina Genome Analyzer IIx (Homo sapiens).

Data pre-processing

Raw fastq read sequences were demultiplexed, and cutadapt was used to remove adapter sequences from DNA high-throughput sequencing data. Afterwards, the expressed genes were aligned to a reference genome using Bowtie [19], an ultrafast genome aligning tool that was used for the mapping of the raw reads to enesembl () hg19 genome with maximum of 2 mismatches in every reads and those who had only 1 match were reserved.

Screening and annotation of HSF1-binding sites

Model-based analysis of ChIP-seq (MACS) was a computational algorithm that identified genome-wide locations of transcription/chromatin factor binding or histone modification from ChIP-seq data [20]. The HSF1-binding sites (also known as peaks) were respectively obtained via MACS2 in the OS cells heat shocked for 10 and 20 min compared with the control samples, and recorded as Peek-10 and Peek-20, respectively. The threshold was P<10−5 and the default ultrasonic interrupt length was 300 bp. ChIPseeker was an R package for annotating ChIP-seq data analysis [21]. Here, based on the ChIPseeker package, their related genes, distance to the closest transcription start sites (TSS), and genome features were assigned to the peaks. Moreover, some genes of Peek-10 and Peek-20 were separately selected out when HSF1-binding sites were located in the promoter regions, and their overlaps were found.

Functional and pathway enrichment analysis of overlaps

The Database for Annotation, Visualization and Integrated Discovery (DAVID) [22] () is a widely used web-based tool for functional and pathway enrichment analysis. In this study it was used to perform Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis of the above overlaps. Moreover, the GO terms and KEGG pathways were selected out with the criteria of P value <0.05.

Construction of the miRNA-gene network

TargetScan is a widely used database that predicts biological targets of miRNAs by searching for the presence of conserved 8mer, 7mer, and 6mer sites that match the seed region of each miRNA [23,24]. The miRNA, which might regulate the overlaps, as well as the miRNA-gene pairs, were screened out with TargetScan. Furthermore, the miRNA-gene-regulated network was constructed based on the above miRNA-gene pairs. However, the selected miRNA-gene pairs might be overwhelming, and the miRNA-gene regulated network too complicated, so that it was not conducive to analyze and the accuracy was not high. Hence, module analysis of the network was conducted and the modules were visualized using Cytoscape software [25].

Results

HSF1-binding sites

A total of 19 626 and 18 015 HSF1-binding sites were obtained in the OS cells heat shocked for 10 and 20 min, respectively, compared with the control samples (recorded as Peek-10 and Peek-20). Their distributions in the genome-wide are shown in Figure 1, showing 10.017% Peek-10 and 7.921% Peek-20 located in the promoter regions. In other words, 1880 genes of Peek-10 and 1283 genes of Peek-20 were located in the promoter regions. Besides, 889 overlaps of them were found, and 40 of the overlaps are shown in Table 1.

Figure 1

The genome-wide distributions of Peek-10 and Peek-20.

Table 1

40 out of the overlapped genes in the promoter regions.

Gene
SGMS2	SOD1	DNAJC7	SPIN3
CGGBP1	GNAS	SPOPL	ACER2
FTL	NCF4	FAM109B	C1GALT1C1
HNRNPH2	ACHE	SMAD1	IL32
JAK3	ADRA2C	OR13D1	KLK6
MEFV	ALDH3B1	GPX8	USP21
MLH1	OPRD1	ZNF761	FAM13A
PGK1	PPIC	TPST2	KAZN
PLOD1	PRIM2	CCDC147	NARG2
KRT5	RPL24	LPPR5	DMAP1

Enriched GO terms and KEGG pathways of overlaps

A total of 122 GO terms were found to be enriched in the overlapped genes, and the top 10 most significant GO terms are exhibited in Figure 2. The overlaps were enriched in 3 KEGG pathways, namely p53 signaling pathway, methane metabolism and viral myocarditis.

Figure 2

The top 10 most significant GO terms of the overlaps.

The miRNA-gene network

The miRNA-gene regulated network of the overlaps was constructed, including 13 657 miRNA-gene-regulated pairs. The network was divided into 4 modules, and they were recorded as Module 1, Module 2, Module 3, and Module 4. Modules and their attributes are shown in Figure 3 and Table 2. SGMS2 and CGGBP1 were included in Module 1 and were regulated by more miRNAs.

Figure 3

Modules of the miRNA-gene network.

Table 2

The 4 modules of the miRNA-gene network and some their attributes.

Cluster	Score	Nodes	Edges
Module 1	8.59	62	262
Module 2	6.625	33	106
Module 3	4.556	37	82
Module 4	4	29	56

Discussion

Heat shock transcription factor 1 (HSF1) is localized to the cytoplasm as an inactive monomer, and heat shock can induce the activation of HSF1 [26]. The active HSF1 leads to synthesis and cytoprotection of heat shock proteins (HSPs) [27]. However, an increase in HSPs might render hyperthermia less effective [28,29]. In this study, the HSF1 binding sites in promoter regions were analyzed after 10 and 20 min of hyperthermia using the ChIP-Seq approach. After the functional enrichment analysis, more overlapped genes were enriched in processes of regulation of programmed cell death and regulation of cell death (Figure 2). Programmed cell death, especially the expressions of some receptors and ligands in the process, play critical roles in the progression and resistance to oncolytic adenovirus of OS [30-32]. Regulations the cell death processes of tumor cells and normal cells has become an important way to increase efficacy and reduce resistance to traditional treatments in OS [33-35]. Studies reported that hyperthermia induced apoptosis and enhanced cell death in human OS cells [36,37]. Therefore, processes of programmed cell death and cell death might be involved in the mechanism of hyperthermia in OS, and their rational regulation might contribute to a better curative effect.

The overlaps were found enriched in 3 KEGG pathways, p53 signaling pathway, methane metabolism, and viral myocarditis. The P53 signaling pathway was involved in many tumors, and OS is no exception. It was reported that the p53 signaling pathway participates in cell proliferation, metastasis, and angiogenesis of OS [38-40]. Furthermore, the p53 pathway is closely associated with cell cycle and apoptosis, which are involved in the enriched GO terms used in this study [41-43]. Matsumoto et al. [44] found that nitric oxide is an initiator of intercellular signal transduction after hyperthermia in mutant p53 cells of human glioblastoma. Few researchers have explored the functions of methane metabolism and viral myocarditis in OS. Our results indicate that the p53 signaling pathway, methane metabolism, and viral myocarditis might play critical roles in gene changes of OS regulated by HSF1. However, further in-depth research is needed, especially on the relationships between OS and pathways of methane metabolism and viral myocarditis.

After TargetScan analysis of the overlaps, the miRNA-gene-regulated network was constructed, including 13 657 miRNA-gene regulated pairs that might be unreliable. Consequently, the network was divided into 4 modules (Figure 3). Module 1 contained the most nodes and pairs, and SGMS2 and CGGBP1 were regulated by more miRNAs than other nodes. Therefore, SGMS2 and CGGBP1 are important nodes in the miRNA-gene-regulated network of the overlapped DEGs. SGMS2 encoded sphingomyelin synthase 2 (SGMS2), which is a risk factor of liver steatosis and atherosclerosis [45,46]. Hailemariam et al. [47] found that SGMS2 is a modulator of NFκB activation. Moreover, NFκB signaling pathways is closely related to the proliferation and progression mechanism of many cancers, and targeting it has been used in cancer prevention and therapy [48-50]. CGGBP1 encoded CGG triplet repeat-binding protein 1 (CGGBP1), which is a nuclear and midbody protein regulating abscission [51]. CGGBP1 expression is important for cell cycle progression in various cancer cell lines [52]. Cell cycle directly affects cell proliferation and apoptosis, which is important in the progression and prognosis of OS. Although few studies have verified that SGMS2 and CGGBP1 are associated with effects of hyperthermia on OS, the present study indicates that they might be key genes of OS treated with heat shock.

The GSE60984 data set was provided by Janus et al. [18], who mainly researched potential genes regulated by TNF-α cytokine and heat shock based on ChIP-Seq analysis, and explored functions of NF-κB signaling pathways on OS in the process of hyperthermia. In this article, some different methods were used and some novel findings were obtained. For example, functional and pathways enrichment analysis were performed, and some different biological processes (e.g., programmed cell death and cell death) and more pathways (e.g., p53 signaling pathway, methane metabolism, and viral myocarditis) were found to be involved in the effects of hyperthermia on OS. Moreover, the miRNA-gene network was constructed and analyzed with TargetScan and Cytoscape, which contained more interactions and thus might be more reliable. We also screened out some different biomarkers (e.g., SGMS2 and CGGBP1).

Conclusions

SGMS2 and CGGBP1 might be biomarkers of heat shocked OS, and programmed cell death and cell death might be involved in the mechanism of hyperthermia in OS. The p53 signaling pathway, methane metabolism, and viral myocarditis might play critical roles in gene changes in OS regulated by HSF1.