Main

Glioblastoma is the most frequent and aggressive primary brain tumor.1 Despite optimal treatment with a combination of surgical resection, radiotherapy and chemotherapy, patients with glioblastoma have a median survival only of 12–15 months.2 Although ionizing radiation (IR) is considered as the gold-standard adjuvant treatment for prolonging survival in patients,3 glioblastomas are often characterized by radiation resistance; therefore, universal recurrence remains a therapeutic challenge.4 Comprehensive understanding of the response of glioblastomas to IR and detailed radioresistance analysis of gliomas may help to identify radiosensitizing agents for combination treatment of the disease.

In the past decade, many advances have been made in understanding IR-induced DNA damage response (DDR), such as the activation of checkpoint pathways, repair and cell death. A previous study showed that preferential activation of the DNA damage checkpoint response and increased DNA repair capacity contributed to glioma radioresistance.5 Furthermore, targeting of the DDR signaling network in gliomas was also found to sensitize tumors to radiation therapy and reverse therapeutic resistance.6 For example, inhibition of PARP1, which is actively involved in the single-stranded DNA (ssDNA) excision repair process, can disturb the cellular repair of radiation-induced ssDNA, and therefore, PARP1 inhibitors have been incorporated into radiation therapy for glioblastomas.7 Furthermore, inhibition of WEE1 or MYT kinases has also been shown to abrogate the G2/M checkpoint and increase radiosensitivity.8, 9 Moreover, it has been found that glioma cells underwent apoptotic cell death in response to IR via both p53-dependent and p53-independent pathways and that evasion of apoptosis led to radioresistance in glioblastomas.10, 11 IR has also been reported to activate NF-κB, which transactivates its target genes, including cox-2, bcl-2, bcl-xL, XIAP and survivin,12 and the expression of these antiapoptotic genes disrupts apoptosis signaling, thereby mediating radioresistance. Vellanki et al.13 demonstrated that therapeutic targeting of XIAP increased the radiosensitivity of glioblastomas by promoting apoptosis. Furthermore, other mechanisms and pathways, including aberrant p21 regulation,14 the Notch signaling pathway,15 the Wnt pathway,16 radiation-induced Akt activation,17 as well as abnormal p53 function,18 have been correlated with radioresistance in human glioblastoma cells.

Although multiple mechanisms have been proposed for radioresistance in glioblastoma cells, the analysis of molecular signaling events is still not comprehensive. To date, advances in high-throughput sequencing methodology have provided a large amount of information regarding gene expression at the transcriptome level, as well as the underlying molecular events in response to irradiation.

In this work, we used the RNA-seq technique to investigate irradiation-responsive genes and the dynamic DDR in the glioma cell line U251 MG. We compared the genome-wide expression between the irradiation group and the corresponding control at 6, 12, 24 and 48 h after irradiation. Functional categories of differentially expressed genes (DEGs) were also analyzed based on gene ontology. Finally, we focused on the apoptosis pathway activated by irradiation and investigated the significance of proapoptosis and antiapoptosis genes in glioblastoma radioresistance. Together, these should provide prime information for research on radiotherapy of glioblastomas.

Results

Characteristics of human glioma cells treated with γ-irradiation

U251 MG cells were treated with a series of γ-irradiation doses from 0 to 30 Gy, followed by the survival assay. On the basis of the results shown in Figures 1a and b, 7 Gy, the half lethal dose at 24 h, was selected as the transition dose for further studies. After treatment with 7 Gy, the cell number remained almost unchanged until 96 h after γ-irradiation (Figure 1c). Furthermore, irradiated cells seemed to lose the proliferative activity at 48 h, with observable morphological changes, such as increased granularity (Figure 1d). Therefore, the time points before 48 h were chosen to capture the gene network controlling the growth arrest and survival.

Figure 1
figure 1

Irradiation doses characterized by cell growth arrest. (a) Cell survival rates detected by the MTT assay at 24 h after treatment with 5, 7, 15, 20 and 30 Gy. (b) The cell survival analysis at 24, 48, 72 and 96 h at the dose of 7 Gy. (c) Survival cell numbers in (b). (d) Morphological changes in U251 cells at 48 h after irradiation with 7 Gy. All experiments were repeated at least three times. **P<0.01, compared to the relevant control

Global changes in gene expression in response to irradiation

Genes with twofold changes or greater (P<0.0005) at 6, 12, 24 and 48 h were defined as DEGs, as shown in the Supplementary Information (Supplementary Dataset 1). Ten DEGs were confirmed by quantitative real-time PCR (Q-RT-PCR) for independent validation and all showed a good similarity between RNA-seq and Q-RT-PCR results. We presented six genes in Figure 2a, indicating the good quality of our data. The dissimilarities among different samples on a global scale are also demonstrated by the scatter plot and principal component analysis (Supplementary Figures S1a and b). The results of PCA analysis are consistent with the number of DEGs that show altered expression in each treatment, relative to the control (Supplementary Dataset 1). Taken together, these observations indicate that there existed much more changed genes in sample at 48 h after irradiation (R48 h) compared with control (C48 h).

Figure 2
figure 2

Changes of the gene expression in U251 at different times after irradiation. (a) Q-RT-PCR versus RNA-seq analyses of the expression for representative six genes (CDKN1A, MCM2, RAD51, CDC25A, HMGB1 and HMGB2). (b) The total number and upregulated/downregulated number of DEGs (≥twofold change, P<0.0005) at 6, 12, 24 and 48 h after irradiation. (c) Venn diagrams showing overlap of irradiation-induced DEGs across different time points. Upregulated and downregulated genes were analyzed separately and have been shown with the number of genes specifically or commonly responsive at different time points. ETU, early transiently upregulated; ECU, early continually upregulated; LU, late upregulated; ETD, early transiently downregulated; ECD, early continually downregulated; LD, late downregulated

As shown in Figure 2b, the number of DEGs induced by irradiation increased with time (180, 365, 698 and 1460 at 6, 12, 24 and 48 h, respectively), and the majority of them were either upregulated or downregulated after 24 h. Venn diagrams were then used to illustrate these radiation-affected genes (Figure 2c). Totally, 1656 genes were determined as radiation-affected genes, with 546 genes upregulated (33%) and 1110 genes downregulated (67%). Among 546 upregulated genes, 81 were upregulated at 6 h, 131 at 12 h, 144 at 24 h and 418 at 48 h; Among 1110 downregulated genes, 99 were downregulated at 6 h, 234 at 12 h, 554 at 24 h and 1042 genes at 48 h. Meanwhile, it was in particular that most of radiation-affected genes were activated exclusively at the 48 h time point, 319 out of 546 upregulated genes (58.42%) compared with only 16 (2.93%), 36 (6.59%) and 43 (7.88%) at 6, 12 and 24 h, similarly, 529 out of 1110 downregulated genes (47.66%) compared with only 7 (0.63%), 10 (0.9%) and 47 (4.23%) at 6, 12 and 24 h, respectively. To better distinguish dynamic regulation and response under radiation stress, 6, 12 and 24 h were chosen as early responsive time points, with 48 h as the late responsive time point. Genes with transcriptional expression changes at any of the time points of 6, 12 and 24 h were designated as early and transiently upregulated or downregulated genes (ETU or ETD, respectively marked as yellow or pink). Genes with expression changes at all time points among 6, 12, 24 and 48 h were grouped as early and continually upregulated or downregulated genes (ECU or ECD, respectively, marked as gray or orange). Furthermore, genes that showed altered expression only at 48 h were considered as late upregulated or downregulated genes (LU or LD, respectively, marked as light or dark green). According to this, it was calculated that among early radiation-affected genes (297), 158 (53.2%) or 139 (46.8%) genes were upregulated (ETU+ECU) or downregulated (ETD+ECD), respectively; among late radiation-affected genes (848), 319 (31.6%) or 529 (62.4%) genes were upregulated (LU) or downregulated (LD), respectively.

Overall, these observations indicate that the time points before 24 and 48 h can be considered as the early and late phases, respectively, of the continuous dynamic process and that both upregulation and downregulation have a role in early radiation responses, whereas downregulation may constitute a major part of the late radiation response in glioma cells under the test conditions.

Dynamic cellular damage responses induced by irradiation

The gene ontology analysis was conducted to analyze these early transiently, early continually and late genes separately (Table 1). The ETU genes were most enriched for genes involved in the induction of apoptosis, p53 signaling pathway and response to radiation, indicating the early stress responses of cells and the activation of related pathways upon radiation exposure. As for the ECU genes, they were most enriched in the p53 signaling pathway, positive regulation of apoptosis and cell cycle arrest, indicating a mechanism of radiation-induced cell growth arrest. During the late response, many genes involved in cell–cell adhesion and various metabolic processes were enriched (LU). Interestingly, many genes involved in the negative regulation of apoptosis were LU, suggesting the counterbalance between antiapoptosis- and proapoptosis-related genes in response to irradiation.

Table 1 Gene ontology analysis of the affected biological processes/pathways based on sets of statistically significant upregulated/downregulated genes (P<0.05)

Most radiation-downregulated genes were ECD genes and LD genes, with only 68 (6.13% of all downregulated genes) being ETD genes that were mostly enriched in the regulation of transcription. ECD genes were highly enriched for genes involved in cell cycle phase, DNA replication and repair, suggesting the regulation of cell cycle and disorders in the DNA metabolic process. During the late response, the most significant events were RNA splicing, indicating radiation-induced alterations in the post-transcriptional regulation of gene expression.

Taken together, we observed that genes involved in cell cycle progression, DNA replication and DNA repair were highly enriched at early time points after irradiation, whereas irradiation-induced proapoptosis-related genes were enriched at the early stage and antiapoptosis-related genes at the later stage.

Irradiation induces G2/M DNA damage checkpoint

It has been demonstrated that DNA damage induced by IR caused cell cycle arrest in proliferating mammalian cells.19, 20 Here, our flow cytometry analysis showed that the cell populations in the G2/M and S phases were significantly higher and lower, respectively, at 24, 48 and 72 h, compared with 0 h control (Figures 3a and b), indicative of cell cycle arrest at the G2 phase. Consistent with previous report,21 we show that irradiation induced the downregulation of many genes related to DNA replication (Supplementary Table S1), thereby inhibiting the process of DNA replication and cell cycle progression. Next, genes associated with the G2/M checkpoint were selectively analyzed. CDKN1A (p21) and GADD45A, two downstream target genes of p53 in the G2 checkpoint, were found to be increased at the transcriptional level at 6 h after irradiation, whereas the expression levels of CHEK1, WEE1, E2F1, E2F2, CDK1, CCNB1 and CDC25C were decreased in a time-dependent manner (Figure 3c). Previous studies have reported that increased p21 expression led to the repression of cyclin B1 and Cdc2 promoters and that increased GADD45A expression inhibits Cdc2 activity, thereby mediating G2/M arrest.22, 23 The G2/M checkpoint protects cell viability by allowing time for DNA repair.24, 25, 26, 27 Interestingly, according to our analysis (Supplementary Table S2), genes involved in the DNA repair pathways, including recombinational repair, base excision repair, mismatch repair and nucleotide excision repair, showed systematic repression, indicating the inefficiency of DNA repair during a certain period after irradiation.

Figure 3
figure 3

Analysis of cell cycle arrest induced by γ-irradiation. (a) The representative distribution and (b) the percentage of cells in G0/G1, S and G2/M phases at 12, 24, 48, 72 and 96 h after irradiation. Histograms of the cell percentage data are also shown (n=3). **P<0.01, ***P<0.001, compared with the control group (0 h). (c) RNA-seq data showing transcriptional expression changes of some key factors in the G2/M checkpoint at 6, 12, 24 and 48 h

Irradiation induces the counteraction of proapoptosis and antiapoptosis

The inactivation or inefficiency of DNA repair may accelerate cell death, and it is known that cell death by apoptosis has a pivotal role in γ-irradiation-induced cytotoxicity.28, 29 In this study, the activities of caspase cascades in irradiated cell lysates were investigated. We observed that the caspase-3 zymogen was cleaved into active fragments (17/19 kDa) after irradiation (Figure 4a). The enzymatic activities of caspase-3/7 and caspase-8 were also higher than the control after irradiation (Figure 4b). However, annexin V-FITC and PI staining, followed by the flow cytometric analysis, indicate that irradiation did not exert significant apoptotic effects on cells with time (Figure 4c), suggesting the existence of antiapoptosis mechanisms in response to irradiation.

Figure 4
figure 4

Apoptotic effects induced by γ-irradiation. (a) Western blotting of pro-caspase-3 and cleaved caspase-3 at 12, 24, 48 and 72 h after irradiation. (b) Proteolytic activities of caspase-8 and caspase-3/7 induced by irradiation with time. (c) The absence of apoptosis analyzed by double staining with annexin V-FITC/PI, followed by flow cytometry. Representative images are shown (n=3). (d) The heat map illustrating the dynamic regulation of apoptosis genes by irradiation, including genes involved in proapoptosis, antiapoptosis, or both positive and negative regulation of apoptosis, or genes related to apoptosis in an unclear way. The green-to-red gradient bar represents log2 values of fold-changes in the gene expression induced by irradiation

On the basis of the gene ontology analysis, we identified 83 genes that associated with apoptosis, of which 43 and 27 genes may function as positive and negative regulators, respectively, and the other 13 genes may regulate apoptosis either in both ways or in an unclear manner, according to the reported information (Supplementary Table S3). Nevertheless, the expression changes shown by the heat map revealed that radiation tended to induce proapoptosis genes (ETU and ECU) at early phases, whereas antiapoptosis genes (LU) were induced at late phases (Figure 4d).

Among proapoptosis genes, p53-dependent target genes, including TP53I3, BBC3, AEN, CYFIP2 and PHLDA3, were induced by irradiation at early time points. In contrast, death receptor genes (FAS and TNFRSF10B), as well as Bcl-2 family proapoptosis genes (BAX and MAGED1, as inhibitors of the IAP family members), were upregulated early or late after irradiation. Furthermore, antiapoptotic genes, TRAF2, MAP3K14, and members of the IAP family (such as BCL2 and BIRC3) were LD (Supplementary Figure S2).

Although the altered expression of genes described above could bring about apoptosis, genes involved in antiapoptosis were significantly induced at a late time point by irradiation. For example, the upregulation of PTGS2, NOTCH1,and BNIP3 may have important roles in radioresistance or antiapoptosis of glioma cells.30, 31, 32, 33, 34, 35

Downregulation of genes related to nuclear DNA fragmentation after irradiation

Table 2 summarizes a selective list of potential target genes whose upregulation or downregulation may promote radioresistance of glioma cells by inhibiting apoptotic cell death after irradiation, according to the reports. We found that three genes, HMGB1, HMGB2 and TOP2A were repressed which were involved in the regulation of apoptotic DNA fragmentation.36 HMGB1 and HMGB2 have been reported to activate the DFF40/CAD nuclease activity during apoptosis,37, 38, 39 whereas TOP2A (encoding topoII) has been found to have a role in the formation of DNA fragments at advanced nuclear apoptotic stages.40 On the basis of these findings, repression of these genes may result in the inhibition of apoptosis downstream of the caspase cascade. That is to say, although there have been upregulated apoptosis signals in the cytoplasm after irradiation, negative feedback of nucleic apoptosis signals might interrupt signal transmission downstream of caspase cascades to the nucleus and trigger the antiapoptosis effects. The significantly downregulated expression of these kinds of genes may reflect the cellular antiapoptosis mechanism induced by irradiation in the nuclear apoptotic stage (Supplementary Figure S2a).

Table 2 Putative target genes inhibiting apoptotic death in glioma cells after irradiation

Targeting of the TRAIL/death receptor 5 pathway increases the radiosensitivity of glioma cells

Our RNA-seq data revealed that TNFRSF10B (the gene encoding Death Receptor 5, DR5) was induced by irradiation (Supplementary Table S3), which was further confirmed by RT-PCR and western blotting (Figures 5a and b). As DR5 is known to be primarily expressed in malignant glioma cells,41, 42 we next investigated whether combination treatment with tumor-related apoptosis-induced ligand (TRAIL), a DR5 ligand, and irradiation destroys the balance between pro- and antiapoptotic factors by activating the TRAIL/DR5 apoptotic pathway. We found that after the cells were treated with TRAIL (20 ng/ml) after irradiation, the number of cells that survived was significantly lower than that for either treatment alone (Figure 5c). The FACS assay results also reflected that the early apoptotic population was increased significantly in a time-dependent manner after the combination treatment (Figure 5d; 24 h, 34.88%±3.13%; 48 h, 45.69%±2.72%), compared with the irradiation alone group (24 h, 3.4%±1.05%; 48 h, 5.73%±0.97%). Markedly increased catalytic activities of caspase-8 and caspase-3/7 were also observed in the early apoptotic population (Figure 5e). Interestingly, we found that radiation-induced downregulation of HMGB1, HMGB2 and TOP2A was significantly attenuated after combination treatment (Figure 5f), suggesting that these genes (involved in DNA fragmentation) may have important roles in the apoptotic process of caspase cascades transmission to the nucleus. Together, these results indicate that the activation of the TRAIL/DR5 signaling pathway in response to irradiation is efficacious for killing radioresistant glioma cells.

Figure 5
figure 5

Increased radiosensitivity by activating the TRAIL/DR5 pathway. (a) RT-PCR and (b) western blotting of the upregulated TRAIL receptor, DR5 (encoded by TNFRSF10B), at 6, 12, 24, 48 and 72 h after irradiation. (c) The number of surviving cells (mean±S.D.) after the indicated treatment (n=3). *P<0.05 and **P<0.01, compared with the radiation group. (d) Apoptosis detected by the FACS assay of the control, irradiation or TRAIL alone treatment, and the combined treatment at 24 and 48 h (n=3). (e) The relative enzymatic activities (mean±S.D.) of caspase-8 and caspase-3/7 at 24 or 48 h after combination treatment of TRAIL and irradiation, compared with the irradiation alone group (n=3; ***P<0.001). (f) The significantly recovered gene expression of HMGB1, HMGB2 and TOP2A detected by Q-RT-PCR after the combination treatment, compared with radiation treatment alone. Data presented are mean±S.D. values of at least three independent experiments. **P<0.01 and ***P<0.001

Discussion

DDR induced by IR was validated to be dose- and time-dependent. Understanding these complex cellular responses requires the detailed analysis of global molecular expression profiling. Here, by using high-throughput RNA sequencing, the dynamic gene expression network in glioma cells induced by irradiation was explored. We identified 1656 genes that displayed significant expression changes at or before 48 h after irradiation. The results also revealed that genes related to DNA repair, DNA replication and cell cycle arrest were transcriptionally modulated and highly enriched at early time points, whereas proapoptosis-related genes were enriched at the early stage and antiapoptosis-related genes were enriched at the later stage after irradiation.

First, our studies demonstrated irradiation-induced repression of numerous genes associated with DNA replication (Supplementary Table S1) and G2 phase arrest. The G2/M checkpoint which may be explained by transcriptional expression changes of some crucial factors involved (Figure 3) was reported to have a protective role against DNA damage-induced cytotoxicity by allowing time for damage repair, the abrogation of which may increase radiosensitivity.43, 44 In our study, U87 also sustained G2/M arrest (data not shown). Collectively, it reflected that cell cycle arrest and increased DNA repair ability may contribute to the radioresistance of glioma cells.

Second, we detected apoptotic effects induced by irradiation. Analysis of apoptosis genes subsequently revealed that the progress of apoptosis was accompanied by changes in both proapoptosis and antiapoptosis gene expression, consistent with the observation in other irradiated glioblastoma cell lines.44 Interestingly, the dynamic gene expression pattern, as shown in Figure 4d, indicated that proapoptosis genes showed early upregulation, and were followed by the late upregulation of antiapoptosis genes. Meanwhile, based on the high-throughput RNA sequencing, we explored the relation of the dynamic gene expression network to radioresisitance in different glioma cells, such as U87 at the dose of 10 Gy (which has worse radiosensitivity than U251 cell line). As shown in Supplementary Figure S2c, most of the proapoptosis genes showed early upregulation (ECU), followed by the similar early upregulation of antiapoptosis genes (ECU). These changes meant that the activation of antiapoptosis genes in U87 cell line was earlier than that in U251 cell line, which might contribute to the stronger radioresistance of U87 glioma cells.

In addition, our data indicated that p53 signaling pathway genes were mainly enriched at the early stage, although they were involved in wide range of columns, including ECU, ETU and ECD (Table 1), and associated with pro- or antiapoptosis responses to DNA damage (Supplementary Figure S2). It seemed like that the radioresistance of glioma cells was mainly due to the regulation at the later stage. Therefore, in spite of the fact that a significant status of p53 regulation in the radioresistance mechanism,8, 9, 10, 11, 14, 18, 24, 45 more studies are needed.

Furthermore, we analyzed factors whose gene expression alterations may contribute to the radioresistance of glioma cells with a focus on the upregulated antiapoptosis genes and the downregulated proapoptosis genes. First, among the upregulated antiapoptosis genes, PTGS2, encoding COX-2, has been reported to increase the radiosensitivity of glioblastoma cells by regulating Ku expression30, 31 and activating the PI3K and PKA pathways.46 Although BNIP3 is a pro-cell death member of the Bcl-2 family, it was found to be localized in the nucleus of glioblastoma cells and repress the gene expression of apoptosis-inducing factor (AIF), thereby preventing cell death.34, 35 In addition, NOTCH1 could activate gene expression programs via the translocation of the intracellular NOTCH domain (NICD).33, 47 Thus, the upregulated gene expression of these three molecules may inhibit the apoptotic effect caused by irradiation and thereby increases the radioresistance of glioma cells. Second, among the downregulated proapoptosis genes, p73 (encoded by TP73), with functions similar to those of p53, was reported to participate in the apoptotic response to DNA damage and had an important role in cellular radiosensitivity.45, 48, 49 Moreover, according to the fact that the downregulation of HMGB1, HMGB2 and TOP2A transcriptional levels was attenuated significantly after combination treatment of TRAIL and irradiation (Figure 5f), we further studied the changes of these three proteins. As shown in Supplementary Figure S3a, HMGB1 indeed decreased at 48 h in a dose-dependent manner (Supplementary Figure S3b), and no significant changes for the other two proteins (data not shown). HMGB1 restored after 15 Gy, which was consistent with the significant activation of caspase-3 at 48 h (Supplementary Figure S3c). Although these results indicate the corresponding signaling cascade from outer space to inter space of the nucleus (Supplementary Figure S3d), it is still unknown whether this type of regulatory network is a key mechanism of the apoptosis-resistance and cell survival in glioma cells after irradiation. We propose that these three genes may be potential targets for increasing radiosensitivity.

On the other hand, p53-targeting genes (PUMA, AEN, BAX and MAGED1) and death receptor genes (FAS and TNFRSF10B) involved in glioma cell death progress were upregulated. The enhancement of these proapoptosis signaling may also increase its radiosensitivity. In particular, the TRAIL/DR5 signaling pathway has been confirmed in glioma cells.50 We found that irradiation could induce significant upregulation of DR5 (encoded by TNFRSF10B) and the combination treatment with TRAIL and irradiation obviously reduced the number of surviving cells and activated caspase cascades, compared with that of irradiation treatment alone (Figures 5c–e). This meant that the activation of TRAIL/DR5 signaling could increase the radiosensitivity and cell-killing efficacy of glioma cells, and the combination of radiation and TRAIL led to additive effects and low synergy in apoptosis induction only, consistent with the results reported by Kuijlen et al.51

Meanwhile, it is noticeable that many genes involved in RNA splicing (based on gene ontology, Table 1) were significantly downregulated during late phase responses. It was reported that γ-irradiation-induced stress affected the splicing of many genes by repressing and redistributing splicing factors and other components of the transcription machinery, which might be used by cells as a survival mechanism to adapt to the stress environment by producing specific mRNA variants.52, 53, 54, 55, 56

In summary, we first identified the radiation-induced dynamic gene expression network involved in various DDRs at a dose that induced cell growth arrest. In particular, apoptotic response was counteracted by dynamic changes in proapoptosis and antiapoptosis gene expression, as well as by the complex interactions between signaling molecules in the cytoplasm and nucleus. These provide important references and sources for exploring radiosensitizing targets for glioblastoma therapy. Studies on the dynamic network of molecules involved in cell death at different irradiation doses are under way to clarify the regulatory mechanisms in glioma cells after irradiation and identify the candidate targets.

Materials and Methods

Cell culture, irradiation and cell proliferation assay

The human glioma cell line U251 MG and U87 was purchased from Cell Center of Peiking Union Medical College (Beijing, China) and cultured in MEM supplemented with 10% fetal bovine serum (Hyclone, Beijing, China), 100 units/ml penicillin-100 μg/ml streptomycin (Beijing Solarbio Science & Technology Co., Beijing, China). The culture was maintained at 37 °C in a humidified atmosphere containing 5% CO2. Cells were irradiated using γ-ray of Co-60 source under atmospheric pressure and ambient temperature (Peking University, Beijing, China). The dose rate of 1 Gy/min was used. Cells were seeded in T-25 flask (Corning, NY, USA) with a density of 2.4 × 105 or 4.3 × 105 cells per flask. For combination treatment, cells treated with 20 ng/ml TRAIL (PeproTech, Rocky Hill, NJ, USA) immediately after irradiation. Cell proliferation and diameter were assessed by counting live cells after enzyme digestion using Scepter hand-held automated cell counter (Milipore, Billerica, MA, USA).

MTT (dimethyl thiazolyl diphenyl tetrazolium) assays

MTT assay was used to evaluate cell survival rate after irradiation. U251 cells were seeded at 5.0 × 103 cells per well of 96-well plates with at least three replicate wells, and then were treated with different irradiation dose at 24 h or irradiated with 7 Gy at the time points of 24, 48, 72 and 96 h. Twenty microliters MTT (Sigma, St. Louis, MO, USA) reagent was added to each well, and then cells were incubated for 4 h. After that, the supernatant was removed and 100 μl dimethyl sulfoxide (DMSO) was added. OD value of each sample was measured at a wavelength of 570 nm. The cells without irradiation were used as control at different time point. All experiments were repeated at least three times and finally the survival rate was calculated as following:

Cell survival rate (%)=OD(the experimental group)/OD (control) × 100%

RT-PCR and quantitative real-time PCR (Q-RT-PCR)

Total RNA from cells treated with or without irradiation was isolated using Trizol reagent (Invitrogen, Carlsbad, CA, USA). RNA integrity was detected by analyzing 18S and 28S bands on a 1% agarose gel. Purity of all RNA samples were determined by the OD260/OD280 rations using a Spectrophotometer (Eppendorf, Hamburg, Germany). Two micrograms RNA from different samples was reverse transcribed to cDNA for gene expression detection by reverse transcriptase PCR (RT-PCR) or Q-RT-PCR. The relatively quantitative PCR reaction was set up in a total volume of 15 μl containing 3 μl of diluted cDNA, 7.5 μl of 2 × SYBR Green Real-time PCR Master mix (Toyobo Co. Ltd., Osaka, Japan), 0.25 μl 10 μmol/l of each primer and 4 μl ddH2O. Q-RT-PCR runs were performed on MX3000P instrument (Stratagene, La Jolla, CA, USA). GAPDH levels were used to for normalization and non-reversed transcribed RNA was used to correct for the presence of genomic DNA. The relative expression of genes was calculated as follows: fold-change=2 ^ (Treatment (Tc target−Tc GAPDH)−control (Tc target−Tc GAPDH)). Primer pairs for 13 genes are detailed in Supplementary Dataset 2.

Illumina sequencing and bioinformatics analysis

Total mRNA of eight samples (including samples at 6, 12, 24 and 48 h after radiation and their corresponding controls) were sent to whole transcriptome sequencing by Illumina HiSeq 2000. The gene expression level was calculated using RPKM (reads per kilobase per million reads). If there is more than one transcript for a gene, the longest one is used to calculate RPKM for the whole genes. DEG were identified by combination of fold-change (≥2) and FDR (false discovery rate; ≤0.001) (Supplementary Dataset 1). The Gene Ontology enrichment analysis was performed using DAVID Bioinformatics Resources 6.7, NIAIS/NIH (http://david.abcc.ncifcrf.gov/). The heat map illustrating dynamic changes of apoptosis-related genes (Figure 4d) was performed using Mev. The scatter plot and PCA analysis (Supplementary Figure S1) were performed via R project.

Analysis of cell cycle progression

Cells were collected after irradiation at different time point (0, 12, 24, 48, 72 and 96 h), and washed with PBS twice, then fixed in 1 ml 70% cold ethanol in test tubes and incubated at −20 °C overnight. After incubation, cells were centrifuged at 1000 r.p.m. for 5 min and the cell pellets were resuspended in 500 μl propidium iodide (PI)/Triton X-100 staining solution containing 50 μg/ml PI (Roche, Penzberg, Germany), 0.1% (v/v) Triton X-100 (Sigma) and 100 μg/ml RNase (Sigma). Then cells were incubated on ice for 30 min and given assay by FC 500 MCL. Cell cycle distribution was calculated from 10 000 cells with ModFit LT software (Becton Dickinson, San Jose, CA, USA).

Annexin V-FITC/PI FACS apoptotic assay

According to the instruction of FITC-Annexin V Apoptosis Detection Kit (Pharmingen BD, San Diego, CA, USA), cells were trypsinized and washed twice with cold PBS, then resuspended in 200 μl 1 × binding buffer. Hundred microliters of cell suspension was transferred to a 5 ml culture tube and incubated with 5 μl of FITC-Annexin V and 10 μl of PI (10 μg/ml). Gently vortex the cells and incubate for 15 min at RT in the dark. Five hundred microliters 1 × binding buffer was added to each tube and the cells were analyzed with flow cytometry FC 500 MCL within 1 h.

Western blotting

Cells were washed twice in cold PBS, and lysed in ice-cold lysis buffer containing 150 mM NaCl, 1.0% Nonidet-P40 and 50 mM Tris-Cl (pH 8.0). Hundred micrograms of the whole-cell protein were electrophoresed by SDS-PAGE and then blotted to an Immobilon-P 0.22 μm membrane (Millipore, Billerica, MA, USA), followed by blocked with 5% nonfat milk in 1% TBST and incubated with primary antibodies overnight. The secondary antibody was conjugated to horseradish peroxidase (CWBIC, Beijing, China), and antigen–antibody complexes were detected by the chemiluminescence (Millipore). Both the primary and secondary antibodies were diluted in 5% nonfat milk in TBST. The primary antibodies were anti-caspase-3 (Abcam, Cambridge, MA, USA), anti-DR5 (Santa Cruz Biotechnology, Santa Cruz, CA, USA) and anti-β-actin (Sigma).

Caspase activity assays

The caspase-3/7 and caspase-8 activation was performed using Apo-One homogenous caspase-3/7 assay and caspase-8 assay (Promega, Madison, WI, USA), respectively, according to the manufacturer’s instructions. Briefly, 1 × 104 cells (treated with or without irradiation) were collected at different time point (0, 12, 24, 48, 72 and 96 h) and lysed in the manufacturer-provided homogeneous caspase-3/7 or caspase-8 reagent. The lysates were incubated at room temperature for 1.5 h before reading in a fluorometer at 485/530 nm. The relative caspase activity was given to evaluate the fold-changes of samples at 6, 12, 24, 48, 72 and 96 h (compared with sample at 0 h).

Statistical analyses (excluding Illumina data)

All experiments were performed at least three times and each time was done in triplicate. The results are shown as the mean values±standard deviation (S.D.), and statistical significance was evaluated by Student’s t-test and ANOVA assay; P values were considered to be statistically significant when less than 0.05.