Experimental validation of Geant4 nuclear interaction models in dose calculations of therapeutic carbon ion beams

Abstract

Background: The choice of nuclear interaction models in Monte Carlo simulations affects the dose calculation accuracy for light ion beam therapy.

Purpose: This study aimed to evaluate the dose calculation accuracy and simulation time of three GATE-RTiON/Geant4 physics lists for therapeutic carbon ion beams, assessing their suitability for independent dose calculation in patient-specific quality assurance (PSQA).

Methods: The normalized beam models for physics lists QGSP_BIC_HP_EMZ, QGSP_INCLXX_HP_EMZ, and Shielding_EMZ were validated against measurements regarding the accuracy of range, spot size and reference dose. Normalized transversal dose profiles ( $D/D_{max}$ ) and field size factor (FSF) were compared with measurements. The accuracy of simulated target dose in 103 fields (various energies, field sizes, depths, and dose gradient ${\nabla }_D$ complexity) of energy-modulated scanned beams was evaluated at 3181 positions. The median of global dose difference $med\left({\Delta }_D\right)$ was calculated at different depth ranges.

Results: The three physics lists with validated beam models showed similar accuracy in $D/D_{max}$ and FSF in the Bragg peak region and proximal depths, while QGSP_INCLXX_HP agreed most closely for $D/D_{max}$ in the fragmentation tail. Accounting for ${\nabla }_D$ -related uncertainty, $med\left({\Delta }_D\right)$ remained within ±1.1% for QGSP_INCLXX_HP, while exhibiting an overall increasing trend with depth for QGSP_BIC_HP (up to 2.3%) and a decreasing trend for Shielding (down to -4.1%), respectively. By tuning the number-of-primaries/monitor unit conversion ( $k_{N/\mathrm{MU}}$ ) as a function of energy, $med\left({\Delta }_D\right)$ of QGSP_BIC_HP was reduced to within ±1.3%, at the cost of reduced accuracy in the simulated reference dose. The simulation time of Shielding was 1.8 times that of QGSP_BIC_HP and 1.5 times that of QGSP_INCLXX_HP.

Conclusions: QGSP_INCLXX_HP demonstrated high dosimetric accuracy in the target region of energy-modulated fields. QGSP_BIC_HP and Shielding showed physics model-related inaccuracies in simulated target dose. Additional $k_{N/\mathrm{MU}}$ tuning improved their target dose calculation accuracy with a trade-off of reference dose accuracy. The computationally efficient QGSP_INCLXX_HP and QGSP_BIC_HP are viable candidates for dose calculation applications of carbon ion beam therapy, such as in silico PSQA.

Keywords: GATE‐RTiON; Geant4; carbon ion radiotherapy; dosimetric accuracy.

FIGURE 1

The measured and simulated IRPD of HBL in water, normalized to the measured dose at 14 mm depth. Each curve is plotted till $1.5R_{80}$, with the nominal energy (MeV/u) marked above the Bragg peak. The uncertainty of simulated dose is represented by the error band. BIC: QGSP_BIC_HP, INCLXX: QGSP_INCLXX_HP, BERT: QGSP_BERT_HP. HBL; horizontal beamline, IRPD; integrated radial profiles as a function of depth.

FIGURE 2

The relative difference of simulated and measured dose at reference condition $D_{ref}$ as a function of nominal energy for HBL (left) and VBL (right). BIC: QGSP_BIC_HP, INCLXX: QGSP_INCLXX_HP, BERT: QGSP_BERT_HP. The subscript ref indicates the normalized beam models, and N/MU indicates additional N/MU conversion tuning. The gray error band represents the relative standard uncertainty of measured $D_{ref}$, $u_{D_{ref}}^{meas}=$ 2.7%. The gray lines represent ±3%. HBL; horizontal beamline, VBL; vertical beamline.

FIGURE 3

The normalized transversal dose profile ($D/D_{max}$) of 402.8 MeV/u carbon ion beam plotted over the lateral distance off the central axis (mm). (a) to (f) represent profiles at various depths in water. ${\sigma }_c$ represents core width. The uncertainty of measured $D/D_{max}$ (black) is shown in the upper axes with shaded region. The uncertainties of $\mathrm{\Delta }(D/D_{max})$ for QGSP_BIC_HP (cyan), QGSP_INCLXX_HP (yellow), and Shielding (red) are shown in the bottom axes with shaded region of the line color. The uncertainty representations are smaller than the symbols in many cases.

FIGURE 4

The normalized transversal dose profile ($D/D_{max}$) of 200.0 MeV/u carbon ion beam plotted over the lateral distance off the central axis (mm). (a) to (e) represent profiles at various depths in water. ${\sigma }_c$ represents core width. The uncertainty of measured $D/D_{max}$ (black) is shown in the upper axes with shaded region. The uncertainties of $\mathrm{\Delta }(D/D_{max})$ for QGSP_BIC_HP (cyan), QGSP_INCLXX_HP (yellow), and Shielding (red) are shown in the bottom axes with shaded region of the corresponding color. The uncertainty representations are smaller than the symbols in many cases.

FIGURE 5

The distribution of ${\mathrm{\Delta }}_D$ for type I fields comparing the three configurations with normalized beam models (left) and the two with tuned N/MU conversion (right). The planned measurement positions ranged from 34.0 to 248.9 mm. The number of points in each depth group is marked beneath the boxes. The light gray and gray horizontal lines indicate ±3% and ±7%, respectively. The lines with markers indicate the $med({\mathrm{\Delta }}_D)$ of the depth group. Both the notches and the error band indicate the $C{I}_{95}$, that is, uncertainty of $med({\mathrm{\Delta }}_D)$. The statistics of ${\mathrm{\Delta }}_D$ is reported in Table S5 in the Supporting Information.

FIGURE 6

The distribution of ${\mathrm{\Delta }}_D$ (left) and $\widehat{{\mathrm{\Delta }}_D}$ (right) for type II fields. The planned measurement positions ranged from 21.9 to 265.9 mm. The number of points in each depth group is marked beneath the boxes. The light gray and gray horizontal lines indicate ±3% and ±7%, respectively. The lines with markers indicate the $med({\mathrm{\Delta }}_D)$ and $med(\widehat{{\mathrm{\Delta }}_D})$ of the depth group. The notches indicate the $C{I}_{95}$ of ${\mathrm{\Delta }}_D$ and $\widehat{{\mathrm{\Delta }}_D}$. The error bands with half width of $|med({\mathrm{\Delta }}_D)-med(\widehat{{\mathrm{\Delta }}_D})|$ indicated the uncertainty of $med({\mathrm{\Delta }}_D)$ of the corresponding bin. The statistics of ${\mathrm{\Delta }}_D$ and $\widehat{{\mathrm{\Delta }}_D}$ are reported in Tables S6 and S7.

FIGURE 7

The distribution of ${\mathrm{\Delta }}_D$ (left) and $\widehat{{\mathrm{\Delta }}_D}$ (right) for type III fields. The planned measurement positions ranged from 11.9 to 274.7 mm. The number of points in each depth group is marked beneath the boxes. The light gray and gray horizontal lines indicate ±3% and ±7%, respectively. The lines with markers indicate the $med({\mathrm{\Delta }}_D)$ and $med(\widehat{{\mathrm{\Delta }}_D})$ of the depth group. The notches indicate the $C{I}_{95}$ of ${\mathrm{\Delta }}_D$ and $\widehat{{\mathrm{\Delta }}_D}$. The error bands with half width of $|med({\mathrm{\Delta }}_D)-med(\widehat{{\mathrm{\Delta }}_D})|$ indicated the uncertainty of $med({\mathrm{\Delta }}_D)$ of the corresponding bin. The statistics of ${\mathrm{\Delta }}_D$ and $\widehat{{\mathrm{\Delta }}_D}$ are reported in Tables S8 and S9.

FIGURE 8

The predicted ($\widehat{t}$) and actual ($t$) normalized simulation time of the monoenergetic pencil beam training plans with and without range shifter, and the verification plans. Simulations were performed for QGSP_BIC_HP_EMZ with ${\mathrm{BIC}}_{ref}$ beam model. The gray line represents $\widehat{t}=t$.