Transcriptomic profiling and pathway analysis of cultured human lung microvascular endothelial cells following ionizing radiation exposure

The vascular system is sensitive to radiation injury, and vascular damage is believed to play a key role in delayed tissue injury such as pulmonary fibrosis. However, the response of endothelial cells to radiation is not completely understood. We examined the response of primary human lung microvascular endothelial cells (HLMVEC) to 10 Gy (1.15 Gy/min) X-irradiation. HLMVEC underwent senescence (80–85%) with no significant necrosis or apoptosis. Targeted RT-qPCR showed increased expression of genes CDKN1A and MDM2 (10–120 min). Western blotting showed upregulation of p2/waf1, MDM2, ATM, and Akt phosphorylation (15 min–72 h). Low levels of apoptosis at 24–72 h were identified using nuclear morphology. To identify novel pathway regulation, RNA-seq was performed on mRNA using time points from 2 to 24 h post-irradiation. Gene ontology and pathway analysis revealed increased cell cycle inhibition, DNA damage response, pro- and anti- apoptosis, and pro-senescence gene expression. Based on published literature on inflammation and endothelial-to-mesenchymal transition (EndMT) pathway genes, we identified increased expression of pro-inflammatory genes and EndMT-associated genes by 24 h. Together our data reveal a time course of integrated gene expression and protein activation leading from early DNA damage response and cell cycle arrest to senescence, pro-inflammatory gene expression, and endothelial-to-mesenchymal transition.


Results
Accelerated senescence observed in endothelial cells after 10 Gy X-irradiation. We previously reported that the primary response of human pulmonary artery endothelial cells (HPAEC) to 10 Gy X-irradiation was accelerated senescence 24,30 . However, it was postulated that ECs from large vessels have different responses to radiation compared with microvascular ECs 14 . Therefore, we investigated the effects of radiation on human lung microvascular endothelial cells (HLMVEC) exposed to 10 Gy X-irradiation. Accelerated senescence, necrosis, and apoptosis were measured. Senescence was significantly increased 24-72 h (Fig. 1a-f). Senescent cells displayed characteristic "fried egg" morphology with expression of senescence-associated beta galactosidase ( Fig. 1a,b). Cells did not display significant necrosis, as indicated by lactate dehydrogenase release in the medium in response to radiation from 24, 48, and 72 h (Fig. 1c). We were unable to detect significant apoptosis using western blotting for cleaved caspase-3 (Fig. 1d). We next performed nuclear morphology analysis to detect low levels of apoptosis occurring over the time course (Fig. 1e,f). Nuclear blebbing is consistent with late apoptotic events. We observed a trend toward 3.6-11.9% apoptosis at 24-72 h post-irradiation, respectively, although these levels did not reach significance due to variation between experiments (Fig. 1f). Our laboratory and others demonstrated the regulation of a variety of signaling pathways in vitro response to radiation, including pathways regulating growth regulation, DNA repair, and cell death 24,30,47 . We investigated alterations in HLMVEC mRNA changes in an early time course (10-120 min) following 10 Gy X-ray irradiation exposure (Fig. 2). The most rapid change was in IGF1R mRNA, which increased significantly at all time points tested with maximum expression (~ 5-6.5-fold, p < 0.05) from 10 min through 2 h post-irradiation. www.nature.com/scientificreports/ (d) HLMVEC were grown to 70% confluence and X-irradiated at 10 Gy. Cell lysates were prepared at the indicated time points, and western blots were performed for cleaved caspase 3 as an indication of apoptotic signaling. A positive control is provided for cleaved caspase 3 (pos control). (e) Cells stained with DAPI to visualize nuclear morphology to identify apoptotic nuclei. (f) Nuclei were scored from random fields to determine percentage of apoptotic nuclei; graph shows average of percent apoptosis ± SEM; NS = not significant. www.nature.com/scientificreports/ Cyclin dependent kinase inhibitor 1A (CDKN1A) increased by fourfold, and murine double minute 2 (MDM2), an E3 ubiquitin ligase regulator of p53, increased ~ threefold, both at 2 h (both p < 0.05). Phosphatidylinositol 3-kinase (PI3K) proteins integrate cell growth in response to stress. The PIK3CD isoform of PI3K was increased at 20 min post-irradiation (p < 0.05). We also observed trends toward increased expression in several stress response and survival genes: homeodomain interacting protein kinase-2 (HIPK2), mammalian target of rapamycin (MTOR), Ras association domain family member (RASSF5), and ribosomal protein S6 kinase B1 (RPS6KB1) 48,49 . Interestingly, we observed no significant changes gene expression for genes encoding proteins associated with DNA repair, stress response, or cell growth regulation: the DNA damage response serine/threonine kinase ataxia-telangiectasia mutated (ATM) gene; the transcription factor tumor protein p53 (TP53), which Figure 2. Gene expression changes in irradiated HLMVEC. HLMVEC were grown to 70% confluence and exposed to 10 Gy X-ray irradiation. RNA was obtained at the indicated time points. mRNA levels in irradiated HLMVEC were assessed by RT-qPCR at indicated time points post-irradiation. Graph represents means, ± SEM from n = 3 independent experiments. p < 0.05 is indicated by *, p < 0.01 is indicated by **. www.nature.com/scientificreports/ regulates senescence-associated genes; and sirtuin 1 (SIRT1), which is activated by oxidative stress and DNA damage to regulate survival. Western blotting was performed to examine the regulation of selected proteins by radiation in a time course from 15 min through 72 h (Fig. 3). Consistent with the RT-qPCR findings, we observed an upregulation of p21/ waf1, (gene CDKN1A) and MDM2 within 2 h post-irradiation that was sustained through 72 h. Phosphorylated Akt, often a surrogate of PI3K activation, was increased at 2 h post-irradiation, consistent with the time point for upregulation of PI3KCD mRNA. Total and phosphorylated ATM were increased within 15 min, consistent with the role of ATM in DNA damage response. Increase in phosphorylated ATM was sustained through 4 h, but the increase in total ATM was maintained through 72 h. The level of p53 and IGF1R displayed trends of reduced levels at early time points, that returned to basal levels ~ 2 h post-irradiation. We unable to observe significant increases in IGF1R phosphorylation, but this may be because cells were not synchronized prior to irradiation.

Scientific Reports
Genome-wide transcriptional responses to 10 Gy X-irradiation. Targeted gene expression analysis provides information regarding known genes with predicted functions in biological processes. To discover previously unidentified gene regulation, and to expand known pathways in primary HLMVEC in response to ionizing radiation, we used comprehensive transcriptome profiling by RNA sequencing (RNA-seq). Gene expression profiles from sham-irradiated (control) HLMVEC were compared with irradiated HLMVEC at 2, 4, 6, 8, and 24 h. Comparative differential expression analysis identified differential gene expression of 9,151 transcripts over all time points, with a q-value < 0.05 (Supplemental Table S1). In order to interpret the functional relationship of candidate DEG in response to ionizing radiation, we filtered for coding transcripts that were differentially expressed in excess of 1.5-fold as compared to the control group. Based on these criteria, we observed 3,581 unique significant DEGs in the irradiated cells irradiated vs the control cells, irrespective of time point (Supplemental Table S1).
The differences in gene expression patterns post-irradiation compared to controls were very robust at 2, 8 and 24 h, indicating significant gene expression activity at these time points (Fig. 4a). Because of this we focused on these three time points for our initial analysis to gain an overall picture of changes in gene expression, and later considered all of the time points in targeted biological pathway analysis. Two Venn diagrams of the data from 2, 8 and 24 h illustrate the overlap of the differentially expressed genes up-and down-regulated in response to 10 Gy irradiation (Fig. 4b). At 2 and 8 h post-irradiation, downregulated genes outnumbered upregulated genes, with 442 and 732 genes downregulated at 2 and 8 h respectively and 322 and 521 upregulated, respectively. At 24 h post-irradiation the number of DEGs upregulated vs downregulated were reversed with 1,651 genes upregulated and 1,271 genes downregulated. mRNA expression levels of selected genes were analyzed using real-time qPCR. We selected genes (ACTA2, BBC3, BCL2A1, CCNE2, E2F1, ESPL1, FAS, IGFBP3, ORC1, SESMN1, and VIM) that were highly regulated in our sequencing data and represented in the six biological pathways we identified as affected by 10 Gy X-irradiation (cell cycle, apoptosis, DNA damage, inflammatory response, senescence, and EndMT). The results from RT-PCR correlated well in the direction of regulation with the sequencing-derived mRNA expression levels. Overall, the RNA-seq results were in good agreement with the qPCR, although there were a few instances in which the results diverged (Supplemental Fig. S1). The 2 h time point provided most of the differences, with expression measured using qPCR higher than RNA-seq. In many of these cases, the qPCR expression is also more variable than the RNA-seq results. The other three times points do not show significant differences suggesting the issue was in the samples or execution of the qPCR at one time point.
Gene ontology and pathway analysis post-irradiation. To identify enrichment of pertinent gene ontology (GO) term clusters, we analyzed 2, 8 and 24 h gene sets using DAVID (v 6.8) 40,41 . Analyses were focused on GO terms and KEGG pathways relevant to cellular response pathways and to the endothelium. The biological processes (BP) graphs show clusters with enrichment scores greater than 1.6 and the y-axis is labeled with the term having the highest fold enrichment in that cluster (Fig. 5). For the KEGG pathway graphs we included all of the relevant to cell biology; these pathways were visualized using DAVID (Supplementary Figs. S2-S4). GOrilla (v 4.1) 43,44 was used to identify biological processes enriched (p ≤ 10 -5 ) in our data, with a focus on categories relevant to cell biology with cluster enrichment scores greater than 1.6, labeling the graph with terms having the highest fold enrichment in the cluster. The pathways with the most enrichment for downregulated pathways are shown in Fig. S5, and the pathways with the most enrichment for upregulation are shown in Fig. S6. DAVID analysis showed that at 2 h post-irradiation, the largest changes (both up-and down-regulation) in BP terms were cell cycle, apoptosis, and DNA damage terms, as well as the related terms, cell proliferation and transcription regulation (Fig. 5a). Consistent with these, we observed enriched KEGG pathways for p53, TGF-β, and TNF signaling, as well as cell cycle at 2 h (Fig. 5b). These pathways, as well as FoxO, and Hippo, suggested the regulation of genes associated with cellular senescence and apoptosis. GO terms related to hypoxia, cell migration, and mesenchymal cell development were also increased. Remarkably, analysis of down-regulated GO term clusters enriched showed that all major down-regulations at 2 h and the majority of the clusters at 8 and 24 h involve cell cycle terms (Supplementary Figs. S5). The GOrilla GO term analysis agreed with the DAVID GO clustering and KEGG pathway analysis, with the most enriched terms associated with cell cycle events (spindle checkpoint, sister chromatid segregation, and metaphase/anaphase transition) and signal transduction by p53 class mediator, and transcription regulation. GOrilla analysis also revealed enrichment in terms related to the TGF-β signaling pathway, positive regulation of cell migration, and response to hypoxia.
At 8 h post-irradiation, DAVID analyses showed that the most enriched BP cluster was response to hypoxia, followed by terms related to regulation of proliferation, motility, and migration, suggesting a continuation of the response seen at 2 h post-irradiation (Fig. 5c). Clusters of gene regulation terms for apoptosis, cell cycle, and DNA damage were present but to a lesser degree. At  www.nature.com/scientificreports/ Figure 3. Changes in protein levels and phosphorylation in irradiated HLMVEC. HLMVEC were grown to 70% confluence and exposed to 10 Gy X-ray irradiation. Protein lysates were prepared at the indicated time points and used for western blotting for the indicated proteins and phospho-proteins. www.nature.com/scientificreports/ followed by the cell cycle and tumor necrosis factor (TNF) and FoxO signaling pathways (Fig. 5d). GOrilla analysis identified enrichment of terms associated with hypoxia response, RNA processing, kinase activity, and cell cycle phase transition. The clusters, pathways, and terms enriched at 8 h suggest that the cells responded to www.nature.com/scientificreports/ radiation by arresting the cell cycle, possibly through p53 signaling, that may be coordinated with DNA repair and the induction of accelerated senescence through FoxO and TNF signaling pathways. At 24 h post-irradiation, there was an increase in the numbers of regulated genes (2,922 genes > 1.5-fold [q ≤ 0.05]) compared with prior time points, and with higher enrichment scores than at 8 h. In addition to hypoxia, cell proliferation and cell cycle terms that were observed at other times, DAVID analysis of genes regulated at 24 h indicated that cell motility and migration, extracellular matrix (ECM) organization, cytoskeleton organization, and protein secretion processes were more prominent. The increased alterations in genes affecting structure, movement, and the exterior of the cell suggested a new phase of cellular response at 24 h (Fig. 5e). KEGG analysis showed that at 24 h, p53 signaling, DNA replication, cell cycle, and homologous recombination were most enriched, however TNF and TGF-ß pathways were also significant (Fig. 5f). GOrilla analysis also Figure 5. GO term cluster and KEGG pathway enrichment analyses of differentially expressed genes in HLMVEC following 10 Gy X-irradiation. HLMVEC were grown to 70% confluence and exposed to 10 Gy X-ray irradiation. RNA was prepared and used for RNA-seq. www.nature.com/scientificreports/ showed regulation of cell migration and adhesion, DNA damage response, replication and repair, focal adhesion assembly, and NF-kappaB transcription factor activity. Metascape (https:// metas cape. org 46 ) was used to create an image of clustered GO terms present in the 1000 genes with the lowest q-value regardless of time point (Fig. 6). This image allows us to visualize the important gene functions that are regulated in HLMVEC in response to 10 Gy irradiation.
Focused analysis of RNA-seq after 10 Gy X-irradiation. To further analyze the DEGs, we used g:Profiler, web based software that performs functional enrichment analysis, and a literature search to identify pathways in the 3000 genes with the lowest q-value (≤ 1.18E-14). We identified four significant pathways activated following 10 Gy X-irradiation: cell cycle regulation; apoptosis; DNA damage response; and cellular senescence (Supplemental Table S2). Based on literature reports, we also found significant regulation of genes related to two additional processes: inflammation and endothelial-to-mesenchymal transition (EndMT) 50 . Focused heat maps were constructed using ClustVis that included genes regulated within each pathway to determine overall pathway activation or suppression over the experimental time course 45 . Cell cycle regulation. Radiation has a broad effect on the cell cycle in a wide variety of cell types 21 . Analysis of HLMVEC gene expression changes following radiation exposure identified 54 genes with roles in cell cycle regulation (Fig. 7a). Nine cyclin genes (CCNA1, CCNA2, CCNB1, CCNB2, CCND1, CCND2, CCNE1, CCNE2, and CCNH) displayed distinct expression patterns in response to ionizing radiation. CCNA2, that promotes transition through S/G 2 and G 2 /M, decreased 1.7-fold at 24 h. CCNB1 and CCNB2, that regulate progression through the G 2 /M transition, decreased at 2-6 h post-irradiation, but returned to near basal levels at 24 h. In contrast CCND1 and CCND2, that regulate G 1 phase transition, were elevated over full the time course, upregulated 1.9-and 4.5-fold at 24 h post-irradiation. CCNE1 and CCNE2, that regulate G 1 /S phase transition, were increased at 2 h (1.5-and 1.8-fold respectively), but then were downregulated by 24 h (2.4-and 3.9-fold respectively). Three cyclin-dependent kinases (CDKs) CDK1 and CDK2 were downregulated 1.6-and 1.5-fold at 24 h. Cyclin-dependent kinase inhibitors CDKN1A and CDKN2B, were upregulated over the entire time course, and increased 4.2-and 4.7-fold, respectively, at 24 h. CDKN1B, CDKN2C, and CDKN2D, that control G 1 progression, were downregulated over the time course, decreasing 1.6-,2.6-, and 2.5-fold, respectively, by 24 h postirradiation.
E2F genes encode transcription activators that control transcription of cell cycle genes 51 . In particular, G 1 -S phase transition requires E2F transcription factor activation. E2F1 and E2F2 transcriptional activators were Heatmap of biological pathway gene expression in HLMVEC following radiation exposure. RNAseq was used to identify gene expression changes in HLMVEC following exposure to 10 Gy X-ray irradiation. ClustVis (https:// biit. cs. ut. ee/ clust vis) was used to generate heat maps. Rows are centered; unit variance scaling is applied to rows. Rows and columns are clustered using correlation distance and average linkage. Set of 3000 genes with lowest q-value was submitted to g:Profiler (https:// biit. cs. ut. ee/ gprofi ler/) to identify pathways. www.nature.com/scientificreports/ initially upregulated (1.4-and 1.7-fold, respectively) and then downregulated (2.9-fold and 3.8-fold respectively) at 24 h. E2F8, that regulates G 1 to S phase progression, was downregulated 3.2-fold by 24 h. Interestingly, the repressive E2F family member E2F4, which suppresses cyclin D expression, was downregulated 1.4-fold, correlating with the increased expression of cyclin D gene expression. Two ORC subunit genes, ORC1 and ORC6, were downregulated 2.7-and 2.1-fold, respectively, at 24 h postirradiation. Expression of MCM complex genes 3-7 also decreased through 24 h in a similar pattern. Additionally, our further investigation of the data showed that MCM2, MCM8, and MCM10 that also participate in DNA replication, although not identified in the program analysis, were also suppressed following X-irradiation.

Regulation of apoptosis. Pathway analysis showed that both pro-and anti-apoptotic proteins involved
in the initiation of the extrinsic and intrinsic apoptotic pathways were regulated, with an early wave of apoptotic gene regulation at 2 h, and a second wave of regulation at 24 h post-irradiation (Fig. 7b). The extrinsic pathway-related death receptors TNFRSF10B, TNFRSF21 and TNFRSF4 were increased 1.7-, 3.7-, and 10.7-fold respectively, at 24 h. However, the anti-apoptotic TNFRSF10C, a decoy death receptor that does not contain a death domain, had early and sustained expression, and was increased 2-to 3.3-fold from 4 to 24 h. Genes encoding TNF receptor ligands were also upregulated inclugin TNFSF4, TNFSF9, TNFSF10, TNFSF12, TNFSF13, TNFSF15, and TNFSF18 (3.6-, 4.0-, 4.2-, 1.5-, 1.5-, 44.0-, and 2.3-fold respectively). LTB, encoding a cytokine that binds a TNF receptor 3 (TNFRSF3) was upregulated 26.5-fold. Genes encoding anti-apoptotic proteins acting to counteract extrinsic apoptotic pathway signaling were also increased, including TNF Receptor Associated Factor 1 (TRAF1), a negative regulator of apoptosis signaling by TNFR superfamily members, increased 7.2-fold by 24 h.
The Bcl2 family of proteins are pro-or anti-apoptotic, depending upon their effects on cytochrome c release from the mitochondrion. Genes encoding pro-apoptotic members of the Bcl-2 family were upregulated, including: BBC3, BCL2L11, BCL2 Interacting Protein 3 Like (BNIP3L), and BOK (5.2-, 6.9-, 3.2-and twofold, respectively) at 24 h. However, at the same time, the anti-apoptotic BCL2A1 gene was upregulated 3.5-fold. The www.nature.com/scientificreports/ transcript for BCL2L1 was upregulated 2-and 2.1-fold at 8 h and 24 h, respectively, although alternative splicing can result in either a pro-or anti-apoptotic activity for this transcript. We also observed mixed regulation of proteins required for downstream apoptotic signaling. Apoptotic protease activating factor 1 (APAF1) required for caspase activation, was increased 1.6-fold at 24 h. The gene for Jun proto-oncogene, AP-1 transcription factor subunit (JUN), also involved in caspase activation, was downregulated at 2 h, but increased 2.2-fold at 24 h. Cathepsins, lysosomal proteolytic enzymes that can directly or indirectly activate caspases, were also upregulated by radiation 52 : CTSK, CTSO, and CTSS were upregulated 3-, 2.2-, and 4.1-fold, respectively, at 24 h. Caspases 1, 3, 6, 7, and 10 were upregulated at 24 h post-irradiation (4.6-, 1.9-, 1.6-, 2.6-, and 1.5-fold, respectively). Caspase 2 was upregulated 1.3-fold at 2 h post-irradiation, but was subsequently downregulated. Our findings of mixed regulation of apoptotic genes is consistent with our finding of only very low levels of apoptosis in HLMVEC in response to radiation.
Regulation of DNA damage response. GO analysis identified 48 genes associated with DNA damage that were regulated > 1.5-fold (Fig. 7c). Cyclins and their related kinases and kinase inhibitors make up 25% of the DNA damage response genes; these were described in Cell cycle regulation section, above. The data set also includes genes encoding DNA damage response and repair enzymes.
Radiation increased mRNA levels for genes in the excision repair pathway including DNA damage binding protein 1 (DDB1), its partner DDB2, and excision repair cross-complementing rodent repair deficiency (ERCC6), all critical for DNA excision repair, upregulated 1.6-, 2.2-, and 3.3-fold, respectively, throughout the time course. Growth arrest and DNA damage inducible 45A (GADD45A), an early DNA damage response gene with multiple functions in cell cycle and survival, was increased 3.4-fold at 24 h.
Interestingly, genes encoding proteins for homologous recombination repair and long-patch excision repair were downregulated. Homologous recombination repair genes that were downregulated included radiation 51 (RAD51) and RAD52, both reduced 1.5-fold at 24 h. Breast cancer 1 (BRCA1) and BRCA2, that also contribute homologous recombination repair, were downregulated 1.5-and 1.8-fold at 24 h. RAD1 and RAD9A, involved in long-patch base excision repair, were downregulated 2.5-and 1.6-fold. Additionally, 9 members of the Fanconi anemia complementation group, involved in intra-strand crosslinking repair, were downregulated from 1.3-to 2.8-fold at later time points. These data suggest that HMLVEC utilized only specific pathways for DNA repair.

Radiation-induced inflammatory response.
Inflammation is a critical in vivo response to radiation exposure 13,53 . The normal function of the immune-inflammatory response is to prevent infection, instigate removal of dead or damaged cells, and initiate normal tissue repair 13 . However, chronic activation of inflammation post-irradiation can lead to an amplification of initial radiation damage to areas outside of the radiation field, leading to fibrosis and tissue necrosis 13,53 . We identified inflammation-related 52 genes with altered expression. Of these, 32 genes were upregulated > two-fold, and 17 were upregulated > four-fold at 24 h (Fig. 7d). These genes included transcription factors, interleukins (ILs), and chemokine receptors and ligands. Tumor necrosis factor receptors and ligands were discussed in the section on apoptosis (see above).
Transcription factors ETS1 and NFATC2 that regulate cytokine expression were upregulated 3.3-and 7.2-fold at 24 h. We identified 19 interleukin (IL) genes, targets of ETS1 and NFATC2, that were significantly regulated from 2 to 6 h post-irradiation. Interleukin 1A (IL1A), a potent pro-inflammatory cytokine, was upregulated 8.5fold at 24 h. T cell regulatory cytokine LIF was upregulated 2.6-fold at 2 h and 4.2-fold at 24 h. CSF1, macrophage colony-stimulating factor, was also upregulated 3.9-fold at 24 h. Chemoattractants for immune cells were also increased. CXCR4, that mediates migration of leukocytes, was increased 8.9-fold at 24 h. CXCL12, that mediates migration of lymphocytes and macrophages, was upregulated early, but was steadily reduced to 9.5-fold decrease at 24 h. CCL2, a monocyte, T cell, and dendritic cell attractant, was upregulated sixfold at 24 h. Senescence related gene expression induced by radiation. GO pathway analysis identified 50 genes involved in senescence (Fig. 7e). As for apoptosis where both pro-and anti-apoptotic genes were regulated, our data set for senescence also showed both pro-and anti-senescence gene regulation. Senescence pathway genes overlap with cell cycle regulation genes and DNA damage genes, previously discussed (see above sections). Interestingly, inflammation-related genes involved in senescence are some of the most upregulated in our data set; the inflammatory genes regulated by radiation were discussed in the prior section.
Insulin-like growth factor (IGF) signaling has been demonstrated to play a major role in age-related senescence as well as accelerated senescence 54 . Although immediate activation of the IGF1 receptor is associated with the induction of senescence in cultured cells following biological stress, the IGF-binding proteins (IGFBPs) that affect the half-lives of IGFs and alter their interactions with cell surface receptors are associated with both proand anti-senescent activities 24,54 . The anti-senescence protein IGFBP1 55 was downregulated 19-fold. IGFBP3, that can mediate proliferation, apoptosis, or senescence 56,57 , was also highly downregulated (> 15-fold) at 24 h. In contrast, IGFBP6, important in the progression of senescence, was upregulated 1.7-fold at 24 h. IGFBP7, that is secreted by senescent cells and can induce senescence in neighboring cells 58 , was also upregulated 1.7-fold at 24 h.

Radiation regulation of endothelial-to-mesenchymal transition. EndMT is the process in which
ECs lose their normal characteristics and exhibit a mesenchymal-like phenotype, including loss of cell-cell junctions, increased migration, and increased invasive capacity 50,59 . EndMT is not currently described as a distinct pathway in GO terms, so we conducted a literature search to identify genes associated with EndMT 60,61 . We found 64 genes in our data set that were regulated > 1.5-fold in at least one time point (Fig. 7f). Forty-four of the genes were upregulated at 24 h, 28 of them > twofold. Twenty-two genes were downregulated at 24 h, seven > twofold. www.nature.com/scientificreports/ Consistent with mesenchymal transition, EC markers were suppressed and mesenchymal (MC) markers were upregulated within 24 h post-irradiation. EC markers FLT1 was downregulated 2.8-fold at 24 h. Other markers, PECAM1, an intercellular junction adhesion molecule and CDH5, a cadherin that controls cohesion and organization of intercellular junctions, were upregulated 2.2-fold at 24 h suggesting radiation effects cell-cell adhesion. Several genes encoding cell surface proteins associated with recruitment of inflammatory cells, that were associated with EndMT in other systems, were also increased: ICAM1 and SELE (2.1-and 6.2-fold respectively). The mesenchymal cell surface markers ACTA2, VIM, CDH2 and CDH11 62 were upregulated 5.4-, 6.1-, 2.3-and 4.6-fold, respectively, at 24 h. Surprisingly, although we observed significant upregulation of mesenchymal markers, we also observed simultaneous upregulation of endothelial markers PECAM1 and CDH5 at 24 h post-irradiation (as stated above).
Cytokines that contribute to EndMT and intracellular downstream signaling molecules were upregulated following HLMVEC irradiation. BMP2 was upregulated 11.7-fold at 24 h post-irradiation. Interestingly, BMP4, shown to inhibit epithelial-to-mesenchymal transition 63

Discussion
Radiation injury to the vasculature is a key factor in delayed radiation damage to tissues. Damage to the endothelium can result in reduced blood flow and lower oxygen perfusion of the tissues 9,16,17 . Additionally, increased permeability of the microvascular endothelial barrier results in increased fluid extravasation, tissue edema, and increased inflammation of the tissues 10,15 . Thus, microvascular damage is considered to be a critical component for persistent radiation tissue damage and for the expansion of injury outside of the original field of radiation exposure. Here we examined an early time course of gene expression changes in cultured primary HLMVEC to gain insight into the pathways activated by acute radiation exposure. An overall summary of our GO analysis shows the regulation of KEGG and Wiki pathways, including cell cycle, apoptosis, DNA damage, and senescence (Fig. 7). We perused the literature for genes that have roles in EndMT and the inflammatory response, and interrogated our RNA-seq data for relevant gene regulation of those pathways. The time course of gene expression shows the coordination of these complex, yet interdependent responses to X-irradiation.
Regulation of genes regulating the cell cycle was notable within 2 h following irradiation. The majority of the genes identified in the cell cycle pathway were cyclin-related genes, origin recognition complex (ORC) subunit genes, and minichromosome maintenance complex components (MCM). Cyclins, cell division cycle genes, cyclin dependent kinases, and cyclin dependent kinase inhibitors cooperate to regulate the cell cycle progression and transitions. ORC protein subunits are required for assembly of the pre-replication complex for initiation of DNA replication during S phase of the cell cycle and serves as a platform at the origin of replication for the assembly of initiation factors and MCM proteins. MCMs control the cell cycle by regulating DNA replication 68 . The data show a complex pattern of regulation of cell cycle controlling genes, and overall suggest a downregulation of genes required for DNA replication and for S, G 2 , and M phase progression, and at the same time an increase in expression of genes encoding proteins involved in early G 1 progression.
DNA damage response for the repair of single-and double-stranded breaks is critical for cell survival following radiation exposure. Multiple pathways of DNA repair have been recognized, including non-homologous end joining, homologous recombination, base excision repair, nucleotide excision repair, mismatch repair, the Fanconi anemia pathway, which corrects DNA intra-strand crosslinks, and DNA demethylating enzymes 47,69 . DNA repair activities are closely coordinated with cell cycle regulation, which is hypothesized to allow the time for DNA repair enzyme activity and to prevent the proliferation of cells with mutations or chromosomal aberrations 70 . We observed the activation of DNA damage response protein ATM within 15-30 min post-irradiation, with other DNA damage response genes upregulated 2-24 h post-irradiation, consistent with rapid and sustained DNA damage response.
Following the interruption of the cell cycle in the presence of DNA damage, cells may undergo a variety of responses including necrosis, apoptosis, or senescence. We did not observe significant necrosis in our biological assays, and we detected only low levels of apoptosis. However, we observed the activation of a significant number of pro-and anti-apoptotic genes with early regulation of both pathways at 2 h, followed by a second wave of activation of genes in both pathways at 24 h including genes for both the intrinsic and extrinsic apoptotic pathways 71 . The extrinsic pathway is usually triggered by activation of the death receptor family, including Fas receptor, and tumor necrosis factors activated by soluble or membrane-bound ligands 72,73 . In contrast, in the intrinsic apoptotic pathway various stimuli trigger the release of cytochrome c from the mitochondria, a process that can be regulated by pro-and anti-apoptotic members of the Bcl-2 protein family. In both apoptotic pathways, caspase proteases are activated, leading to DNA cleavage and disintegration of cellular structures 74  www.nature.com/scientificreports/ gene regulation that we observed for pro-and anti-apoptotic genes in both the intrinsic and extrinsic apoptosis pathways suggests that cells may be receiving multiple conflicting signals for survival and programmed death. This is in agreement with our findings of low levels of apoptosis from 24 to 72 h post-irradiation. Ultimately, the integrated response to gene expression changes in the apoptosis pathway appears to result in HLMVEC cell survival.
Our results indicate several groups of senescence-related genes are highly regulated in response to 10 Gy X-irradiation. Normal cells undergo an ageing process that induces irreversible replicative arrest 54,75 . Genotoxic stress and other stressors can cause the cell to permanently exit the cell cycle, resulting in accelerated senescence 54,75,76 . An increasing body of evidence suggests that induction of cellular senescence by ionizing radiation may be a driver of the pathogenesis of RILI 77 . The senescent phenotype is characterized by a variety of physiological alterations including changes in the production of extracellular matrix proteins, cell-cell connectivity, resistance to apoptosis, and the secretion of pro-inflammatory cytokines and factors 78 . Senescence signaling pathways overlap with cell cycle regulation and DNA damage response, since response to genotoxic stress, such as is induced by radiation exposure, can lead to DNA damage response that induces a pause of the cell cycle as an initial step in DNA repair, and failed DNA repair can result in senescence 54 . Insulin-like growth factor (IGF) signaling has been demonstrated to play a major role in age-related senescence as well as accelerated senescence 54 , although we did not observe IGFR activation in the HLMVEC.
The use of GO and KEGG for the analysis of RNA-seq data is limited by the molecular functions, biological processes, and cellular components that have been annotated. Because EndMT is not currently identified in GO terms, we used published findings to identify EndMT associated genes. EndMT process plays a role in organ development and in adult wound healing 79 . Dysregulation of EndMT is associated with the secretion of abnormal ECM in fibrotic organ diseases and in tumor microenvironments where it is believed to contribute to cancer progression and metastasis 59,79,80 . Signaling transduction cascades that contribute to EndMT include those induced by TGF-β family members, Notch, and Wnt ligands, oxidative stress, and inflammation 79 . MMPs are important in a variety of developmental processes including mediation of cell-cell adhesion, tissue remodeling, cell migration, and proliferation 81 . Our identification of EndMT in the biological response to radiation in HLMVEC may indicate that these cells can play a role in fibrotic remodeling following radiation exposure. Further investigation of EndMT may help to determine whether this pathway can be inhibited or modulated to reduce delayed tissue injuries following radiation exposure.
Together, our data indicate a complex and integrated regulation of biological processes leading mostly to cellular senescence, a low levels of apoptosis, increased inflammation, and EndMT. Further investigation of these pathways may lead to the identification of targets for prevention, mitigation and treatment of radiation injury to normal tissues.