A unified multi‐activation (UMA) model of cell survival curves over the entire dose range for calculating equivalent doses in stereotactic body radiation therapy (SBRT), high dose rate brachytherapy (HDRB), and stereotactic radiosurgery (SRS)

Purpose Application of linear‐quadratic (LQ) model to large fractional dose treatments is inconsistent with observed cell survival curves having a straight portion at high doses. We have proposed a unified multi‐activation (UMA) model to fit cell survival curves over the entire dose range that allows us to calculate EQD2 for hypofractionated SBRT, SRT, SRS, and HDRB. Methods A unified formula of cell survival S=n/eDDo+n‐1 using only the extrapolation number of n and the dose slope of Do was derived. Coefficient of determination, R2, relative residuals, r, and relative experimental errors, e, normalized to survival fraction at each dose point, were calculated to quantify the goodness in modeling of a survival curve. Analytical solutions for α and β, the coefficients respectively describe the linear and quadratic parts of the survival curve, as well as the α/β ratio for the LQ model and EQD2 at any fractional doses were derived for tumor cells undertaking any fractionated radiation therapy. Results Our proposed model fits survival curves of in‐vivo and in‐vitro tumor cells with R2 > 0.97 and r < e. The predicted α, β, and α/β ratio are significantly different from their values in the LQ model. Average EQD2 of 20‐Gy SRS of glioblastomas and melanomas metastatic to the brain, 10‐Gy × 5 SBRT of the lung cancer, and 7‐Gy × 5 HDRB of endometrial and cervical carcinomas are 36.7 (24.3–48.5), 114.1 (86.6–173.1),, and 45.5 (35–52.6) Gy, different from the LQ model estimates of 50.0, 90.0, and 49.6 Gy, respectively. Conclusion Our UMA model validated through many tumor cell lines can fit cell survival curves over the entire dose range within their experimental errors. The unified formula theoretically indicates a common mechanism of cell inactivation and can estimate EQD2 at all dose levels.


INTRODUCTION
Linear-quadratic model (LQ) that can fit the initial shoulder of cell survival curves at the low-dose domain has been extended for calculation of equivalent dose in 2-Gy fractions (EQD2) to the tumor and organs at risk (OARs) receiving large fractional doses during hypofractionated stereotactic body radiation therapy (SBRT), high dose rate brachytherapy (HDRB), and intracranial stereotactic radiation therapy (SRT) or stereotactic radiosurgery (SRS). Such an extrapolation by using the same α/β ratio derived at 2-Gy fractions might be inappropriate to predict cell survival at 7-20-Gy fractions, which are most likely located within the straight or almost straight portion of observed cell survival curves. 1,2 Correct estimation of radiation responses of human tumors and organs at risk (OARs) based on the intrinsic radiosensitivity of the cell lines measured in-vivo or invitro is required in design of any radiation therapy schemes particularly for hypofractionated radiation therapy with a high dose per fraction. Gurrero and Li 2 have extended the LQ model pertinent to stereotactic radiotherapy by modifying the β parameter into a product of the β and a function of time and repairing-rate constants. This modified LQ model has improved the fitting of cell survival curves at high doses as does by a four-parameter-based lethal and potential lethal (LPL) damage repairing model. [3][4][5] However, models with more parameters and complicated functions are more difficult to use than the simple LQ and multitarget (MT) models each using only two parameters. Garcia et al. 6 have derived the α and β parameters by fitting cell survival curves for different dose ranges. However, a prediction from such a model is restricted to the individual dose ranges. Park et al. 7 have proposed a universal survival curve (USC) model to calculate equivalent dose by applying an LQ model at the low-dose domain, a MT model at the high-dose domain, and an analytical formula for determination of the transition dose point between the two dose domains. Astrahan 8 has then added a linear tail function at the low dose domain and Kehwar et al. 9 have modified the transition dose formula using D o , n, and α to improve the USC model. A review of the radiobiology for SBRT by Garau et al. 10 has outlined some differences between the LQ and USC models in the calculation of BED and EQD2. More importantly, application of the USC uses different descriptions of radiosensitivity of tumor or tissue cells when their fractional doses change across the transition dose during SBRT, SRT, and HDRB. In addition, there is no evidence of cell responses suddenly switching at the "transition dose point". A better understanding of cell responses to avoid a systematic discrepancy between any model predictions and experimental observations at all dose levels is urgently needed for evaluation of hypofractionated radiotherapy 11,12 in light of recent technological advances.
This article attempts to resolve the radiobiological discrepancies (or catastrophes) between observed survival curves and model predictions, at either high-dose domain by the LQ model or at low-dose domain by the MT model, by introducing a UMA model to fit all cell survival curves over the entire dose range. A unified formula using only two parameters is empirically derived and then numerically validated with many observed survival curves that we could find and redraw from published data that have either repeated measurements or estimated experimental errors. Such a formula using doseindependent parameters allows us to estimate EQD2 and the α and β parameters for the LQ model at any doses per fraction. Figure 1 illustrates three repeatedly measured x-ray survival fractions of a human melanoma cell line by Weichselbaum et al. 13 which were refit by the MT model using the linear least square regression of the natural logarithm of the cell survival fraction S ¼ ne ÀD=D o for doses > 3 Gy with D o = 1.35 Gy and n = 5.4. The MT model overestimated cell survival at the initial shoulder (low doses) that results a poor R 2 = 0.04 and a large residual r = 1.22. The initial shoulder was better fit by the second-order polynomial trend line of the natural logarithm of the LQ model 14 S ¼ e ÀαDÀβD 2 for the low dose data with results of α = 0.13 Gy −1 and β = 0.06 Gy −2 , respectively. However, the LQ model underestimated the survival at high doses ≥ 7 Gy. The first author has proposed a cell survival formula of S ¼ a= e bD þ c À Á , where a, b, and c are dose-independent parameters. We now determine the formula through an empirical approach with no need of a theoretical derivation to be presented in another report.

2.A. A UMA Model derived from an observed cell survival curve
Prior to irradiation at dose of D ¼ 0 Gy, we have S ¼ 1 ¼ a e bD þc ¼ a 1þc , thus, a ¼ 1 þ c. At high doses with e bD >> c, the survival curve is almost straight as S≈ne À D Do ¼ a=e bD so that a = n, b = 1/D o and c ¼ n À 1. Now, we have obtained a unified formula of the survival curve using only two parameters as: The parameters of n and D o are similar to that of the MT model but iteratively determined through first setting of lnðS i ÞÀlnðS iÀ1 Þ and n ¼ S iÀ1 e D iÀ1 =D o at the two high dose points of D i and D i-1 , then adjusting the n and D o with the linear regression using least squares of Ln(n/S i + 1-n) verses D i /D o at all measured dose points for having the lowest relative residual (r) of modeled responses of S(D i ) from the observed S i as r ¼ 1 Þ 2 are the residual sum of squares and the corrected sum of squares of survival 15 , respectively. The relative residual of r is compared with the relative experimental errors of e = 1 I ∑ i¼I i¼1 ΔS i =S i , where ΔS i is one standard deviation of the multiple survival measurements at a dose point of i. I is the total number of dose points in the survival curve. The relative residuals are the second index of good fitness in addition to the R 2 since R 2 is insensitive to the discrepancy at the high doses with very low survival fractions as shown in Fig. 1 between the LQ and UMA models sharing the same R 2 but different r. Typically, a good fitting of a survival curve should have R 2 in a range of 0.9-1.00 and the relative residual smaller than that of the relative experimental errors, that is, r < e. In congruence with survival curves in literatures, mean values and one standard deviations of the survival fractions at individual dose points are redrawn and used in our curve fitting. Survival fraction normalized to 1 with no error bar at the zero dose is used for all CSCs.

2.B. Analytical Predictions of EQD2 for any radiotherapy schemes
There is only a single term of e D=D o that contains the independent variable of dose D in Eq. (1) for modeling the cell survival, thus, the biologically effective dose (BED) can be the physical dose D and there is no need for BED calculation in our UMA model. To achieve the same ending survival S end by a course of D-Gy N fractions as that treated with 2-Gy N 2 fractions, that is, Since n and D o are independent of the fractional dose of D, Eq. (2) allows us to calculate the EQD2 of the tumor or tissue cells receiving any fractional doses. If D changes with different sessions, the EQD2 at the end of treatment is just the summation of individual fractional EQD2. Another advantage of the proposed UMA model is to deal with tumor cells under hypoxic condition by simply changing the D o and n for the hypoxic tumor cells in Eq. (2).
The LQ model has also assumed S end ¼ S D Þso that Ln(S) (the natural logarithm of survival fraction) will linearly correlate with the BED that compose factors of the total physical dose D total and 1 þ D= α=β ð Þ ð Þ : To achieve the same ending survival (or BED if there is no change of α with the fractional dose) as that of 2-Gy fractions, one could determine the equivalent dose in 2-Gy per fraction: ÞÞ. The EQD2 estimated by the LQ model would have the same value or dose as Eq. (2) if both the α and α/β ratio in the LQ model were dose independent.
2.C. α and α/β ratio varied with D LQ model defines Ln S ð Þ ¼ ÀαD À βD 2 and by replacing S with our Eq. (1), we can determine the α and β parameters by the first and second derivatives of Ln(S) as: Clearly, α, β and α/β ratio change significantly with the fractional dose of D and our UMA determines the EQD2 by Eq. (2) without using the dose dependent α, β and α/β ratio.

2.D. Interpretation of the new UMA model
There are a number of radiobiological models and the choice of a model particularly for the prediction with extrapolations from low dose of~2 Gy to high doses of 10 Gy has been a major concern in clinical implementation. 16 The proposed UMA model is conceptually different from other models (including the latest USC model) by fitting both the shoulder and straight portions without any mechanism changes. In fact, recent studies have demonstrated that radiation induced DNA double-strand breaks (DSB) have the same nonhomologous end-joining (NHEJ) as the mainstay for sublethal damage repairing (SLDR) after a low-dose fraction and potentially lethal damage repairing (PLDR) after a high-dose fraction. 17,18 If both SLDR and PLDR depend mostly on NHEJ of DSB of DNA, it is reasonable to explore the UMA model for a simple and fundamental quantification of radiosensitivity of cells. For which, we have adapted the mean inactivation dose (MID) previously introduced for fractionated radiotherapy and brachytherapy by Fertil et al. 19 and we have obtained an analytical solution of The MID is determined by the D o and n of the cells receiving the radiation. One can compare the MID with the mean life span of a radioactive source which is an integration of "radioactive nuclei or excited atoms" over the time from 0 to ∞ with a result of T = 1/λ = 1.44 T 1/2 . Where λ and T 1/2 are decay constant and half-life, respectively. The life span T is a single parameter to describe the activity change with time and the MID may be used as the single parameter to describe the cell response to the radiation dose. Certainly, MID is more complicated than the life span and the n may serve as the number of activated cell death pathways of a live cell by the radiation dose. In fact, if n = 1, the survival curve is exactly given by S ¼ e ÀD=D o and the MID = D o . The number of active pathways, n, may change with the type of radiation as shown by Fig. 2 for the change of D o from 0.80 to 1.49 Gy and n from 1.0 to 3.8 to the same breast cancer cell line irradiated with α particles to γ-rays, respectively. Eq. (6) tells us that the change of MID of the breast cancer from 0.80 Gy for α-particle irradiation to 2.7 Gy for the γ-ray irradiation is magnified by the changes of D o and n. Thus, the activation number n plays an important role in MID or cell killings from the ionization radiation. We have named the new model as the unified multi-activation (UMA) model instead of unified multi-target model since the cell survival corresponds to the inactivation of nactivated pathways by the specific radiation and our modeled n and D o are differently from that of the MT model.

3.A.. Survival curves of human cells fitted with the UMA model
Our UMA model has provided the best fit among the LQ, MT and UMA models for in-vivo and in-vitro cell lines we have found in literature search with presentation of error bars or repeated measurements as shown in Fig. 1 with three measurements on individual dose points. The UMA model has also the best fit of three skin fibroblast cell lines from courtesy of Weichselbaum et al. 21 (not presented here). Survival curves of some malignant tumor cell lines of interest in intracranial SRS are presented in Fig. 3(a) for in-vivo hypoxic and oxic as well as in-vitro plated Bell cells 22 and plated MelH cells 13 two metastatic melanoma lines and two plated U373MG cellsa glioblastoma line 6,23 and in Fig. 3(b) for four metastatic uveal melanomas of OMM-1, OMM 2-2, OMM2-3, and OMM2-6 cells. 24 Apparently, invivo hypoxic and aerobic Bell cells have high D o of 6.5 and 5. 5 Gy, respectively. The extrapolation number n varies significantly with the cell lines. All thin lines are LQ model fitting curves with R 2 > 0.98 but systematically differing from thick lines of UMA model curves at high doses. Our α/β ratios should agree with that of the original publications but results of 24.6 and 20.5 Gy for the OMM-1 and OMM 2-3 cell lines are respectively doubled the original values of 14.1 AE 2.6 and 10.3 AE 1.6 Gy 24 which is due to the exclusion of the highest dose point for the LQ model fitting of the curves in the original report. Thus, α/β ratio of LQ is extremely sensitive to the dose range.
Results of tumors treated with extracranial SBRT are grouped in Fig. 4(a) for survival curves of squamous cell cancer (SCC) cell lines: 24-hour delay and no delay plated SW1573a lung SCC, 25 DaFua head and neck SCC 26 and SKXa base of tongue SCC 26 and in Fig. 4(b) for invivo lung cancer lines other than SCC including HX147large-cell carcinomas, HX149Mvariant small cell carcinomas, HX144 -Adenocarcinomas, and HC12classical small cell lung cancer (SCLC). 27 The invert shoulder having a negative curvature of SKX cells in Fig. 4(a) was fitted well with D o = 2.55 Gy and n = 0.2. Figure 5(a) represents in-vivo and in-vitro survival curves of HT29 cellsa colorectal cancer cell line, under hypoxic, aerobic, and no delay plate conditions. 22 Fig. 5(b) shows two prostate cell lines of DU145 and CP3 cells 6 with significantly different D o = 2.20 and 1.08 Gy and n = 2.8 and 20.0, respectively. Figures 6(a) and 6(b) address the tumor cell lines for HDRB with some irregular survival curves of the cervical adenocarcinomas irradiated under extreme hypoxia with or without pre-irradiation of the contact medium 28 and some in-vitro survival curves of the endometrial carcinomas, 29 respectively. For all of the tested cell lines, our new model has R 2 > 0.97 and relative residuals r < e. In other words, the UMA model has described all cell survival curves within the experimental uncertainties. The survival curves of a SKX cell line in Fig. 4(a) and two metastatic endometrial carcinomas (EC) of UM-EC-2 and UT-EC-2 in Fig. 6   A big advantage of UMA model is that the parameters of n and D o can be determined from preclinical measurements of cell survival curves within experimental dose ranges while α and α/β ratios determined with LQ models change significantly with dose ranges. There are negative α/β ratios for the tumors with a negative curvature (β < 0 and n < 1) for which hyperfractionated 1.2-Gy twice-daily radiation (BID) irradiation is desired. EQD2 for 20 Gy SRS of the primary glioblastomas or metastatic melanomas, range from 24.3 to 48.5 Gy with an average of 36.7 Gy, much lower than 50 Gy from the LQ model. EQD2 of 10 Gy × 5 SBRT of the lung carcinomas, except for the large cell lung carcinomas with 86.6 Gy, ranging from 101 to 173 Gy are much higher than 90 Gy from the LQ model (EQD2 = 83.3 Gy if α/β = 10 Gy). EQD2 for 3 Gy × 18 of H&N SCC, 5 Gy × 5 of colorectal adenocarcinomas and 2.66 Gy × 16 of the breast cancer are all higher than that from the LQ model. EQD2 for 7 Gy × 5 HDRB of cervical and endometrial carcinomas, except for the mixed cell population survival curve of NHIK3025, varies from 40 to 52.6 Gy. EQD2 for any tumors with mixed cell populations is calculated by using the low EQD2 of the most radioresistant cell population since the tumor cells would be dominated by the Using LQ model of Ln(S) = αBED, EQD2 calculation should not only use clinical α/β ratio but also α value that varies with the fractional dose.

3.C. Fractional dose selection
The benefits of a desired fractional dose are demonstrated by Figs. 7(a) and 7(b), presenting the ratios of S 10Gy /(S 2Gy ) 5 for the hypofractionated 10-Gy fractions and (S 1.2Gy ) 2/1.2 / S 2Gy for the hyperfractionated 1.2-Gy fractions in the D o and n plane, respectively. The lower the isocontour value, the higher tumor control is, by delivering the same physical dose in comparison with the 2-Gy fractions. Fig. 7(a) illustrates that tumor control in SBRT/SRT at 10-Gy fractions rapidly decreases with the value of D o but has a peak value across the value of n at ∂ S 10Gy = S 2Gy À Á 5 h i =∂n ¼ 0 for a given D o (>2 Gy). This is important since x-ray survival curves of  Fig. 7(b), the tumor control rate in hyperfractionated radiotherapy of 1.2 Gy per fraction rapidly decreases with the value of n and changes less with D o , in favor of smaller n and smaller D o . Both Figs 7(a) and 7(b) have the highest value of 1 located along the line of n = 1. Thus, the α-particle irradiation of the breast cancer with n = 1 has no difference among different fractional doses and it may be an advantage for taking continuously internal irradiation in theranostics using Bi-213 alphaparticle radiation. Certainly, the normal tissue complication should be considered in order to increase the therapeutic ratios. Due to large variations of dose and dose distribution to OARs in different treatment plans and delivery techniques, normal tissue complication requires further studies. The values of D o and n along the radiosensitivity (or MID) of the tumor and normal tissue cells are essential in optimal selection of fractional doses.

DISCUSSION
This study has validated a newly proposed UMA model that can fit available survival curves of human tumor cell lines over the entire dose range with R 2 > 0.97 and accurate predictions within the experimental uncertainties. A unified formula using only two dose-independent parameters and a single term in a linear form of the physical dose has enabled us to reliably estimate the EDQ2 based on experimentally measured survival curves of tumor cell lines. To compare with the LQ model, we have fitted the same cell survival data by using the second order polynomial with zero intercept for LQ formula of Ln S ð Þ ¼ ÀαD À βD 2 : We have obtained very good fitting with R 2 > 0.98 for all curves but the bending of a EQD2 for α-particle irradiation is determined by using α-particle parameters in the numerator and γ-ray parameters in the denominator of Eq. (2).
Medical Physics, 48 (4), April 2021 fitted curves beyond the measured dose range is of great concern. For example, LQ-modeling of OMM2-2 cell survival curve in Fig. 3(b) has results of R 2 = 1.00, α = 0.497 Gy −1 , β = 0.0827 Gy −2 , α/β = 6.01 Gy and EDQ2 of 65 Gy in 20 Gy SRS which significantly differs from our model result of 37.4 Gy. This is due to the low α/β ratio determined directly from the survival curve in a low-dose range. Similarly, the LQ modeling of the lung cancer cell survival curves in Fig. 4 gives results of R 2 > 0.98, α ∈ (0.07, 0.15) Gy −1 , β ∈ (0.016 0.069) Gy −2 and relatively low values of α/β ∈ (1.5, 7.2) Gy. By fitting the DU145 and CP3 cell survival curves in Fig. 5  in SRT of hundreds of acoustic neuromas 34 but the early 5 Gy × 5 scheme might have decreased EDQ2 to the benign tumors. Using the UMA model for tumor and normal tissue cell lines could provide us a theoretical analysis of the SRS/ SRT of brain tumors. In the last decade, high local control rate of early stage lung cancer treated with SBRT have been experienced by many institutions. 35,36 High EDQ2 from our UMA model and other modified LQ models 33 might provide some explanations. EQD2 for 7 Gy × 5 HDRB of the cervical and endometrial cancer and EQD2 for 7.4 Gy × 5 SBRT of the prostate cancer decrease with D o and increase slightly with n similar to the contours in Fig. 7(a). Our model calculated EQD2 for hyperfractionated radiotherapy of the tumors with inverted shoulders mainly depend on the extrapolated number n, either higher or lower than that of LQ model prediction for cells with n = 0.2 or n = 0.9, respectively. This indicates that the benefit for the selection of continuous hyperfractionated accelerated radiotherapy (CHAR) regimen depends on the characteristics of the cell survival curves. The UMA modelprediction in Fig. 7(b) tells us which tumor cell line may benefit the most from a CHART regimen. Thus, the proposed UMA model is useful in the design and evaluation of any new fractionated radiation therapy schemes.
Our UMA model, as any new models, is not matured yet. 37 The model has not included many factors, such as dose rate effect, cell repopulation and regroup, tumor regrowth with treatment delay, and synergistic effects for combination with other treatment modalities. Certainly, one can add more parameters and functions to improve the fitting of survival curves and address the radiation damage repairing, cell repopulation, and redistribution during a radiotherapy course. In fact, many cell survival curves have been measured with the effects of damage repairing, under hypoxic/aerobic conditions, contact, blood circulation changes such as pre-irradiation of the medium or capillary blood vessel damages as well as altering cell phases or multiple populations of cells during in-vivo and in-vitro observations. If there are data from split dose experiments to simulate the current SBRT or HDRB procedures, we can check the temporal effects by remodeling the survival curves with different dose schemes. A generalized LQ model 5 adds a modification function to the β parameter of LQ model to deal with dose rate or damage repairing for fractionated treatments but cannot explain the α and β changes with dose level in modeling the cell survival curves in comparison with the UMA modeling of the measured cell survival curves using the same dose rate but different dose levels. If the model works for observed cell survival curves in various conditions particularly for solid tumors under hypoxic, contacted, and/or pre-irradiation conditions during the course of radiation therapy, the prediction should be applicable to clinics with similar situations. Thus, the personalized radiation therapy scheme could be achieved with patientspecific tumor and tissue cell survival curves.
It is important to consider some mixed radiations in clinical situations such as proton therapy with tumor cells irradiated with shooting through as low LET irradiation and at the Bragg peaks as high LET irradiation. Most recently, Pfuhl et al. 38 have proposed a local effect model for prediction of cell survival irradiated by mixed radiation based on the assumption that the same spatial DNA double-strand break (DSB) distribution in the cell nucleus leads to the same effects, independent of the radiation quality. Our unified formula applicable for all types of ionization radiation over the entire dose range indicates a possible common mechanism of cell killing among various ionization radiations with their determinable MID = nLn (n)D o /(n-1)a single parameter for the description of the cell response to the radiation dose.
Better understanding of the molecular pathways of radiation-induced cell death with multi-activations is even more important in the design of combination therapy with the increasing use of novel agents in chemo or immune therapy. 12 A change of biological model from a low dose of 2-Gy fractions to a high dose of 10-Gy fractions as does by the USC model has brought considerations for the possibility of different cell killing mechanisms in radiation therapy. Clinically observed double median survival for oligometastatic diseases 39 and more than double the response rate to immunotherapy 40 based on EQD2 predictions of the LQ model have led to the hypothesis that the biology of tumor response to irradiation is different when a high dose per fraction is given. Our new model has indicated that the improved local control seen in SBRT might be the results of an EQD2 much higher than the EQD2 estimated with the LQ model and there is no need of different mechanisms among different dose ranges. We have also found that some cell survival curves for combination of radiotherapy with chemotherapy or immunotherapy as well as hyperthermia 23,25 could mostly be fitted well with the unified formula (to be presented with our other report). Such a unified formula representing the same mechanism through the entire dose range or combination of multiple treatment agents (including radiation) allows us to describe the dose response without changing of the mechanisms or dose prescription.

CONCLUSIONS
The feasibility of the proposed UMA model has been validated through the fitting of survival curves of many tumor cell lines available to us. The capability of modeling survival curves of in-vivo and in-vitro human cell lines over the entire dose range within their experimental errors provides us a new way to calculate EQD2 of SRS, SBRT, HDRB and even hyperfractionated radiation therapy courses. In comparison with the current LQ model estimations, this study has found EQD2 that is lower for intracranial SRS but higher for SBRT of lung cancer using parameters extracted from some preclinical cell survival curves. Most importantly, the unified formula has resolved the catastrophe of the traditional LQ and MT models and it theoretically indicates a common mechanism of cell killings from ionization radiation and possibly from other agents at all dose levels.