Vol 28, No 3 (2023)
Research paper
Published online: 2023-07-06

open access

Page views 1045
Article views/downloads 340
Get Citation

Connect on Social Media

Connect on Social Media

research paper

Reports of Practical Oncology and Radiotherapy

2023, Volume 28, Number 3, pages: 416–428

DOI: 10.5603/RPOR.a2023.0044

Submitted: 30.12.2022

Accepted: 06.07.2023

© 2023 Greater Poland Cancer Centre.

Published by Via Medica.

All rights reserved.

e-ISSN 2083–4640

ISSN 1507–1367

Optimal values of the Electron Monte Carlo dose engine parameters

Jacek Wendykier1Ewa Wojtyna2Barbara Bekman1Adam Bekman3Bożena Woźniak3Beata Niewiadomska3Joanna Prażmowska3Piotr Wendykier4Marcin Bieniasiewicz5Aleksandra Grządziel1
1Radiotherapy Planning Department, Maria Skłodowska-Curie National Research Institute of Oncology, Gliwice Branch, Gliwice, Poland
2Department of Medical Physics, St. John of Dukla Cancer Center, Lublin, Poland
3Medical Physics Department, Maria Skłodowska-Curie National Research Institute of Oncology, Gliwice Branch, Gliwice, Poland
4Independent Scholar, Warsaw, Poland
5Radiotherapy Center, Department of Medical Physics, Multidisciplinary Hospital, Gorzów Wielkopolski, Poland

Address for correspondence: Ewa Wojtyna, St. John of Dukla Cancer Center, Department of Medical Physics, Lublin, Poland, e-mail: ewaa.wojtyna@gmail.com

This article is available in open access under Creative Common Attribution-Non-Commercial-No Derivatives 4.0 International (CC BY-NC-ND 4.0) license, allowing to download articles and share them with others as long as they credit the authors and the publisher, but without permission to change them in any way or use them commercially

Abstract
Background: The aim of this study was to indicate the most favorable in terms of to the time of calculation and the uncertainty of determining the dose distribution values of the parameters for the Electron Monte Carlo (eMC) algorithm in the Eclipse treatment planning system.
Materials and methods: Using the eMC algorithm and the variability of the values of its individual parameters, calculations of the electron dose distribution in the full-scattering virtual water phantom were performed, obtaining percentage depth doses, beam profiles, absolute dose values in points and calculation times. The reference data included water tank measurements such as relative dose distributions and absolute point doses.
Results: For 63 sets of calculation data created from selected values of the parameters for the eMC algorithm, calculation times were analyzed and the absolute calculated and measured doses were compared. Performing a statistical analysis made it possible to determine whether the differences in the values of deviations between the actual dose and the calculated dose in individual regions of the percentage depth dose curve and the beam profile are statistically significant between the analyzed sets of parameters.
Conclusions: Taking into account obtained results from the analysis of the discrepancy between the distribution of the calculated and measured dose, the correspondence of the absolute value of the calculated and measured dose and the duration of the calculation of the dose distribution, the optimal set of parameters was indicated for the eMC algorithm which allows obtaining the dose distribution and the number of monitor units in an acceptable time.
Key words: Electron Monte Carlo algorithm; calculation time; uncertainty of dose determination
Rep Pract Oncol Radiother 2023;28(3):416–428

Introduction

The uncertainty of the calculated dose distribution by the computerized treatment planning system affects the actual dose that a patient will receive during the therapeutic session. Compliance between the measured (reference) value of the dose and its calculated value is extremely important in order to carry out radiation therapy in both safe and effective [1, 2] way. The use of increasingly accurate algorithms in computerized treatment planning systems plays an important role in determining the three-dimensional dose distribution for different anatomical structures of a patient and different beam geometries based on the data obtained from the measurements. Their task is to reproduce the dose distribution in a precise way with regard to the reference data [3]. A breakthrough solution for calculating the dose distribution, conditioned by increasing speed and efficiency of computers, turned out to be the Monte Carlo computational method (MC), which allows regardless of the complexity of geometry the interaction of particles that make up the therapeutic beam to be simulated [4].

Monte Carlo is a commonly used method in radiotherapy, allowing for precise reproduction of the interaction of radiation with matter. The use of this technique in recent decades has increased almost exponentially [4]. In this method, a large number of simulated particles is used to statistically obtain dose distributions from a predetermined computational uncertainty. It is worth noting that Monte Carlo is the only simulation method that allows the implementation of all interactions of particles occurring in a heterogeneous medium such as the human body [5, 6]. Simulations for different particles (electrons, photons, etc.), energy ranges, their centers of motion and types of interactions are possible due to the variety of MC simulation codes [5, 7].

The Electron Monte Carlo (eMC) algorithm available in the Eclipse treatment planning system (Varian Medical Systems, Palo Alto, CA) is an efficient implementation of the MC method used to calculate the electron beam dose distribution [8]. This algorithm allows the user to select the values of seven parameters such as: statistical uncertainty (floating-point number from 0.2 to 8.0), size of the calculation grid (0.10, 0.15, 0.20, 0.25, 0.30, 0.40, 0.50 cm), starting number of random generator, number of particle history, dose threshold for uncertainty (10, 20, 30, 40, 50, 60 70, 80%), smoothing method (no smoothing, 2-D Median and 3-D Gaussian), smoothing options and smoothing level (low, medium and strong). Only four of them have a significant impact on both the calculation time and the uncertainty of determining the dose distribution: statistical uncertainty, the size of the calculation grid, the smoothing method and the level of smoothing [9].

The effect of the dose threshold for uncertainty on the time of calculation was not in the area of interest of this study. In this work, it was assumed that in the eMC algorithm the threshold can be selected to the clinical situation. If the only interest is to cover, for example, isodose 80%, and there are no critical organs nearby, then one can take a threshold of 80% and then the volume of tissue with a dose greater than the threshold will have a dose determined with a certain uncertainty. If the clinical situation is such that a critical organ is close to the irradiated volume, then a smaller threshold value can be selected to achieve a reliable estimate of the dose absorbed in the critical organ under consideration.

These parameters of the eMC algorithm have been described in detail in the instructions prepared by the manufacturer [9] and in publications [6, 10].

The aim of the study was to find optimal values for the parameters of the Electron Monte Carlo algorithm in the Eclipse treatment planning system for an electron beam with an energy of 6 MeV as a function of calculation time and computational uncertainty determination of the dose distribution.

Materials and methods

Calculations of the dose distribution in the computerized treatment planning system using the Electron Monte Carlo v. 15.6 algorithm, implemented in the Eclipse treatment planning system, and measurements made during system verification, were carried out for the TrueBeam linear accelerator from Varian Medical Systems.

Measurements

The reference data were acquired during the verification process of the treatment planning system before the accelerator was approved for clinical work. The measurements, including both relative dose distributions and absolute point doses, were made in the full-scattering water phantom PTW MP3 (PTW Freiburg). The selected geometry is SSD equal to 100.0 cm and all original applicators with square fields: 6 × 6, 10 × 10, 15 × 15, 20 × 20, 25 × 25 cm2.

The transverse and the longitudinal beam profiles were obtained using a Semiflex ionization chamber (PTW Freiburg) with an active volume of 0.125cm3 (directed parallel to the beam axis). The depth of maximum dose, the dose at depth of 80% and 50% relative to the normalized maximum dose were chosen which correspond to values R100 = 1.3 cm, R80 = 1.9 cm and R50 = 2.4 cm. The resolution of acquired data was 2.0 mm in the plateau and the umbra area and 1.0 mm in the penumbra region. The percentage depth dose measurements with a step of 1.0 mm were performed using a plane parallel ionization chamber (Markus type) recommended by international protocols as the most suitable for its spatial resolution and presence of guarding ring. The absolute dose values were calculated from measurement with a plane parallel chamber, which is a reference tool for absolute dosimetry of electron beams, especially for low energies such as 6 MeV. The measuring depth was zref, recommended by Report 398 [11].

Calculations

The dosimetric model of the accelerator head has been prepared based on only the mandatory data. These data include the air profile for a 40 × 40 cm2 field without an applicator and a percentage depth dose for the same field in the water. For each applicator, a characterization of dose change with depth and the dose rate at the depth of the maximum dose are necessary. Calculations of the electron dose distribution were made in a water phantom created in the treatment planning system with dimensions of 50 × 50 × 50 cm3. To achieve speedup, calculations were performed on the graphics cards and were divided into several parallel processes. Percentage depth doses and beam profiles were extracted from exported three-dimensional dose distributions from the phantom (RTDOSE files in DICOM format [12]). The dose values in points and the calculation time were obtained directly from the Eclipse treatment planning system. The entire phantom was covered by a calculation grid. The effect of the external size of the calculation grid on the calculation time has not been studied.

Due to the wide range of statistical uncertainty (SU) values, preliminary studies were conducted. In order to estimate the effect of statistical uncertainty on the calculation time and quality of the output, calculations were carried out for SU = 0.2, 0.25, 0.3, 0.4, 0.5, 0.75, 1.0, 2.0 and 3.0%. The limit of 3% is due to the fact that this is the recommended tolerance for profiles in the therapeutic field. The calculations were carried out for a grid of 0.25 cm, applicators: 6 × 6, 15 × 15 and 25 × 25 cm2 and SSD = 100.0 cm at 80% uncertainty limit and without smoothing the output. Eight computation threads were used for the calculations.

Based on preliminary results (Tab. 1), three statistical uncertainty values (1%, 2%, 3%) were selected for further studies due to the shortest calculation times, relevant especially for the largest applicator. The times given are averaged values from two measurements. The selected statistical uncertainty values were then checked in sets with other parameters affecting the calculation time and the uncertainty of the dose distribution: calculation grid size (0.10 cm, 0.25 cm, 0.50 cm), smoothing method (none, 2-D Median, 3-D Gaussian), smoothing level (low, medium, high). Table 2 shows that 63 data sets were obtained from calculations described above. These sets were then analyzed both in terms of the duration of the dose distribution calculations in the Eclipse computerized treatment planning system, the consistency of the calculated and measured absolute dose values, as well as the dose distributions (relative values) determined experimentally with the results of computer calculations, indicating in the end the optimal set of parameters of the eMC algorithm for groups of parameters selected above.

Table 1. Calculation time as a function of statistical uncertainty values. The calculation was carried out using 8 threads and for the grid of 0.25 cm without smoothing with a threshold of 80%. The time measurement was repeated twice

SU

<t>

t(25x25)/t(6x6)

t(25x25)/t(15x15)

t(15x15)/t(6x6)

(%)

[s]

6 ×6

15 ×15

25 ×25

0.20

342

2220

6177

18.1

2.8

6.5

0.25

236

1456

4073

17.3

2.8

6.2

0.30

178

1036

2849

16.0

2.8

5.8

0.40

116

590

1575

13.6

2.7

5.1

0.50

77

392

1063

13.9

2.7

5.1

0.75

55

190

486

8.8

2.6

3.5

1.00

43

136

323

7.5

2.4

3.2

2.00

34

58

124

3.6

2.1

1.7

3.00

33

48

82

2.5

1.7

1.5

Table 2. The time of calculation of the dose distribution by the computerized treatment planning system for individual sets of parameters of the Electron Monte Carlo (eMC) algorithm. In the last column in the fields with the value “1”, the absolute value of the difference between the maximum calculated dose Dmax,cal and the reference Dmax,ref is within the specified uncertainty interval for each applicator

No.

Statistical uncertainty [%]

Grid size [cm]

Smoothing method

Level of smoothing

Time [s]

|Dmax,cal - Dmax,ref|

1.

1

0.10

5749

1

2.

1

0.25

192

1

3.

1

0.50

73

1

4.

2

0.10

2496

1

5.

2

0.25

106

0

6.

2

0.50

63

1

7.

3

0.10

2031

1

8.

3

0.25

62

0

9.

3

0.50

62

0

10.

1

0.10

2-D Median

Low

5840

1

11.

1

0.25

2-D Median

Low

220

0

12.

1

0.50

2-D Median

Low

80

0

13.

2

0.10

2-D Median

Low

2667

1

14.

2

0.25

2-D Median

Low

107

1

15.

2

0.50

2-D Median

Low

67

0

16.

3

0.10

2-D Median

Low

2225

1

17.

3

0.25

2-D Median

Low

96

0

18.

3

0.50

2-D Median

Low

62

0

19.

1

0.10

2-D Median

Medium

0

20.

1

0.25

2-D Median

Medium

197

0

21.

1

0.50

2-D Median

Medium

69

0

22.

2

0.10

2-D Median

Medium

4927

0

23.

2

0.25

2-D Median

Medium

113

0

24.

2

0.50

2-D Median

Medium

57

0

25.

3

0.10

2-D Median

Medium

0

26.

3

0.25

2-D Median

Medium

108

0

27.

3

0.50

2-D Median

Medium

56

0

28.

1

0.10

2-D Median

Strong

0

29.

1

0.25

2-D Median

Strong

196

0

30.

1

0.50

2-D Median

Strong

66

0

31.

2

0.10

2-D Median

Strong

0

32.

2

0.25

2-D Median

Strong

146

0

33.

2

0.50

2-D Median

Strong

54

0

34.

3

0.10

2-D Median

Strong

3163

0

35.

3

0.25

2-D Median

Strong

119

0

36.

3

0.50

2-D Median

Strong

52

0

37.

1

0.10

3-D Gaussian

Low

6214

1

38.

1

0.25

3-D Gaussian

Low

183

1

39.

1

0.50

3-D Gaussian

Low

67

0

40.

2

0.10

3-D Gaussian

Low

2772

1

41.

2

0.25

3-D Gaussian

Low

93

1

42.

2

0.50

3-D Gaussian

Low

67

0

43.

3

0.10

3-D Gaussian

Low

2142

1

44.

3

0.25

3-D Gaussian

Low

88

0

45.

3

0.50

3-D Gaussian

Low

55

1

46.

1

0.10

3-D Gaussian

Medium

0

47.

1

0.25

3-D Gaussian

Medium

188

0

48.

1

0.50

3-D Gaussian

Medium

69

0

49.

2

0.10

3-D Gaussian

Medium

4864

0

50.

2

0.25

3-D Gaussian

Medium

99

0

51.

2

0.50

3-D Gaussian

Medium

58

0

52.

3

0.10

3-D Gaussian

Medium

3525

0

53.

3

0.25

3-D Gaussian

Medium

78

0

54.

3

0.50

3-D Gaussian

Medium

55

0

55.

1

0.10

3-D Gaussian

Strong

6038

0

56.

1

0.25

3-D Gaussian

Strong

216

0

57.

1

0.50

3-D Gaussian

Strong

73

0

58.

2

0.10

3-D Gaussian

Strong

2786

0

59.

2

0.25

3-D Gaussian

Strong

98

0

60.

2

0.50

3-D Gaussian

Strong

62

0

61.

3

0.10

3-D Gaussian

Strong

2140

0

62.

3

0.25

3-D Gaussian

Strong

97

0

63.

3

0.50

3-D Gaussian

Strong

60

0

Calculations were carried out for all five applicators simultaneously and time was measured for such calculations. Due to the comparison of the profile with the depth of R50, a threshold of statistical dose uncertainty of 50% was chosen.

Comparison of measured dose distributions and those calculated in the treatment planning system was performed according to the recommendations published by the Polish Society of Medical Physics [13] and based on international recommendations [14]. The Alfard package [15] was the main software used for the processing of the measurement data.

Wendykier-1.png
Figure 1. Graphical representation of the calculation time as a function of statistical uncertainty values, error bars denote 2 standard deviations (SD)
Wendykier-2.png
Figure 2. Beam profile as the statistical uncertainty (SU) function for 6 MeV at dmax for 10 × 10 cm2 applicator (normalization in the beam axis intentionally done so that the profiles do not obscure each other and noise becomes visible)

For the purpose of this work, the measurements and calculations were normalized to 100% in the central axis of beam profiles and up to 100% at the maximum point in the case of depth doses. As recommended, the measurement data were divided into areas differentiated according to the size of the dose and its gradient regions with literature designations δ1, δ2, δ3, δ4. These regions are presented on the percentage depth dose curve (Fig. 3) and electron beam profile (Fig. 4). Alfard software made it possible to find boundaries between these areas, based on a gradient value of 0.5%/mm [16, 17].

Wendykier-3.png
Figure 3. The percent depth dose (PDD) curves of the electron beam dose with marked areas where deviations between calculations and measurements are compared. A tenfold increased gradient for the PDDs are also shown. Graphs for energies of 6 MeV (top panel) and 22 MeV (bottom panel) are presented
Wendykier-4.png
Figure 4. Beam profile with marked regions where deviations between calculations and measurements are compared using different tolerance values

Profiles from the treatment planning system as well as from measurements have been symmetrically transformed relative to the beam axis. The differences between the measurement and calculated graphs for specific intervals were computed in the last step, and the resulting data were used for analysis. The above method was described in more detail, on the example of determining the boundaries of the region of a large and small gradient of the photon beam profile, in the works [16, 17].

The available literature does not provide a way to divide the electron beam percent depth dose (PDD) into the appropriate regions. Therefore, in this work, a division into three regions was proposed based on the value of the dose and its gradient. The first region stretched from the surface of the phantom to a point with a gradient value of 0.5%/mm, but deeper than the maximum, and it was defined as a high dose-small gradient (δ1). The second region was the region of a large gradient (δ2) and included a dose gradient greater than 0.5%/mm, located deeper than the maximum dose. The third low-dose region was the deepest and was defined by a gradient smaller than the assumed limit value (δ4). Schematically, these areas are depicted in Figure 3. The dashed vertical lines represent the boundaries of the regions.

In the case of low energies, the gradient in the first part of the slope is greater than the accepted limit value (top panel Fig. 3), but in order to maintain consistently one tolerance value in a given region, it was decided to go for the above solution and mark the region as high dose and small gradient. Energies other than 6 MeV were not in the area of interest for this study, but it should be noted that for high energies the situation is different. In the first part of the PDD curve, the gradient is definitely smaller, as can be shown for a beam with the energy of 22 MeV (bottom panel Fig. 3). Another argument for adopting such a definition of this area is Table 3. It shows that the tolerance of 10.0% is recommended for a large gradient, which seems inappropriate for a region that is important from a therapeutic point of view.

Table 3. Recommended deviation values for different regions of the electron beam [14]

Region designation

Type and location of the region

Deviation value for simple geometry, homogeneous environment

δ1

High dose, small dose gradient

Central axis of the radiation beam

2.0%

δ2

High dose, large dose gradient

Dose drop area in the central axis, penumbra region

10.0%

δ3

High dose, small dose gradient

Outside the central axis of the radiation beam

3.0%

δ4

Low dose, small dose gradient

Region beyond the beam boundary

3.0%

Both the differences between the profiles and the percentage depth doses were calculated with a resolution of 1.0 mm. The value of the deviation of the δ between the point of the target dose and the measured dose is expressed as a percentage and is determined according to the formula (1):

154655.png (1)

where:

Dc dose calculated at a given point of the phantom,

Dm dose measured at the same point in the phantom,

D dose:

  • Dr reference dose in the beam axis for the profile,
  • Dmax maximum dose for PDD.

For specific regions: δ1, δ2, δ3, δ4, the mean deviation δmean was determined in the Alfard software, due to the fact that each region contains many points and then the average confidence limit Δmean [14], also expressed as a percentage (%), which is described by the formula (2):

154666.png (2)

where:

δmean mean value of the deviation,

σ standard deviation of the mean.

For individual sets of parameters of the eMC algorithm, the values of the mean confidence limits Δ were obtained. To assess the discrepancy between the actual and the calculated values in individual regions, the criteria set out in TRS 430, which are defined in Table 3, were used [13, 14].

In order to find the sets in which the deviations between the calculated and measured dose distributions are the smallest, the average of the obtained error values (mean confidence limits Δmean) for both beam profiles and PDDs was calculated, and then the obtained values were added up. For 7 sets with the least errors, a statistical analysis was performed, which checked whether there were significant statistical differences in error values between them. The statistical hypotheses were verified using the Jamovi package version 1.1.9.0 [18]. First, the analysis of basic descriptive statistics was performed along with the Shapiro-Wolf normality test, followed by the one-way analysis of variance. The threshold value was α = 0.05.

For some of the variables analyzed, the result of the Shapiro-Wolf test turned out to be statistically significant, which means that these distributions deviate from the normal distribution. However, in none of the cases analyzed, the skew value exceeded the contractual absolute value of 2, therefore the distributions were considered sufficiently close to normal distributions [19]. Due to this, parametric statistical tests were used. In order to test the assumption of the equality of variance, a Levene test was performed. In the case when the result of this test was statistically significant, i.e. there were differences between the variances in the studied groups, the result of the univariate analysis of variance was read from the Welch correction (the variances are not equal), while otherwise the Fisher correction was used (the variances are equal). Based on the results obtained from the univariate analysis of variance, it was determined whether the compared sets of eMC parameters differed statistically significantly in a given region (δ1, δ2, δ3, δ4) of the beam profile and the percentage depth dose. Then, in order to detect cases that differed statistically significantly, post-hoc tests were used with the Games-Howell correction (due to uneven variances) or with the Tukey correction (due to equal variances).

In order to analyze the absolute dose values, 10 simulations of calculations were performed in a computerized treatment planning system for several sets of parameters of the eMC algorithm to determine the standard deviation of the maximum dose Dmax (SD). Based on the results obtained, an absolute dose uncertainty (2SD) of 0.02 Gy was assumed. Then, the absolute value of the difference between the maximum calculated dose Dmax,c and the reference Dmax,ref for all cases was calculated, after which it was checked whether it falls within the assumptions of this 2SD uncertainty range.

Results

Preliminary tests

The results of the preliminary tests are presented in Table 1 and in Figure 1. In Table 1 also includes calculation time quotients for different applicators at the same SU. For the selected four SU values, beam profiles are presented (Fig. 2).

Fundamental research

In order to indicate the optimal set of parameters of the eMC algorithm, an analysis was carried out for 63 sets of computational data (Tab. 2), taking into account the calculation time, the discrepancies obtained between the calculated and measured dose distributions and the results obtained from comparisons of the absolute value of the calculated and reference doses.

The results of the time of the calculation of the dose distribution by the computerized treatment planning system for individual cases are presented in Table 2. The last column shows aggregated data for the given dataset. Number one means that for each applicator, the absolute value of the dose difference calculated and measured was less than 0.02 Gy, which is 2% of the set dose at dmax. Number zero indicates those cases in which for at least one applicator this difference exceeded the set value.

In Table 4 the results of the preliminary analysis are presented. Based on these results, seven sets of parameters for the eMC algorithm were selected, in which the deviations between the distribution of the calculated and measured doses were the smallest. For the indicated sets, a statistical analysis was carried out, based on which it was determined whether the obtained differences in deviation values are statistically significant between sets. The analysis was performed for individual areas of electron beam profiles and from the percentage depth doses.

Table 4. Sets of electron Monte Carlo (eMC) algorithm parameters with the smallest deviation values along with the corresponding time of calculation of the dose distribution by the computerized treatment planning system. Determinations: the mean value of all obtained mean confidence limits for electron beam profiles; the mean value of all the mean confidence limits obtained for the percentage depth doses

No.

Statistical uncertainty (%)

Size nets [cm]

Smoothing method

Level Honing

μprof

μPDD

μprof + μPDD

Time

1.

1

0.25

3-D Gaussian

Low

2.36

1.70

4.06

3 mins 3 s

2.

1

0.25

2-D Median

Low

2.65

1.65

4.30

3 mins 40 s

3.

2

0.25

3-D Gaussian

Low

2.44

2.05

4.49

1 mins 33 s

4.

1

0.25

lack

2.62

1.90

4.52

3 mins 12 s

5.

1

0.25

3-D Gaussian

Medium

1.96

2.73

4.69

3 mins 8 s

6.

1

0.10

3-D Gaussian

Low

2.16

2.59

4.75

1 h 43 mins 34 s

7.

2

0.10

3-D Gaussian

Low

2.15

2.61

4.76

46 mins 12 s

In order to verify the hypothesis of a significant impact of individual sets of eMC parameters on the differences in deviation values between the actual dose and the calculated dose in the area of δ4(l), δ2(l), δ3, δ2(r), δ4(r) of the beam profile, the univariate analysis of variance was performed and, on its basis, it was found that the deviations between the actual dose and the calculated dose in the compared sets differ statistically significantly in the area of δ4(l), δ2(l), δ2(r), δ4(r) of electron beam profiles, i.e. in the so-called shadow and penumbra regions (l left side of the analyzed profile, r right side). However, in the most important region, which is the therapeutic region δ3, the differences between the analyzed sets are not statistically significant. In order to find groups that differ from each other, multiple post-hoc comparisons were carried out. On their basis, a statistically significant difference was observed in individual areas between the eMC parameter sets. Based on the results obtained from the statistical analysis, it was found that the best set which differs statistically significantly from each other in the area of δ4(l), δ2(l), δ2(r), δ4(r) is a set for statistical uncertainty of 1%, size of the calculation grid: 0.25 cm, smoothing methods: 3-D Gaussian, level of smoothing: mean. In each of the analyzed regions, the mean value of the deviation between the actual dose and the calculated dose for this case is the smallest.

In the case of the percentage depth dose, a one-way analysis of variance was also performed to test another hypothesis of a significant influence of individual sets of eMC parameters on the differences in deviation values between the reference dose and calculated dose in the area of δ1, δ2, δ4. Analogous calculations were made as in the case of beam profiles and, on their basis, it was found that the deviations between the actual dose and the calculated dose in the compared sets were statistically significantly different in all analyzed regions of the PDDs. Based on the results obtained from the multiple post-hoc comparisons and the average error values obtained for the compared sets, which differ statistically significantly, it was concluded that the set for statistical uncertainty of 1%, the size of the calculation grid: 0.25 cm, the method of smoothing: 3-D Gaussian, the mean level of smoothing, is the least favorable for region δ1 and δ4 (scattered photons). Significant statistical differences between the sets were also observed in the area of δ2 percentage depth doses, in which cases for statistical uncertainty of 1%, calculation grid size: 0.10 cm, smoothing methods: 3-D Gaussian, level of smoothing: low and statistical uncertainty of 2%, calculation grid size: 0.10 cm, smoothing methods: 3-D Gaussian, level of smoothing: low have the highest average error value.

Discussion

The calculation time can be roughly described by a power function, dependent on SU, as presented on Figure 1. In the case of the smallest applicator, the increase in calculation time caused by the change in the SU value from 3.0 to 0.2% is eleven times. For the largest applicator, the same change causes the calculation to take 76 times longer. Correlation of preliminary data from Table 1 with Figure 2 allows a conclusion that the best results (the lowest noise) are obtained for very long times. However, in routine clinical work, such times would be unacceptable, especially for large applicators. Therefore, for further studies, SU values were adopted, which gave acceptable results and computation time was not too long.

In a simplified mental model it can be assumed that voxels with a high dose will only be in the area covered by the applicator. The “layer thickness” of voxels of a certain dosage in the direction of depth in the phantom is approximately constant for given beam energy and does not depend on the size of the applicator. In such a model, the calculation time will be proportional to the surface area of the applicator. The quotient of the applicator surface areas is 252/62 = 17.4, 252/152 = 2.8 and 152/62 = 6.3. Such a relationship is approximately satisfied only for the smallest SU values and tends to approach unity with an increase in the value of statistical uncertainty.

From the obtained results, it can be concluded that it is not possible to clearly indicate the set of parameters that would be fully satisfactory. Changing one of the parameters may improve one of the characteristics of the beam, but at the same time it results in deterioration of another.

For the duration of the calculation (Tab. 2) under constant conditions of the other parameters for the eMC algorithm, it was found that the statistical uncertainty and the size of the calculation grid have an impact on the calculation speed. As both the statistical uncertainty value and the mesh size increase, the time it takes for the computerized treatment planning system to calculate the dose distribution is reduced. Changing the level, as well as the method of smoothing, does not significantly affect the duration of the calculation. In the case of the smallest grid (0.10 cm) and the smallest analyzed statistical uncertainty value, the calculation lasted about 1 h 40 min., regardless of the other parameters. Changing the size of the grid to 0.25 cm caused the reduction of the calculation time to less than 4 min., which is a significant change. If the grid was maintained at 0.10 cm and the statistical uncertainty changed to 2.0%, the calculation time was about 40 min. In two cases, the calculations took more than an hour, but this was due to other computational processes sharing resources. For the statistical uncertainty of 3.0% at 0.10 cm grid, the calculation times are also of the order of thirty minutes. In all other cases, the calculation time is in the order of a few minutes. In several cases in Table 2 information about the calculation time is missing (rows 19, 25, 28, 31, 46) because the system indicated errors and stopped the calculation after more than 2 hours of calculation. The error was repeatable. The reason for the behaviour of the Eclipse system remained unknown. For such cases the calculations were performed individually for single applicators to receive point dose and dose distribution.

For profiles (Tab. 4), the only set with an average level of smoothing out of the selected seven turned out to be the best suited for the measurements. In that case, the average of all obtained mean confidence limits Δmean for electron beam profiles was below 2.00%. The highest value of the average confidence limits was obtained by the set with a low level of smoothing, but using the 2D Median method (2.65%). Just behind it, with a value equal to 2.62%, was the case without smoothing. In the case of smoothing profiles, a small value of the calculation grid had a positive effect, regardless of statistical uncertainty, and gave a result of 2.16 and 2.15%, respectively (the last two items in Tab. 4).

The winning set for the profiles simultaneously reached the highest value of the mean confidence limits Δmean for the percentage depth doses (2.73%). Not much worse are the sets ranked 2nd and 3rd in the profile comparisons, with 2.61% and 2.59%, respectively. The set of parameters best reflecting the dose decrease with depth reaching an average confidence limit of Δmean of 1.65% and simultaneously was the set with the worst model fit for the measurements for the profiles.

It has been observed that an increase in the level of smoothing causes an increase in the discrepancy between the maximum dose measured at the point and the maximum dose calculated by the computerized treatment planning system (Tab. 2). For the medium and strong levels of smoothing, regardless of the smoothing method used, not a single set was observed for which the calculated differences would be less than those set for each of the applicators. Another factor that affects these discrepancies is the calculation grid. The smaller the grid size, the better the fit. Out of the sixteen sets in which there was a consensus for all applicators, nine had a grid of 0.10 cm. There were three cases in which the grid was equal to 0.50 cm. Therefore, the best results were obtained for unsmoothed data sets (for which the level of smoothing does not matter) or smoothed ones with a low level.

Conclusions

The above analysis of matched profiles and dose decreases confirms the opposite effect of specific parameters on the obtained results. Improving the fit of profiles with the help of smoothing, clearly causes a deterioration in the matching of percentages depth doses. In addition, smoothing results in greater discrepancies between absolute dose values in points.

Taking into account the duration of the calculations, the correspondence of the absolute dose calculated with the measured one within the determined uncertainty limits and the results obtained from the statistical analysis of relative doses, the optimal values of the parameters of the eMC algorithm are as follows: statistical uncertainty: 2.0%, calculation grid size: 0.25 cm, 3-D Gaussian smoothing method with low level of smoothing. Such a set of parameters allows the distribution of the dose and the number of monitor units to be obtained in a clinically acceptable time of several minutes. Based on the obtained results, it can also be concluded that statistical uncertainty, the size of the calculation grid, the method of smoothing and the level of smoothing have a significant impact on the dose distribution, and the first two factors also mostly affect time of calculation.

When analyzing the results, a suggestion arises that the process of planning the distribution of the dose with electron beams should be carried out in two stages. In the first one, the values of the parameters for the Electron Monte Carlo algorithm described in the above paragraph can be used to quickly obtain a satisfactory result, not fundamentally deviating from the one confirmed by measurements, without waiting too long for the results of individual geometry changes (necessary changes of arm angle, collimator angle, applicator dimension, frame shape, use of bolus, its various thicknesses, etc.), which are an integral part of the treatment planning process. At this stage, the planning process can be completed, bearing in mind the optimal compliance of the data calculated with the measurements.

However, purists creating treatment plans may then feel unsatisfied with the noise visible on the profile with a low level of smoothing and a relatively high SU value. Therefore, following the example of other treatment planning systems, a second step can be taken, the effect of which would be the distribution with the least possible level of noise and, at the same time, the most convergent with the measurements. To obtain such a final distribution, the SU value should be significantly reduced, preferably to the minimum value and then the results obtained would be more similar to the measurements, both in terms of the three-dimensional dose distribution and the designated number of monitor units. This may involve several dozen minutes of calculations, but with the final geometry of the plan already in place, only this single calculation remains to be performed.

Acknowledge

To Tadeusz Jedynak for providing Alfard software.

Ethical permission

Ethical approval was not necessary for the preparation of this article.

Conflicts of interest

The authors have no conflict of interest to declare.

Funding

This publication was prepared without any external source of funding.

References

  1. Mijnheer BJ, Battermann JJ, Wambersie A. What degree of accuracy is required and can be achieved in photon and neutron therapy? Radiother Oncol. 1987; 8(3): 237–252, doi: 10.1016/s0167-8140(87)80247-5, indexed in Pubmed: 3107087.
  2. Dutreix A. When and how can we improve precision in radiotherapy? Radiother Oncol. 1984; 2(4): 275–292, doi: 10.1016/s0167-8140(84)80070-5, indexed in Pubmed: 6522641.
  3. Malicki J, Ślosarek K. Planowanie leczenia i dozymetria w radioterapii. 2nd ed. Via Medica, Gdańsk 2016.
  4. Seco J, Verhaegen F. Monte Carlo Techniques in Radiation Therapy. 1st ed. CRC Press 2016.
  5. Matuszak N. [Monte Carlo jako jedna z metod symulacyjnych w radioterapii]. Lett Oncol Sci. 2019; 16(2): 15–22, doi: 10.21641/los.2019.17.2.91.
  6. Tertel J, Wulff J, Karle H, et al. Verification of a commercial implementation of the Macro-Monte-Carlo electron dose calculation algorithm using the virtual accelerator approach. Z Med Phys. 2010; 20(1): 51–60, doi: 10.1016/j.zemedi.2009.11.001, indexed in Pubmed: 20211424.
  7. Konefał A. [Symulacje metodą Monte Carlo za pomocą oprogramowania GEANT4]. Postępy Fizyki. 2006; 57(6): 242–251.
  8. Chamberland E, Beaulieu L, Lachance B. Evaluation of an electron Monte Carlo dose calculation algorithm for treatment planning. J Appl Clin Med Phys. 2015; 16(3): 4636, doi: 10.1120/jacmp.v16i3.4636, indexed in Pubmed: 26103470.
  9. Eclipse Photon and Electron Algorithms Reference Guide, P1020505-003-C, Varian Medical Systems (accessed: July 2018).
  10. Zhang A, Wen N, Nurushev T, et al. Comprehensive evaluation and clinical implementation of commercially available Monte Carlo dose calculation algorithm. J Appl Clin Med Phys. 2013; 14(2): 4062, doi: 10.1120/jacmp.v14i2.4062, indexed in Pubmed: 23470937.
  11. IAEA Technical Reports Series no. 398.. Absorbed Dose Determination in External Beam Radiotherapy: An International Code of Practice for Dosimetry based on Standards of Absorbed Dose to Water. International Atomic Energy Agency, Vienna 2006.
  12. DICOM PS3.3 2023b — Information Object Definitions. http://dicom.nema.org/medical/dicom/current/output/chtml/part03/sect_C.8.8.3.html (17.02.2022).
  13. Dybek M, Winiecki J, Iwanicki T, et al. [Kontrola Systemów Planowania Leczenia 3D W Radioterapii Wiązkami Zewnętrznymi Fotonów I Elektronów]. Pol J Med Phys Eng. 2014; 20(1): 1–32, doi: 10.2478/pjmpe-2014-0001.
  14. IAEA Technical Reports Series no. 430. . Commissioning and Quality Assurance of Computerized Planning Systems for Radiation Treatment of Cancer. International Atomic Energy Agency, Vienna 2004.
  15. ALFARD. http://alfard.eu5.net/alfard/index.html (08.01.2022).
  16. Wendykier J, Bieniasiewicz M, Grządziel A, et al. Determination of boundaries between ranges of high and low gradient of beam profile. Rep Pract Oncol Radiother. 2016; 21(3): 168–173, doi: 10.1016/j.rpor.2015.12.007, indexed in Pubmed: 27601946.
  17. Wendykier J, Grządziel A, Bekman B, et al. Using beam profile inflection point in process of treatment planning system verification. Rep Pract Oncol Radiother. 2021; 26(4): 553–562, doi: 10.5603/RPOR.a2021.0049, indexed in Pubmed: 34434571.
  18. Jamovi. Open statistical software for the desktop and cloud. https://www.jamovi.org (17.10.2020).
  19. George D, Mallery P. IBM SPSS statistics 23 step by step: a simple Guide and Reference. Routledge, New York 2016.