Ethical Approval

All applicable international, national, and/or institutional guidelines for the care and use of animals were followed. All procedures performed in studies involving animals were in accordance with the ethical standards of the institution or practice at which the studies were conducted.

Animals

Eighteen baboons were bred by the Centre National de la Recherche Scientifique (Rousset sur Arc, France) for the purpose of biomedical research. In the non-human primate facility of the French Army Biomedical Research Institute, the baboons were placed in individual cages at 21°C, with a relative humidity of 55% and a 12:12 h light-dark schedule. The animals received fresh fruit and solid food twice a day and had access to water ad libitum. The male baboons were an average age of 8.1 years (±3.3 years) and weighed 23.7 ± 5.2 kg. The experiment was approved by the French Army Animal Ethics Committee (no. 2010/12.0). All baboons were treated in compliance with the European legislation related to animal care and protection to minimize pain and suffering.

Irradiations

The animals were anesthetized with a combination of tiletamine and zolazepam (6 mg/kg intramuscularly, Zoletilt 100; Virbac, Carros, France) before irradiation. Then, the baboons were placed in restraint chairs, sitting orthogonally, front to a horizontal and homogeneous field of gamma rays delivered by a ⁶⁰Co source (IRDI 4000; Alsthom, Levallois, France) for either TBI or PBI. To attain different patterns of PBI, a 20-cm-thick lead screen was used to shield different parts of the body. Of four baboons that received TBI, two animals received 5 Gy and two animals received 2.5 Gy. Of 16 baboons that received PBI, eight different exposure patterns were simulated and two baboons were exposed per pattern, corresponding to an equivalent TBI dose of 2.5 or 5 Gy. Two dose rates were used (8 cGy/min for 5 Gy TBI and 5 Gy 50% PBI, and 32 cGy/min for all other exposure scenarios) because the ⁶⁰Co source was changed during this study. Moreover, to achieve the same homogeneous radiation field, whatever the dose rate, all baboons were irradiated at the same distance from the source. Consequently, irradiations lasted between 8 and 62 min. The midline tissue (right anterior iliac crest) dose in air was measured using an ionization chamber. Delivered doses were controlled by alumina powder thermoluminescent dosimeters placed on different cutaneous areas (thorax, thoracic and lumbar vertebrae, head, tibia, femur, femoral and head).

Blood Collection, Determination of H-ARS Severity Classifications and Other Clinical Symptoms

Changes in blood cell counts were used to determine H-ARS severity classifications (0 = unexposed; 1–4 = low-fatal) in accordance with METREPOL. The H-ARS classification was based on changes in differential blood counts taken up to 22 times over a course of 7 to 203 days postirradiation. Often, changes in lymphocyte counts or thrombocytes over time indicated H-ARS differing from each other so that intermediates between, e.g., H2–3 had to be defined. Six baboons experienced further clinical symptoms including transient (2 h) moderate vomiting (4 animals), erythema and hair loss (3 animals) and a transient weight loss of 4–10% (3 animals).

RNA Extraction and Quality Control

Whole blood samples (2 × 2.5 ml) were processed using the PAXgene Blood RNA system (BD Diagnostics, PreAnalytiX, Hombrechtikon, Switzerland). In brief, blood was drawn into a PAXgene Blood RNA tube at the French Army Biomedical Research Institute. The tube was gently inverted (10 times), stored at room temperature overnight and then stored at –20°C. After all samples were collected, the PAXgene tubes were sent to Germany for further processing. After thawing, washing and centrifugation, cells in the supernatant were lysed (proteinase K) followed by addition of lysis/binding solution taken from the mirVana Kit (Life Technologies, Darmstadt, Germany). With the mirVana kit total RNA, including small RNA species, was isolated by combining a phenol-chloroform RNA precipitation with further processing using a silica membrane. After several washing procedures DNA residuals became digested on the membrane (RNAse free DNAse Set; QIAGEN, Hilden, Germany). RNA was eluted in a collection tube and frozen at –80°C. Quality and quantity of isolated total RNA were measured spectrophotometrically (NanoDrop; PeqLab Biotechnologie, Erlangen, Germany). RNA integrity was assessed using the Agilent 2100 Bioanalyzer (Life Science Group, Penzberg, Germany) and DNA contamination was controlled by conventional PCR using an actin primer pair. We used only RNA specimens with a ratio of A260/A280 nm ≥2.0 (Nanodrop) and RNA integrity number (RIN) ≥7.5 for whole genome microarray (IMGM Laboratories, Martinsried, Germany) or RIN ≥7.3 for qRT-PCR analyses.

Screening for mRNA Species: Whole Genome Microarray

Whole genome screening for differentially expressed genes (protein coding mRNAs) was performed using the Agilent Oligo Microarray GE 83 60 K (Agilent Technologies, Waldbronn, Germany) combined with a one-color-based hybridization protocol of GeneSpring GX12 software for data analysis. We analyzed gene expression by quantile-normalized log2-transformed probe signals as an outcome. We used the nonparametric Mann-Whitney test to compare gene expression across different H-ARS classification groups with the H0 group as the reference (control). Only those gene transcripts that had a call “present” in at least 50% of RNA specimens were included in the analysis of gene expression and only genes with Mann-Whitney P values ≤0.05 and a ≥2-fold gene expression difference among compared groups were considered to represent a candidate gene for validation in phase II using qRT-PCR. Due to the explorative nature of this study and the low sample size we did not correct for multiple comparisons on screening (phase I) of the study, but considered this within our bioinformatic approach as well as the validation (phase II) of our study, where the numbers of hypotheses tested in parallel are reduced from approximately 20,000 in phase I to less than 100 genes in phase II. Gene expression data used in our study were deposited in the NCBI's Gene Expression Omnibus (GEO accession no. GSE77254;  http://1.usa.gov/1OQEw5e).

Validation of mRNA Candidate Genes via qRT-PCR

For validating the candidate genes (mRNA) from whole genome screening, only the remaining RNA samples were used. Also, the methodology was changed from microarrays to the use of a custom low-density array [(LDA), high-throughput qRT-PCR platform] and TaqMan^™ chemistry. An RNA aliquot (typically 1 µg) of RNA sample was reverse transcribed using a two-step PCR protocol (high-capacity kit). cDNA (400 µl; 1 µg RNA equivalent) was mixed with 400 µl 2 × RT-PCR Master Mix and pipetted into the eight fill ports of the LDA. Cards were centrifuged twice (1,200 rpm, 1 min, Multifuge 3S-R; Heraeust, Hanau, Germany), sealed and transferred into the 7900 qRT-PCR instrument. The qRT-PCR was run for 2 h using the qRT-PCR protocol for a 384-well LDA format. All measurements were performed in duplicate. All technical procedures for qRT-PCR were performed in accordance with standard operating procedures implemented in our laboratory in 2008 when the Bundeswehr Institute of Radiobiology became certified according to DIN EN ISO 9001/2008. All chemicals for qRT-PCR using TaqMan chemistry were provided by Life Technologies.

Prior to this study, we had established the upper limit of the linear-dynamic range of our qRT-PCR using replicate measurements on a separate sample. The upper limit occurred at threshold cycles (C_t) ≤30 in accordance to previously published work. C_t values were normalized relative to the diluted 18S rRNA (0.01 ng/10 µl) measured in an aliquot of the RNA samples using a 96-well format TaqMan qRT-PCR platform. We have found that this approach to normalization was more robust compared to the use of the internal control (GAPDH and 18S rRNA) spotted on the LDA.

Screening and Validation of miRNAs

The screening and validation of miRNAs was performed using the same commercially available miRNA LDA (type A and B). In detail, a 384-well qRT-PCR platform (LDA) was used that provided the simultaneous detection of 380 different miRNAs. Two different LDAs (type A and B) were combined so that the detection of 667 miRNA species was possible. Aliquots from each RNA sample (in general 2 µg total RNA/LDA type A/B) were reversely transcribed without preamplification over 3 h using “Megaplex pools without preamplification l for microRNA expression analysis protocol.” Using different sets of primers, two kinds of cDNAs suitable for each of both LDAs were created. In a second step, the whole-template cDNA and 450 µl 2× RT-PCR master mix were adjusted to a total volume of 900 µl by adding nuclease-free water, and aliquots of 100 µl were pipetted into each fill port of a 384-well human LDA. Cards were centrifuged twice (12,000 rpm, 1 min; Multifuge3S-R, Heraeus, Germany), sealed, transferred into the 7900 qRT-PCR instrument and a specific qRT-PCR protocol was run over 2 h using the 384-well LDA format.

We used a previously established upper limit of the linear-dynamic range of our qRT-PCR, threshold cycles (CT) <30. Normalization was performed using the median miRNA expression on each LDA separately, because this proved to be the more robust and slightly more precise method compared to a normalization approach using a housekeeping miRNA species provided on the LDA (data not shown). The median miRNA expression was subtracted from the CT value of each of the spotted genes, following the ΔCT-quantitative approach for normalization purposes. Only twofold ratios (≥2/≤0.5) were considered to represent differentially expressed genes.

All technical procedures for qRT-PCR were performed in accordance with standard operating procedures implemented in our laboratory in 2008 when the Bundeswehr Institute of Radiobiology became certified according to DIN EN ISO 9001/2008. All chemicals for qRT-PCR using TaqMan chemistry were provided by Life Technologies.

Due to the explorative nature of this study and the small sample size we did not correct for multiple comparisons on the screening step (phase I) of the study. We did multiple comparison adjustment (Bonferroni correction) in the validation step (phase II) of our study.

Bioinformatics

All genes associated with P values ≤0.05 and a ≥2-fold gene expression difference (up- or downregulated) relative to the reference underwent gene set enrichment analyses using PANTHER pathway version 10 ( http://www.pantherdb.org/). Using PANTHER, genes with similar biological function are grouped based on their annotation (reference list was the current homosapiens GO database). For these P values, we corrected for multiple testing by employing the Bonferroni algorithm or the false discovery rate (FDR).

For investigation of miRNA-mRNA interactions, a data set was downloaded encompassing 22.9 million predicted miRNA-mRNA interactions from the TargetScan 7.0 database ( http://www.targetscan.org). Visualization of these interactions and calculated cumulative weighted context++ scores, which indicate how effectively a given miRNA may target mRNAs, were visualized in Cytoscape.

For identification of potential target genes, respectively interacting proteins, the Search Tool for the Retrieval of Interacting Genes/Proteins [STRING version 11.0;  https://string-db.org/ (23)] was utilized. The interactions include direct (physical) and indirect (functional) associations and stem from computational predictions, from knowledge transfer between organisms, and from interactions aggregated from other (primary) databases. Gene enrichment analysis for gene ontologies and pathways were visualized for these potential targets. Multiple test correction with an FDR and P < 0.05 was performed.

Statistical Analysis

For statistical analysis of candidate genes, only quantitative qRT-PCR gene expression results were used. Descriptive statistics (n, mean, standard deviation, minimum and maximum) and P values for group comparisons (t test or Kruskal-Wallis test, where applicable) were calculated. Logistic regression analysis using H-ARS categories as binary outcome variable (H0 served as the reference) was run for each of the variables (genes) separately (univariate) or combined (multivariate) in some of our studies, considering parsimony and significantly associated variables only. Odds ratios, 95% confidence intervals (CI) and corresponding P values (Wald Chi-squared) were calculated. We also determined the area under a receiver operating characteristic (ROC) curve providing a reasonable indication of overall diagnostic accuracy. ROC areas of 1.0 indicate complete agreement among the predictive model and the known classifications. Based on the probability function of the ROC curves reflecting true positives (TP), true negatives (TN), false positives (FP) and false negatives (FN), we calculated positive predictive values [(PPV) = TP × 100/(TP + FP)] and negative predictive values [(NPV) = TN × 100/(TN + FN)]. PPV and NPV are best thought of as the clinical relevance of a test. PPV is the accuracy of the test among cases that test positive and indicates how seriously to take a positive result. An 85–100% range for PPV and NPV was considered a “good” test in this context. Samples identified in this range in certain of our publications were superimposed in the ROC curves for better visualization of correct predictions of the logistic regression model in predicting a certain fraction of measurements. All calculations were performed using SAS version 9.2 (Cary, NC).

Conserving High-Quality RNA

Blood samples from irradiated baboons were taken before (day 0) and at days 1, 2, 7, 28 and 75–106 postirradiation (Fig. 2). Peripheral whole blood was immediately drained into PAXgene blood RNA tubes (QIAGEN, PreAnalytiX GmbH, Hilden, Germany). At that time, alternative RNA blood tube systems existed. PAXgene tubes proved to be superior, because of the ease in use (invert filled tubes 10 times, store at room temperature over several hours and freeze afterwards) and remarkable RNA quality (RIN between 8–9) and quantity (approximately 6 µg total RNA/2.5 ml whole blood). Two PAXgene tubes were taken per time point and in the laboratories only one tube was used for RNA extraction at the same time. This robust workflow helped us to overcome unforeseen methodological problems. For example, in one instance the robot system did not perform the DNA digestion step during RNA isolation, but no error was indicated by the machine. Once recognized, this had to be repeated manually, but total RNA during an additional RNA isolation step was lost. In another instance, usually-reliable columns for RNA isolation malfunctioned in that the membrane in the column moved unrecognized from a horizontal into a vertical position. This resulted in total RNA absent of any mRNA and small RNA species, as we found later when performing the microarray and qRT-PCR analysis.

Two-Phase Study Design: Considerations and Rationale

Screening approach (phase I). Newly available technologies such as microarrays, RNAseq or small RNAseq next generation sequencing (NGS) enable examination of radiation-induced whole genome gene expression changes at the same time. We utilized them for screening purposes (phase I), but always found it important to validate them (phase II) by using another technique to avoid false positive results and prove that the method worked (Fig. 3A). Quantitative real-time PCR (qRT-PCR) using expensive chemistry like TaqMan primer-probe assays (compared to SYBR^™ Green) is still considered the gold standard for quantifying gene expression and, based to our experiences, preferable. Over time, we learned that validation of microarray as well as NGS results demanded an extensive search for the appropriate TaqMan assay. Usually, we selected them with a preference for an exon-exon junction primer probe design (designed to span an exon-intron boundary) and thus, avoided amplification of residual DNA copy numbers. Selecting “inventoried” assays guaranteed the applicability of the ΔΔC_t approach for quantifying gene expression differences and no extensive additional experiments to prove that were required. The “best coverage” inventoried assays target exon regions, representative for high baseline gene expression gene regions observed in normal tissues. However, recent NGS validation studies indicated that targeting those exon regions, where radiation-induced changes occurred (we called them “radiation-responsive exons”), must also be considered to avoid false positive or false negative results. Those radiation-induced “radiation-responsive exons” are not necessarily regions with highest baseline gene expression (Fig. 3B). Impressive correlations of gene expression measurements among both technologies were generated and can be further improved when considering “radiation-responsive exons” (Fig. 3B).

FIG. 3

Methodological considerations are shown with respect to the two-phase study design. Using microarrays, RNA from irradiated and nonirradiated baboon samples was hybridized and gene fluorescence intensity per gene (dot) compared (panel A). A significant fold change as shown was used as a filter to identify candidate genes, which were forwarded for validation using qRT-PCR. Using next generation sequencing (NGS) in the context of another project, we identified radiation-induced “exon driver regions” and exon regions predominantly contributing to the baseline gene expression level (panel B, left side). Panel B (right side) reflects the agreement of differential gene expression measured in microarrays (baboons) and qRT-PCR. The low-density arrays (LDA) were loaded with three technical replicates (indicated by different symbols and colors) and expected overlap of replicate measurements is shown between 15–30 raw C_t values (threshold cycle), but not afterwards where technical replicate measurements vary over a wide range of raw C_t values (panel C). This procedure allowed defining of the linear-dynamic range as well as the range where measurements are reproducible and convert into quantitative values.

Validation approach (phase II). In study scenarios where dozens or several hundred candidate genes and larger sample numbers are analyzed, formats other than the use of 96-well plates must be utilized. So-called customized low-density arrays (LDA, Thermo Fisher Scientific^™ Inc., Waltham, MA; 384 qRT-PCRs are run simultaneously, Fig. 3C) or recent developments such as the 12k open array (OA; Thermo Fisher Scientific) format (12,000 qRT-PCRs are run simultaneously) can be applied. In our hands, LDA proved to be reliable, but 12k OA must be run in triplicate due to technical limitations inherent to this technology.

Normalization and control of qRT-PCR data is essential. Using similar amounts of RNA from each sample (typically 1 µg total RNA) for cDNA synthesis and per qRT-PCR reaction (typically 30 ng cDNA) reduced methodological-driven gene expression changes. When using, e.g., high-abundance 18S rRNA for normalization purposes, the expected C_t values (threshold cycles) showed similar values provided the same amount of RNA was introduced per qRT-PCR reaction. For instance, when using 30 ng cDNA per qRT-PCR reaction (and holding other qRT-PCR parameter constant, e.g., the baseline), then the 18S rRNA C_t values always ranged between 18–20. This strategy was used to control our qRT-PCR and for normalization purposes. Of note, 18S rRNA had to be diluted accordingly and cDNA synthesis had to be performed with random hexamers (the structural 18S rRNA includes no poly-A-tail). This caused additional burden, but it gained control over the validation using qRT-PCR. In the LDAs, primer and probes for each gene are filled and lyophilized in each well, thus ensuring single qRT-PCR reactions and identification of up to 384 different genes per LDA. A specific primer probe design filled into LDA wells identifies 18S rRNA which can be chosen for normalization purposes. However, C_t values due to 18S rRNA high abundance occur early (approximately 6–7 raw C_t values). They lie outside the linear-dynamic range of the methodology (typically ranging between 15–30 raw C_t values, Fig. 3C) and render these LDA data uninformative. Indeed, 18S rRNA has to be examined in a separate step from the LDA if chosen for normalization purposes. Alternatively, we used the median expression on each LDA for normalization purposes or an additional internal control spotted on the LDA (like GAPDH).

To confirm the reproducibility of the LDAs and check the linear-dynamic range before starting the experiment, we ran one exemplary sample in triplicate (Fig. 3C). As expected, measurements of genes within the linear-dynamic range of the technique produced similar results (stacked over one another), but C_t values of replicate measurements outside of it (around C_t = 30) were disparate (Fig. 3C).

We found that with qRT-PCR, approximately one third to two thirds of assays did not show amplification plots. One reason may have been that the gene copy number was too low. This issue can be solved to a certain extent by loading more cDNA, but two times more cDNA translates in a shift of one C_t value only. We noted this problem when using custom LDAs or commercially available miRNA LDA type A/B covering 667 distinct miRNA species. It is an unavoidable issue and must be considered when selecting an appropriate number of candidate genes. A better selection of appropriate candidate genes considering radiation-responsive exon might lead to further improvements of detectable genes.

Sample distribution on screening (phase I) and validation phase (II). Screening in phase I with microarrays (or NGS) using up to five samples per group and independent validation in phase II using the remaining samples (Fig. 3A) proved an efficient strategy for several reasons.

Retrospective Dosimetry: A Challenging Approach

In all of our baboon studies we did not perform retrospective dosimetry. However, in this section we address the challenges related to this approach.

Inter-individual differences in radiosensitivity related to survival after high-dose radiation exposures are important to consider and understand. For instance, with a whole-body exposure of 3–4 Gy, approximately 50% of humans will die within 60 days without treatment [LD_50/60 (7)]. Knowing the dose provides a 50% chance to predict the clinical outcome, which is no better than tossing a coin.

The relationship of radiation exposure with ARS depends on many factors including radiation exposure characteristics (e.g., radiation quality, fractionation, dose rate, partial/total-body irradiation) and inherent biological processes (e.g., radiosensitivity, cell cycle dependency, oxygenation, Fig. 4) (8) as well as other aspects such as age at exposure or pre-exposure health conditions.

FIG. 4

Box plots representing dose estimates using the dicentric chromosomal assay performed on radiation therapy patients originating from a medical database are plotted against the H-ARS severity which was generated based on changes in blood cell counts stored at a database comprising clinical data of irradiated, radiation-accident patients following the METREPOL (MEdical TREatment ProtocOLs) approach (panel A). Radiation exposure characteristics must be considered when estimating health effects such as the H-ARS (panel B). It is hypothesized that biological changes, such as clinical signs and symptoms, gene and protein expression changes, have the potential to integrate those radiation exposure characteristics (as well as biological aspects, e.g., radiosensitivity, which is not shown here) and provide a simplified approach for effect prediction, but not dose estimation.

That does not render absorbed dose uninformative, but reflects the complexity of the concept when using exposure as a predictor of clinical effects. It is a reminder that dose estimation is not the primary aim. Dose estimates are a “surrogate” for risk estimation (e.g., late effects) or acute health effects. The relationship of dose with acute effects requires further delineation. Recently, our group examined the relationship of dose to different H-ARS severity degrees using an archive comprising radiation accident case histories (9). All individuals received whole-body exposures (single exposures predominated) and individuals were drug-treated, which challenges a comparison with the LD_50/60 (not considering treatment) mentioned above. Exposures below 1 Gy and doses above 5 Gy roughly corresponded with no- or a low-grade H-ARS and severe grade H-ARS, respectively, and were consistent with medical expectations (Fig. 4A). However, whole-body doses between 1–5 Gy corresponded poorly to different H-ARS degrees of severity, making an individual recommendation based only on dose essentially impossible [Fig. 4A (9, 10)]. We interpreted this dose range of 1–5 Gy as an exposure where inter-individual radiosensitivity drives the clinical outcome. At lower or higher exposures, all humans respond similarly. Knowing the dose in these situations allows for a good estimate for acute health effects.

We recently introduced the concept of radiation-related bioindicators for ARS (effect) prediction [Fig. 4B, (10, 11)]. In our model we considered bioindicators to be factors that integrate multiple radiation exposure characteristics as well as cell- and molecular-based processes to improve clinical prediction in persons with, e.g., ARS (Fig. 4B). For instance, altered peripheral blood cell counts are caused by cell death in irradiated radiosensitive bone marrow stem cells. If partial-, and not whole-body irradiation occurred with either neutrons or gamma rays, bone marrow stem cell death would be reduced. The number of surviving bone marrow stem cells would directly affect peripheral blood cell counts and the subsequent events of H-ARS severity, thus integrating differences in radiation exposure pattern. Likewise, biological factors, such as cell type or individual radiosensitivity, would be included (integrated) by measuring bioindicators of effect. Different groups have already shown that changes in blood cell counts (10, 12, 13) and blood protein markers (4, 5) can be used for effect prediction. In a similar fashion, gene expression levels can be exploited for the same purpose, as explained below.

Predicting H-ARS Based on Radiation-Induced Gene Expression Changes

ARS categorization using METREPOL (MEdical TREatment ProtocOLs). As outlined above, it was not our purview to undertake retrospective dosimetry on 18 baboons. Instead, our goal was to predict acute clinical effects, namely H-ARS.

Approaches like METREPOL are used to categorize H-ARS into four severity degrees based on blood cell count changes in the weeks that follow radiation exposure. These are: no H-ARS (H0, BCC in the normal range referenced in METREPOL), low (H1), medium (H2), severe (H3) and fatal (H4) H-ARS (14). H-ARS severity degrees 2–3 are associated with recovery, while no recovery is expected in H4. These severity degrees are associated with treatment decisions outlined in METREPOL. Thus, predicting H-ARS severity will provide important insight for clinical diagnosis and treatment decisions. In our experience, categorization of H-ARS degrees after irradiation was often challenging, because, e.g., lymphocyte counts over time (e.g., H-ARS 1) indicated a different H-ARS severity compared to thrombocytes (e.g., H-ARS 2) or granulocytes (e.g., H-ARS 2, Fig. 5). In these cases, we aggregated our findings to an intermediate severity category, namely H-ARS 1–2. Since H-ARS develops with a delay over time, it is categorized based on the entire clinical follow-up of BCC spanning about 60 days. Some groups intended to examine correlations of biomarkers with H-ARS categories obtained at the same early time points after irradiation. Since H-ARS develops with a delay over time, it is categorized based on the entire clinical follow-up of BCC spanning about 60 days. Correlating biomarkers with H-ARS categories obtained at the same early time points after irradiation is an approach without clinical significance, since the later-developing H-ARS severity degree is not predicted (personal communications). This approach describes an actual status of the H-ARS, for early time points after irradiation, but the delayed manifestation of the H-ARS at later time periods is not addressed with this approach. Therefore, from the clinical perspective, a H-ARS prediction of the delayed H-ARS based on early biological changes appears preferable.

FIG. 5

Changes in blood cell counts over time are shown. They correspond to METREPOL (MEdical TREatment ProtocOLs) H-ARS (hematological acute radiation syndrome) severity categories (H1–3) as outlined. The highest H-ARS 4 severity category according to SEARCH develops at lethal doses and not at sublethal exposures applied to the baboons in our studies.

Rather than relying solely on the less common METREPOL ARS score, we also examined radiation-induced gene expression changes that predicted pancytopenia. Pancytopenia is a well-known clinical end point and does not require further introductions among clinicians. However, examinations in baboons for radiation-induced mRNA or miRNA species (most promising was miR-574-3p) associated with pancytopenia did not result in genes any better than those related to different H-ARS scores and are, therefore, not further described here (15).

First generation gene expression signatures. With the decrease in granulocytes and thrombocytes in the peripheral blood over time, H-ARS is characterized mainly by immune suppression and hemorrhage. We hypothesized that the depletion of BCC would precede changes in gene expression causally or timely related to their later decline and, therefore, could serve as an indication of the eventual H-ARS severity score.

Using preirradiation day 0 (reference, H0) and day 1 and 2 postirradiation baboon blood samples, we looked for early radiation-induced gene expression changes and tried to predict the final H-ARS severity degree based on the BCC changes over the entire follow-up time (Fig. 5). The first three postirradiation days are considered to represent the diagnostic timeframe when clinicians require guidance for medical management decisions such as continuing diagnostics, hospitalization and whether to begin intensive treatment regimens. Given the low number of animals and different clinical follow-ups, we identified three H-ARS severity groups, namely H0 (blood samples before exposure also used as the reference), H-ARS 1–3 and H-ARS 2–3 severity. H-ARS 4 severity (observed at lethal doses) was not expected given the sublethal exposures (Fig. 5). With H-ARS 1–3, we identified irradiated animals requiring different degrees of medical support and this in turn helped to identify H-ARS 0 (nonirradiated) animals, not requiring medical support. Identifying H-ARS 0 is important in situations where individuals believe they have been exposed, but were not. H-ARS 2–3 disease requires hospitalization and intensive medical support. After confirmation with qRT-PCR, 22 genes were associated with H-ARS 1–3 and seven genes with H-ARS 2–3 (16). Each of these 29 radiation-induced genes detected in the peripheral blood of baboons within the first two days after irradiation could independently predict H-ARS. Most notable were genes such as POU2AF1 and WNT3 that were constantly (over two days) about 30–40-fold downregulated and indicated a more severe H-ARS 2–3. These two genes were not seen in the nonirradiated (H-ARS 0) animal samples. We used the criteria of fold change, FC ≥ |2|, P ≤ 0.05 as a cut-off value to define candidates that were constantly up- or downregulated over both days to establish a simple and robust diagnostic tool. Because transcriptional changes are supposed to change constantly over time, with a whole genome screening approach we reasonably expected to identify genes not appreciated for biodosimetric purposes or H-ARS effect prediction. To the best of our knowledge, at that time we were the first to examine the applicability of gene expression changes for H-ARS effect prediction. Therefore, we were not surprised to find that both genes had not been identified as radioresponsive at that time. We also examined FDXR, a well-known and established radioresponsive gene in the field of biodosimetry (17, 18). Strangely, FDXR in irradiated baboons appeared constantly downregulated as opposed to human in vitro peripheral blood models (19). At this time of our studies we questioned if the in vivo baboon results reflect gene expression changes also found in irradiated humans (in vivo) and we therefore had to expand our research to include leukemia patients who had received a whole-body radiation exposure (results discussed in the next section).

Validation of baboon data and further studies. Additionally, we wondered about the origin of altered gene expression. Where did it come from? Was it related to irradiated peripheral BCC? Did it reflect a response of irradiated body areas releasing mRNA and miRNA species into the blood similar to the rationale of “liquid biopsy” where tumors release their molecular “fingerprint” in the peripheral blood (20, 21)? Thus, in the next steps we aimed to: 1. Compare baboon and human in vivo data for interspecies validation purposes; 2. Perform corresponding baboon and human in vitro radiation experiments to better understand where radiation-induced altered gene expression originates; and 3. Construct dose-response curves to reconfirm the hypothesized radiation responsiveness of our genes and to support the association of the gene's altered expression with H-ARS categories. For instance, for WNT3 marginal changes in gene expression at low doses and pronounced changes at higher doses would be expected to be corresponding to a predicted H-ARS 2–3 severity degree occurring usually at high-dose and high-dose-rate exposures. Due to limited access to blood of patients with severe diseases like leukemia, we only examined nine genes (19). We chose six promising baboon candidate genes (WNT3, POU2AF1, CCR7, ARG2, CD177, WLS), as well as three genes commonly used in ex vivo whole blood experiments (FDXR, PCNA, DDB2). To cover a wide dose range, we used 24-h blood samples from patients receiving diagnostic computed tomography (0.004–0.018 Gy), radiotherapy for prostate cancer (0.25–0.3 Gy) or radiotherapy for leukemia (TBI, 3–4 Gy). Additionally, peripheral whole blood samples of baboons and healthy human donors were cultivated ex vivo over 24 h and irradiated at 0–4 Gy. Results can be summarized as follows. 1. The six baboon candidate genes were confirmed in leukemia patients. 2. FDXR was downregulated in baboons and upregulated in humans after radiation exposure. FDXR measurements underscore the importance of independent assessments even when candidates from animal models have striking gene sequence homology to humans. 3. WNT3 and POU2AF1 gene expression changes in both species were similar after in vivo WBI (peripheral whole blood as well as the whole body were irradiated) and in vitro experiments irradiating the peripheral blood only (19). This indicates that peripheral blood cells represent the radiation-responsive targets. CCR7, ARG2, CD177 and WLS showed different patterns between these in vivo and in vitro irradiation experiments and, therefore, appeared to be altered due to radiation-responsive targets other than the whole blood cells (19). 4. Linear dose-response relationships of FDXR, WNT3 and POU2AF1 over a wide dose range (0.001–5 Gy for POU2AF1) underscore that these genes are associated with radiation exposure. This, however, also indicates diagnostic limitations, because only pronounced (10-fold and higher) downregulation of, e.g., WNT3 or POU2AF1 will correspond to a severe H-ARS. Other, already established biodosimetry genes, namely FDXR and DDB2, proved to be strongly radioresponsive, indicating that an exposure occurred. In the absence of a FDXR or DDB2 response, the implication may be as a means to identify unexposed individuals.

Narrowed gene set applicable for H-ARS diagnosis. These intermediate results narrowed the potential set to four genes (FDXR, DDB2, WNT3 and POU2AF1) to be applicable as a diagnostic tool for H-ARS prediction. However, questions remained. What are baseline values for these promising validated candidate genes and could inter-individual variation interfere with the diagnostic applicability in humans? We pursued these questions by analyzing peripheral blood samples of 200 healthy male and female donors (22). Interestingly, for these genes, gender and age contributed marginally (less than twofold) to the inter-individual variance expressed as the fold change calculated between the lowest (reference) and the highest gene expression measurements (C_t values). FCs ranging between 10 and 17 were observed for most of the ten genes examined, but fold change for WNT3, POU2AF1, FDXR and DDB2 were 37.1, 17.0, 9.9 and 9.4, respectively. When plotting, e.g., WNT3 or POU2AF1 normalized gene expression values of 200 healthy nonirradiated humans versus normalized gene expression values of diseased irradiated baboons, the complete separation of diseased and healthy groups remained. Assuming a similar response in humans and baboons, we concluded that WNT3 and POU2AF1 baseline variance measured in humans did not impact the complete separation of H-ARS categories measured in baboons (22).

The diagnostic window and origin of gene expression changes. In our most recently published work, we elucidated the diagnostic window (onset and duration) and the origin of these four genes on healthy human peripheral blood samples irradiated in vitro: How early after irradiation can we use these genes for diagnostic purposes and which cell populations in the peripheral blood contribute to the results? We performed these experiments in vitro with peripheral blood irradiated at 0.5 Gy (related to a H-ARS 1 severity), 2 Gy or 4 Gy [related to a H-ARS 2–4 degree severity, (23)]. The earliest gene expression changes were detected 2–4 h after irradiation (FDXR), followed by 4 h (DDB2 and POU2AF1) and 8 h (WNT3) indicating that the diagnostic window onset for prediction of the subsequent H-ARS requiring intervention may be available as early as 2–8 h after exposure and might be predictive up to three days after exposure (23). Further examinations revealed T-lymphocytes as the predominant cell population for radiation-induced upregulation of FDXR/DDB2 (74.8%/80.5%) and B-lymphocytes for downregulated POU2AF1/WNT3 [97.1%/83.8%, (24)]. Regarding the onset and length of the diagnostic window, we are fairly certain of our in vitro data. But due to methodological difficulties the in vitro model becomes unstable over time, for which qRT-PCR can compensate only to a certain degree at three days after irradiation. We suggest that the 3-day postirradiation end of the diagnostic window must be validated in another model. To that end we are currently analyzing irradiated rhesus macaque NHP for whom blood samples collected 3 days after irradiation are available.

Regarding the origin of radiation-induced gene expression changes, it is interesting to have identified immune competent T- and B-lymphocytes as the predominant radiation cell targets. Differences in radiosensitivity of peripheral blood cell populations and decreased blood cell counts are known (25–27), however, in our in vitro experiments radiation-induced changes in T- and B-lymphocyte blood cell counts were marginal and mostly not significant (24). This may be another limitation of our in vitro model requiring further effort for improvement, e.g., separate cultivation of T or B cells when proceeding with mechanistic work on WNT-catenin signaling and cell death in B cells. Nevertheless, in our examinations, transcriptomic changes preceded cell death processes that occurred later after irradiation and led to onset of H-ARS.

H-ARS prediction, features. A diagnostic tool for H-ARS prediction should address features like early and high-throughput diagnosis. To address the high-throughput applicability of gene expression measurements for early H-ARS diagnostics, we recently employed targeted NGS in vitro on 1,000 blood samples from healthy human donors (28). All 1,000 samples were processed within 30 h and classification of H-ARS severity categories was based on FDXR, DDB2, WNT3 and POU2AF1. Radiation-induced gene expression changes and classification of H-ARS severity categories had an overall agreement ranging between 90–97% (28). Another group identified a gene set including DDB2 for dose estimation, considering medical interventions for treatment of comorbidities in patients, indicating their applicability even in diseased individuals (29). However, much more research in this regard is required.

The diagnostic tool. A combination of four genes allows an early and high-throughput diagnosis of H-ARS severity including the identification of unexposed individuals (Fig. 6). With this tool, there are three goals:

FIG. 6

Description of the diagnostic tool based on radiation-induced gene expression changes and predicted clinical outcomes of the hematological acute radiation syndrome (H-ARS severity degree). Gene annotations and detection assays using qRT-PCR are summarized in the table below.

Merging existing H-ARS categories in those three goals appears appropriate since it addresses urgent clinical questions (prioritizing hospitalization as well as restricted medical resources) and considers low number of irradiated baboons. For the tool, we introduced pairwise redundant genes (FDXR/DDB2 and WNT3/POU2AF1) showing the same association with H-ARS severity degrees in order to increase the robustness of our diagnostic tool.

For identification of H-ARS 0 (unexposed individuals) we proposed an algorithm where all four radiation-induced gene expression changes relative to unexposed do not exceed a fold change of 2. This is the fold change considered to adjust for methodological variance.

For identification of H-ARS 1, we proposed an algorithm characterized by ≥2-fold upregulation of FDXR and DDB2, while no gene expression changes of WNT3 and POU2AF1 below the methodological variance (FC 2 downregulated) are expected.

For identification of H-ARS 2–4, we proposed an algorithm characterized by >2-fold upregulation of FDXR and DDB2 and a pronounced downregulation of WNT3 and POU2AF1. Based on our currently available data, we expect a 10-fold downregulation of WNT3 or POU2AF1 to predict H-ARS 2–4 with a positive predictive value (PPV) of approximately 100% (22). A >2 and <10-fold downregulation of WNT3/POU2AF1 according to our current knowledge would reduce the PPV from 100% to approximately 90%, but these data are based on 17 irradiated baboons and 200 nonirradiated human healthy donor samples only. Additional data will be generated soon from 63 irradiated rhesus macaques that we are currently investigating.

With a low number of genes (FDXR/DDB2 and WNT3/POU2AF1), clinically-relevant decisions about hospitalization and based on H-ARS severity category could be made without complex gene signatures of tens or hundreds of genes required (30, 31). While our strategy does not rule out the use of signatures, it suggests that a parsimonious approach appears to provide robust predictions. Using downregulated genes (WNT3/POU2AF1) for diagnostic purposes adds challenges, but the precision gained in H-ARS prediction outweighed the disadvantage.

These data provide the matrix for establishing the diagnostic tool in other laboratories worldwide. But each laboratory must customize it in their way with regard to, for example, the amount of input cDNA and normalization strategies. This is why we are not providing meaningless raw C_t values, which apply to our setting, because settings will differ in other laboratories.

Considerations and next steps in NHP and human models. Although validated in human patient models and using extensive in vitro analysis, we acknowledge the low baboon and human sample numbers examined so far. Additional analysis for confirmation is required. We have been entrusted with pre- and postirradiation rhesus macaque blood samples (mentioned above) for which NGS analyses are underway in 63 animals.

Upcoming analyses include blood samples from 92 chemotherapy-related female breast cancer patients. Chemotherapy and radiotherapy are both targeting cell proliferation and induce DNA damage (32, 33). We will assess H-ARS severity (based on BCC changes according to METREPOL) and examine our candidate genes. After these examinations in two different human patient groups (leukemia and breast cancer patients) and two NHP species (baboons and rhesus macaque), a final judgment of our diagnostic tool can be made.

Predicting H-ARS is an important step in diagnosis of the ARS. But other organ systems are also affected by ARS such as in the gastrointestinal or the neurovascular system. Those are not covered by our diagnostic tool. However, the bone marrow is considered the most radiosensitive and when identifying individuals who will develop H-ARS, we automatically hospitalize individuals who may also present with the gastrointestinal or neurovascular syndrome and who will require medical support.

As outlined above, changes in gene expression used for diagnostic purposes of H-ARS could be a timely predictor or causally related to BCC decline. The canonical WNT signal transduction pathway via catenin stabilization and consecutive upregulation of genes such as c-myc are already known to impact proliferation (apoptosis) of gastrointestinal crypt stem cells or pre-B-lymphocytes (34–36). Whether this can be observed in matured B cells is unclear and might provide a new therapeutic target to fight H-ARS immune deficiency after radiation exposure.

It is unlikely that these genes are expressed only after ionizing radiation. However, in a radiological or nuclear scenario they can certainly be considered and further validated by other measures such as dicentric chromosomal analysis, known to increase almost exclusively after ionizing radiation exposure (37, 38).

For future developments, we are currently working on high-throughput qRT-PCR platforms, e.g., automatic RNA isolation combined with 12,000 qRT-PCR measurements performed simultaneously within 2 h. These complex workflows will be established in specialized laboratories. In parallel, we are also working on point-of-care diagnostics, where the workflow from the blood drop to RNA-isolation, cDNA synthesis and qRT-PCR will be performed automatically on a microfluidic card in hospitals, where potentially exposed individuals will arrive after a radiologic or nuclear event.

When considering the COVID-19 pandemic one lesson we should learn is to be prepared for unexpected events. This applies to radiological or nuclear events as well.

Rationale

Several years ago, we examined modes of cell death such as apoptosis, micronucleation, mitotic catastrophe or necrosis in ten cell lines of different origins (39). Interestingly, the cell death mode shifted with increasing dose but mostly depended on the cell type reflected by its transcriptomic activity. For instance, radiation exposure in mouse fibroblasts (L929 cells) almost exclusively increased the micronucleus frequency only, while apoptosis frequency predominated in cells of blood origin (e.g., human leukemia, HL-60 cell lines). We assumed there was an “intrinsic” program dictating the mode of cell death in different cell lines. Intrinsic radiosensitivity is a well-accepted mechanism explaining clonogenic survival (or death as the mechanism for ARS development) and is reflected by the D_o term of the D_o/D_q model (7). It is furthermore well established that cells irradiated in different cell cycle phases differ in their radiosensitivity (survival). Its importance is reflected by its contribution to the 4Rs (repair, repopulation, redistribution, reoxygenation) in radiobiology (7). Extrapolating from these well-accepted concepts and findings, where differences in irradiated cell activity might dictate their fate after irradiation, we assumed that the pre-exposure transcriptomic status somehow reflects the intrinsic (D_o) or DNA repair (D_q)-related radiosensitivity of irradiated individuals.

Radiosensitivity Predicted in Baboons and Further Studies in Rhesus Macaques

We examined preirradiation gene expression (mRNA and miRNA species) in baboons (n = 17) that developed H-ARS without pancytopenia (n = 12) or a more aggravated H-ARS with pancytopenia (n = 5) after irradiation (40). Interestingly, two of the five baboons that developed the more aggravated H-ARS (with pancytopenia) had received only one half the dose (2.5 Gy) of that received by most of the animals (5 Gy) that subsequently developed a less aggravated H-ARS. Thus, it appears that dose may not be a major predictor of H-ARS, with or without pancytopenia, and that other factors related to inherent radiosensitivity might play a role.

Measurements on mRNA species. From microarray analysis in our baboon (male) study, only 20 mRNA species met the selection criteria (≥2-fold difference among groups, P ≤ 0.05) and were forwarded for validation. Only 12 mRNAs were able to be examined using qRT-PCR and four mRNAs (DCAF12, HBE1, SLCO4C1, TSPO2) were successfully validated (40), since fold change showed the same direction (2–5-fold different between both groups), but P values were only borderline significant. Of note, in this study, radiosensitivity was referring to the prediction of a more aggravated H-ARS. We are currently conducting another study on a total of 142 rhesus macaque NHP (male, n=72; female, n=70) to examine radiosensitivity in preirradiation samples. However, here, our clinical end point is survival status due to lethal radiation exposure (LD_66/60). In an additional analysis, we examined the 20 baboon mRNA candidates (microarray data) in the male rhesus macaques (NGS data). We found SECTM1 (baboon FC 0.2 vs. rhesus macaques FC 0.4) and PATE2 (FC 2.7 vs. 3.6) showed a similar fold change (≥|2|-fold difference among groups, P ≤ 0.05) in both species. The other candidate genes from the baboon study were not differentially (≥|2|-fold difference) expressed in the rhesus macaque study. Nevertheless, another five genes showed a tendency of gene expression in the same direction (PQLC3 baboon FC 0.3 vs. rhesus macaques FC 0.7; PLSCR1 0.4 vs. 0.6, SLCO4C1 0.4 vs. 0.6, XAF1 0.4 vs. 0.6; TSPO2 2.6 vs. 1.5).

In the rhesus macaque study in males, in particular SLC22A4, EPX, IGF2BP1 appeared promising, but in the baboon microarray data, SLC22A4 was not differentially expressed and EPX and IGF2BP1 appeared inversely regulated when comparing the two interspecies data sets (upregulated in baboons and downregulated in rhesus macaques) after irradiation (41).

Measurements on miRNA species. We examined 667 miRNA species in baboons. Six of nine selected candidate miRNAs remained significantly deregulated during validation. In particular, miR-425-5p showed nearly complete discrimination between H-ARS groups with and without pancytopenia. This miRNA has been previously associated with chemo- and radiosensitivity (42, 43), cell proliferation (44, 45) and cell death (42). Target gene searches of miR-425-5p identified known and new potential mRNAs and associated biological processes linked with radiosensitivity (40). We further detected transport-related genes (heterogeneous nuclear ribonucleoprotein (hnRNP)-family, SYNCRIP) as well as mRNAs linked to the Wnt-β-catenin signaling pathway. Interestingly, in previous work we had identified HNRNPA1, belonging to the hnRNP-family, which appeared to be involved in atherosclerotic processes in Mayak workers with internal plutonium exposure (46). WNT3 was identified as a robust bioindicator that predicted more severe H-ARS as described above.

In an additional analysis, we examined the nine selected candidate miRNA species from the screening phase of the baboon model (microarray data) in the male rhesus macaques (NGS data). These examinations revealed only mir-192, which showed a downregulation in baboons with fold change of 0.4 as well as in male rhesus macaques with fold change of 0.5.

Interpretation

In the male baboon study, few validated candidate mRNAs were found, other than in the larger rhesus macaque study (41). This may reflect the 4-times-lower sample size (17 preirradiation peripheral blood samples) in baboons versus 72 male rhesus macaques. However, even the validated genes in rhesus macaques sometimes showed similar gene expression values by survival status, rendering discrimination of clinical outcome difficult and suggesting that larger sample sizes may be needed. For instance, in the rhesus macaque study, the most promising gene (MBOAT4 in female macaques) discriminated survival status among 74% (29/39) of NHP samples correctly (NPV and PPV of 85–100%). No prediction could be made for the remaining 26% of irradiated female NHP. That differs with respect to genes such as WNT3 or POU2AF1 where complete separation was possible in prediction groups of different H-ARS categories in baboons (see above). Clearly, predicting radiosensitivity in nonirradiated animals is more challenging than identifying genes altered during the radiation-induced biological process giving rise to severe acute health effects, namely H-ARS.

Our recently published study of 142 rhesus macaques suggested that different gene sets were involved in radiosensitivity prediction for males and females (41), but the studies had many other differences besides species and end points (H-ARS severity vs. survival status). Comparisons of data sets of baboons and rhesus macaques previously identified a low number of mRNA and miRNA species similarly regulated, so common ground was limited. In future examinations on a second rhesus macaque group (n = 63) we attempt to replicate of our current results.

It is difficult to conduct similar studies in humans. By taking advantage of preirradiation blood samples from 92 chemotherapy treated female breast cancer patients (mentioned above), we aim to examine chemosensitivity in this cohort to compare it with our radiosensitivity results. We hypothesize that a common mechanism of cell death of both treatments might share similar features regarding chemo- and radiosensitivity.

Rationale

For meaningful prediction of acute radiation health effects, several characteristics related to radiation exposure must be considered, including exposure pattern of the body [Fig. 2 (8, 11)]. For instance, during the Georgian accident, local skin exposures of up to 150 Gy led to acute injuries (ulcerations), but individuals survived (47). Such would not be the case with a whole-body dose of the same magnitude.

Dictated by low baboon sample numbers (4–5 animals per exposure group), we restricted the exposure pattern to three body areas, namely upper-body (legs shielded), left hemibody or total-body irradiation [Fig. 2 (48)]. The equivalent whole-body dose was either 2.5 or 5 Gy for each exposure. Differentially expressed gene sets relative to the preirradiation sample were constructed per exposure pattern group and compared with each other at days 1, 2, 7, 28 and 75–106 postirradiation (Fig. 2). As a proof-of-concept, we expected to identify differentially expressed genes associated with the exposed body area. Along similar lines, we also aimed to assess gene expression changes corresponding to the percentage of exposed body area (TBI is equivalent to 100% exposed body area). This allowed simultaneous evaluation of all 17 samples with linear regression analysis and an increase in statistical power over categories using non-parametric tests. We examined body exposure patterns associated with mRNA species of the whole genome (transcriptional level) and 667 miRNA species (post-transcriptional level) and validated our findings using qRT-PCR.

Examinations in Baboons

Results can be summarized as follows:

1. Over all time points, we identified 55 miRNAs and approximately 20 times more mRNA species associated with exposed body areas (48, 49). Robust fold changes (5–10 times for 22 miRNAs and 10–40 times for the 14 most promising mRNAs) and either complete or almost complete separation of exposure groups was observed.

2. However, for miRNAs, a different gene set was required for each time point and for mRNA species, a few genes (0–10, depending on the comparison over time) could be used for discrimination of exposure groups over a larger time range, namely 2–75 days postirradiation. Bioinformatic analysis confirmed a shift of biological processes occurring over time after irradiation with almost no genes in common over the time points examined. Furthermore, a wave-like pattern in the number of differentially expressed mRNAs or miRNAs over time in all group comparisons was seen. In the comparison of TBI versus PBI, similar features emerged: While most differentially expressed miRNAs were exclusively observed on day 7, a corresponding peak in the number of deregulated mRNAs was found on day 28 for mRNAs (Fig. 7). This pattern in time could indicate a causal miRNA-mRNA association which is well established in the literature, especially considering cancer pathways (50–53). However, bioinformatic reanalysis showed no statistically significant negative correlation between miRNAs and mRNAs differentially deregulated.

FIG. 7

The number of deregulated (up- and downregulated) miRNAs and mRNAs examined in baboons over time after irradiation are shown. The highest number of deregulated miRNAs and mRNAs was found on days 7 and 28, respectively.

We were curious about the cellular reason for these observed gene expression changes. Our analysis of the number of deregulated miRNAs over time indicated a potential origin from the proportion of irradiated lymphocytes, which would explain the unusual pattern (peak at day 7 postirradiation, corresponding to lowered lymphocyte counts observed in partial- relative to total-body exposure at the same time). It can be speculated, that some of the miRNAs originated from the proportion of peripheral lymphocytes irradiated, as we found for FDXR, DDB2, WNT3 and POU2AF1 (see above). However, due to the low-dose rate and extended exposure times, large volumes of circulating blood cells will enter (and exit) the local radiation area, become irradiated and, following this rationale, no gene expression changes related to the local exposure pattern would be expected. Since we observed strong linear associations of miRNA expression with increased exposed body area (in percentages), it is likely that locally exposed bone marrow stem cells were the actual targets of ionizing radiation. With increasing exposed body area more bone marrow stem cells will be irradiated, which is then reflected in the peripheral blood.

3. Several of the miRNA and mRNA species showed gene expression changes significantly associated with the percentage of the exposed body area. We interpret this as a finding in support of gene expression measurements providing hints on the exposure pattern of the body. However, these strong associations were usually observed for one point in time only and generally were absent before or after this time point.

For a few genes, these correlations were observed in preirradiation samples, rendering these genes and their applicability uninformative. This was a humbling reminder that false positive associations can survive conservative statistical efforts to reduce chance findings due to multiple comparisons.

Our analysis strategy certainly reduced the chance of false positives by using Bonferroni corrections, but these analyses suffered from low sample sizes and, therefore, must be validated in other settings. False positives are one concern, but to overlook associations due to small sample sizes (false negatives) is another limitation and inherent to our baboon study. Our ongoing analyses in two other rhesus macaque cohorts of 63 and 142 animals will mitigate some small sample size concerns.

In summary, these findings point to the activation of several transient biological processes, which have, in most instances, returned to baseline during the follow-up period. Due to the transient gene expression changes, a different set of candidate mRNAs and miRNAs appeared to be required at each day after irradiation, but further validation in other species and larger studies is urgently required.

Rationale

Are individuals really “back to normal” after surviving life-threatening acute health effects? Recently, a series of published articles identified delayed effects of acute radiation exposure (DEARE) due to degenerative and inflammatory conditions in multiple organs [see below (54–56)]. Reports of chronic radiation syndrome (CRS) observed after chronic protracted doses exceeding the natural background 100-fold or after ARS argues against it (57). Fatigue, one among other symptoms of CRS, is reported to last over months or even years after cancer therapy (58). It was reported in patients after treatment of tumor entities such as breast cancer (59–61), prostate cancer (62) or lung cancers (62, 63). In addition, sublethal partial-body irradiation of baboons bears exposure features with radiotherapy cancer patients who may be at risk of second cancers. Finally, understanding long-lasting expression changes of genes identified as being of clinical significance (e.g., predictor of H-ARS), would improve our knowledge of underlying mechanisms and further describe the diagnostic window.

With our baboon model, using the two-phase study design and searching the whole genome for deregulated mRNAs and 667 miRNAs before, days 7, 28 and 75–106 postirradiation, we aimed to describe persistence of gene expression [Fig. 2 (64)]. Samples from days 1 and 2 were not used for this analysis, since they are more reflective an acute response. Preirradiation samples were used as the reference for calculating differential gene expression at the indicated time points. In other words, do we detect gene expression changes and even persistent gene expression changes which hypothetically might suggest molecular mechanisms related to acute (e.g., fatigue) or late health effects (tumorigenesis)? Since neither outcome of acute (e.g., fatigue) or late health effects was part of our baboon study, our results provide only a fragile basis for further studies examining a possible relationship of gene expression changes with other health effects.

Examinations in Baboons

Results can be summarized as follows:

The pattern, in which we observed strong and persistent alterations of miRNA gene expression over time and almost no corresponding transcriptomic response (mRNA), can be partly explained by a weakness in our validation procedure, as described in Methods. Otherwise, our findings align with previously reported studies indicating that miRNA changes perform better than mRNA changes in the classification and diagnosis of tumors (65).

The baboon study employed complicated body exposure patterns to mimic several possible real exposure scenarios. However, it was the goal of this research to examine persistent changes in gene expression, which was independently found irrespective of partial- and whole-body exposure patterns. Such exposure scenarios cannot be performed in healthy humans, making our findings particularly interesting.

Our results indicate a marked radiation-related modification of the post-transcriptome, observed in nonirradiated cell generations more than three months after exposure. Because our measurements were performed on peripheral whole blood, presumably our results reflect either altered blood stem cells or indicate a selection for less radiosensitive stem cells that reconstitute the peripheral blood after radiation exposure. That would explain the rather high fold changes observed in several of the miRNA species and their persistence in nonirradiated cellular progeny. For the gastrointestinal tract, two functionally distinct pools of intestinal stem cells are proposed (34). Actively cycling, radiosensitive, WNT pathway modulated stem cells vis a vis slowly cycling, radioresistant, WNT pathway refractory and injury-inducible stem cells. After radiation injury, a shift in the composition of these differently radiosensitive stem cells after irradiation may occur. A similar mechanism in the hematopoietic system could explain our results.

While the possible clinical implications of our findings were not examined explicitly, we can offer some inferences, below. Interestingly, many of the miRNA species we identified are already known to be linked with radiosensitivity (e.g., miR-22, miR-29c, miR-195, miR-212) and most of these miRNAs are thought to be associated with risk of leukemia and other cancer types (64). Tumorigenesis represents a multi-step process. We suggest that this shift in the composition of cells might make the exposed individuals more prone to developing a tumor during follow-up, and that this might represent the first step toward tumorigenesis. While it is uncertain which animal might get a “second hit” and finally develop a tumor, radiation-exposed individuals are undoubtedly at higher risk for developing secondary tumors and the changes observed might represent the first step of tumor initialization. This hypothesis is very much in line with established concepts such that several populations of adult stem cells have been identified as the cell of origin for cancer (36, 66). The capacity of stem cells for self-renewal via oncogenes or loss of tumor suppressor genes makes them prone to hijack these functions and to transform these cells (67, 68). The question remains whether these persistent gene expression alterations are associated with risk of primary or secondary tumors, and further research is needed to address this hypothesis. Therefore, the concept of persistent radiation-related gene expression changes may provide an avenue for early prediction and diagnosis of cancer, which, again, must be shown in future studies.

Several of our candidate genes are known to impact immunological and inflammatory mechanisms (e.g., miR-130b, miR-29c, miR-212) and this corresponds to a reported persistent altered state of the immune system (based on gene expression measurements) up to 30 days after irradiation in NHP (69). Immunological alterations have been implicated as a mechanism for persistent fatigue after treatment (70, 71). In this regard, miR-212 is of interest due to the persistently 48–78-fold downregulation over the entire time period and its known effect on radiosensitivity, its immune modulating effect and its association with cancer. The development of DEARE due to degenerative and inflammatory conditions affecting multiple organs (54), cerebrovascular injury (55) or lung injury (56) are other known health effects in line with our findings.

We are currently preparing for validation of these data in two larger independent rhesus macaque NHP cohorts.