Dose-Dependent Shift in Relative Contribution of Homologous Recombination to DNA Repair after Low-LET Ionizing Radiation Exposure: Empirical Evidence and Numerical Simulation

1. Introduction

Ionizing radiation (IR) affects DNA structure by inducing various types of damage, of which DNA double-strand breaks (DSB) are considered to be the most deleterious ones [1,2]. DSBs are considered to be critical DNA lesions because their misrepair can lead to severe mutations, oncogenesis, cell death or senescence [3]. In this regards, selection of an optimal repair pathway is crucial for the cell in terms of a final outcome [4]. Non-homologous end joining (NHEJ) and homologous recombination (HR) are two major DSB repair pathways in mammalian and human cells. NHEJ fuses the two broken ends, with little regard for homology, leading to deletions and other rearrangements [5]. In contrast, HR typically copies the missing information from the sister chromatid into the break site, resulting in exact reconstitution of the original sequence [5]. In the last decades, at least three alternative pathways of DSB repair were suggested to be operational following the ionizing radiation exposure. These distinct pathways, the single-strand annealing (SSA), microhomology-mediated end-joining (MMEJ), and alternative end-joining (Alt-EJ), are distinguished based on the amount of DNA sequence complementarity used to align DNA ends [6].

Revealing the mechanistic nature of DSB repair pathways and evaluating their manifestation in response to different radiation doses aims to improve our understanding the deleterious effects of ionizing radiation and our ability to predict them [7,8]. Development of our predictive capabilities impacts the accuracy of radiotherapy, radiation-induced cancer risk estimations, as well as relevant for health issues of people living in areas with high background radiation and for future manned space exploration where astronauts can be exposed to complex radiation fields.

One of the initial events in the complex process of DNA repair from DSBs is the phosphorylation of the core histone H2AX by kinases of the PIKK family (phosphatidylinositol 3 kinase related kinases) ATM, ATR, DNA-PKcs in the flanking regions of DSBs of chromatin [9]. Dynamic microstructures containing thousands of copies of the phosphorylated histone H2AX (γ-H2AX) called “foci” in the literature are easily visualized using immunocytochemical staining [10-12]. Analysis of γH2AX foci provides information on the number of DNA repair sites from DSBs, their distribution over the nuclear volume, and on the kinetics of repair. However, it is equally important to evaluate not only changes in the total number of DSBs, but also the proportion of DSBs repaired by the HR correct pathway. For this purpose, the analysis of the key HR protein Rad51 foci is most often used [13].

Computational and mathematical modeling associated with cancer radiotherapy and radiation risk assessments has been undertaken for the decades and proceeds along with the new biophysics phenomena that have been identified [14-18]. Along with the modeling of radiation risk itself, there are a series of models designed to understand kinetics of radiation-induced DNA damage repair and associated secondary cancer risk. Combining the numerical simulation with experimental methods allows one to identify new mechanistic properties of the DNA repair process, which are hardly accessible to experimental studies but could be important for detailed understanding of the variety of outcomes induced by ionizing radiation.

One of the intriguing tasks being potentially solved via combining the experimental measurements with simulation techniques is the activation of specific DSB pathways following the low-dose radiation exposure. Currently, there are many works devoted to changes in the number of DNA damage foci in mammalian cells irradiated at low doses [19-21]. However, the question about DSB repair efficiency after low-dose radiation exposure is one of the most discussed and controversial topics in modern radiation biology.

The current study is aimed at the experimental probing of relative contribution of the HR pathway to DSB repair, induced by X-ray exposure in the dose range from 20 to 1000 mGy, with the subsequent numerical evaluation of γ-H2AX foci remaining considered as the main DSB repair marker.

2. Materials and Methods

2.1. Experimental study

2.1.1. Cell culture

The studies were performed with a human skin fibroblast culture (Cell Applications, United States, Cat.no. 106K-05a). Cells were cultured at 37°C and 5% CO₂ in a standard DMEM culture medium with a high glucose content (4.5 g/L) (Thermo Fisher Scientific, United States) supplemented with 2 mmol/L L-glutamine (Thermo Fisher Scientific, United States), 100 U/mL penicillin, 100 μg/mL streptomycin (Thermo Fisher Scientific, United States) and 10% fetal bovine serum (Thermo Fisher Scientific, United States). Fibroblasts were then plated on 4 cm² coverslips in 35 mm Petri dishes (SPL Lifesciences, South Korea) with a density of 10⁴ cells/cm². 3-5 passages of cells in the phase of exponential growth (cell population density ~ 60–70%) were used for the experiments.

2.1.2. Irradiation

The cells were irradiated using a RUB RUST-M1 X-irradiator (Russia), equipped with two X-ray emitters, at a dose rate of 40 mGy/min, voltage of 100 kV, current of 0.8 mA, 1.5 mm Al filter, and a temperature of 4°C (LAB ARMOR BEADS thermal granules were used for cooling). The error in the exposure dose was estimated to be within 15%. After irradiation, cells were incubated under the standard conditions of a CO₂ incubator (37°C, 5% CO₂) for 0.25–24 hours.

2.1.3. Immunocytochemistry

Cells were fixed on coverslips in 4% paraformaldehyde in PBS (pH 7.4) for 15 min at room temperature followed by two PBS rinses (pH 7.4) and permeabilization in 0.3% Triton-X100 (in PBS, pH 7.4) supplemented with 2% bovine serum albumin (BSA) to block non-specific antibody binding. For immunocytochemical staining, the slides were incubated for 1 hour at room temperature with a primary rabbit monoclonal antibody against γH2AX (dilution 1 : 200, clone EP854(2)Y, Merck-Millipore, United States), Rad51 (dilution 1 : 200, АВЕ257, Merck-Millipore, United States) or CENPF (dilution dilution 1:200, ab5, Abcam, USA) in PBS (рН 7.4), containing 1% BSA. After rinses with PBS, cells were incubated for 1 h at room temperature with the goat anti-mouse Alexa Fluor 488 conjugated secondary antibodies IgG (H + L) (Life Technologies, United States), dilution 1 : 600 and goat anti-rabbit rhodamine conjugated antibodies (Merck-Millipore, United States), dilution 1 : 400 in PBS (pH 7.4) with 1% BSA. ProLong Gold medium (Life Technologies, United States) was used with DAPI for DNA counter-staining and for the prevention of photo fading. Cells were viewed and imaged using a Nikon Eclipse Ni-U microscope (Nikon, Japan) equipped with a high definition ProgRes MFcool camera (Jenoptik AG, Germany). The filter sets UV-2E/C (340–380 nm excitation and 435–485 nm emission), B-2E/C (465–495 nm excitation and 515–555 nm emission) and Y-2E/C (540–580 nm excitation and 600–660 nm emission) were used. At least 200 cells were analyzed for each data point. γ-H2AX and Rad51 foci were counted by manual scoring and using DARFI (https://github.com/varnivey/darfi).

2.1.4. Statistical Analysis

Statistical and mathematical analyses of the data were conducted using the Statistica 8.0 software (StatSoft, Tulsa, OK, USA). The results are presented as means of three independent experiments ± standard error (SE).

2.2. Evaluating the percent contribution of HR to DSB repair

In order to evaluate the percent contribution of HR to the repair process, the quantitative model of radiation-induced DSB repair was used [22]. The model consists of five main parts. The first one evaluates the initial yield of radiation-induced DSBs. The other parts are the quantitative models of NHEJ, HR, SSA, and Alt-EJ repair systems. To simulate the processing of DNA lesions by repair enzymes the mass-action chemical kinetics approach is used.

The kinetics of DSB induction and repair is simulates as

$$\frac{d{N}_0}{dt}=\alpha (L)\frac{dD}{dt}{N}_{ir}-V_{\mathrm{NHEJ}}-V_{\mathrm{HR}}-V_{\mathrm{SSA}}-V_{\mathrm{micro}-\mathrm{SSA}}-V_{\mathrm{Alt}-\mathrm{EJ}}$$

(1)

where N₀ = N_ncDSB + N_cDSB; V_NHEJ, V_HR, V_SSA, V_micro-SSA, and V_Alt-EJ are the terms characterizing elimination of DSBs by the NHEJ, HR, SSA, micro-SSA, and Alt-EJ repair pathways respectively. Detalization of these terms is given in Eqs. A1 of Appendix A. In Eq. (1), N_ir is the fraction of irreparable DSBs which corresponds to the level of γ-H2AX foci remained in the cell 24 h postirradiation. The rate of initial DSB induction is calculated as$d{N}_0/dt=\alpha (L)dD/dt$, similarly as it is done in [23-25]. Here, D is the dose of ionizing radiation (Gy) and α(L) is the slop coefficient of linear dose dependence which describes DSB induction per unit of dose (Gy^–1 per cell).

Enzymatic interactions occurred in the course of repair are simulated as follows. Within the NHEJ pathway, the stage of Ku binding to a DSB is represented by the following kinetic equation

$$\left[\mathrm{DSB}\right] + [\mathrm{Ku}]\underset{K_{-1}}{\overset{K_1}{\rightleftarrows }}[\mathrm{DSB} · \mathrm{Ku}],$$

(2)

where quantities in brackets denote time-dependent intracellular concentrations of repair complexes; K values with an appropriate subscript used to represent the dimensional reaction rate-constants. Here [DSB] is the number of DSBs that undergo binding by Ku, $[\mathrm{DSB} · \mathrm{Ku}]$ is the level of resulting intermediate complex.

The stage of recruitment of the DNA-dependent protein kinase catalytic subunit (DNA-PKcs) and Artemis to a DSB site is described as

$$[\mathrm{DSB} · \mathrm{Ku}] + [\mathrm{DNA}-\mathrm{PKcs}]\underset{K_{-2}}{\overset{K_2}{\rightleftarrows }}[\mathrm{DSB} · \mathrm{DNA}-\mathrm{PK}/\mathrm{Art}],$$

(3)

In this kinetic equation [DNA-PK] denotes a complex of Ku and DNA-PKcs. Art indicates that the mentioned intermediate complex is formed in the presence of Artemis.

The autophosphorylation of DNA-PKcs is represented by

$$[\mathrm{DSB} · \mathrm{DNA}-\mathrm{PK}/\mathrm{Art}]\stackrel{K_3}{\to }[\mathrm{DSB} · {\mathrm{DNA}-\mathrm{PK}/\mathrm{Art}}^{\mathrm{P}}],$$

(4)

where the superscript P defines the phosphorylated DNA-PK product.

The subsequent end-bridging process is described as a junction of two ${[\mathrm{DSB}}_{\mathrm{NHEJ}} · {\mathrm{DNA}-\mathrm{PK}/\mathrm{Art}}^{\mathrm{P}}]$ constructs formed at the previous repair stage:

$$[\mathrm{DSB} · {\mathrm{DNA}-\mathrm{PK}/\mathrm{Art}}^{\mathrm{P}}] + [\mathrm{DSB} · {\mathrm{DNA}-\mathrm{PK}/\mathrm{Art}}^{\mathrm{P}}] \underset{K_{-4}}{\overset{K_4}{\rightleftarrows }} \left[\mathrm{Bridge}\right],$$

(5)

Here the [Bridge] intermediate complex characterizes the final product of the reaction.

The further assembling of the NHEJ repair complex promotes the recruitment of LigIV with its associated factors XRCC4 and XLF, and subsequent involvement of the polynucleotide kinase phosphatase (PNKP) with a break site:

$$\left[\mathrm{Bridge}\right] + [\mathrm{LigIV}/\mathrm{XRCC}4/\mathrm{XLF}] \underset{K_{-5}}{\overset{K_5}{\rightleftarrows }} [\mathrm{Bridge} · \mathrm{LigIV}/\mathrm{XRCC}4/\mathrm{XLF}],$$

(6)

$$[\mathrm{Bridge} · \mathrm{LigIV}/\mathrm{XRCC}4/\mathrm{XLF}] + [\mathrm{PNKP}] \underset{K_{-6}}{\overset{K_6}{\rightleftarrows }} [\mathrm{Bridge} · \mathrm{LigIV}/\mathrm{XRCC}4/\mathrm{XLF} · \mathrm{PNKP}].$$

(7)

The final step of the NHEJ pathway implying recruitment of a polymerase can be denoted as

$$\begin{array}{l}[\mathrm{Bridge}·\mathrm{LigIV}/\mathrm{XRCC}4/\mathrm{XLF}·\mathrm{PNKP}]+[\mathrm{Pol}]\underset{K_{-7}}{\overset{K_7}{\rightleftarrows }}[\mathrm{Bridge}·\mathrm{LigIV}/\mathrm{XRCC}4/\mathrm{XLF}·\mathrm{PNKP}·\mathrm{Pol}]\stackrel{K_8}{\to }\\ \stackrel{K_8}{\to }[\mathrm{dsDNA}]+[\mathrm{LigIV}/\mathrm{XRCC}4/\mathrm{XLF}]+[\mathrm{Pol}]+[\mathrm{PNKP}]+ [\mathrm{DNA}-\mathrm{PKcs}]+[\mathrm{Ku}]. \end{array}$$

(8)

In this consideration, after gap filling and ligation are finalized, it is accepted that the repair complexes dissociate, leaving the recovered double-strand DNA (dsDNA).

The first stages of HR associated with action of MRN, its co-factors (CtlP, ExoI, Dna2) and ATM are described as follows in the model:

$$\left[\mathrm{DSB}\right] + [\mathrm{MRN}/\mathrm{CtIP}/\mathrm{ExoI}/\mathrm{Dna}2]\underset{P_{-1}}{\overset{P_1}{\rightleftarrows }}[\mathrm{DSB}·\mathrm{MRN}/\mathrm{CtIP}/\mathrm{ExoI}/\mathrm{Dna}2],$$

(9)

$$[\mathrm{ATM}]\stackrel{P_2}{\to }{[\mathrm{ATM}}^{\mathrm{P}}],$$

(10)

$$\begin{array}{l}[\mathrm{DSB}·\mathrm{MRN}/\mathrm{CtIP}/\mathrm{ExoI}/\mathrm{Dna}2]+[{\mathrm{ATM}}^{\mathrm{P}}]\underset{P_{-3}}{\overset{P_3}{\rightleftarrows }}[\mathrm{DSB}·\mathrm{MRN}/\mathrm{CtIP}/\mathrm{ExoI}/\mathrm{Dna}2·{\mathrm{ATM}}^{\mathrm{P}}]\stackrel{P_4}{\to }\\ \stackrel{P_4}{\to }[\mathrm{ssDNA}]+[\mathrm{MRN}/\mathrm{CtIP}/\mathrm{ExoI}/\mathrm{Dna}2]+[{\mathrm{ATM}}^{\mathrm{P}}]\end{array}$$

(11)

where MRN (Mre11-Rad50-Nbs1) complex interacting with other repair factors is considered as a single complex and superscript P defines the phosphorylated species.

An involvement of the replication protein A (RPA) eliminating secondary structure of DNA and protecting single-stranded regions from other enzymatic activities is denoted by

$$\left[\mathrm{ssDNA}\right] + [\mathrm{RPA}]\underset{P_{-5}}{\overset{P_5}{\rightleftarrows }}[\mathrm{ssDNA} · \mathrm{RPA}],$$

(12)

The next HR step associated with formation of Rad51 filament is described with the following kinetic equation:

$$\begin{array}{l}[\mathrm{ssDNA} · \mathrm{RPA}] +\text{} [\mathrm{Rad}51/\mathrm{Rad}51\mathrm{par}/\mathrm{BRCA}2]\underset{P_{-6}}{\overset{P_6}{\rightleftarrows }}[\mathrm{ssDNA} · \mathrm{RPA} · \mathrm{Rad}51/\mathrm{Rad}51\mathrm{par}/\mathrm{BRCA}2]\stackrel{P_7}{\to }\\ \stackrel{P_7}{\to }[\mathrm{Rad}51 \mathrm{filament}] + \left[\mathrm{RPA}\right] \end{array}$$

(13)

where Rad51par abbreviation denotes five biologically important Rad51 paralogs (Rad51B, Rad51C, Rad51D, XRCC2, and XRCC3), and [Rad51 filament] defines the $[\mathrm{ssDNA} · \mathrm{Rad}51/\mathrm{Rad}51\mathrm{par}/\mathrm{BRCA}2]$ complex.

Formation of a displacement loop (D-loop) and two ways of its subsequent resolution are represented as follows:

$$\begin{array}{l}[\mathrm{Rad}51 \mathrm{filament}] +\text{} [{\mathrm{DNA}}_{\mathrm{inc}}]\underset{P_{-8}}{\overset{P_8}{\rightleftarrows }}[\mathrm{Rad}51 \mathrm{filament} · {\mathrm{DNA}}_{\mathrm{inc}}]\stackrel{P_9}{\to }\\ \stackrel{P_9}{\to }[\mathrm{D}-\mathrm{loop}] + [\mathrm{Rad}51/\mathrm{Rad}51\mathrm{par}/\mathrm{BRCA}2]\end{array}$$

(14)

$$[\mathrm{D}-\mathrm{loop}] \stackrel{P_{10}}{\to }\left[\mathrm{dHJ}\right]\stackrel{P_{11}}{\to } \left[\mathrm{dsDNA}\right] + {[\mathrm{DNA}}_{\mathrm{inc}}],$$

(15)

$$[\mathrm{D}-\mathrm{loop}] \stackrel{P_{12}}{\to } \left[\mathrm{dsDNA}\right] + {[\mathrm{DNA}}_{\mathrm{inc}}],$$

(16)

where [DNA_inc] defines the incoming DNA duplex, $[\mathrm{Rad}51 \mathrm{filament} · {\mathrm{DNA}}_{\mathrm{inc}}]$ complex is assumed to contain Rad54 protein, and dHJ denotes the double Holiday junction.

The SSA pathway is given as the following set of kinetics equations. The first step assuming interaction with Rad52 is denoted by

$$[\mathrm{ssDNA} · \mathrm{RPA}] + [\mathrm{Rad}52]\underset{Q_{-1}}{\overset{Q_1}{\rightleftarrows }}[\mathrm{ssDNA} · \mathrm{RPA} · \mathrm{Rad}52],$$

(17)

where the $[\mathrm{ssDNA} · \mathrm{RPA}]$ complex is the same as in Eq. (12).

The junction between Rad52 heptamer rings and each ssDNA termini allowing formation of a flapped structure is represented as

$$[\mathrm{ssDNA} · \mathrm{RPA} · \mathrm{Rad}52] + [\mathrm{ssDNA} · \mathrm{RPA} · \mathrm{Rad}52] \stackrel{Q_2}{\to } [\mathrm{Flap}] + [\mathrm{RPA}].$$

(18)

The cutting of the flapped ends by the ERCC1/XPF endonuclease and final ligation of a damaged site with LigIII complex are simulated by

$$[\mathrm{Flap}]+[\mathrm{ERCC}1/\mathrm{XPF}]\underset{Q_{-3}}{\overset{Q_3}{\rightleftarrows }}[\mathrm{Flap}·\mathrm{ERCC}1/\mathrm{XPF}]\stackrel{Q_4}{\to }[{\mathrm{dsDNA}}_{\mathrm{nicks}}]+[\mathrm{Rad}52]+[\mathrm{ERCC}1/\mathrm{XPF}], \left(19\right)$$

(19)

$$[{\mathrm{dsDNA}}_{\mathrm{nicks}}] + [\mathrm{LigIII}] \underset{Q_{-5}}{\overset{Q_5}{\rightleftarrows }} [{\mathrm{dsDNA}}_{\mathrm{nicks}}· \mathrm{LigIII}] \stackrel{Q_6}{\to } [\mathrm{dsDNA}] + [\mathrm{LigIII}].$$

(20)

Here, Rad52 and ERCC1/XPF assumed to dissociate from the processing site.

The simulation of the alternative end-joining pathways is based on the hypothesis suggesting two different Ku-independent repair mechanisms [26]. First of them is represented by MMEJ which is sometimes considered as an independent end-joining mechanism distinct from the other pathways. Meanwhile, experimental evidence suggests that this type of repair represents a specific subclass of the SSA pathway called as micro-SSA [27,28]. On the basis of these considerations, in our model, rejoining via MMEJ is simulated as the additional part of DSBs which follow the SSA pathway described above.

To simulate the Alt-EJ pathway, an additional mechanistic explanation is proposed. According to the recent hypothesis, Alt-EJ requires activity of MRX complex and possibly exhibits the same initiation steps as HR [26]. Therefore, the initial stages of Alt-EJ can be described by Eqs. (9)-(11). After the production of ssDNA, activity of PARP1 recruited to the DSB site is characterized by Eq. (21).

$$\left[\mathrm{ssDNA}\right] + [\mathrm{PARP}1] \underset{R_{-1}}{\overset{R_1}{\rightleftarrows }} [\mathrm{ssDNA} · \mathrm{PARP}1] .$$

(21)

The kinetics of microhomology production via a polymerase activity is simulated with

$$[\mathrm{ssDNA} · \mathrm{PARP}1] + [\mathrm{Pol}] \stackrel{R_2}{\to } [\mathrm{ssDNA} · \mathrm{Pol}] +[\mathrm{PARP}1],$$

(22)

$$[\mathrm{ssDNA} · \mathrm{Pol}]\stackrel{R_3}{\to } [\mathrm{MicroHomol}] +[\mathrm{Pol}].$$

(23)

In Eq. (23), [MicroHomol] denotes the yield of microhomology produced.

The final step of Alt-EJ is believed to rely on LigI activity [26]. In the model, this stage is represented as follows

$$[\mathrm{MicroHomol}]+[\mathrm{LigI}]\underset{R_{-4}}{\overset{R_4}{\rightleftarrows }}[\mathrm{MicroHomol} · \mathrm{LigI}]\stackrel{R{}_5}{\to }[\mathrm{dsDNA}]+[\mathrm{LigI}].$$

(24)

Kinetics of induction of γ-H2AX foci is simulated by summing up all active forms of DNA-PKcs and ATM, which are considered in the model, similarly as it was done in [23]

$$V_{\mathsf{\gamma }-\mathrm{H}2\mathrm{AX}+}=\frac{K_9 [\mathrm{Sum}] [\mathrm{H}2\mathrm{AX}]}{K_{10}+[\mathrm{Sum}]},$$

(25)

where [H2AX] is the level of histone variant H2AX and

$$\begin{array}{l} [\mathrm{Sum}]=[\mathrm{DSB} · {\mathrm{DNA}-\mathrm{PK}/\mathrm{Art}}^{\mathrm{P}}]+[\mathrm{Bridge}]+[\mathrm{Bridge} · \mathrm{LigIV}/\mathrm{XRCC}4/\mathrm{XLF}]+\\ + [\mathrm{Bridge} · \mathrm{LigIV}/\mathrm{XRCC}4/\mathrm{XLF} · \mathrm{PNKP}]+[\mathrm{Bridge}·\mathrm{LigIV}/\mathrm{XRCC}4/\mathrm{XLF}·\mathrm{PNKP}·\mathrm{Pol}]+\\ + [\mathrm{DSB}·\mathrm{MRN}/\mathrm{CtIP}/\mathrm{ExoI}/\mathrm{Dna}2·{\mathrm{ATM}}^{\mathrm{P}}].\end{array}$$

(26)

The mechanism of γ-H2AX foci dephosphorylation is assumed to be proportional to the amount of repaired DNA and its spontaneous decay with the corresponding rate constants K₁₁ and K₁₂

$$V_{\mathsf{\gamma }-\mathrm{H}2\mathrm{AX}-}= K_{11}[\mathrm{dsDNA}]+K_{12}[\mathsf{\gamma }-\mathrm{H}2\mathrm{AX}].$$

(27)

To describe the interactions of repair enzymes with their substrates, the mass-action chemical kinetics approach was used. A dynamic change of intracellular concentrations of main repair complexes was expressed as the following differential equations

$$\frac{dX}{dt}=V_{+}(X_i,N_0)-V_{-}(X_i,N_0),$$

(28)

where $X_i$ $(i=1,\mathrm{...},n)$ is the intracellular level of the i-th repair complex, $t$ is time, the functions $V_{+}$ and $V_{-}$ describe the complex accumulation and degradation, respectively. The dimensionless form of the system of ordinary differential equations (ODE) referred to simulation of NHEJ, HR, SSA and two alternative repair pathways as well as its parameters and initial conditions are presented in Appendix A. In the present study, the solution of the ODE system is performed using the fourth-order Runge–Kutta method. The integration time-step is chosen as 10^–10 s in order to satisfy the fastest reactions being simulated in the model.

The initial yield (N⁰) of DSBs was calculated as

$$\frac{d{N}_0}{dt}=\alpha (L)\frac{dD}{dt},$$

(29)

where D is the radiation dose (Gy) and $\alpha (L)=a\exp (-bL)$ is the slop coefficient of linear dose dependence which describes DSB induction per unit of dose (Gy^–1 per cell). Parameters a and b used in this equation are also presented in Appendix A. The particular details of parameter evaluation as well as simulation of each repair pathway can be found in [22].

In order to simulate DSB rejoining in the asynchronous cell culture, we set a model cell cycle distribution similar to the distribution of HF19 cells observed in the work by [29] where 40% of cells were in G₀/G₁, <10% in S phase and 50% in G₂/M. Assuming that S-phase cells include two equal sub-fractions of early-S and late-S cells, the final distribution was set as following: 45% of cells is in G₀/G₁ and early S phases, and 55% is in late S and G₂/M. Based on the data discussed in [25,27,30,31], the model provides the following pathways for the repair of non-clustered and clustered DSBs. In G₀/G₁ and early S phases non-clustered DSBs are likely to follow rejoining via NHEJ, while the limited number of complex DSBs may undergo PARP1-dependent Alt-EJ. The other pathways are assumed to be unavailable or masked by NHEJ in these phases. In the late S and G₂/M, non-clustered DSBs may follow NHEJ and may be a substrate for SSA and micro-SSA pathways. The clustered DSBs are suggested to proceed through HR, SSA, and Alt-EJ. Finally, in our calculations the models of corresponding repair mechanisms are processed according the following scheme

$$\begin{gathered} {\mathrm{G}}_0/{\mathrm{G}}_1, \mathrm{and} \mathrm{early} \mathrm{S} \left\{\begin{array}{l}0.45 T{N}_{\mathrm{ncDSB}}\to \mathrm{NHEJ},\\ 0.45 T{N}_{\mathrm{cDSB}}\to \mathrm{Alt}-\mathrm{EJ},\end{array}\right. \\ \mathrm{late} \mathrm{S} \mathrm{and} {\mathrm{G}}_2/\mathrm{M} \left\{\begin{array}{l}0.55 T{N}_{\mathrm{ncDSB}}\to \mathrm{NHEJ}, \mathrm{SSA},\\ 0.55 T{N}_{\mathrm{cDSB}}\to \mathrm{HR}, \mathrm{SSA}, \mathrm{Alt}-\mathrm{EJ},\end{array}\right. \end{gathered}$$

(30)

where numbers 0.45 and 0.55 represent the fractions of cells in G₀/G₁ and early S, and in late S and G₂/M phases, respectively, N_ncDSB and N_cDSB denote amounts of non-clustered and clustered DSBs, and T is the total size of the exponentially growing cell culture set to be 1×10⁵ cells, as mentioned in [29]. The model accounts for the very cell cycle-dependence of the repair pathways through the initiating the computation procedure with the levels of N_ncDSB and N_cDSB parameters according to the cell cycle distribution of the cell culture used in the experiment.

The evaluation of the percent contribution of HR to DSB repair was done through the obtaining the calculated time-courses of γ-H2AX and Rad51 foci for the particular dose of X-rays. The HR contribution to PHR was evaluated as the following ratio

$$P_{\mathrm{HR}}=100\times \overline{y_9}/ \overline{x_{14}},$$

(31)

where $\overline{y_9}$ and $\overline{x_{14}}$ are the mean numbers of Rad51 and γ-H2AX foci respectively, counted per 24-hour simulation period.

To obtain a dependence of PHR on X-ray dose, the ODE system was run with a sufficiently small dose step equalling to 0.1 mGy within the range of 0–1000 mGy. This simulation procedure yielded a curve of PHR (D) dependence that is expressed by the formula

$$P_{\mathrm{HR}}(D)=100\times \overline{y_9(D)}/ \overline{x_{14}(D)}.$$

(32)

3. Results

3.1. Experimental results

In Figure 1, γ-H2AX foci yields are shown, as scored at 0.25–24 h after exposure to different doses of X-rays ranged from 20 mGy to 1000 mGy. Foci yields scored in cells exposed to 500 mGy and 1000 mGy tend to peak values at 0.25–1.0 h, while exposure to doses of 20–250 mGy results in maximum levels of foci at 0.25–2.0 h. Closer to 4 h after exposure, the signal of γ-H2AX sharply decreases with time due to completion of the fast DSB repair being a feature of NHEJ. Then, up to 24 h after exposure, the slow DSB repair occurs, which is typically associated with complex DSBs being repaired via HR pathway.

It is noteworthy that 24 h after exposure to 40–80 mGy the level of γ-H2AX foci does not drop to control values, indicating a slowdown in the DSB repair kinetics following low doses of ionizing radiation. On the other hand, 24 h after exposure to 500 mGy and 1000 mGy the observed levels of γ-H2AX foci are slightly below the control values. A possible explanation is that the exposure to low and moderate doses of ionizing radiation affects cell proliferation in an oppositely directed manner. Low doses stimulate cell proliferation, while higher doses lead to its reduction.

One of the sources of DNA DSBs in proliferating cells is the collapse of replication forks in S phase of the cell cycle [32]. The repair of DSBs induced in this process follows through the HR pathway [32,33]. The histone H2AX in S phase is phosphorylated mainly by ATR kinase [34,35]. As a consequence, the percentage of cells in S phase contributes to the average amount of γ-H2AX foci in asynchronous cell populations. This leads to an overestimation of γ-H2AX foci level in cells exposed to low doses and to an underestimation of it in cells irradiated with moderate doses. The observed pattern generally meets our previous findings in human fibroblasts and mesenchymal stromal cells exposed to X-rays [36,37].

To estimate the contribution of HR to the entire process of DNA DSB repair, we analyzed changes in the level of radiation-induced Rad51 foci being the key marker of HR (Figure 2). A statistically significant increase in the number of Rad51 foci was observed 2 h after irradiation and reached its maximum values by 6 h. After 6 h, a decrease in the number of Rad51 foci was indicated, and by 24 h after irradiation, the level of foci dropped down to control values in cells exposed to 160–1000 mGy and remained higher than control in cells irradiated with doses of 40 mGy and 80 mGy.

3.2. Percent contribution of HR to DSB repair

To confirm the model validity for radiation doses in the range of 20–1000 mGy, the calculated time-courses of γ-H2AX and Rad51 foci were compared to the experimental data for normal human skin fibroblasts exposed to X-rays at the doses of 20, 40, 80, 160, 250, 500, and 1000 mGy. The results of comparison are shown in Figure 3 and Figure 4, which depict kinetics of these foci in the range of 0–24 h after irradiation.

Although the exposure to doses of 20–250 mGy results in relatively low absolute levels of γ-H2AX, the entire repair process does not seem to occur faster than in the case of higher considered doses. The proportion of γ-H2AX foci observed after exposure to 20–250 mGy remains elevated for a similarly long period as following the exposure to the doses of 500–1000 mGy. We hypothesize that the observed pattern of DSB repair occurs due to preferential activation of HR pathway, which takes longer to eliminate the damage.

Interestingly, that following the exposure to doses of 160 mGy and higher the levels of γ-H2AX foci remaining 24 h repair after irradiation were lower than the background value. This suggests that DSBs normally induced as a consequence of non-radiation factors can undergo a more effective elimination if cells enter into the process of radiation-induced DSB repair. In could be also expected that in cells that have successfully completed repair of a certain portion of radiation-induced DSBs the background levels of γ-H2AX foci can be suppressed as compared to non-irradiated cells at least for a residual period after exposure.

Overall juxtaposition of simulated and experimentally measured data suggests validity of the considered DSB repair model for relatively low radiation doses, enabling calculation of time-courses of γ-H2AX and Rad51 foci within the full dose range of interest (0–1000 mGy) to obtain a continuous dependence of P_HR on D.

The results yielding the percent contribution of HR to DSB repair are shown in Figure 5. The general pattern of HR involvement in the DSB repair is characterized by near-exponentially decreasing function, which is supposed to depend on a cell line, type of ionizing radiation exposure and likely on other factors. The dependence of the HR on different factors requires more detailed examination.

3.3. Dose-dependent changes in the S/G2-phase сell fractions

DNA repair by homologous recombination occurs predominantly in the S/G2 phases of the cell cycle. Therefore, for the correct interpretation of the obtained results, it was important to estimate the changes in the proportion of cells in the S/G2 phases in the irradiated cell populations. For this purpose, we used immunohistochemical analysis of a protein-marker of cells in S/G2 phases - Centromere protein F (CENPF). This protein, being a component of the nuclear matrix during G2 phase, has been used as a marker of S/G2 cells [38,39]. Its synthesis commences in early S phase and ceases in M phase, with a peak at G2 phase [39].

From the results presented in Figure 6, it can be seen that irradiation at low doses (20-80 mGy) leads to a slight increase in the proportion of CENPF positive cells. However, these changes are not statistically significant (p>0.05). Whereas irradiation at doses of 250, 500 and 1000 mGy leads to a significant decrease in the proportion of CENPF positive cells by 2.15 (p=0.004), 3.69 (p=0.0003) and 6.93 (p=0.0002) times, respectively (Figure 6). In general, after irradiation at doses of 160-1000 mGy, the pattern of change in the proportions of S/G2-phase fibroblasts is characterized by near-exponentially decreasing function. The obtained results are in good agreement with dose-dependent changes in the relative contribution of HR in irradiated cells.

4. Discussion

Our study reveals distinct patterns of γ-H2AX and Rad51 foci dynamics following the application of low and moderate doses of X-rays. We identified several parameters of foci dynamics, which demonstrate different regularities in response to low and moderate doses, including periods of reaching peak levels of γ-H2AX foci, numbers of γ-H2AX and Rad51 foci remaining 24 h post-irradiation as compared to control levels. All these characteristics of DNA DSB repair demonstrate a pronounced dose-dependent shift. Within the dose steps selected for the analysis, the dose-dependent shift associated with reaching peak levels of γ-H2AX foci is observed between 250 mGy and 500 mGy. The doses triggering the levels of remaining γ-H2AX and Rad51 foci appeared to be within the ranges of 80—500 mGy and of 80—160 mGy respectively.

Juxtaposition of our results with other findings [40] suggests another confirmation of a dose-dependent shift in the activity of DNA repair systems, at least with regards to the repair of radiation-induced DSBs. On the one hand, the number of DSBs increases linearly with radiation dose and the entire cell response depends on the correctness of DNA repair. Our results show that relatively low doses of low-LET ionizing radiation lead to a higher contribution of the error-free repair via HR pathway. Since HR takes longer than NHEJ, low-dose mediated activation of DNA repair mechanisms may not only protect DNA from the immediate damage, but also result in prolonged adaptive responses, protecting the cell from future oxidative insult (such as high-dose radiation), as is also confirmed in [40-46]. The observed shift to error-free HR allows restoration of the genome with maximum fidelity.

The mechanistic nature of the relative contribution of HR at low doses remains under debates. Nevertheless, our results support one of the possible models of DSB repair pathway engagement in responses to low doses of low-LET ionizing radiation [47]. According to these findings, error-prone pathways including NHEJ, Alt-EJ and SSA are partly or completely suppressed and likely only operate when HR fails to process the damage. Increase in the dose of low-LET ionizing radiation leads to the suppression of HR by mechanisms that remain to be identified, while the contribution of NHEJ increases and becomes the first choice. The decrease in the fraction of cells in the S/G2 phases of the cell cycle after irradiation at doses of 250–1000 mGy, shown in our work, suggests that the simplest mechanistic mechanism for reducing the contribution of HR with increasing radiation dose is G1/S arrest of the cell cycle. The dose-dependent G1/S cell cycle arrest leads to a decrease in the faction of cells in which HR repair is possible. Our study was limited to the radiation dose of 1 Gy, while the model in [47] postulates that DNA end resection remains active showing signs of suppression only above 20 Gy. Therefore, increased engagement of error-prone Alt-EJ and SSA under conditions of persisting resection and partially suppressed HR can be identified here only through the simulation approach.

Considering plenty of well-established facts on a strong compensatory power of the DBS repair in eukaryotes, it is unlikely to assume that the radiation effects related to DSB removal follow the linear dose dependence. Even the mechanistic nature of the assembly of DSB repair super-complexes in response to the appearance of radiation-induced lesions suggests a non-linear kinetics of this process, which may blur the resulting outcome within a certain dose range. It may suggest that genotoxicity and late risks associated with the quality of DSB repair after exposure to low-LET radiation can also undergo some radiation dose threshold.