Advertisement

Comparison of projection- and image-based methods for proton stopping power estimation using dual energy CT

Open AccessPublished:September 04, 2017DOI:https://doi.org/10.1016/j.phro.2017.08.001

      Abstract

      Background and Purpose

      Several strategies for estimating stopping power ratio (SPR) from dual-energy CT (DECT) have been proposed to improve accuracy of proton dose calculations. DECT methods can mainly be categorized into projection-based methods, where material decomposition is performed prior to image reconstruction, and image-based methods, where decomposition takes place after image reconstruction. With the advent of photon-counting and dual-layer technology, projection-based methods could be considered for SPR estimation. In this simulation-based study we compared the SPR accuracy of one projection- and three image-based DECT methods.

      Materials and Methods

      X-ray CT projections of the female ICRP phantom were simulated using two different X-ray spectra with a realistic detector response and noise levels. ICRP slices at four different locations were selected. Reference SPR-maps were computed at 200 MeV. The SPR comparison was based on percentage deviation inside ROIs and relative range errors calculated with Radon transform of difference maps.

      Results

      SPR root-mean-square errors (RMSE) over the selected ROIs were 0.54% for the projection-based method and 0.68%, 0.61% and 0.70% for the three image-based methods. The RMSE for the relative range errors were slightly smaller for the projection-based approach, but close to zero for all decomposition domains as positive and negative errors averaged out over the slice.

      Conclusions

      SPR estimations with the projection-based method produced slightly better results (though not statistically significant) than the three image-based methods used in this simulation-based study, therefore, with the advent of technological developments, projection-based methods could be considered for SPR estimation if projection data is available.

      Keywords

      1. Introduction

      Treatment planning dose calculations of proton therapy require conversion of each voxel of the computed tomography (CT) scan into proton stopping power ratios (SPRs) relative to water [
      • Yang M.
      • Virshup G.
      • Clayton J.
      • Zhu X.R.
      • Mohan R.
      • Dong L.
      Theoretical variance analysis of single- and dual-energy computed tomography methods for calculating proton stopping power ratios of biological tissues.
      ]. However, since photons and protons interact differently inside the human body, there is no direct correspondence between CT numbers and SPRs. The uncertainty of the conversion from CT numbers to SPR may lead to proton range errors up to 3.5% [
      • Yang M.
      • Zhu X.R.
      • Park P.C.
      • Titt U.
      • Mohan R.
      • Virshup G.
      • et al.
      Comprehensive analysis of proton range uncertainties related to patient stopping-power-ratio estimation using the stoichiometric calibration.
      ]. It is therefore important to estimate the SPR accurately to reduce treatment margins and to exploit the full potential of proton therapy.
      The use of dual-energy CT (DECT) was investigated [
      • Alvarez R.E.
      • Macovski A.
      Energy-selective reconstructions in X-ray computerized tomography.
      ,
      • Rutherford R.A.
      • Pullan B.R.
      • Isherwood I.
      Measurement of effective atomic number and electron density using an EMI scanner.
      ] already a few years after the introduction of the CT scanner [
      • Hounsfield G.N.
      Computerized transverse axial scanning (tomography). Part I. Description of system.
      ]. Both Alvarez and Macovski [
      • Alvarez R.E.
      • Macovski A.
      Energy-selective reconstructions in X-ray computerized tomography.
      ] and Rutherford et al. [
      • Rutherford R.A.
      • Pullan B.R.
      • Isherwood I.
      Measurement of effective atomic number and electron density using an EMI scanner.
      ] observed that the linear attenuation coefficient, μ, could be separated into terms describing the two main interaction processes for X-ray photons in the energy range used in diagnostic CT: Photoelectric absorption and Compton scattering. Material decomposition can therefore be performed when having attenuation measurements at two different energies [
      • Heismann B.J.
      • Schmidt B.
      • Flohr T.
      Spectral computed tomography.
      ]. The decomposition schemes proposed by Alvarez and Macovski [
      • Alvarez R.E.
      • Macovski A.
      Energy-selective reconstructions in X-ray computerized tomography.
      ] and Rutherford et al. [
      • Rutherford R.A.
      • Pullan B.R.
      • Isherwood I.
      Measurement of effective atomic number and electron density using an EMI scanner.
      ] differed both in the suggested basis functions and in the sense that Alvarez and Macovski based their estimation scheme on the projection data directly obtained from the CT scan and then reconstructed these basis images, whereas Rutherford et al. relied on reconstructed CT data of the linear attenuation for their estimations. These two approaches can be categorized as projection-based and image-based methods, respectively. The main difference in the approaches is therefore the order of the decomposition and the reconstruction step. In projection-based methods, material decomposition is performed prior to image reconstruction, and for image-based methods the decomposition takes place after image reconstruction.
      Both projection-based and image-based approaches can be further divided into two main types: Type (1) Model-based: prior knowledge of source energy spectra and detector response is required to solve analytically the system of equations; Type (2) Calibration-based: X-ray measurements of materials with known properties are performed to obtain parameters that characterize the energy spectra and the detector response.
      Several algorithms based on DECT have been proposed to predict the SPR or the intermediate variables required to compute the SPR, e.g., the relative electron density (RED), the effective atomic number (Zeff) or the logarithm of the mean excitation energy [
      • Yang M.
      • Virshup G.
      • Clayton J.
      • Zhu X.R.
      • Mohan R.
      • Dong L.
      Theoretical variance analysis of single- and dual-energy computed tomography methods for calculating proton stopping power ratios of biological tissues.
      ,
      • Saito M.
      Potential of dual-energy subtraction for converting CT numbers to electron density based on a single linear relationship.
      ,
      • Hünemohr N.
      • Krauss B.
      • Tremmel C.
      • Ackermann B.
      • Jäkel O.
      • Greilich S.
      Experimental verification of ion stopping power prediction from dual energy CT data in tissue surrogates.
      ,
      • Hünemohr N.
      • Paganetti H.
      • Greilich S.
      • Jäkel O.
      • Seco J.
      Tissue decomposition from dual energy CT data for MC based dose calculation in particle therapy.
      ,
      • Landry G.
      • Seco J.
      • Gaudreault M.
      • Verhaegen F.
      Deriving effective atomic numbers from DECT based on a parameterization of the ratio of high and low linear attenuation coefficients.
      ,
      • Bourque A.E.
      • Carrier J.F.
      • Bouchard H.
      A stoichiometric calibration method for dual energy computed tomography.
      ,
      • Farace P.
      Comment: experimental verification of ion stopping power prediction from dual energy CT data in tissue surrogates.
      ,
      • van Abbema J.K.
      • van Goethem M.J.
      • Greuter M.J.W.
      • van der Schaaf A.
      • Brandenburg S.
      • van der Graaf E.R.
      Relative electron density determination using a physics based parameterization of photon interactions in medical DECT.
      ,
      • Garcia L.I.R.
      • Azorin J.F.P.
      • Almansa J.F.
      A new method to measure electron density and effective atomic number using dual-energy CT images.
      ,
      • Han D.
      • Siebers J.V.
      • Williamson J.F.
      A linear, separable two-parameter model for dual energy CT imaging of proton stopping power computation.
      ,
      • Möhler C.
      • Wohlfahrt P.
      • Richter C.
      • Greilich S.
      Range prediction for tissue mixtures based on dual-energy CT.
      ,
      • Taasti V.T.
      • Petersen J.B.B.
      • Muren L.P.
      • Thygesen J.
      • Hansen D.C.
      A robust empirical parametrization of proton stopping power using dual energy CT.
      ]. At present, no consensus is reached within the community to prefer one method over the others, nor a decomposition domain (the projection- or the image-domain) for the SPR estimation. Nearly all DECT methods are image-based [
      • van Elmpt W.
      • Landry G.
      • Das M.
      • Verhaegen F.
      Dual energy CT in radiotherapy: current applications and future outlook.
      ], possibly because image-based methods have the clear advantage that they can be performed on all DECT scanners, including conventional SECT scanners by acquiring two consecutive images at different kVp settings. In contrast, for projection-based methods the low and the high energy projections need to be acquired at the same angle or be accurately interpolated. Nevertheless, with the advent of dual-layer systems and photon-counting detector technology, projection-based methods are becoming more attractive since the projections from the different energy spectra are obtained without angular separation and there cannot be any motion between the two (or more) projections per source position. Another reason is that projection-based methods require access to dual-energy raw-data, which are not accessible in most commercial CT systems without the vendor agreement or reverse engineering.
      It has been shown that DECT has a potential to improve the accuracy of the SPR estimation compared to conventional single energy CT (SECT) [
      • Yang M.
      • Virshup G.
      • Clayton J.
      • Zhu X.R.
      • Mohan R.
      • Dong L.
      Theoretical variance analysis of single- and dual-energy computed tomography methods for calculating proton stopping power ratios of biological tissues.
      ,
      • Hansen D.C.
      • Seco J.
      • Srensen T.S.
      • Petersen J.B.B.
      • Wildberger J.E.
      • Verhaegen F.
      • et al.
      A simulation study on proton computed tomography (CT) stopping power accuracy using dual energy CT scans as benchmark.
      ,
      • Hünemohr N.
      • Krauss B.
      • Dinkel J.
      • Gillmann C.
      • Ackermann B.
      • Jäkel O.
      • et al.
      Ion range estimation by using dual energy computed tomography.
      ,
      • Zhu J.
      • Penfold S.N.
      Dosimetric comparison of stopping power calibration with dual-energy CT and single-energy CT in proton therapy treatment planning.
      ], however, these studies were based on image-based methods. Tremblay et al. [
      • Tremblay J.-É.
      • Bedwani S.
      • Bouchard H.
      A theoretical comparison of tissue parameter extraction methods for dual energy computed tomography.
      ] concluded that combining the output of projection-based methods with the image-based stoichiometric calibration proposed by Bourque et al. [
      • Bourque A.E.
      • Carrier J.F.
      • Bouchard H.
      A stoichiometric calibration method for dual energy computed tomography.
      ] gave better accuracy for RED and Zeff than a projection-based method alone. The aim of the present study was therefore to investigate the impact on the SPR accuracy when using a projection-based method. We compared the performances of the SPR estimation and the corresponding proton range errors of one projection- and three image-based DECT methods, through a simulation-based CT acquisition framework with the measured detector response and dual-energy spectra of a commercial DECT scanner and realistic noise levels.

      2. Methods and materials

      2.1 SPR estimation methods

      Four DECT-based SPR estimation methods were investigated, including one projection-based and three image-based approaches. The projection-based method relied on the two-material decomposition proposed by Alvarez and Macovski [
      • Alvarez R.E.
      • Macovski A.
      Energy-selective reconstructions in X-ray computerized tomography.
      ] to estimate the RED and the conversion from RED to SPR proposed by Kanematsu et al. [
      • Kanematsu N.
      • Inaniwa T.
      • Koba Y.
      Relationship between electron density and effective densities of body tissues for stopping, scattering, and nuclear interactions of proton and ion beams.
      ]. We will refer to this method as “AMK” for Alvarez and Macovski, and Kanematsu. The first image-based method used a parametrization for the SPR proposed by Taasti et al. [
      • Taasti V.T.
      • Petersen J.B.B.
      • Muren L.P.
      • Thygesen J.
      • Hansen D.C.
      A robust empirical parametrization of proton stopping power using dual energy CT.
      ]. This image-based method will be referred to as “SPP” for Stopping Power Parametrization. The second image-based approach was a combination of the methods proposed by Saito [
      • Saito M.
      Potential of dual-energy subtraction for converting CT numbers to electron density based on a single linear relationship.
      ] to compute the RED and by Kanematsu et al. [
      • Kanematsu N.
      • Inaniwa T.
      • Koba Y.
      Relationship between electron density and effective densities of body tissues for stopping, scattering, and nuclear interactions of proton and ion beams.
      ] to derive the SPR. We will refer to this method as “SK” for Saito and Kanematsu. The third image-based method is a two-material decomposition method with different basis materials for soft and bone tissues, suggested by Han et al. [
      • Han D.
      • Siebers J.V.
      • Williamson J.F.
      A linear, separable two-parameter model for dual energy CT imaging of proton stopping power computation.
      ]. The four SPR estimation methods are described in Appendix A.

      2.2 Virtual patient

      The Adult Female (AF) reference computational phantom of the International Commission on Radiological Protection (ICRP) [

      ICRP. Adult Reference Computational Phantoms. ICRP Publication 110. Ann. ICRP 39 (2).

      ] was used as a virtual patient to evaluate the performance of the SPR methods. This phantom, with voxel dimensions of 1.775 × 1.775 × 4.84 mm3, represented an average female subject structured with 140 organs made of 53 standard human tissues. The mass density and elemental weights of each tissue were provided in the ICRP Publication 110 [

      ICRP. Adult Reference Computational Phantoms. ICRP Publication 110. Ann. ICRP 39 (2).

      ]. To evaluate a large number of tissue densities and compositions, four anatomical regions were selected: head, sternum, breast and pelvis (Fig. 1).
      Figure thumbnail gr1
      Fig. 1Reference SPR for the four slices used in the comparison. Placement of ROIs are shown by blue circles, three to four ROIs are placed in each slice. (For interpretation of the references to colour in this figure caption, the reader is referred to the web version of this article.)

      2.3 Calculation of reference SPR

      Reference SPR values of the ICRP phantom were computed from the chemical composition and density pixel-by-pixel using the Bethe equation without correction terms, as described by Schneider et al. [
      • Schneider U.
      • Pedroni E.
      • Lomax A.
      The calibration of CT Hounsfield units for radiotherapy treatment planning.
      ]. A proton beam of 200 MeV kinetic energy was considered and the mean ionization energy of water was set to 78 eV [

      Sigmund P, Schinner A, Paul H. Errata and Addenda for ICRU Report 73. Stopping of ions heavier than helium. J ICRU 2009;5(1).

      ]. The mean ionization energies of the tissues were calculated using the Bragg additivity rule [
      • Schneider U.
      • Pedroni E.
      • Lomax A.
      The calibration of CT Hounsfield units for radiotherapy treatment planning.
      ] and the I-values given in Tables 2.8 and 2.11 in the ICRU49 report [

      ICRU. Stopping Power and Ranges for Protons and Alpha Particles. ICRU Report No. 49.

      ].

      2.4 CT imaging and CT reconstruction

      Virtual CT projections of the Gammex and the ICRP phantom were generated in GATE [
      • Jan S.
      • Santin G.
      • Strul D.
      • Staelens S.
      • Assié K.
      • Autret D.
      • et al.
      GATE: a simulation toolkit for PET and SPECT.
      ] v7.2 (based on Geant4 v10.1, physics-list: emlivermore) using the Fixed Forced Detection Actor. This module deterministically computes digitally reconstructed radiographs using the reconstruction toolkit (RTK) v1.3.0 [
      • Rit S.
      • Vila Oliva M.
      • Brousmiche S.
      • Labarbe R.
      • Sarrut D.
      • Sharp G.C.
      The Reconstruction Toolkit (RTK), an open-source cone-beam CT reconstruction toolkit based on the Insight Toolkit (ITK).
      ] and the Geant4 database of X-rays cross sections. Fan-beam projections with 600 views of 2052 pixels subsequently re-binned to 1026 pixels of 1 mm2 were acquired. The thickness of the detector row was 1 mm for the Gammex acquisitions and 2 mm for the ICRP phantom. The source-to-isocenter distance was 626 mm and the source-to-detector distance was 1026 mm. For the simulations of the CT projections, the measured dual-energy spectra (low energy (LE): 100 kVp, high energy (HE): 150 kVp + 0.6 mm Sn) and the measured detector response for the SOMATOM Force dual-source CT scanner were used, kindly provided by Siemens Healthcare (Forchheim, Germany). To have a fair comparison between the projection-based and the image-based methods, energy spectra with high mean energies were chosen to reduce the beam-hardening for the image-based methods. To represent a realistic scenario, Poisson distributed noise was applied to the projections. For each slice, a total central dose of 20 mGy was delivered, with an equal dose split between the two energy spectra, thus a central dose of 10 mGy per energy spectrum; the calculation of the number of photons needed to deliver a given dose is described in Appendix C. No bowtie-filter was simulated.
      Prior to CT reconstruction, the LE and the HE sinograms of the Gammex and the AF phantom were corrected for beam-hardening using the method described in Ref. [
      • Brooks R.A.
      • Di Chiro G.
      Beam hardening in X-ray reconstructive tomography.
      ], since this was found to be mandatory for the image-based methods. A look-up table linearizing the input projection values for water [
      • Brooks R.A.
      • Di Chiro G.
      Beam hardening in X-ray reconstructive tomography.
      ] was created using the energy spectra and the detector response and increasing water thicknesses.
      CT images of the Gammex phantom and the AF slices were reconstructed using the LE and the HE sinograms, for the image-based methods; and soft tissue and cortical bone density images of the AF slices, for the projection-based method. All reconstructions were performed using the filtered-backprojection (FBP) reconstruction of RTK; for the reconstructions of the AF slices preserving the initial pixel dimension of 1.775 × 1.775 mm2 and for the Gammex the pixel dimension of 1 × 1 mm2.

      2.5 SPR comparison

      The performance of the SPR estimation was evaluated based on different criteria: relative SPR differences over defined ROIs, each covering a single tissue, and relative range errors computed over a whole slice. The placement of the ROIs can be seen in Fig. 1, and their reference SPRs are listed in Table B.5 in Appendix B. For each ROI, the mean error was calculated and its precision was assessed using the standard error of the mean (SEM). To have a direct comparison of the three methods, the root-mean-square error (RMSE) over the mean relative SPR difference for the thirteen ROIs was taken; this gave a combined error measure for the thirteen defined ROIs.
      To estimate the influence of the SPR deviations on the range calculation, we calculated the range errors along the proton beam path through the entire slice using the Radon transform implemented in MATLAB (The MathWorks Inc., Natick, MA). The Radon transform was computed on the absolute SPR difference images, for angles in the interval from 0–179° in steps of 1°. The SPR difference for pixels outside the body outline was set to zero. The Radon transform was also performed on the reference SPR maps to calculate the water-equivalent path length (WEPL). To exclude beam paths entirely outside the body, an exclusion penalty was applied to pixels in the reference SPR sinogram with a WEPL of zero. The results for the range errors were given relative to reference WEPL maps, to take into account that the beam paths in different directions and slices were not of equal length.

      3. Results

      The RMSE over the mean relative error for each of the thirteen ROIs was 0.54% for the AMK method, 0.68% for the SPP method, 0.61% for the SK method, and 0.70% for Han’s method (Fig. 2). The bias (given by the mean error) for the head, sternum and breast slices was smallest for the SPP method, however, for the pelvis slice (ROIs 4a-d) this method gave the largest errors (Table 1). Considering all slices, the SK method had the largest bias. For all methods, the relative errors for ROIs in the pelvis slice (ROIs 4a-d) were larger than for the other slices.
      Figure thumbnail gr2
      Fig. 2Relative SPR errors for each of the thirteen ROIs (for the placement see ). The center squares show the mean of the relative SPR error over the ROI, and the error-bars show the standard error of the mean (SEM).
      Table 1Comparison of relative SPR deviations over the defined ROIs. The upper part of the table gives the results averaged over the individual slices, while the lower part gives the results averaged over all four slices together. (PB: Projection-based, IB: Image-based).
      Per sliceHeadSternumBreastPelvis
      PB (AMK)Mean0.28%0.32%0.29%−0.45%
      Uns. mean0.28%0.32%0.64%0.66%
      IB (SPP)Mean−0.06%0.03%0.17%−0.99%
      Uns. mean0.46%0.30%0.29%0.99%
      IB (SK)Mean−0.17%−0.06%−0.09%−0.85%
      Uns. mean0.20%0.20%0.40%0.85%
      IB (Han)Mean0.01%0.26%0.79%−0.60%
      Uns. mean0.34%0.26%0.79%0.86%
      All slicesRMSEMeanUns. meanσ
      PB (AMK)0.54%0.07%0.49%0.56%
      IB (SPP)0.68%−0.27%0.55%0.65%
      IB (SK)0.61%−0.33%0.44%0.53%
      IB (Han)0.70%0.06%0.59%0.73%
      The AMK method gave the lowest RMSEs for the relative range errors. For the AMK and SPP methods, the positive and negative range errors were nearly averaged out over the slice, giving mean relative range errors close to 0%; only for the head slice estimated by the SPP method the 0% range error was not within the interquartile range (IQ) (Fig. 3). In contrast, for the image-based SK method, the range error distribution had a negative bias for all slices, and Han’s method gave a slight positive bias for some slices. Even though the mean relative range errors for the SPP method were fairly low, this method produced the widest error distributions (Fig. 3 and Table 2).
      Figure thumbnail gr3
      Fig. 3Box-plot of relative range error distribution, as calculated by the Radon transform for each anatomical slice using each of the four SPR estimation methods (PB: Projection-based, IB: Image-based). The boxes show the interquartile (IQ) range, i.e. from the 25%-percentile (Q1) to the 75%-percentile (Q3), and the black dots represent the median. The whiskers go from Q1-2.5×IQ up to Q3+2.5×IQ. The outliers are not shown in the figure, but the percentages (excluding the sinogram entries corresponding to a WEPL of zero) are given in the right side of the figure.
      Table 2Comparison of the relative range errors over each image slice. The results are given as the signed mean errors (μ±σ) and the root-mean-square errors (RMSE) in percentage. The range errors are taken relative to the reference water equivalent path length (WEPL); the mean reference WEPL for each slice is given in parentheses beneath the slice name. (PB: Projection-based, IB: Image-based).
      Head (122.8 mm)Sternum (162.5 mm)Breast (162.1 mm)Pelvis (181.7 mm)
      Methodμ±σ (%)RMSE (%)μ±σ (%)RMSE (%)μ±σ (%)RMSE (%)μ±σ (%)RMSE (%)
      PB (AMK)−0.24 ± 0.810.84−0.01 ± 0.640.64−0.04 ± 0.580.58−0.14 ± 0.800.82
      IB (SPP)−0.41 ± 0.860.950.01 ± 0.810.810.04 ± 0.720.720.03 ± 1.101.10
      IB (SK)−0.51 ± 0.800.95−0.28 ± 0.630.69−0.33 ± 0.570.66−0.44 ± 0.800.91
      IB (Han)−0.10 ± 0.820.820.36 ± 0.750.830.35 ± 0.650.740.47 ± 1.191.28

      4. Discussion

      In this study, we have compared one projection-based and three image-based SPR methods. We found that the projection-based method produced slightly better results than the three image-based methods used in this study. The differences were not found to be statistically significant. The image-based SPP method gave mean errors comparable to the results for the projection-based method, but larger standard deviations. In contrast, the image-based SK method resulted in slightly biased results. It has been shown earlier that the SPR accuracy using Saito’s method in combination with an estimation of Zeff is dependent on the calibration and evaluation materials being very similar [
      • Hansen D.C.
      • Seco J.
      • Srensen T.S.
      • Petersen J.B.B.
      • Wildberger J.E.
      • Verhaegen F.
      • et al.
      A simulation study on proton computed tomography (CT) stopping power accuracy using dual energy CT scans as benchmark.
      ].
      In our simulations, clinical polychromatic energy spectra, a realistic detector response and noise levels corresponding to a 20 mGy central dose were used. The material basis for the projection-based approach was optimized to give the best results. To have a direct comparison between projection-based and image-based calculations, we used an image-based method (the SK method) where the conversion from RED to SPR were performed in the same way as for the projection-based method, using the polyline fit proposed by Kanematsu et al. [
      • Kanematsu N.
      • Inaniwa T.
      • Koba Y.
      Relationship between electron density and effective densities of body tissues for stopping, scattering, and nuclear interactions of proton and ion beams.
      ], and we also compared to a two-material decomposition image-based method, though not using the exact same materials as for the projection-based method, as we chose to use the optimal basis materials in each case. Further, we also used one of the latest image-based methods [
      • Taasti V.T.
      • Petersen J.B.B.
      • Muren L.P.
      • Thygesen J.
      • Hansen D.C.
      A robust empirical parametrization of proton stopping power using dual energy CT.
      ], as this method was proven to give smaller SPR errors and to be more robust to noise than two other existing DECT methods [
      • Han D.
      • Siebers J.V.
      • Williamson J.F.
      A linear, separable two-parameter model for dual energy CT imaging of proton stopping power computation.
      ,
      • Hansen D.C.
      • Seco J.
      • Srensen T.S.
      • Petersen J.B.B.
      • Wildberger J.E.
      • Verhaegen F.
      • et al.
      A simulation study on proton computed tomography (CT) stopping power accuracy using dual energy CT scans as benchmark.
      ]. For these reasons, we consider that the comparison was carried out under fair conditions.
      The image-based SK method showed a net negative bias on the range accuracy that was not observed for the image-based SPP method. Range errors seem to correlate well with SPR errors obtained for some tissues, such as breast tissue (ROI 3c), adipose tissue (ROI 4d) or pelvis spongiosa (ROI 4c), which were found in large proportions in the selected slices. This suggests that accuracy errors add up along the proton beam path. The SPP method resulted in noisier SPR images and, thus, wider error distributions. In terms of RMSE, the projection-based approach produced better results for the SPR and the range with respect to the three image-based methods. Therefore, based on these results, we recommend increasing the consideration of projection-based approaches for proton SPR estimation. However, the inherent limitations of each method should be taken into account, and therefore the most well-suited SPR estimation method for a specific facility might depend on the DECT technology available at the facility. It can therefore not be generally stated which method will produce the most accurate results.
      Beam-hardening artifacts occur because the energy dependency of the linear attenuation coefficients is neglected. To overcome this approximation, image-based methods require a beam-hardening correction. In contrast, material decomposition in projection-based approaches is performed prior to reconstruction, and as the energy dependency is factored out in this decomposition process (see Eq. (A.1)) the basis-images should not be affected by beam-hardening. But this is not true unless a good choice of the basis-material decomposition is performed, such that the assumption that μ can be separated into an energy-dependent and an energy-independent part holds true. In the energy range of clinical CT scanners (i.e. from 40 to 150 keV) and for materials with Z<15 this assumption is acceptable [
      • Jackson D.F.
      • Hawkes D.J.
      X-ray attenuation coefficients of elements and mixtures.
      ]. In this study, the implemented strategy to correct for beam-hardening in the image-based method works well for water-equivalent tissues but might not be accurate enough for low- and high-density tissues. A more sophisticated beam-hardening correction might therefore improve the results for the image-based methods. If a better beam-hardening correction had been applied, a LE spectrum with a lower mean energy might be more relevant than the 100 kVp spectrum used in this study.
      The three image-based methods used in this study were calibration-based (type 2), and a calibration was needed to find the effective energies for the LE and HE spectrum. This calibration procedure is sensitive to the calibration phantom being comparable to the object for which the SPR is to be estimated; i.e. the calibration phantom should be of a similar size and consist of materials with a composition comparable to human tissues, such that the phantom hardens the X-ray energy spectrum in the same way as the investigated anatomical site. For this study, only a single calibration phantom with a diameter of 16.5 cm was used for all four anatomical slices. This may be one cause for the larger relative SPR deviations in the pelvis slice for the image-based methods (Table 1).
      For model-based (type 1) image-based methods, the energy spectrum must also be properly hardened to reproduce the beam-hardening of the scanned object. This issue is taken into account by additionally filtering the X-ray source for instance by water [
      • Yang M.
      • Virshup G.
      • Clayton J.
      • Zhu X.R.
      • Mohan R.
      • Dong L.
      Theoretical variance analysis of single- and dual-energy computed tomography methods for calculating proton stopping power ratios of biological tissues.
      ,
      • Bazalova M.
      • Carrier J.-F.
      • Beaulieu L.
      • Verhaegen F.
      Dual-energy CT-based material extraction for tissue segmentation in Monte Carlo dose calculations.
      ]. Here, the amount of water may also be dependent on the anatomical region being investigated, since the performance of these methods will depend on the accuracy of this prior knowledge on the energy spectrum and the detector response [
      • Yang M.
      Dual energy computed tomography for proton therapy treatment planning.
      ]. Therefore, the input parameters must be tuned with care for both types of image-based methods. However, for calibration-based methods (type 2), this is more straight-forward since it only requires using a well-suited calibration phantom of the proper size.
      For the projection-based method applied in this study, the system of Eqs. ((A.4), (A.5)) was solved numerically assuming that the energy spectra and the detector response were well-known (type 1). The same was assumed for the beam-hardening correction of image-based methods. It was shown that projection-based methods are sensitive to the amount of noise [

      Brendel B, Bergner F, Brown K, Koehler T. Penalized likelihood decomposition for dual layer spectral CT. In: Proc 4th intl mtg on image formation in X-ray CT; 2016. p. 41–4.

      ] and to alterations of the energy spectrum [
      • Tremblay J.-É.
      • Bedwani S.
      • Bouchard H.
      A theoretical comparison of tissue parameter extraction methods for dual energy computed tomography.
      ]. Tremblay et al. found that for projection-based methods the accuracy of the RED was more degraded when the spectrum was altered than when using an image-based method [
      • Tremblay J.-É.
      • Bedwani S.
      • Bouchard H.
      A theoretical comparison of tissue parameter extraction methods for dual energy computed tomography.
      ]. The projection-based results in the present work are therefore dependent on accurate knowledge of the energy spectrum and detector response. An alternative would be to approximate the system of Eqs. ((A.4), (A.5)) by a polynomial expansion of two variables [
      • Alvarez R.E.
      • Macovski A.
      Energy-selective reconstructions in X-ray computerized tomography.
      ,
      • Chuang K.-S.
      • Huang H.K.
      Comparison of four dual energy image decomposition methods.
      ]. The unknown coefficients could be determined experimentally through a calibration procedure with increasing thicknesses of two well-known materials (type 2). However, there would also be disadvantages with this approach, including sensitivity to the order of the polynomial, to the use of cross-terms in the polynomial, to the exact thicknesses of the interposed slabs of two materials and to the number of points of the calibration curve. It should also be noted that type 1 approaches always rely on a pre-calibration of the model, e.g., the one presented in [
      • Vilches-Freixas G.
      • Létang J.M.
      • Brousmiche S.
      • Romero E.
      • Oliva M.V.
      • Kellner D.
      • et al.
      Technical note: procedure for the calibration and validation of kilo-voltage cone-beam CT models.
      ].
      The basis functions used for the projection-based method were found to give the best SPR results, after testing several different basis functions. Other basis functions have been suggested in the literature, e.g. Han et al. [
      • Han D.
      • Siebers J.V.
      • Williamson J.F.
      A linear, separable two-parameter model for dual energy CT imaging of proton stopping power computation.
      ] used two different sets of basis functions for soft and bone tissues, respectively, where the separation was based on the ratio of the CT numbers which is not available for projection-based methods. This method was tested in image-space in this study. Tremblay et al. [
      • Tremblay J.-É.
      • Bedwani S.
      • Bouchard H.
      A theoretical comparison of tissue parameter extraction methods for dual energy computed tomography.
      ] found that in most cases a two-material decomposition rather than decomposition into physical processes performed better.
      Another ongoing discussion within the community is how to determine the reference SPR values. In this study, a 200 MeV proton energy was used to determine the reference SPR, as for proton energies in the range between 80 and 300 MeV, the variation of SPR with proton energy is negligible (<1%) [
      • Arbor N.
      • Dauvergne D.
      • Dedes G.
      • Létang J.M.
      • Parodi K.
      • Quiñones C.T.
      • et al.
      Monte Carlo comparison of X-ray and proton CT for range calculations of proton therapy beams.
      ]. In addition, Yang [
      • Yang M.
      Dual energy computed tomography for proton therapy treatment planning.
      ] conducted a study with a 175 MeV initial kinetic energy beam and concluded that neglecting the SPR dependence with energy introduced a 0.5% range error.

      5. Conclusion

      In this simulation-based study, the dual-energy decomposition domains – the projection-domain (decomposition prior to image reconstruction) and the image-domain (decomposition after image reconstruction) – were compared in terms of the performance of extracting the SPR from DECT and the resulting WEPL. In terms of SPR accuracy and range errors, the observed differences between projection- and image-based methods were not statistically significant. If properly calibrated, the two SPR estimation approaches can therefore be assumed to give similar results. With the advent of technological developments, projection-based methods are becoming more attractive. Most of the existing methods to compute the SPR are image-based but, based on these results, projection-based approaches should also be considered for proton SPR determination.

      Disclosure of conflicts of interest

      The authors have no relevant conflicts of interest to disclose.

      Acknowledgments

      This work was partially supported by grant ANR-13-IS03-0002-01 (DEXTER project) from the French National Research Agency (ANR), and by the Danish Cancer Society (Grant No. A5992).

      Appendix A. SPR estimation methods

      A.1 Projection-based method (AMK)

      The method used to compute the SPR from dual-energy projection data is an adaptation of the two-step procedure proposed by Farace [
      • Farace P.
      Comment: experimental verification of ion stopping power prediction from dual energy CT data in tissue surrogates.
      ] in the image-domain. This method was chosen for its simplicity and because it produced comparable results with respect to Hnemohr’s approach [
      • Hünemohr N.
      • Krauss B.
      • Tremmel C.
      • Ackermann B.
      • Jäkel O.
      • Greilich S.
      Experimental verification of ion stopping power prediction from dual energy CT data in tissue surrogates.
      ]. First, RED was derived from the reconstructed mass density images of soft tissue (ST) and cortical bone (CB) [

      P. Linstrom, W. Mallard, NIST Standard Reference Database Number 69, [Online] Available:http://webbook.nist.gov, National Institute of Standards and Technology, Gaithersburg, MD. (retrieved December 6, 2016).

      ], issued from the implementation of the two-material decomposition method in the projection domain proposed by Alvarez and Macovski [
      • Alvarez R.E.
      • Macovski A.
      Energy-selective reconstructions in X-ray computerized tomography.
      ]. Second, RED was converted into SPR using the poly-line relations proposed by Kanematsu et al. [
      • Kanematsu N.
      • Inaniwa T.
      • Koba Y.
      Relationship between electron density and effective densities of body tissues for stopping, scattering, and nuclear interactions of proton and ion beams.
      ].
      The basis material decomposition (BMD) method assumed that the linear attenuation coefficient of the scanned object at a location x can be expressed as a linear combination of two energy-dependent basis functions of two materials with energy-independent coefficients:
      μ(E,x)=μρST(E)·aST(x)+μρCB(E)·aCB(x)
      (A.1)


      where (μ/ρ) represented the mass attenuation coefficient of ST and CB, and a is the mass fraction of ST and CB in the volume at location x – in this BMD these mass fractions represented the energy-independent coefficients.
      The projection of the mass fractions aST and aCB was expressed as follows:
      AST(u,θ)=L(u,θ)aST()d
      (A.2)


      ACB(u,θ)=L(u,θ)aCB()d
      (A.3)


      where L(u,θ) was the line-segment between the source and a detector pixel located at position u for a given projection angle θ.
      When performing a dual-energy acquisition, two sinograms of the same object at two different energy spectra were available:
      ILE(AST,ACB)=EminEmaxSLE(E)exp-ASTμρST(E)-ACBμρCB(E)dE
      (A.4)


      IHE(AST,ACB)=EminEmaxSHE(E)exp-ASTμρST(E)-ACBμρCB(E)dE
      (A.5)


      where ILE and IHE were the measured intensities for the LE and the HE spectrum, respectively, for a given projection angle; SLE and SHE were the normalized energy spectra weighted by the detector response.
      The unknowns AST and ACB were determined by solving the system of Eqs. ((A.4), (A.5)) using the implementation of Nelder and Mead [
      • Nelder J.A.
      • Mead R.
      A simplex method for function minimization.
      ] in RTK [
      • Rit S.
      • Vila Oliva M.
      • Brousmiche S.
      • Labarbe R.
      • Sarrut D.
      • Sharp G.C.
      The Reconstruction Toolkit (RTK), an open-source cone-beam CT reconstruction toolkit based on the Insight Toolkit (ITK).
      ], assuming that the energy spectra and the detector response were known (type 1 decomposition).
      Then, the mass fraction per volume of soft tissue (aST) and cortical bone (aCB) was determined by reconstruction of AST and ACB (see Section 2.4). Afterwards, the RED was derived for each pixel by:
      RED(x)=aST(x)iwiZiAiST+aCB(x)iwiZiAiCBρWiwiZiAiW
      (A.6)


      where Zi was the atomic number, Ai was the atomic mass and wi was the elemental weight fraction for element i of the tabulated compounds ST, CB and water (represented with the index W) [

      P. Linstrom, W. Mallard, NIST Standard Reference Database Number 69, [Online] Available:http://webbook.nist.gov, National Institute of Standards and Technology, Gaithersburg, MD. (retrieved December 6, 2016).

      ]; ρW was the mass density for water.
      To reproduce the piece-wise linear relations suggested by Kanematsu et al. [
      • Kanematsu N.
      • Inaniwa T.
      • Koba Y.
      Relationship between electron density and effective densities of body tissues for stopping, scattering, and nuclear interactions of proton and ion beams.
      ], the SPR and RED were calculated for 92 ICRU46 tissues of known chemical composition and mass density [

      ICRU, Photon, electron, proton and neutron interaction data for body tissues, ICRU Report No. 46.

      ] applying the Bethe equation (see Section 2.3). The fitting lines are presented in Fig. A.1.
      Figure thumbnail gr4
      Fig. A.1Kanematsu fitting lines to convert from RED to SPR/RED.

      A.2 Image-based method – SPR parametrization (SPP)

      The details of the image-based SPP method were described by Taasti et al. [
      • Taasti V.T.
      • Petersen J.B.B.
      • Muren L.P.
      • Thygesen J.
      • Hansen D.C.
      A robust empirical parametrization of proton stopping power using dual energy CT.
      ], but the main principles are described here.
      CT scans of a calibration phantom were required for the characterization of the X-ray energy spectra (type 2 decomposition). This calibration process was performed on simulated CT images of the Gammex 467 electron density calibration phantom (Gammex, Middleton, WI). The dimensions of the phantom were down-scaled to half the size of the real Gammex phantom which has a diameter of 33 cm; this was done to obtain the same beam hardening effect on this calibration phantom as for the ICRP image slices. Average CT numbers over ROIs inside the inserts in the reconstructed CT images were used.
      Based on the known density and elemental composition of the inserts of the Gammex phantom, the effective energies of the used energy spectra were found using the CT numbers (Hj) for the phantom inserts. The effective energy was defined as the energy which maximized the coefficient of determination, R2, for the linear fit:
      μ(Eeff,j)μW(Eeff,j)=HjAjt+Bjt
      (A.7)


      Here, subscript j refers to the energy spectrum (j=LE,HE), and A and B are fitting parameters. The linear attenuation coefficients, μ(E), for the Gammex inserts were calculated based on XCOM data [

      Berger MJ, Hubbell JH, Seltzer SM, Chang J, Coursey JS, Sukumar R, et al., XCOM: photon cross section database (version 1.5), [Online] Available: http://physics.nist.gov/xcom, National Institute of Standards and Technology, Gaithersburg, MD; 2010.

      ]. When the effective energies were found, the Gammex inserts were divided into two categories, soft and bone tissues based on their CT numbers in the LE CT image. Then Eq. (A.7) was refitted for each tissue group, to find two sets of fitting parameters, Ajt and Bjt, where superscript tsoft,bone indicates the tissue group, using the effective energies Eeff,j for each energy spectrum j. In this study we used a separation point between the tissue groups set to HLE=150 HU.
      The attenuation ratios, u, for the 92 reference human tissues [

      ICRU, Photon, electron, proton and neutron interaction data for body tissues, ICRU Report No. 46.

      ] were calculated at the effective energies based on their density and elemental composition:
      uμ(Eeff)μW(Eeff)=ρiwiμρi(Eeff)iwiμρi(Eeff)W
      (A.8)


      These attenuation ratios were fitted to the SPR for the tissues based on the following equations:
      SPRsoftest=(1+x1)uHE-x1uLE+x2uLE2+x3uHE2+x4uLE3+uHE3
      (A.9)


      SPRboneest=(1+x5)uHE-x5uLE+x6uLEuHE+x7uLE2-uHE2+x8uLE3+uHE3
      (A.10)


      where the xi’s are fitting parameters. The fitting parameters used in this study can be found in Table A.1. When these expressions were used to estimate the SPR for the AF phantom, the attenuation ratios were calculated using the fitting parameters found together with the effective energies, ujt=Hj/Ajt+Bjt. The same separation between soft and bone tissue was used for the SPR estimation as for the calculation of the attenuation ratios.
      Table A.1Calibration parameters used in the SPR parametrization; for the energy spectra characterization (Eq. A.7)), and for the SPR estimation (Eq. (A.9), (A.10)). The effective energies of the low energy (LE) and high energy (HE) spectra are given in parentheses.
      Energy spectra characterizationSPR fitting parameters
      LE (64 keV)HE (96 keV)Soft tissuesBone tissues
      Asoft988.8991.3x13.161x50.8251
      Abone971.8984.8x21.176x60.03853
      Bsoft1.0061.007x3−1.136x70.1150
      Bbone0.98031.004x4−0.01883x8−0.008910

      A.3 Image-based method – Saito and Kanematsu’s (SK) method

      This image-based DECT method is a combination of two methods originally presented by Saito [
      • Saito M.
      Potential of dual-energy subtraction for converting CT numbers to electron density based on a single linear relationship.
      ] (calculation of the relative electron density, RED), and by Kanematsu et al. [
      • Kanematsu N.
      • Inaniwa T.
      • Koba Y.
      Relationship between electron density and effective densities of body tissues for stopping, scattering, and nuclear interactions of proton and ion beams.
      ] (conversion from RED to SPR).
      In Saito’s method the RED was estimated as
      RED=a(1+α)HHE-αHLE1000+b
      (A.11)


      The constants a, b and α were found by making calibration fits to the theoretical ρe values for the 92 reference human tissues and their CT numbers calculated from Eqs. ((A.7), (A.8)). To take the low RED values for lung tissue into account the constants were found by minimizing the relative deviations REDtheo-REDest/REDtheo. The parameters used in this study can be seen in Table A.2. The CT numbers for the 92 reference human tissues used in the calibration were calculated from Eq. (A.8) and the effective energies given in Table A.1.
      Table A.2Calibration parameters for Saito’s method, found by minimization of Eq. A.11. These parameters were used for calculation of RED in the SK method applied in this study.
      abα
      1.00851.00910.5202
      The RED estimates were converted into SPR using the same calibration curve as described in Section A.1 for the Kanematsu method (see Fig. A.1).

      A.4 Image-based method – Han’s (Han) method

      The image-based two-material decomposition method used in this study was proposed by Han et al. [
      • Han D.
      • Siebers J.V.
      • Williamson J.F.
      A linear, separable two-parameter model for dual energy CT imaging of proton stopping power computation.
      ], and was also explained in a compressed way in Appendix A.2 in Ref. [
      • Taasti V.T.
      • Petersen J.B.B.
      • Muren L.P.
      • Thygesen J.
      • Hansen D.C.
      A robust empirical parametrization of proton stopping power using dual energy CT.
      ].
      The basis materials used in this method were water and polystyrene or water and a CaCl2 23% aqueous solution. An unknown material was assigned to one of the two basis material sets according to its ratio of reduced CT numbers, ξ=uHE/uLE; materials with ξ0.97 were categorized as soft tissues and assigned to the water-polystyrene group, while materials with ξ<0.97 were categorized as bone tissues and assigned to the water-CaCl2-group.
      The reduced CT numbers for this method were calculated as uj=Hj/Aj+Bj, however unlike in the SPP method, described in Appendix A.2, only one set of fitting parameters, Aj and Bj, was used to characterize all tissues. These parameters are listed in Table A.3.
      Table A.3Effective energies and CT number fitting parameters, from Eq. A.7 used in Han’s method.
      Eeff (keV)AjBj
      100 kVp64992.11.003
      Sn150 kVp96990.21.007
      The “concentrations” of the two basis materials were denoted c1 and c2, and they were found by solving the linear least square problem for c1 and c2:
      u(Eeff,LE)u(Eeff,HE)=u1(Eeff,LE)u2(Eeff,LE)u1(Eeff,HE)u2(Eeff,HE)c1c2
      (A.12)


      where the subscripts 1 and 2 denotes the two basis materials.
      The electron density and the logarithm to the mean excitation energy were calculated as
      ED=c1ED1+c2ED2
      (A.13)


      and
      I=fI·expc1ED1lnI1+c2ED2lnI2c1ED1+c2ED2fI·expIBVM
      (A.14)


      where ED is the electron density. The RED is then be found by dividing by the electron density of water. fI is a correction factor for the I values. It was found as a linear fit between the parameter c1c1+c2 and the ratio ItheoIBVM, fI=aItheoIBVM+b. Here Itheo is the theoretical I-value for the reference human tissue calculated based on the Bragg additivity rule. There is a fit for each of the two tissue groups, and the calibration parameters are listed in Table A.4. If the correction factor, fI, was negative, when estimating for an unknown material, it was set to 1 to avoid complex numbers when calculating the SPR estimates.
      Table A.4Fitting parameters (a and b) for the correction factor fI=aItheoIBVM+b used in Han’s method, see Eq. A.14.
      Soft tissuesBone tissues
      asoft=0.0837bsoft=0.8996abone=-0.0599bbone=1.021

      Appendix B. SPR reference values for the thirteen ROIs

      Table B.5Reference data for the thirteen defined ROIs used for comparison of the relative SPR differences. N is the number of pixels in the ROI.
      ROI nameAF Material-IDTissue typeSPRrefN
      Head1a49Adipose tissue0.97275
      1b32Brain1.051111
      1c8Cranium, spongiosa1.20327
      Sternum2a3Humeri, upper half, spongiosa1.157195
      2b29Muscle tissue1.046147
      2c29Muscle tissue1.04675
      Breast3a50Lung tissue (compressed lungs)0.384195
      3b28Blood1.055195
      3c48Breast (mammary gland)1.04047
      Pelvis4a29Muscle tissue1.046111
      4b49Adipose tissue0.972147
      4c14Pelvis, spongiosa1.10047
      4d9Femora, upper half, spongiosa1.05375

      Appendix C. Calculation of CT dose

      For the CT imaging setup, the level of Poisson noise depends on the number of primary photons per simulation, Nprim, required to deliver a central dose of 20 mGy, D=20mGy, which was determined analytically assuming an homogeneous water medium:
      Nprim=DAbeamES(E)e-μW(E)rμen,W(E)ρWEdE
      (C.1)


      where Abeam is the area covered by the beam at the isocenter, S is the detected energy spectrum with unity area: ES(E)dE=1. μen,W(E)/ρW and μW are the energy-dependent mass energy absorption coefficient and the linear attenuation coefficient of water, respectively, both taken from the NIST database [

      Hubbell JH, Seltzer SM. Tables of X-ray mass attenuation coefficients and mass energy-absorption coefficients (version 1.4), [Online] Available: http://physics.nist.gov/xaamdi. National Institute of Standards and Technology, Gaithersburg, MD; 2004.

      ], and r is the radius of the phantom.

      References

        • Yang M.
        • Virshup G.
        • Clayton J.
        • Zhu X.R.
        • Mohan R.
        • Dong L.
        Theoretical variance analysis of single- and dual-energy computed tomography methods for calculating proton stopping power ratios of biological tissues.
        Phys Med Biol. 2010; 55: 1343-1362
        • Yang M.
        • Zhu X.R.
        • Park P.C.
        • Titt U.
        • Mohan R.
        • Virshup G.
        • et al.
        Comprehensive analysis of proton range uncertainties related to patient stopping-power-ratio estimation using the stoichiometric calibration.
        Phys Med Biol. 2012; 57: 4095-4115
        • Alvarez R.E.
        • Macovski A.
        Energy-selective reconstructions in X-ray computerized tomography.
        Phys Med Biol. 1976; 21: 733-744https://doi.org/10.1088/0031-9155/21/5/002
        • Rutherford R.A.
        • Pullan B.R.
        • Isherwood I.
        Measurement of effective atomic number and electron density using an EMI scanner.
        Neuroradiology. 1976; 11: 15-21
        • Hounsfield G.N.
        Computerized transverse axial scanning (tomography). Part I. Description of system.
        Br J Radiol. 1973; 46: 1016-1022
        • Heismann B.J.
        • Schmidt B.
        • Flohr T.
        Spectral computed tomography.
        SPIE PRESS, Bellingham, Washington, USA2012
        • Saito M.
        Potential of dual-energy subtraction for converting CT numbers to electron density based on a single linear relationship.
        Med Phys. 2012; 39: 2021-2030https://doi.org/10.1118/1.3694111
        • Hünemohr N.
        • Krauss B.
        • Tremmel C.
        • Ackermann B.
        • Jäkel O.
        • Greilich S.
        Experimental verification of ion stopping power prediction from dual energy CT data in tissue surrogates.
        Phys Med Biol. 2014; 59: 83-96
        • Hünemohr N.
        • Paganetti H.
        • Greilich S.
        • Jäkel O.
        • Seco J.
        Tissue decomposition from dual energy CT data for MC based dose calculation in particle therapy.
        Med Phys. 2014; 41: 061714https://doi.org/10.1118/1.4875976
        • Landry G.
        • Seco J.
        • Gaudreault M.
        • Verhaegen F.
        Deriving effective atomic numbers from DECT based on a parameterization of the ratio of high and low linear attenuation coefficients.
        Phys Med Biol. 2013; 58: 6851-6866https://doi.org/10.1088/0031-9155/58/19/6851
        • Bourque A.E.
        • Carrier J.F.
        • Bouchard H.
        A stoichiometric calibration method for dual energy computed tomography.
        Phys Med Biol. 2014; 59: 2059-2088https://doi.org/10.1088/0031-9155/59/8/2059
        • Farace P.
        Comment: experimental verification of ion stopping power prediction from dual energy CT data in tissue surrogates.
        Phys Med Biol. 2014; 59: 7081-7084https://doi.org/10.1088/0031-9155/59/1/83
        • van Abbema J.K.
        • van Goethem M.J.
        • Greuter M.J.W.
        • van der Schaaf A.
        • Brandenburg S.
        • van der Graaf E.R.
        Relative electron density determination using a physics based parameterization of photon interactions in medical DECT.
        Phys Med Biol. 2015; 60: 3825-3846
        • Garcia L.I.R.
        • Azorin J.F.P.
        • Almansa J.F.
        A new method to measure electron density and effective atomic number using dual-energy CT images.
        Phys Med Biol. 2016; 61: 265-279
        • Han D.
        • Siebers J.V.
        • Williamson J.F.
        A linear, separable two-parameter model for dual energy CT imaging of proton stopping power computation.
        Med Phys. 2016; 43: 600-612
        • Möhler C.
        • Wohlfahrt P.
        • Richter C.
        • Greilich S.
        Range prediction for tissue mixtures based on dual-energy CT.
        Phys Med Biol. 2016; 61: N268-N275
        • Taasti V.T.
        • Petersen J.B.B.
        • Muren L.P.
        • Thygesen J.
        • Hansen D.C.
        A robust empirical parametrization of proton stopping power using dual energy CT.
        Med Phys. 2016; 43: 5547-5560https://doi.org/10.1118/1.4962934
        • van Elmpt W.
        • Landry G.
        • Das M.
        • Verhaegen F.
        Dual energy CT in radiotherapy: current applications and future outlook.
        Radiother Oncol. 2016; 119: 137-144https://doi.org/10.1016/j.radonc.2016.02.026
        • Hansen D.C.
        • Seco J.
        • Srensen T.S.
        • Petersen J.B.B.
        • Wildberger J.E.
        • Verhaegen F.
        • et al.
        A simulation study on proton computed tomography (CT) stopping power accuracy using dual energy CT scans as benchmark.
        Acta Oncol. 2015; 54: 1638-1642
        • Hünemohr N.
        • Krauss B.
        • Dinkel J.
        • Gillmann C.
        • Ackermann B.
        • Jäkel O.
        • et al.
        Ion range estimation by using dual energy computed tomography.
        Z Med Phys. 2013; 23: 300-313
        • Zhu J.
        • Penfold S.N.
        Dosimetric comparison of stopping power calibration with dual-energy CT and single-energy CT in proton therapy treatment planning.
        Med Phys. 2016; 43: 2845-2854
        • Tremblay J.-É.
        • Bedwani S.
        • Bouchard H.
        A theoretical comparison of tissue parameter extraction methods for dual energy computed tomography.
        Med Phys. 2014; 41: 081905https://doi.org/10.1118/1.4886055
        • Kanematsu N.
        • Inaniwa T.
        • Koba Y.
        Relationship between electron density and effective densities of body tissues for stopping, scattering, and nuclear interactions of proton and ion beams.
        Med Phys. 2012; 39: 1016-1020
      1. ICRP. Adult Reference Computational Phantoms. ICRP Publication 110. Ann. ICRP 39 (2).

        • Schneider U.
        • Pedroni E.
        • Lomax A.
        The calibration of CT Hounsfield units for radiotherapy treatment planning.
        Phys Med Biol. 1996; 41: 111-124
      2. Sigmund P, Schinner A, Paul H. Errata and Addenda for ICRU Report 73. Stopping of ions heavier than helium. J ICRU 2009;5(1).

      3. ICRU. Stopping Power and Ranges for Protons and Alpha Particles. ICRU Report No. 49.

        • Jan S.
        • Santin G.
        • Strul D.
        • Staelens S.
        • Assié K.
        • Autret D.
        • et al.
        GATE: a simulation toolkit for PET and SPECT.
        Phys Med Biol. 2004; 49: 4543-4561
        • Rit S.
        • Vila Oliva M.
        • Brousmiche S.
        • Labarbe R.
        • Sarrut D.
        • Sharp G.C.
        The Reconstruction Toolkit (RTK), an open-source cone-beam CT reconstruction toolkit based on the Insight Toolkit (ITK).
        J Phys Conf Ser. 2014; 489: 012079https://doi.org/10.1088/1742-6596/489/1/012079
        • Brooks R.A.
        • Di Chiro G.
        Beam hardening in X-ray reconstructive tomography.
        Phys Med Biol. 1976; 21: 390-398https://doi.org/10.1088/0031-9155/21/3/004
        • Jackson D.F.
        • Hawkes D.J.
        X-ray attenuation coefficients of elements and mixtures.
        Phys Rep. 1981; 70: 169-233
        • Bazalova M.
        • Carrier J.-F.
        • Beaulieu L.
        • Verhaegen F.
        Dual-energy CT-based material extraction for tissue segmentation in Monte Carlo dose calculations.
        Phys Med Biol. 2008; 53: 2439-2456https://doi.org/10.1088/0031-9155/53/9/015
        • Yang M.
        Dual energy computed tomography for proton therapy treatment planning.
        ([Ph.D. dissertation]) The University of Texas, 2011
      4. Brendel B, Bergner F, Brown K, Koehler T. Penalized likelihood decomposition for dual layer spectral CT. In: Proc 4th intl mtg on image formation in X-ray CT; 2016. p. 41–4.

        • Chuang K.-S.
        • Huang H.K.
        Comparison of four dual energy image decomposition methods.
        Phys Med Biol. 1988; 33: 455-466https://doi.org/10.1088/0031-9155/33/4/005
        • Vilches-Freixas G.
        • Létang J.M.
        • Brousmiche S.
        • Romero E.
        • Oliva M.V.
        • Kellner D.
        • et al.
        Technical note: procedure for the calibration and validation of kilo-voltage cone-beam CT models.
        Med Phys. 2016; 43: 5199-5204https://doi.org/10.1118/1.4961400
        • Arbor N.
        • Dauvergne D.
        • Dedes G.
        • Létang J.M.
        • Parodi K.
        • Quiñones C.T.
        • et al.
        Monte Carlo comparison of X-ray and proton CT for range calculations of proton therapy beams.
        Phys Med Biol. 2015; 60: 7585
      5. P. Linstrom, W. Mallard, NIST Standard Reference Database Number 69, [Online] Available:http://webbook.nist.gov, National Institute of Standards and Technology, Gaithersburg, MD. (retrieved December 6, 2016).

        • Nelder J.A.
        • Mead R.
        A simplex method for function minimization.
        Comput J. 1965; 7: 308-313
      6. ICRU, Photon, electron, proton and neutron interaction data for body tissues, ICRU Report No. 46.

      7. Berger MJ, Hubbell JH, Seltzer SM, Chang J, Coursey JS, Sukumar R, et al., XCOM: photon cross section database (version 1.5), [Online] Available: http://physics.nist.gov/xcom, National Institute of Standards and Technology, Gaithersburg, MD; 2010.

      8. Hubbell JH, Seltzer SM. Tables of X-ray mass attenuation coefficients and mass energy-absorption coefficients (version 1.4), [Online] Available: http://physics.nist.gov/xaamdi. National Institute of Standards and Technology, Gaithersburg, MD; 2004.