Clonogenic Survival RBE Calculations in Carbon Ion Therapy: The Importance of the Absolute Values of α and β in the Photon Dose-Response Curve and a Strategy to Mitigate Their Anticorrelation

: The computation of the relative biological effectiveness (RBE) is a fundamental step in the planning of cancer radiotherapy treatments with accelerated ions. Numerical parameters derived analyzing the dose response of the chosen cell line after irradiation to photons (i


Introduction
Cancer radiation therapy with accelerated ions was proposed as an alternative to conventional radiotherapy with X-rays due to its advantageous physical characteristics [1]. Unlike X-rays, the ion energy deposition increases with decreasing ion energy. Therefore, the energy of the ion beams can be precisely tuned so that the particles stop within the tumor volume, thus minimizing the radiation dose to the surrounding healthy tissues [2] and reducing the risk of malignant secondary effects [3]. However, ions are well known to have a denser pattern of energy deposition per unit of path (i.e., a higher linear energy transfer, LET). This higher LET is associated with markedly different biological effects [4], especially for treatments with carbon ions [5,6].
In order to account for this different relative biological effectiveness (RBE), computational models of radiation-induced effects are used during the planning of ion therapy treatments. Among all the developed models, the mixed beam model [7], the first version of the local effect model (LEM I [8]), and the modified microdosimetric kinetic model (modified MKM [9]) are currently implemented in clinics [10][11][12][13][14]. The modeled endpoint is clonogenic survival [15,16] due to its relevance for tumor control calculations [17][18][19]. Although differences are present between the models and their formalism [20,21], the key idea is to describe and predict the cell survival after high LET irradiation by knowing the response of the cell to sparsely ionizing radiation such as photons. A subset of ion-exposed data is generally used for the determination of the model parameters.
In view of the upcoming proton and carbon ion therapy center of Mayo Clinic in Jacksonville (Florida, United States of America) [22], we recently developed a model named Mayo Clinic Florida microdosimetric kinetic model (MCF MKM) [23,24]. The MCF MKM is based on the experience gathered with previous microdosimetric models [9,[25][26][27][28] and the hypothesis that the DNA damage clustering over giant loops of chromatin (subnuclear structures containing approximately 2 Mbp) has a leading role in determining the cell survival after a radiation insult [29]. The latter hypothesis is shared with different RBE models based on the amorphous track structure [30,31] instead of microdosimetry. Thanks to new strategies to determine a priori the model parameters based on the results of morphometric and karyotypic cell measurements [23,24], the RBE calculations with the MCF MKM can be performed knowing the in vitro clonogenic survival after exposure to photons only. It follows that an accurate estimation of the surviving dose response after the photon exposure is of utmost importance.
Since the in vitro clonogenic survival data (i.e., the surviving fraction S as a function of the absorbed dose D) are discrete, an interpolation is generally performed by means of the linear-quadratic model (Equation (1)) [32].
where α and β are exposure-and cell-specific fitting parameters. However, α and β are anticorrelated during the LQM fit [33]. Therefore, different couples of LQM terms can be obtained by fitting the same in vitro survival curve using different minimization algorithms and constraints or by independent users. A preliminary study for conference proceedings [34] showed that the RBE 10% of two human cell lines calculated with the MCF MKM could be significantly different if dissimilar published sets of LQM terms (obtained by different authors fitting the same photon survival curve) were used as input for the calculations. It must be noted that these different LQM sets (α and β for the photon reference radiation, named α ref and β ref here and in the following) could lead to markedly dissimilar α ref /β ref ratios for the same in vitro survival curve. For instance, α ref /β ref values of 3.5 and 6.9 Gy were obtained for photon-exposed human brain glioblastoma cells (A-172 cell line) [34]. Since the α ref /β ref ratio is used as a relevant parameter for RBE calculations with models other than the MCF MKM (i.e., proton phenomenological models [35]) and for fractionation studies [36], it is very important to assess it accurately.
As an example, the biologically effective dose (BED) for a treatment composed of 20 fractions of 2 Gy was calculated with Equation (2) [37] for human brain glioblastoma cells (A-172 cell line) using the two published LQM fits reported for the same photonirradiated survival curve (α ref /β ref respectively equal 3.5 Gy and 6.9 Gy [34]).
n is the number of fractions, d is the absorbed dose per fraction, α and β are the LQM terms.
In the case of accelerated ions, an extension of the BED formalism was proposed to account for their different RBE with respect to photons [38]. Despite being based on the same photon survival curve, the BED calculated using the results of the two different LQM fits was significantly different; namely, 62.9 Gy and 51.6 Gy for the LQM fits with α ref /β ref equal to 3.5 Gy and 6.9 Gy, respectively.
A fundamental parameter in the RBE calculations with MKMs is α 0 , namely α in the limit of LET → 0. Negative values α 0 were reported in the case of cell lines with a low α ref /β ref [24,39]. However, the physical and biological meaning of these negative values of α 0 (which might lead to the computation of unreasonable negative RBE values) is unclear, if any. It is possible that the anticorrelation between α ref and β ref is one of the causes behind the computation of these negative α 0 values.
Finally, previous studies with protons suggest that the α ref /β ref after photon irradiation can be used as an indication of the radiosensitivity of the cell line for ion irradiations, with higher RBE values expected for cell lines with a lower α ref /β ref [35]. Similarly, in the case of carbon ions, higher RBE values were observed at low doses for cell lines with a lower α ref /β ref [6,23]. However, at higher doses (i.e., for a surviving fraction of 10% or 50%), the RBE values appear to not significantly depend on the α ref /β ref [6]. In addition, despite the ion RBE values calculated with the clinically implemented LEM I [8] for a fixed α ref /β ref were reported to be nearly independent of the absolute values of α ref and β ref [40], the knowledge of the absolute values of α ref and β ref appears to be necessary [41] for 12 C ion RBE calculations with the latest version of the LEM (LEM IV [30]) and the modified MKM [9].
Considering all the above, this article aims to: (

Relative Biological Effectiveness
The relative biological effectiveness (RBE) is defined as the ratio between the absorbed dose needed to obtain the same biological effect for the reference radiation (D ref ) and the radiation under investigation (D) (Equation (3)).
In this article, the reference radiation is photons, while the radiation under investigation is carbon ions. The biological endpoint is the in vitro clonogenic survival, as described with the LQM (Equation (1)). The RBE was evaluated as a function of the surviving fraction S (RBE S ) using Equation (4) [28]. The microdosimetric kinetic model (MKM [25]) of clonogenic survival was developed to explain the enhanced killing effectiveness of ions with respect to photons by means of microdosimetry [42]. The cell nucleus was assumed to be composed of subnuclear radiation structures (radius~0.2-0.4 µm) named domains where the accumulation of lethal and sublethal radiation damages leads to cell death. The effect of different radiation types within these domains is assessed by means of microdosimetric density distributions of stochastic quantities such as the lineal energy.
In the MKM [25], the linear term of clonogenic survival (α) represents the effect of lethal lesions, and it is calculated with Equation (5): where α 0 is the cell-specific value of α in the limit of y → 0 where y is the lineal energy, β ref is the quadratic term of the LQM model for the reference photon exposure, y D is the dose-mean lineal energy assessed in a domain with density ρ and radius r d . The dose-mean lineal energy y D is computed as in Equation (6) y D = y d(y) dy (6) where d(y) is the dose probability density of the lineal energy y [42].
On the other hand, the quadratic term of clonogenic survival (β) is assumed to be radiation-independent and equal to that for the photon exposure (Equation (7) [25]), i.e.,: Once α ref and β ref are determined with Equation (1) by fitting the clonogenic survival curve after photon exposure, the two MKM parameters (α 0 and r d ) are generally assessed by fitting the results of ion-exposed in vitro experiments with Equation (5) [43].
Corrections to Equations (5) and (7) were afterward introduced to improve the agreement with the in vitro data by explicitly considering the overkill effects at high linear energy transfer (LET) and the radiation dependence of β [9,27,44,45]. However, this requires the determination of additional parameters, such as, for instance, the saturation parameter y 0 in the modified MKM [9]. Nonetheless, it is worth mentioning that the results of these newer MKMs [25,27] were found to be equivalent to that of the original MKM [25] for LET <~30 keV/µm [46], where Equations (5) and (7) are still valid.
where c(y) (Equation (10), adapted from [27]) is a correction factor accounting for the non-Poisson distribution of lethal lesions.
y is the lineal energy, d(y) is the dose probability density of the lineal energy, α 0 and β 0 are the LQM terms in the limit of y → 0, R n is the mean radius of the cell nucleus, r d is the mean radius of the subnuclear domains, and ρ is the density (=1 g/cm 3 ).

MCF MKM Parameters
Both MCF MKM parameters (R n and r d ) are assessed a priori without using clonogenic survival in vitro data [23,24]. As first choice, R n can be determined by means of morphometric measurements for (unattached) cells presenting a spherical nucleus. However, this information is not always available. If the cross-sectional area A of the cell nucleus is reported, an estimate of the mean radius of the cell nucleus can be obtained with Equation (11) [23,24].
If also A is unavailable, an approximate value of R n can be obtained using an empirical correlation between R n and the mean DNA content of an asynchronized population Γ [Gbp] (Equation (12) [23,24]) The mean DNA content of an asynchronized population Γ is computed with Equation (13) [23,24] Γ = γ p ξ (13) where γ is the normal DNA content for one set of chromosomes (3050 Mbp for humans, 2750 Mbp for rats, 2700 Mbp for Chinese hamsters, and 2650 Mbp for mice [47]), p is the ploidy number, and ξ is a factor accounting for the cell cycle distribution of the irradiated population (ξ = 4/3 for asynchronized cells [23]). The ploidy number p is equal to 2 for healthy cell lines and for diploid cancer cells. Though mutated chromosomes might contain an abnormal amount of DNA, an approximated value of p for aneuploid cancer cells can be calculated with Equation (14) [23,24].
where x is the mean number of chromosomes in the aneuploid cell population and x n is the number of chromosomes in a normal set (23 for humans, 21 for rats, 11 for Chinese hamsters, and 20 for mice [47]). In the MCF MKM, it is hypothesized that the subnuclear domains are giant loops of chromatin, subnuclear structures containing approximately 2 Mbp of DNA [29,48,49] where the accumulation of DNA damage leads to cell inactivation. Therefore, the mean radius of the subnuclear domains r d is computed with Equation (15) [23,24] where R n is the mean radius of the cell nucleus, λ is the mean amount of DNA in a chromatin substructure (2 Mbp), and Γ is the mean DNA content of the irradiated cell population. Finally, since R n and r d are known, and the correction factor c(y) (Equation (10)) is~1 for the reference photon exposure [23], then The following quantities were calculated and compared as a function of the LET for 12 C ions with energy between 1 and 1000 MeV/n: α/α ref , β/β ref , RBE 50% , RBE 10% , and RBE 1% . The ratio α/α ref is known as low dose RBE (RBE α ) and can be derived from Equation (4) in the limit of S → 1. The reference radiation chosen for the RBE calculations was 6 MV X-rays (y D,re f~2 .3 keV/µm [23]).  [33,50] (see Figure 1). The biological endpoints, the reference radiation, and the numerical values of R n , r d , and y D,re f were the same of Section 2.3.

MKM-Based Fit of the Photon Survival Curve
As mentioned before, the anticorrelation between α and β [54] during the LQM fit of the photon clonogenic survival data could produce very different sets of α ref and β ref . Since no physical nor biological constraints are generally applied during the LQM fits, null and negative values of both α ref and β ref were reported in literature [33]. This could immediately result in negative or infinite RBE values at low doses (RBE α ).
Additionally, if α re f < β re f y D,re f ρ π r 2 d then negative α 0 are calculated with Equation (16).
Therefore, if the dose-mean lineal energy for the radiation under investigation (y D ) is lower than that for the reference radiation y D,re f , a negative α value is computed and thus a negative RBE α .
To prevent the occurrence of this effect and to possibly mitigate the anticorrelation between the two LQM terms during the fit of the in vitro clonogenic survival curve for the photon exposure, we propose that the latter is fitted with Equation (18) under the condition of α 0 ≥ 0. Alternatively, the in vitro survival curve for the photon exposure can be fitted with the classic expression of the LQM (Equation (1)) in combination with the following numerically equivalent constraint (Equation (19)): As an example, for the synthetic cell line of Sections 2.3 and 2.4 (r d = 0.29 µm, We performed these MKM-based LQM fits using the Curve Fitting Tool of Matlab R2021 (The MathWorks Inc., Natick, MA, USA) based the "NonlinearLeastSquares" method and the "Trust-Region" algorithm 2.6. Effect of Different Fits of the Same Photon Survival Curve on the Calculated RBE for 12

C Ions
Representative in vitro datasets from the Particle Irradiation Data Ensemble (PIDE) version 3.2 [33,53] were chosen to test the effect of different fits of the same in vitro photon survival curve (α ref and β ref ) on the RBE calculated with the MCF MKM for 12 C ions. PIDE is a database of in vitro clonogenic survival data for human and rodent cell lines exposed to photons and ions up to 238 U. The following information is reported for each entry in PIDE: cell line details (name, human or rodent, tumor or healthy, synchronized or asynchronized), ion exposure details (ion energy, ion LET, monoenergetic or spread-out Bragg peak (SOBP)), the reference photons used (i.e., 200 kVp X-rays or 60 Co γ-rays), and the results of the LQM fits (α, β, α ref , and β ref ). When available, two sets of LQM terms are included for each entry in PIDE: the results of the LQM fits performed by the authors of the original publication and the results of LQM fits performed by the PIDE team. As performed in previous studies [23,24,46], we initially discarded entries for synchronized cells, irradiations along SOBPs, and exposures to ions with energy lower than 1 MeV/n. In the second step, we included only PIDE entries with non-negative values of β for the ion exposure. Negative β could be due to the presence of a more radioresistant cell subpopulation [53].
As a relevant example of datasets leading to negative values of α 0 , we selected the cell lines with most entries (for a specific ion) for which α ref was equal to 0 in one of the two LQM sets (i.e., fit by the original author of the publication or by the PIDE team) and greater than 0 in the other photon LQM set: human astrocytoma cells (U-251MG cell line) and human mammary epithelial cells (M/10 cell line). In both cases, these two cell lines were irradiated with 12 C ions. It must be stated that the RBE for the M/10 and the U-251MG cell lines was computed in previous work with the MCF MKM [24] using only the photon LQM terms with α ref = 0. These published RBE results were compared with the novel RBE calculations using the MCF MKM in combination with the photon LQM terms with α ref = 0 and the new LQM terms obtained fitting the in vitro photon survival curve with the MKM-based approach of Section 2.7. The a priori assessed model parameters (r d and R n ) for both cell lines were extracted from [24] and listed in Table 1 together with the model parameters derived analyzing the different LQM fits of the photon in vitro survival curve (α ref , β ref , α 0 , β 0 ). Additionally, in vitro data for human brain glioblastoma cells (A-172 cell line) irradiated with 12 C ions were also included since marked differences were observed between the photon LQM terms (α ref and β ref ) obtained by the original authors of the publication and by the PIDE team [34]. Different RBE 10% values were computed when these two sets of photon-exposed LQM terms were used in combination with the MCF MKM [34]. In this work, we use these computed RBE 10% s as a benchmark for the novel RBE calculation based on the α ref and β ref assessed with the MKM-based fit of Section 2.5. Furthermore, we extend the analysis to other RBE levels (RBE α , RBE 50% , and RBE 1% ). The reference photon radiation was 137 Cs γ-rays for the A-172 and U-251MG cell lines from [55], while 137 Cs or 60 Co γ-rays (not specified) for the M/10 cell line data from [56].

Computer Simulations
The computer simulations were carried out using the Particle and Heavy Ion Transport code System (PHITS [57]) version 3.24, following the methodology described in detail in previous works [23,24]. Therefore, only a summary is given here. The calculations were performed over an infinitesimal layer of water (no slowing down and no nuclear reactions) for monoenergetic 12 C ion beams with energy between 1 and 1000 MeV/n.
The microdosimetric PHITS function [58] implemented in the PHITS [T-SED] tally was used to analytically compute the lineal energy distributions for a spherical liquid water target with radius equal to 0.30 µm. The minimum energy deposition considered in the microdosimetric calculations was 10.9 eV (i.e., one event of one ionization [58]). The results of the PHITS microdosimetric function were validated against experimental measurements with gas detectors [59][60][61]. Furthermore, the PHITS microdosimetric function was previously used to successfully model the response of LiF:Mg,Ti and LiF:Mg,Cu,P thermoluminescent detectors [62][63][64][65], Al 2 O 3 :C optically stimulated luminescent detectors [66], and an optical-fiber-based BaFBr:Eu detector [67]. The unrestricted linear energy transfer (LET) in water was calculated with the PHITS [T-LET] tally. All results of this study are plotted as a function of a generic "unrestricted LET in water" because of the lacking details and the contrasting approaches used in the LET assessment between the authors of the different studies presenting the in vitro data. A previous study reported a good agreement between the PHITS-calculated LET values and corresponding literature results for 12 C ions [68]. More details on the simulations can be found in [23,24] and in the PHITS manual (https://phits.jaea.go.jp/manual/manualE-phits.pdf, accessed on 24 January 2023). These synthetic survival curves were used to determine the MCF MKM parameters α 0 and β 0 using Equations (16) and (17). As described in the methodology, R n and r d were fixed at 4.7 and 0.29 µm. Since the biological weighting function α(y) (Equation (8)) is the most important factor in the RBE calculations with the MCF MKM, Figure 2B

Effect of the Absolute Values of α ref and β ref on the Calculated RBE for 12 C Ions
The

Effect of Different Fits of the Same Photon Survival Curve on the Calculated RBE for 12 C Ions
The effect of the different LQM fits on the MCF MKM calculations for carbon ions is shown in Figures 5-7 [55,56] and by the subsequent reanalysis performed by the PIDE team (open red diamonds) [53] on the same raw data. All the in vitro results are reported without error bars since the PIDE database, including the results of the analyses by original authors and by the PIDE team, currently does not provide uncertainty intervals for the clonogenic survival data. As previously discussed [24], this is likely due to the fact that most published biological papers do not include a systematic uncertainty analysis nor include sufficient information to perform a retrospective uncertainty study.
As shown in Figure 5A,B, the description of the photon survival curve of the A-172 cell line with the MKM-based fit is very similar to that from the original article [55]. On the other hand, the fit by the PIDE team appears to better describe the data point at 6 Gy, but it reproduces the in vitro data at low doses less accurately ( Figure 5B). The LQM terms after the carbon exposures ( Figure 5C  The results of the different LQM fits for the photon-exposed M/10 cell line are plotted in Figure 6A,B. Relevant differences are present between the three fits, especially for that by the PIDE team (dashed red lines). The fit by the PIDE team appears to well reproduce the in vitro data points at 4 and 6 Gy, but it significantly overestimates the in vitro data below 4 Gy ( Figure 6B). As a consequence of this unreliable description of the photon survival curve at relatively low doses, the experimental in vitro RBE 50% and RBE 10% results are not included in Figure 6E,F. Except for the carbon-irradiated data point at~300 keV/µm (there is a minor underestimation), the MCF MKM results agree well with the in vitro data ( Figure 6C,G) when the photon LQM terms from the original article [56] and the MKM-based fit are used as input to the MCF MKM calculations. By contrast, the in silico calculations based on the photon LQM terms obtained by the PIDE team (dashed red lines in Figure 6C,G) overestimate the corresponding in vitro results.
As shown in Figure 7A,B (U-251MG cell line), the LQM fit from the original article and that by the PIDE team are quite similar. These two fits seem to well describe the in vitro data (black diamonds) over the whole dose range (0-8 Gy). On the other hand, the MKM-based fit appears to overestimate the in vitro photon data at 6 and 8 Gy. The LQM terms from these three fits of the photon curve were used as input to the MCF MKM to compute the cell response to carbon ions ( Figure 7C,G). Though a reasonable agreement between the in silico and the in vitro results is present at 20 and 40 keV/µm, the model calculation appears to overestimate the experimental data at~80 and~100 keV/µm. The underestimation is the largest for the MCF MKM calculations based on the photon LQM determined by the PIDE team, while it is the smallest in the case of the MKM-based fit.

Discussion
In order to perform RBE calculations with the MCF MKM, four parameters must be known together with the microdosimetric spectra: a geometrical measure of both the cell nucleus (R n ) and the subnuclear domains (r d ) and the numerical values of the LQM terms in the limit of zero lineal energy (α 0 and β 0 ). While the first two parameters (R n and r d ) can be assessed by means of morphometric measurements of the cell nucleus and estimates of the mean DNA content, α 0 and β 0 are derived from the in vitro LQM terms of the photon dose-response of the cell line under investigation (α ref and β ref ) [23,24]. Since the two LQM terms are anticorrelated during the LQM fit of the photon dose response [33], different sets of α ref and β ref might be obtained from the same in vitro survival curve. This affects the results of the following calculations for photons and the computation of the RBE of ion therapy treatments.
In the first case (photons), different values of the biologically effective dose (BED, Equation (2) [69,70] can be affected by the different LQM fits. This happens because these MKMbased repair models [69,70] rely on the application of a temporal correction factor to β, while α is left unchanged. The anticorrelation between the two LQM terms could therefore result in different sets of LQM terms which lead to a different estimation of the temporal dependence (i.e., fractionation, dose rate, beam interruption) of the radiationinduced effects.
Secondly, since α ref and β ref are used as input to the RBE calculations (i.e., via Equations (16) and (17)  In the framework of implementing a variable RBE in proton therapy practice, many phenomenological RBE-vs-LET models were developed. A numerical comparison between the results of these models revealed large differences between the results of these LET-based proton RBE models [35]. As concluded in a subsequent review of these phenomenological models [71], these differences could be attributed to the significant disagreement between the role of α  Table 1) and therefore negative α 0 values with Equation (16). As a result, the MCF MKM α(y) functions computed with Equation (8) could assume negative values, as shown in Figure 2B. This might lead to the calculation of unreasonable negative RBE α values for sparsely ionizing radiation. Though α 0 = 0 denotes that all lethal lesions arise from the combination of sublethal lesions [72], the biological meaning of negative α 0 values is, at best, unclear.
In this work, we presented and used a set of novel MKM-based constraints (Equations (18) and (19)) for the LQM fit of photon-exposed in vitro data (panels A and B of Figures 5-7). We do not claim that this MKM-based fit for the reference photon radiation is universally valid nor mathematically the best one. On the other hand, we introduced this method as a strategy to mitigate the occurrence of negative α 0 values, to avoid the calculation of negative RBE values for ions whose dose-mean lineal energy is smaller than that of the reference photons, and to prevent the calculations of infinite RBE α for cell lines whose LQM fit of the reference photon exposure resulted in α ref = 0 (Table 1). When the photon LQM terms (α ref and β ref ) obtained with this MKM-based fit were used as input to the MCF MKM calculations, the predicted cell response to carbons (solid blue lines in panels C to G of Figures 5-7) were generally found to be in reasonable agreement with the in vitro data. The results obtained using the other two LQM fits of the photon survival curve (short-dashed red lines and dashed black lines in panels C to G of Figures 5-7) were found to describe the in vitro data in a comparable or worse way (i.e., some large deviations were observed especially in case of the photon fits with α ref = 0).
Finally, as previously discussed, we benchmarked the RBE results obtained with the MCF MKM in combination with the new MKM-based photon LQM fit against in vitro RBE data for carbon ions only. This was performed because carbon ions are the focus of this article and since the most interesting in vitro datasets were found for this ion (see Section 2.6). Nonetheless, this approach can be used in combination with in vitro experiments dealing with other ions. As an example, we applied this method to published in vitro data of proton-irradiated human glioblastoma cells (U-87 cell line [73]; see Figure Table S1 of the Supplementary Materials. In conclusion, this study indicates that the knowledge of both the absolute values of α ref and β ref and their ratio α ref /β ref appears to be necessary for an accurate calculation of the RBE in carbon ion therapy. Care must be taken in fitting the reference in vitro photon survival curve. It is suggested to pay particular attention to the relatively-low dose range since α 0 and β 0 represent the LQM parameters in the limit of very sparsely ionizing radiation and low dose. However, at low doses (i.e., below~1-2 Gy), phenomena such as low-dose hypersensitivity [74], likely due to intercellular signaling (bystander effects) [75], could further complicate the fitting process. Furthermore, it should be remembered that LQM fits over different dose ranges could produce different sets of LQM terms [76]. This is particularly evident at relatively high doses where the in vitro clonogenic survival curves transition to a fixed slope [77]. The MKM-based fit of the in vitro photon survival curve appears to be a useful tool to determine