Regression fitting megavoltage depth dose curves to determine material relative electron density in radiotherapy

The relative electron density (RED) parameter is ubiquitous throughout radiotherapy for clinical dosimetry and treatment planning purposes as it provides a more accurate description of the relevant radiological properties over mass density alone. RED is in practice determined indirectly from calibrated CT Hounsfield Units (HU). While CT images provide useful 3D information, the spectral differences between CT and clinical LINAC beams may impact the validity of the CT-ED calibration, especially in the context of novel tissue-mimicking materials where deviations from biologically typical atomic number to atomic weight ratios 〈Z/A〉 occur and/or high-Z materials are present. A theoretical basis for determining material properties directly in a clinical beam spectrum via an electron-density equivalent pathlength (eEPL) method has been previously established. An experimental implementation of this approach is introduced whereby material-specific measured percentage depth dose curves (PDDs) are regressed to a PDD measured in a reference material (water), providing an inference of 〈Z/A〉, which when combined with the physical density provides a determination of RED. This method is validated over a range of tissue-mimicking materials and compared against the standard CT output, as well as compositional information obtained from the manufacturer's specifications. The measured PDD regression method shows consistent results against both manufacturer-provided and CT-derived values between 0.9 and 1.15 RED. Outside of this soft-tissue range a trend was observed whereby the 〈Z/A〉 determined becomes unrealistic indicating the method is no longer reporting RED alone and the assumptions around the eEPL model are constrained. Within the soft-tissue RED range of validity, the regression method provides a practical and robust characterisation for unknown materials in the clinical setting and may be used to improve on the CT derived RED where high Z material components are suspected.


Introduction
A core challenge of clinical radiation therapy (RT) is measuring and predicting the dose deposition in various tissues and materials, and from different beam modalities.This is exacerbated by the advent of 3D printing, enabling an ever-expanding range of novel materials [1].Both reference dosimetry protocols and treatment planning systems (TPS) are dedicated to ensuring reference conditions or managing deviations therein.However, these have limitations in both theoretical and empirical underpinnings.For instance, while IAEA TRS-398 provides correction factors for non-reference phantom materials, the guidance is limited to 'water equivalent' phantoms, and many other materials that might be encountered in the clinic e.g., biological tissues, are not included [2].
For novel materials, including tissue mimicking materials (TMM) used in anthropomorphic phantoms, the characterisation of attenuation properties forms a key part of clinical end-to-end testing [3,4].ICRU 44 recommends that material composition is reported [5].This has been reinforced in recent reviews [1].A recent trend to obtain material composition through scanning electron microscopes is a welcome step forward towards reducing uncertainties [6,7] however this kind of direct assessment is rarely available in the clinic.Various techniques to infer composition based on radiation attenuation or CT-derived measurements have been developed [3,4,[8][9][10][11][12][13].Early dose correction methods involved influence factors based only on mass density, which have validity where the ratio of atomic number to atomic weight (Z/A) is approximately constant, covering the elements found in the biological tissues commonly encountered in radiation therapy [14,15].An issue arises for some biological tissues which significantly deviate from this ratio e.g., adipose and cortical bone, as well as high Z materials where this ratio trends down with higher Z.
Methods limited to physical density and where the variation in 〈Z/A〉 is omitted can result in significant attenuation differences for high Z materials [14].Use of electron density, and its water-normalised counterpart, RED, better accommodate materials across a range of Z/A [14].Given the predominance of Compton interactions from typical MV therapy beams in common tissues encountered in RT, the use of RED is generally appropriate for dose calculation and comparing the dosimetric properties of materials [1].This assumption generally holds despite the spectral differences between CT and LINAC.However, significant differences in CT-derived RED can result from the presence of high-Z materials [11].This presents a risk when novel materials of unknown composition are introduced, which is common with patient prosthetics, bolus, and the plastics associated with the introduction of 3D printing techniques [9,16].
In characterising non-reference materials for dose prediction purposes, the radiation transport of both primary and secondary particles in the energy range of interest should be considered [17].Methods to characterize materials directly in the clinical beam spectrum have been developed and reported [3,4,7,9,10,11,12,13].Single point measurements can be used to deduce a simple scaling factor RED eff , providing quick and easy attenuation coefficients [11,18].However, this approach has limited applicability beyond the specific assembly material and geometry measured.Others have compared PDDs measured in TMMs to that measured in water or calculated in the TPS to evaluate water equivalence and identify likely dose errors [9,10].These methods are readily implemented in the clinic given the ease of acquiring PDD data, and clinically relevant dose metrics such as the depth of maximum dose (D max ), the depth where the dose is 90% of the maximum dose (R 90 ), direct dosedifferences, as well as other qualitative methods which can be used to estimate the dose error or corrections needed [9,10].PDDs have also been used via a regression method to evaluate the water equivalence of 3D printed materials [19].However, these methods did not result in broader radiological characterisation in terms of material RED.
A theoretical electron equivalent path length (eEPL) method was introduced to model the primary beam Compton interactions up to 10MV, thus factoring the variation in 〈Z/A〉 between materials [14].This provides a framework for deriving RED from a clinical beam under narrow beam geometry.Broad beam geometry can be catered for using the O'Connor theorem, which relates the dose in two media of different physical densities and equal atomic composition irradiated by the same beam.In practice this can be applied as a field-size scaling factor [15].
Beyond single point measurement techniques, empirical implementation of the eEPL method for material characterization has not been reported.A method for characterising a material via regression of measured PDDs based on the work of Karl, Seco and O'Connor is presented [14,15,19].The utilisation of non-linear regression analysis leads to a robust characterisation of attenuation properties for materials directly in the clinical beam spectrum.The method is applied to measured PDDs of commonly available radiotherapy TMMs and the results compared to manufacturerprovided and CT-derived RED values.

Basic material characteristics and arrangement
Slabs of Plastic Water and tissue-mimicking materials: Adipose, Lung (Inhale and Exhale), and Bone (Trabecular and Cortical) were acquired from CIRS (Norfolk, VA USA).Slabs of High Equivalency Solid Water (HE SW), LN 300 Lung, IB3 inner bone and SB cortical bone were acquired from Gammex (Wisconsin, USA).The CIRS tissue mimicking material dimensions were 20 cm wide by 20 cm long.The Plastic Water and Gammex tissue mimicking materials were 30 cm wide by 30 cm long.
The mass of each slab was determined with a set of digital scales and the dimensions were measured with a digital vernier calliper, arriving at a physical density (g/cm 3 ) for each slab.

RED from manufacturers data
The RED values supplied by the manufacturer are noted in Table 1.
In addition, RED values derived from the manufacturer's reported material composition were determined using the following methodology.The relative electron density of a medium is defined as where emed is the electron density of the medium, and eH 2 O is the electron density of water.Now where med is the mass density of the medium and ⟨Z∕A⟩ med is the ratio of effective atomic number to atomic weight of the medium: where for each element i of the medium, f i is the fractional weight of the element and (Z∕A) i is the ratio of the atomic number to the atomic mass [14].
Except for HE and Plastic Water, manufacturers did not specify tolerance or uncertainties with their stated RED values.Clarification sought from the manufacturers was not forthcoming.For the HE and plastic water the stated uncertainties were ± 0.2% (RED) and ± 0.5% (dosimetrically to water), respectively.No coverage factors were specified.For all the other materials, uncertainties were inferred under the assumption that these would meet ICRU 44 specifications [5], i.e., ± 1%.Coverage factors of k = 2 were attributed to these uncertainty values.For compositional analysis, the manufacturers state the materials meet ICRU 44 [5].Here the coverage factors were assumed as k = 2. (2

RED from CT imaging
CT-image derived RED values for the CIRS materials were acquired with a Siemens SOMATOM Confidence CT (Siemens, Munich, Germany) using the "RT Chest" CT protocol (for imaging parameters see Table 2).The CIRS target materials were sandwiched between two slabs of 2 cm thick solid water, one on top and one in between the TMM slabs and the surface of the treatment couch.CT-image derived RED values for the Gammex materials were acquired with a Siemens Biograph mCT (Siemens, Munich, Germany) using the "Pelvis" protocol.The imaging parameters are listed in Table 2.
The DICOM CT images were imported into ImageJ (U. S. National Institutes of Health, Bethesda, Maryland, USA) for CT number sampling.A rectangular region of interest was defined on the central CT slice.Partial volume effect and higher density shell artefacts were avoided by excluding pixels within 1cm of the edge of the slab material volume [20].Mean and standard deviation (SD) CT numbers for each slab were sampled.CT-ED conversion for the CIRS materials was established by the relationship from 23 previous assessments of an RMI phantom (Gammex Inc, Wisconsin, USA).For the Gammex materials a single assessment of the Advanced Electron Density Phantom (Gammex) was used.
An estimate of uncertainty for the CT-derived RED values was established by sampling the variance of Hounsfield Units (HU) obtained from 23 measurements on the RMI phantom [21], with a coverage factor k = 2.The plot of CT-RED was used to identify the linear portions of the relationship, noting a small discontinuity around RED = 1 (near water materials).The second standard deviations of sampled CT-numbers were plotted, and a linear fit determined for upper and lower limits, on each side of the discontinuity.The linear fitting equations were used to interpolate from the CT-number standard deviation to RED standard deviation, and a range was taken between the upper and lower values, whereby relative U (k = 2) could be derived for inter-scan uncertainty.Additionally, each image was sampled using the region of interest (ROI) tool in ImageJ software; resulting in a sample size, mean, and standard deviation, providing an estimate of intra-scan uncertainty.Finally, the tolerance supplied by the manufacturer was assumed as a representation of U (k = 2) for the stated RED.These three components of uncertainty were added in quadrature, as an estimate of total expanded uncertainty, k = 2, inherent to RED derived from CT images for the various materials.

RED from measured PDDs
PDDs for each material were measured over a practical physical range (see Table 1) using a PPC40 plane-parallel ion chamber with an active volume of 0.4 cm 3 (IBA Dosimetry GmbH, Schwarzenbruck, Germany).The chamber was housed in a 20 mm thick slab of water equivalent plastic with a compatible chamber insert (Fig. 1).Accumulated charge was measured under various arrangements of the material slabs to achieve depths between 1 and 225 mm.The 1 mm PMMA chamber window was accounted for.The plastic water chamber slab was placed on an additional 10 cm thick plastic water slab to provide back scatter.Clinical LINACs Elekta Versa HD (Elekta, Stockholm, Sweden) and Varian Truebeam (Varian, Palo Alto, USA) were used to deliver 200 MU from a 6 MV photon beam collimated to 10 cm × 10 cm impinging the slab arrangement for each PDD measurement (see Table 2).An SSD of either 90 cm or 100 cm was maintained using the optical distance indicator cross-checked by the digital readout on the treatment couch for the Elekta and the calibrated front pointer on the Varian.Dose-difference between the measured TMM PDD points and the water reference data at each depth were regressed via least squares fit to determine the summed squared residuals.The Generalized Reduced Gradient (GRG) algorithm, developed by Frontline Systems Inc and implemented as part of the Solver plug-in for Microsoft Excel (Microsoft, Redmond, Washington, United States), was used to optimise over key parameters to minimise the sum of the squared residuals between datasets.The reference data was interpolated using a linear interpolation between datapoints.The coarser measurement dataset was regressed to the higher resolution reference dataset to minimise error from interpolation between neighbouring points.
The optimisation results in three factors, , and , in fitting the TMM PDD (PDD′) to the water PDD (PDD), where = ⟨Z∕A⟩ med of the TMM, is a depth scaling factor; is the dose normalisation scaling factor; is a depth shift, ( 4) to account for relative translations of the dose distributions; and z is the depth determined from the electron-density equivalent pathlength.The field size of the reference water PDD used in Eq. ( 4) was determined by scaling the mass density of the TMM to water.Now, where l e i is the physical ray line length through slab i of the medium and RED i is the relative electron density of slab i. Substituting Eqs. ( 1) and (2) into Eq.( 5), and assuming each slab has elemental composition ⟨Z∕A⟩ med but can vary in density, med i , gives The average RED for each TMM is subsequently determined via An indication of uncertainty for each RED value resulting from the regression method was determined from an assessment of the components of uncertainty in the measured PDD inputs using confidence intervals calculated using Fisher's F distribution [22].A plot of the sum of the square of the residuals (SSR) was generated by varying around its optimised value while still optimising and .A threshold SSR value is then determined based on Fisher's F distribution taking into consideration the degrees of freedom of each measurement, the required confidence, and the minimum SSR.Using the threshold value, the upper and lower limits of the ⟨Z∕A⟩ med confidence interval are determined.Confidence intervals for all material REDs were converted to standard deviations, whereby an expanded U (k = 2) were derived and added in quadrature with other components of uncertainty relating to measurement of the physical densities.

Analysis of variance
The results across the materials from each method were assessed via a one-way analysis of variance (ANOVA).This was performed using the summary data for each-sample size, mean, and standard deviation, with an assumption of normality underlying all distributions around the mean REDs.Where not directly measured, standard deviations were inferred from the total expanded U (k = 2).For the manufacturer's data, the sample size was set arbitrarily to 100 indicating high confidence in the mean value.For the CT data the sample size was set as the number of pixels contained within the ROI on the central slice.For the measured PDD data the sample size was the number of points for each material depth dose scan used in the regression.The desired confidence level was set to 95%.P-values and significance levels via a Tukey-Kramer HSD post-hoc test were reported.

Results
The results required to determine material RED and estimate uncertainties using the PDD method are presented and then summarised with the other methods.
For the measured-PDD method, Fig. 2 presents two examples (a) Gammex LN 300 Lung and (b) SB Cortical Bone fitted to the corresponding water PDD.Here, the resulting depth scaling factor = ⟨Z∕A⟩ med is 0.4724 and 0.8322 respectively (cf.Gammex: 0.5356 and 0.5145 respectively); , the dose normalisation point scaling, is − 0.2% and − 0.1% different from the maximum charge value measured respectively; , the depth shift, is -0.5 mm and 0.2 mm respectively.Differences in and were less than 0.5% and 1 mm for all materials studied.
The SSR difference is shown for a range of ⟨Z∕A⟩ med for the examples of Gammex LN 300 Lung in Fig. 3(a) and SB Cortical Bone in Fig. 3(b).
Figure 4. compares the resultant 〈Z/A〉 based on the manufacturer's stated RED and the PDD fitting method for all 10 materials, including uncertainty bars that indicate 95% confidence interval The compositional analysis and the accompanying uncertainty are reported for the full range of materials in Table 3.
Figure 5 presents a comprehensive summary of the results for all four methods used to determine REDs: the regressionfitted PDD method, measured CT numbers, manufacturerstated values, and those calculated based on the composition.

Discussion
The application of non-linear regression provides a mathematically robust way of comparing the PDDs measured in various materials to that of the reference data (water), and a novel method of estimating material RED via the theoretical formalism of Seco and O'Connor.
Generally, manufacturers don't outline the methods used to determine the RED stated for TMMs.However, where they provide the material composition the RED can be derived via the effective atomic weight ratio and the mass density.In most cases we saw agreement around 0.5% between the stated RED and compositionally derived values, however with the exception of HE solid water, tolerances associated with the stated values were not provided.For comparison with other methods, we used the stated values, though further information from the vendors is desirable.
For the two water equivalent materials, the fitted PDD method agreed with manufacturer-stated or compositionally derived RED.For Gammex High Equivalence solid water, the PDD method RED was 1.01 ± 0.01, which was not statistically significantly different (p value 0.25) from the manufacturer's RED of 1.0000 ± 0.0002.For CIRS Plastic Water, the PDD method RED was 0.99 ± 0.01, The solid red horizontal line is the threshold value, determined by the Fisher distribution.The 〈Z/A〉 at which the parabola intersects the threshold value indicates the confidence interval which was not significantly different (p value 0.30) from the manufacturer's RED of 0.9991 ± 0.0005.At lower densities, the RED stated by the manufacturer for CIRS Adipose was 0.934 ± 0.001, which was not significantly different (p value 0.12) to the value determined using the measured PDD method (0.94 ± 0.03).This was found to be the lower limit of agreement, as differences between these methods for the CIRS and Gammex lung materials (RED approx.0.2-0.5)were statistically significant (p values < 0.00001, for Gammex LN300 lung, and CIRS Inhale and Exhale lung).
For higher densities, the RED stated by the manufacturer for Gammex IB3 inner bone was 1.090 ± 0.001, which was not significantly different to the value determined using the measured PDD method of 1.08 ± 0.01 (p value 0.47).The RED stated by the manufacturer for CIRS Trabecular bone TMM was 1.156 ± 0.001, which was not significantly different to the value determined using the PDD method of 1.16 ± 0.04 (p value 0.99998).This was found to be the upper limit of agreement and statistically significant differences were observed for the CIRS and Gammex cortical bone (RED 1.7-1.8)(p values < 0.00001, for Gammex SB Cortical bone and CIRS Cortical bone).
Overall, the range of RED agreement between the manufacturer stated RED and measured PDD methods was 0.93 and 1.15, which covers the range of soft tissues defined by ICRU 44.Previous work noted a divergence between CT and single-point measurement methods at higher RED, and suggested high-Z additives as the cause.
A parallel was observed with previous work where RED eff derived from single point measurements and that from CT diverged outside the range of 1.7 and 0.5 RED [11].As only the Compton component is considered for the eEPL and the RED eff methods, whereas CT mainly based on photoelec- tric absorption, the presence of high Z additives may explain the discrepancies observed [11].
The range of applicably may relate to a limitation of Seco's initial work which only considers primary beam attenuation.Field size scaling provides a partial correction by factoring in secondary scatter via the consideration of mass density ratios, nevertheless, the results at the more extreme RED values were significantly different.Investigation into the derivation of RED from the eEPL model embodied in the measured PDD method showed that the determined 〈Z/A〉 became physically unreasonable for RED outside this range, as evident in Fig. 4.

Table 3
Atomic composition to 1 D.P. of TMMs (supplied by manufacturers), ordered from lowest to highest mass density value, ρ The effective atomic number to atomic weight ratio of each material is ⟨Z∕A⟩ i as defined by Eq. ( 3).Chemical elements presented in the table are identified by their chemical symbol.Electron densities are given by Eq. ( 2) and RED are given by Eq. (  A distinct pattern was observed where the uncertainty increased for RED below 1.0, with the most pronounced result in LN300 (see Fig. 4).This trend, however, was not reciprocated for RED above 1.0.The primary source of this pattern appears to be the inherent uncertainty in the regression itself.This is demonstrated in the SSR plots for LN300 where the broader minima, indicates a higher uncertainty in the resulting regression (see Fig. 3).Nevertheless, the minima, although broad, remained well-defined.This suggests that despite the increased uncertainty, the estimation based on the regression remains valid.
The inherent limitations of the FS scaling method also contribute to the overall uncertainty in material RED.While the method partially adjusted for LN300's discrepancy by shifting Dmax deeper into the tissue, it failed to provide full correction.Understanding the FS that correlates with the measured data for LN300 may be informative in discerning the underlying causes of this limitation.Tailoring the approach to address the unique properties of LN300 could yield a more precise model of beam behaviour across various materials.
Seco suggested that additional factors could be added to the eEPL model to take into consideration the higher energy incident beams and the Z of the materials [14].This accounts for Pair production and improves the modelling of Compton and Pair production interactions, overall improving the accuracy of correction to material heterogeneities under the primary beam across a wider range of energies.
However, this extension would introduce significant complications relating to computational requirements and sample size for the measured PDD regression method.Secondary scatter contains a spectrum of energies based on their initial interactions in the primary beam.The additional factors proposed would still not attend to the energy spectrum of the secondary scatter.
The agreement between the CT-derived RED values and the manufacturer stated RED was generally good in the mid to low density range, however CIRS cortical bone, CIRS Plastic Water, and IB3 Inner bone all reported statistically significant differences (p values < 0.00001).This result was reflected for the manufacturer-stated RED, which is used to determine the CT to ED curve.None of these materials are considered high Z materials where the spectral properties of CT may arise, and the differences between material manufacturers is unlikely to contribute, as the HU results for these lay well within the uncertainty range determined on the CT calibration curve.However, the result is expected for CIRS Plastic Water where the "Original" Plastic Water was not designed for water equivalence below 150 keV.The reported significance of the differences here may be a function of the arbitrary components of uncertainty attributed to the manufacturer's stated and compositional data.Further information from the manufacturers on these fronts is required to investigate further.For all RED methods, an attempt to infer uncertainties under standard methods was made.However, different uncertainty approaches were required to accommodate the various methodologies.For instance, the method of calculating confidence intervals using Fisher's F distribution is robust and well suited to regression approaches, for the PDD method [22].In contrast, the manufacturer-stated RED had little to no tolerances or uncertainties provided, so for these data the assumption of meeting ICRU 44 tolerances was made, and for the purposes of the ANOVA a high confidence was attributed by way of assuming a larger n value.For the manufacturer's compositional data, the main component of uncertainty arose from the limitations in measuring physical density and slab volumes in the clinic.For the CT, the uncertainty was conservatively generated from an analysis of the sampled RED of the Gammex calibration materials over 23 scanning sessions.The inter-scan uncertainty was the main component, far outweighing the uncertainty associated with the variance observed within the materials sampled in each scan.Taken as a whole, the ANOVA was useful to identify and investigate differences between methods in each material, but further information is required regarding components of uncertainty, especially from the manufacturers, before attributing high confidence to the results reported as significant.
Some limitations to the robustness of the PDD-measured results were due to practical considerations-available slab arrangements limited both the sampling frequency and range along the measured PDD's.These limits influence the power of the regression method to fully capture the features of the measured PDD for comparison to the reference data.Due to its high density and limited range, the build-up region of CIRS cortical bone was not well sampled, reducing confidence in the result.Similarly due to its low-density and limited range (depth), no measurement of the PDD beyond the build-up region was obtained for CIRS inhale lung and only a limited number of samples of CIRS exhale lung were possible.Fortunately, similar materials with adequate sampling were also present in these ED ranges, specifically Gammex SB cortical bone and Gammex LN300 lung.In addition, even when the overall sampling range is well balanced, i.e., either side of Dmax, the regression can be biased by over or under sampling different parts of the curve.Further, the limited possible measurement locations made precise measurement of Dmax impractical which led to the need of an extra normalisation variable used with the solver.This highlights the constraints of the method imparted from limited quantities of material available.
The main parameter optimised by the regression method was the dose difference.Other comparison methods are available e.g., Gamma [23], Delta [24], etc., which may be more suitable and should be explored in the context of regression optimisation methods.Data interpolation may also impact the quality of the PDD regression results.In this study a linear interpolation between nearest neighbours was applied to the higher resolution reference data for matching to the measured PDD data.Other interpolation methods may be more suitable e.g., spline fitting [25], etc., further improving the regression optimisation results.
Overall, the quality of the regression of measured-PDD is impacted by practical challenges relating to sampling in discrete material arrangements, but improvements are possible and may form further work.A suggestion for dealing with new or novel materials in a clinical radiotherapy setting, CT could be used initially to confirm the spatial geometry and provide an estimate RED [26].The RED could be confirmed based on the behaviour of the material in the MV beam and then compared to the TPS.While a typical approach has been to use single points as a comparison, using PDDs regressed to a reference provides a more robust evaluation of the behaviour of the material over a range of geometries providing a result valid over a greater range of applications.

Conclusion
Determination of material properties in the clinical beam remains a challenge.The regression method using measured PDD data showed results consistent with both manufacturer-stated values and CT-derived RED in the range around 0.9 and 1.15.Outside of this soft tissue range a trend was observed whereby the determined 〈Z/A〉 becomes physically unrealistic indicating the method is no longer reporting RED alone and the assumptions in the eEPL model are constrained.Within this range of validity, the regression method provides a practical and robust characterisation for unknown soft tissue materials in the clinical setting and may be used to improve on the CT-derived RED where high Z material components are suspected.

Fig. 1
Fig. 1 PDD measurement arrangement, showing the layout for two different depth measurements of a tissue mimicking material

Fig. 2 aFig. 3 a
Fig. 2 a Scaled PDD for Gammex LN 300 Lung and b scaled PDD for SB Cortical Bone quantifying the fit as the residual percentage differences between the measured and reference PDDs

Fig. 4
Fig.4 Comparison of 〈Z/A〉 results calculated from the manufacturer stated RED (Manu) and the PDD fitting method (PDD) for a range of tissue mimicking materials

Fig. 5
Fig.5 Comparison of results using all four methods to determine RED for a range of tissue mimicking materials: RED calculated from the measured mass density and manufacturer stated composition (Manufacturer Composition), RED stated by the manufac-

Table 1
Material characteristicsThe depth range available for measurements, measured mass density, manufacturer supplied relative electron density (RED), assumed total uncertainty (U), and corresponding O'Connor Field Size (FS) scaling factor used for the PDD fitting method

Table 2
Measurement setup-LINAC and CT parameters