Experimental validation of Monte Carlo based treatment planning system in bone density equivalent media

Accuracy of the TPS Monaco ver. 5.11 calculation algorithm was experimentally verified using a semi-anthropomorphic CIRS Thorax phantom (CIRS Inc., Norfolk, VA, USA) consisting of a body made of water equivalent material (ρ = 1.003 g/cm³), lung equivalent parts (ρ = 0.207 g/cm³), and bone equivalent part (ρ = 1.506 g/cm³) with cylindrical holes for placement of ionization chamber into interchangeable rod inserts having three different densities.² The phantom was scanned using the Somatom Open CT simulator. Acquired CT images were used for the delineation of volumes of interest and subsequent dose calculations. Measurements of absorbed dose were performed at ten measuring positions within the phantom (Figure 1) for 15 different irradiation set-ups (Table 1), using a PTW Farmer-type ionization chamber. All measurements were carried out at the central part of the selected radiation fields, excluding the regions of high dose gradients.

Figure 1

Photo of the semi-anthropomorphic CIRS Thorax phantom with interchangeable rod inserts (left) and its CT image (right). Positions of 10 interchangeable rod inserts are marked with numbers from 1 to 10. Five measuring points are in the water equivalent part of the phantom (grey area), four points are in the lung density equivalent material (black area), and one point is in the bone density equivalent part of the phantom (white area).

Measured doses were compared to the corresponding doses obtained by both calculation options, D_m and D_w. Dose differences δDm$\delta {{D}_{m}}$and δDw$\delta {{D}_{w}}$between measured and calculated values for dose-to-medium and dose-to-water calculation approach, were calculated according to the IAEA methodology^1,2 as:

[1]δDm=100⋅Dm−DmeasDmeas,ref$$\delta {{D}_{m}}=100\cdot \frac{{{D}_{m}}-{{D}_{meas}}}{{{D}_{meas,ref}}}$$

[2]δDw=100⋅Dw−DmeasDmeas,ref$$\delta {{D}_{w}}=100\cdot \frac{{{D}_{w}}-{{D}_{meas}}}{{{D}_{meas,ref}}}$$

where D_meas denotes measured absorbed dose at the selected measuring point, while D_meas,ref stands for the absorbed dose measured at the reference point, which was chosen on the central axis of the beam at the isocenter (Table 1).

Dose differences δDm$\delta {{D}_{m}}$and δDw$\delta {{D}_{w}}$between calculated and measured doses were analysed for both calculation options through the comparison of the respective average values δD¯m${{\overline{\delta D}}_{m}}$and δD¯w${{\overline{\delta D}}_{w}}$

[3]δD¯m=1n∑i=1i=nδDm,i$${{\overline{\delta D}}_{m}}=\frac{1}{n}\sum\limits_{i=1}^{i=n}{\delta {{D}_{m,i}}}$$

[4]δD¯w=1n∑i=1i=nδDw,i$${{\overline{\delta D}}_{w}}=\frac{1}{n}\sum\limits_{i=1}^{i=n}{\delta {{D}_{w,i}}}$$

The index i stands for a particular dose difference for i-th dose measurement and corresponding calculated dose for two different calculation modes in the selected part of the CIRS Thorax phantom (water equivalent part, lung equivalent part, or bone density equivalent part).

Throughout the study, all calculations within Monaco TPS were performed on a 0.2 cm calculation grid, with 0.5% statistical uncertainty per control point.

In the second part of the study, we were aiming to determine differences between D_m and D_w calculation approaches in the Monaco ver. 5.11 TPS in the bone equivalent part of the CIRS Thorax phantom, following the same methodology as described in the preceding section.

Three irradiation geometries (single asymmetric rectangular fields having different gantry angles: 0°, 90°, and 180°) were selected for this part of the study (Table 1, set-ups 6, 7, and 8). For each of those irradiation geometries, two phantom assemblies were used to analyze differences between the two calculation approaches with respect to the measurements performed by PTW 30013 Farmer type ionization chamber in the bone density equivalent (BDE) part of the CIRS Thorax phantom. In the first assembly, referred to as non-standard, the water equivalent insert with the ionization chamber was placed into the BDE part of the phantom (Figure 2A). In this way, the measuring point in the phantom was surrounded by water equivalent material of sufficient thickness to fulfill conditions required by the Bragg-Gray cavity theory for the determination of absorbed dose in terms of dose to water. In the second assembly, referred to as standard, the BDE insert was placed in the BDE part of the phantom (Figure 2B).

Figure 2

CT image of the CIRS Thorax phantom: water equivalent insert inside BDE part of the phantom (A); a BDE insert inside bone density equivalent (BDE) part of the phantom (B) and cross-section of small “water cylinders” of different dimensions delineated inside BDE part of the phantom to find limits for calculating geometry where cavity theory applies (top right).

In the last part of the study, the phantom assembly was additionally virtually modified for the calculation purposes in the Monaco TPS: cylinders of various volumes (constant length and different diameters) were delineated inside the BDE insert on the CT scans (Figure 2, top right). This approach was utilized to obtain the limits above which the differences between D_m and D_w calculation approaches become non-significant and in agreement with experimentally determined absorbed doses. The length of the cylinders was set equal to the length of the cavity volume of the PTW 30013 ionization chamber, while the electron density of such cylinders was set to be equal to the electron density of the water. According to the IAEA TRS-398 Code of practice²², the charge measured by an ionization chamber calibrated in terms of absorbed dose to water is directly proportional to the absorbed dose in water at the point of measurement in the absence of the chamber. By delineating cylinders having the electron density of water inside the BDE part of the phantom, we have tried to simulate the mentioned theoretical situation to different degrees.

To verify the accuracy of dose-to-medium and dose-to-water calculation modes, we have analyzed differences δDm$\delta {{D}_{m}}$and δDw$\delta {{D}_{w}}$between calculated and measured absorbed doses for both calculation modes and different volumes of “water cylinders” smaller than the volume of the PTW 30013 ionization chamber’s cavity volume (0.6 cm3), using Eqs. [1] to [4]. We were aiming to find the volume of “water cylinder,” above which there will exist an agreement between calculated and measured doses without a statistically significant difference between both calculation approaches. Our final challenge was to define an addendum to the existing TPS verification methodology based on the described method and experimental findings from the present work.

The uncertainty of δD¯m${{\overline{\delta D}}_{m}}$was estimated as the combination in quadrature of the statistical uncertainty of δD¯m${{\overline{\delta D}}_{m}}$and the uncertainty of Monte Carlo calculation of 0.5% (1 SD) for D_m, using a coverage factor k = 2 (2 SD). The uncertainty of δD¯w${{\overline{\delta D}}_{w}}$was calculated in the same manner.

We considered that the D_m and D_w calculation modes differed significantly within 95% confidence limits (two standard deviations – 2 SD, i.e., coverage factor k = 2) if the relation

[5]|δD¯m−δD¯w|>uc(k=2)$$\left| {{\overline{\delta D}}_{m}}-{{\overline{\delta D}}_{w}} \right|>{{u}_{c}}\left( k=2 \right)$$

was satisfied. u_c is a combined uncertainty which was determined as the combination in quadrature of the individual uncertainties of δD¯m${{\overline{\delta D}}_{m}}$and δD¯w${{\overline{\delta D}}_{w}}$. This estimation was considered conservative due to the fact that the uncertainties of the terms D_meas/D_meas,ref were included in the compute of the individual uncertainties δD¯m${{\overline{\delta D}}_{m}}$and δD¯w.${{\overline{\delta D}}_{w}}.$

Secondly, we considered that the dose calculations within Monaco TPS were in agreement with the experimentally determined doses if the conditions

[6]|δD¯m|<1%$$\left| {{\overline{\delta D}}_{m}} \right|<1\text{%}$$

[7]|δD¯w|<1%$$\left| {{\overline{\delta D}}_{w}} \right|<1\text{%}$$

were satisfied. At this point we note, that throughout the rest of the paper all combined uncertainties are stated within two standard deviations, i.e., using a coverage factor k = 2.

Differences between calculated and measured absorbed doses for two calculation modes, dose-to-medium D_m and dose-to-water D_w were determined using Eqs. [1] and [2] for all 15 standard irradiation configurations and ten measurement points in the CIRS Thorax semi-anthropomorphic phantom (Table 1). Mean values of percentage dose differences δD¯m${{\overline{\delta D}}_{m}}$and δD¯w${{\overline{\delta D}}_{w}}$calculated by Eqs. [3] and[4] are presented with corresponding uncertainties in terms of two standard deviations in Figure 3, separately for the water equivalent part (five measurement points), lung density equivalent part (four measurement points), and BDE part (one measurement point) of the phantom. Statistical significance of the obtained differences between δD¯m${{\overline{\delta D}}_{m}}$and δD¯w${{\overline{\delta D}}_{w}}$was evaluated using the relations shown in Eqs. [5] to [7].

Figure 3

Mean percentage dose differences δD¯m${{\overline{\delta D}}_{m}}$and δD¯w${{\overline{\delta D}}_{w}}$between calculated and measured doses in different parts of the CIRS Thorax phantom (water, lung, and bone density equivalent materials) for both calculation options built in the Monaco TPS: dose-to-medium D_m and dose-to-water D_w. Error bars represent corresponding combined uncertainties.

Comparison of measured and calculated doses in the water equivalent part of the phantom showed that the mean percentage dose difference for all points was - 0.6% (uc=1.1%)$\left( {{u}_{c}}=1.1\text{%} \right)$for the dose-to-medium calculation mode and - 0.6% (uc=1.1%)$\left( {{u}_{c}}=1.1\text{%} \right)$for the dose-to-water calculation mode (Figure 3). The two calculations were found not to be significantly different within 95% confidence limits since the condition from Eq. [5] was not satisfied: |δD¯m−δD¯w|=0.0%(uc=1.5%).$\left| {{\overline{\delta D}}_{m}}-{{\overline{\delta D}}_{w}} \right|=0.0\text{%}\left( {{u}_{c}}=1.5\text{%} \right).$

Comparison of measured and calculated doses in the lung density equivalent part of the phantom showed that δD¯m=0.1%(uc=1.1%)${{\overline{\delta D}}_{m}}=0.1\text{%}\left( {{u}_{c}}=1.1\text{%} \right)$for the dose-to-medium calculation approach, while δD¯w${{\overline{\delta D}}_{w}}$= 0.0% (uc=1.1%)$\left( {{u}_{c}}=1.1\text{%} \right)$for the dose-to-water mode (Figure 3). Also in this case, the difference between both applied calculation approaches was statistically non-significant within 95% confidence limits: δD¯m−δD¯w=0.1%uc=1.5%. $\left| {{\overline{\delta D}}_{m}}-{{\overline{\delta D}}_{w}} \right|=0.1\text{%}\left( {{u}_{c}}=1.5\text{%} \right).$

In the bone density equivalent part of the CIRS Thorax phantom, the percentage dose differences between the two calculation options were larger than in the previous two cases (Figure 3). Mean difference δD¯m${{\overline{\delta D}}_{m}}$for the dose-to-medium calculation mode was - 2.8% (uc=2.0%),$\left( {{u}_{c}}=2.0\text{%} \right),$while for the dose-to-water calculation approach the mean difference δD¯w${{\overline{\delta D}}_{w}}$was 2.9% (uc=1.8%).$\left( {{u}_{c}}=1.8\text{%} \right).$Consequently and importantly, in the BDE part of the phantom, the absolute differences between the two calculation modes were found to be statistically significant within 95% confidence limits: |δD¯m−δD¯w|=5.7%(uc=2.6%).$\left| {{\overline{\delta D}}_{m}}-{{\overline{\delta D}}_{w}} \right|=5.7\text{%}\left( {{u}_{c}}=2.6\text{%} \right).$

Dose calculations within Monaco TPS were in agreement with experimentally determined doses for water equivalent and lung equivalent parts of the CIRS Thorax phantom, since the conditions from Eqs. [6] and [7] were satisfied. On the contrary, for the BDE part of the phantom, conditions from Eqs. [6] and [7] were not satisfied. Therefore, we can conclude that the dose calculations in Monaco TPS ver. 5.11 were not in agreement with measured absorbed doses for the BDE part of the phantom, regardless of the calculation mode.

The second part of the study was focused on the differences between calculated and measured absorbed doses in the BDE part of the CIRS Thorax phantom. Three simple asymmetric fields with different gantry angles were selected for that purpose utilizing two different phantom assemblies, standard and non-standard, as described in the section Materials and methods and shown in Table 1 (set-ups 6, 7, and 8) and Table 2.

Table 2

Irradiation geometry (field, gantry)	Phantom assembly	δDm[%]$\delta {{D}_{m}}\left[ \text{%} \right]$	δDw[%]$\delta {{D}_{w}}\left[ \text{%} \right]$
(6+8) x 15 cm2Gantry = 0°	standard(2)	- 2.9	2.9
	non-standard(3)	- 0.7	- 0.2
(3+8) x 15 cm2Gantry = 90°	standard(4)	- 3.0	5.1
	non-standard(5)	- 0.7	- 0.1
(4+10) x 15 cm2Gantry = 180°	standard(6)	- 5.7	5.4
	non-standard(7)	0.5	1.3

For non-standard phantom geometry, we did not find statistically significant differences between measured and calculated absorbed doses: δD¯m=−0.3%(uc=1.3%)${{\overline{\delta D}}_{m}}=-0.3\text{%}\,\,\,\,\,\left( {{u}_{c}}=1.3\text{%} \right)$and δD¯w=0.3%${{\overline{\delta D}}_{w}}=0.3\text{%}$(uc=1.3%).$\left( {{u}_{c}}=1.3\text{%} \right).$In this case, the absolute difference |δD¯m−δD¯w|$\left| {{\overline{\delta D}}_{m}}-{{\overline{\delta D}}_{w}} \right|$between both approaches was 0.6% and was statistically non-significant within 95% confidence limits (uc=1.8%).$\left( {{u}_{c}}=1.8\text{%} \right).$

In the standard phantom geometry, however, the differences δD¯m${{\overline{\delta D}}_{m}}$and δD¯w${{\overline{\delta D}}_{w}}$between measured and calculated doses were larger and statistically significant (Table 2). δD¯m=−3.9%(uc=2.1%)${{\overline{\delta D}}_{m}}=-3.9\text{%}\left( {{u}_{c}}=2.1\text{%} \right)$and δD¯w=4.4%(uc=1.9%).${{\overline{\delta D}}_{w}}=4.4\text{%}\left( {{u}_{c}}=1.9\text{%} \right).$

The absolute value of the difference between both approaches was in this case statistically significant: |δD¯m−δD¯w|=8.3%(uc=2.8%).$\left| {{\overline{\delta D}}_{m}}-{{\overline{\delta D}}_{w}} \right|=8.3\text{%}\,\,\,\left( {{u}_{c}}=2.8\text{%} \right).$

As a final point, we investigated differences between calculated and measured doses in the phantom, which was virtually modified for the calculation purposes, as described in the section Materials and methods. Results for five delineated “water cylinders,” including the results for standard geometry (V = 0 cm³), are presented in Table 3. Differences gradually decrease as the volumes of delineated “water cylinders” become larger. The maximal difference was |δD¯m−δD¯w|=8.3%,$\left| {{\overline{\delta D}}_{m}}-{{\overline{\delta D}}_{w}} \right|=8.3\text{%},$for V = 0 cm³ (i.e., BDE plug without delineated “water cylinder”). The smallest difference of 0.1% between δD¯m${{\overline{\delta D}}_{m}}$and δD¯w${{\overline{\delta D}}_{w}}$was found for the largest investigated “water cylinder” of volume 0.573 cm³. This difference was statistically non-significant within 95% confidence limits (uc=1.9%).$\left( {{u}_{c}}=1.9\text{%} \right).$

Table 3

V [cm3]	δD¯m${{\overline{\delta D}}_{m}}$[%]	uc,m[%]	δD¯w${{\overline{\delta D}}_{w}}$[%]	uc,w[%]
0	- 3.9	2.1	4.4	1.9
0.035	- 2.6	1.5	2.5	1.9
0.141	- 1.4	1.3	1.8	1.7
0.279	- 0.3	1.2	1.2	1.5
0.573	0.3	1.4	0.4	1.3

Differences between calculated and measured doses in the water equivalent part of the CIRS Thorax semi-anthropomorphic phantom were within 1% and not significantly different from zero (Eqs. [6] and [7]), regardless of the applied calculation option. The latter is in good agreement with previously published data.^9,16 Similarly, in lung density equivalent material, the calculated mean percentage dose differences were not significantly different than zero for both calculation modes, confirming the results from previously published studies.^3,9,13

The differences between the two calculation approaches, dose-to-medium and dose-to-water, were, however, significant in BDE media (Table 2 and Figure 3). Andreo et al.⁹ have shown that a 10% difference in ICRP bone can be expected for Monaco ver. 5.0 TPS between two calculation modes after conversion of D_m to D_w. Results of the present study confirm those findings as well as the opposite signs of mean percentage dose differences for D_m and D_w reporting modes in the case when Monaco ver. 5.11 TPS has been used. Considerable differences between calculated dose distributions using D_m and D_w calculation approaches have also been reported in clinical studies.^15,23

In the BDE part of the CIRS Thorax semi-anthropomorphic phantom, mean percentage dose differences δD¯m${{\overline{\delta D}}_{m}}$and δD¯w${{\overline{\delta D}}_{w}}$were calculated by applying Eqs. [1] and [2] for two phantom assemblies - standard and non-standard and three selected irradiation geometries, as shown in Table 2. In the case of nonstandard geometry, both δD¯m${{\overline{\delta D}}_{m}}$and δD¯w${{\overline{\delta D}}_{w}}$were within 1%, demonstrating that there is a negligible difference between applied calculation modes.

However, differences between the respective mean values δD¯m=−3.9%${{\overline{\delta D}}_{m}}=-3.9\text{%}$and δD¯w=4.4%${{\overline{\delta D}}_{w}}=4.4\text{%}$were statistically significant if standard geometry was utilized. The latter case was also assumed as our first result in the part of the study where we attempted to find the volume of “water cylinder” delineated in the Monaco ver. 5.11 TPS for which the difference between δD¯m${{\overline{\delta D}}_{m}}$and δD¯w${{\overline{\delta D}}_{w}}$would become non-significant.

For further discussion, analysis, and graphical presentation, the exponential function was selected to fit the data from Table 3. The general form of the fitting function is given as

[8]y=a+b⋅ecx$$y=a+b\cdot {{e}^{cx}}$$

with fitting coefficients, a, b, and c. The dependent variable y denotes average values δD¯m${{\overline{\delta D}}_{m}}$and δD¯w,${{\overline{\delta D}}_{w}},$while x stands for volumes of delineated “water cylinders”. The explicit forms of the exponential fitting functions obtained were

[9]δD¯m=0.397−3.995⋅e−6.274.V$${{\overline{\delta D}}_{m}}=0.397-3.995\cdot {{e}^{-6.274.V}}$$

[10]δD¯w=0.526+3.510⋅e−8.131⋅V$${{\overline{\delta D}}_{w}}=0.526+3.510\cdot {{e}^{-8.131\cdot V}}$$

for D_m and D_w reporting modes, respectively. Both functions from Eqs. [9] and[10] are graphically presented in Figure 4 having residual standard errors of the fit equal to 0.340% and 0.165% (on two degrees of freedom) for D_m and D_w calculation modes, respectively.

Figure 4

Average differences δD¯m${{\overline{\delta D}}_{m}}$and δD¯w${{\overline{\delta D}}_{w}}$between calculated and measured doses in the bone density equivalent (BDE) part of the CIRS Thorax phantom, as a function of the volumes of the simulated “water cylinders” (see Figure 2 and Table 3). δD¯m${{\overline{\delta D}}_{m}}$and δD¯w${{\overline{\delta D}}_{w}}$are presented as individual values/points calculated using Eqs. [1] to[4], and in the form of two analytical functions from Eqs. [9] and [10]. Error bars represent corresponding uncertainties within 95% confidence limits.

Applying Eqs. [9] and[10] for large volumes, we can see that δD¯m${{\overline{\delta D}}_{m}}$and δD¯w${{\overline{\delta D}}_{w}}$converge to the values of the free fitting coefficients a, i.e., δD¯m=am≅${{\overline{\delta D}}_{m}}={{a}_{m}}\cong $0.397% and δD¯w=aw≅${{\overline{\delta D}}_{w}}={{a}_{w}}\cong $0.526%. a_m and a_w denote free fitting coefficients in Eqs. [9] and [10], respectively. Those values are non-significantly different from zero, thus in agreement with experimentally determined absorbed doses. From the latter observations, we can deduct two key facts, which form a basis for the recommended additional procedure to the existing methodology for the verification of the accuracy of the Monte Carlo based TPS we were aiming at. Briefly:

• Differences |δD¯m−δD¯w|$\left| {{\overline{\delta D}}_{m}}-{{\overline{\delta D}}_{w}} \right|$between dose-to-medium and dose-to-water calculation approaches gradually fade away as the volumes of “water cylinders” become larger and closer to the volume of the Farmer chamber;

• |δD¯m|$\left| {{\overline{\delta D}}_{m}} \right|$and |δD¯w|$\left| {{\overline{\delta D}}_{w}} \right|$fall below 1% for volumes of delineated “water cylinders” larger than 0.3 cm³.

Irrespective of the fact that the ionization chamber is calibrated in terms of dose to water, we propose an additional verification test of the accuracy of the Monaco TPS calculation modes for BDE regions considering the mentioned observations:

• One can select three simple irradiation geometries (single fields, different gantry angles) and perform measurements of absorbed doses with the Farmer type ionization chamber in the BDE part of CIRS Thorax semi-anthropomorphic phantom, using a BDE insert (“standard” geometry). The ionization chamber should be positioned at the central part of the radiation field, where the measured signal is sufficiently large.

• Measured doses are compared to the calculated ones using both calculation modes, D_m and D_w applying Eqs. [1] to [4] for the additional four “water cylinders” delineated in the TPS.

• Obtained mean values δDm¯$\overline{\delta {{D}_{m}}}$and δD¯w${{\overline{\delta D}}_{w}}$of the percentage dose differences are fitted by the analytical function from Eq. [8].

Finally, the acceptability of the tested TPS algorithm is based on two conditions, which have to be fulfilled concurrently:

• i)

  Differences |δD¯m−δD¯w|$\left| {{\overline{\delta D}}_{m}}-{{\overline{\delta D}}_{w}} \right|$between dose-to-medium and dose-to-water calculation approaches should fall within 1% for the “water cylinder” of volume 0.6 cm³,i.e.,

[11]|δD¯m(V=0.6cm3)−δD¯w(V=0.6cm3)|<1%$$\left| {{\overline{\delta D}}_{m}}\left( V=0.6\,c{{m}^{3}} \right)-{{\overline{\delta D}}_{w}}\left( V=0.6\,c{{m}^{3}} \right) \right|<1\text{%}$$

Fulfilment of this condition means that both calculation options yield to the same results within statistical uncertainty for large volumes, as expected. Since significant differences do exist for small volumes of delineated “water cylinders,” we have to consider this fact as well. The maximal difference |δD¯m−δD¯w|$\left| {{\overline{\delta D}}_{m}}-{{\overline{\delta D}}_{w}} \right|$can be obtained from the corresponding fitting functions for V = 0 cm³ (in our study, the maximal difference between both calculation options was 7.6%).

• ii)

  Obtained values |δD¯m|$\left| {{\overline{\delta D}}_{m}} \right|$and |δD¯w|$\left| {{\overline{\delta D}}_{w}} \right|$have to fall below 1% (see Eqs [6] and [7]) for large volumes of delineated “water cylinders”. If this condition is fulfilled, one can conclude that TPS dose calculations are in agreement with experimentally determined doses for both calculation modes.

It is important to note that our investigation was limited to the region of charged particle equilibrium (CPE) and for 6 MV photon beam only.