Prioritized efficiency optimization for intensity modulated proton therapy

A high dosimetric quality and short treatment time are major goals in radiotherapy planning. Intensity modulated proton therapy (IMPT) plans obtain dose distributions of great conformity but often result in long delivery times which are typically not incorporated into the optimization process. We present an algorithm to optimize delivery efficiency of IMPT plans while maintaining plan quality, and study the potential trade-offs of these interdependent objectives. The algorithm is based on prioritized optimization, a stepwise approach to implemented objectives. First the quality of the plan is optimized. The second step of the prioritized efficiency optimization (PrEfOpt) routine offers four alternatives for reducing delivery time: minimization of the total spot weight sum (A), maximization of the lowest spot intensity of each energy layer (B), elimination of low-weighted spots (C) or energy layers (D). The trade-off between dosimetric quality (step I) and treatment time (step II) is controlled during the optimization by option-dependent parameters. PrEfOpt was applied to a clinical patient case, and plans for different trade-offs were calculated. Delivery times were simulated for two virtual facilities with constant and variable proton current, i.e. independent and dependent on the optimized spot weight distributions. Delivery times decreased without major degradation of plan quality; absolute time reductions varied with the applied method and facility type. Minimizing the total spot weight sum (A) reduced times by 28% for a similar plan quality at a constant current (changes of minimum dose in the target  <1%). For a variable proton current, eliminating low-weighted spots (C) led to remarkably faster delivery (16%). The implementation of an efficiency-optimization step into the optimization process can yield reduced delivery times with similar plan qualities. A potential clinical application of PrEfOpt is the generation of multiple plans with different trade-offs for a multicriteria optimization setting. Then, the planner can select the preferred compromise between treatment time and quality for each individual patient.

optimization step into the optimization process can yield reduced delivery times with similar plan qualities. A potential clinical application of PrEfOpt is the generation of multiple plans with different trade-offs for a multicriteria optimization setting. Then, the planner can select the preferred compromise between treatment time and quality for each individual patient.
Keywords: prioritized optimization, IMPT treatment planning, delivery efficiency (Some figures may appear in colour only in the online journal)

Introduction
Intensity modulated radiotherapy treatment planning involves optimization of multiple criteria. Dosimetric goals defined in physicians' prescriptions, such as the coverage of the planning target volume (PTV) with a prescribed dose and dose limits to organs at risk (OARs), are often conflicting and require compromises. In addition to those dosimetric goals, short treatment times are a major concern in radiotherapy for various reasons: faster plan deliveries increase patient comfort, cost-effectiveness and the number of patients that can be handled by the clinic. Shorter irradiation times potentially reduce intrafractional dose uncertainties, due to patient movement and organ deformation, and may further have a positive impact on biological response (Paganetti 2005, Bewes et al 2008, Suzuki et al 2011, Mittauer et al 2013, van de Water et al 2015. During the last few years faster irradiation techniques have been developed, volumetric arc therapy being an example of a modality that is being increasingly used due to its short treatment times (Otto 2008).
In intensity modulated radiotherapy with photons (IMXT), as well as protons (IMPT), the degree of modulation often determines the quality of the plan and consequently the irradiation time (RT time). The aim of greater efficiency is correlated with plan quality, such that the major challenge when increasing delivery efficiency is not to compromise dosimetric quality. Treatment planning faces a trade-off problem, not only with respect to conflicting dosimetric goals but also between plan quality and delivery efficiency (Craft et al 2007, Bortfeld and Webb 2009, Mittauer et al 2013, Wilkie et al 2013. Multicriteria optimization (MCO) is a method well suited to investigate the optimization of interdependent objectives. MCO treatment planning is based on the creation of a database of optimal plans, the so-called Pareto front. It allows the user to search interactively by surface navigation for the plan with the best trade-off (Craft et al 2005). Craft et al utilized this approach to study the complexity of IMXT plans as a measure of efficiency versus quality (Craft et al 2007).
Previous studies on correlations between delivery efficiency and quality demonstrated that RT times of IMXT and IMPT plans can often be decreased without compromising plan quality (Coselmon et al 2005, Craft et al 2007, Kang et al 2008, Mittauer et al 2013, Cao et al 2014, van de Water et al 2015. Sometimes plans are optimized to an unnecessary complexity (Bortfeld and Webb 2009). The degeneracy of the solution space allows one to generate similar dose distributions using different fluence maps or spot patterns Webb 2009, Cao et al 2014).
Different methods for decreasing IMPT delivery times have already been published (Kang et al 2008, Cao et al 2014, van de Water et al 2015 and are partly realized in clinical practice. IMPT RT times can be reduced at different stages of the planning process: prior to optim ization by selecting a large spot grid, via post-processing methods, e.g. by eliminating selective spots, or as part of the optimization routine. The initial spot raster, given by the lateral spot distance ∆x and ∆y and the steps in depth ∆z, defines the number of available spots. Larger grids typically reduce application times but concurrently decrease the number of degrees of freedom for the optimization. A trade-off strategy, which keeps the number of spots at a minimum level but still high enough to achieve adequate dosimetric qualities, is required (Hillbrand and Georg 2010). Non-uniform depth scanning decreases the number of required energies, in particular for deep-seated tumors (Kang et al 2008). Recent publications have demonstrated remarkable reductions of IMPT times by eliminating energy layers, especially for facilities with large energy switch times (Cao et al 2014, van de Water et al 2015. IMPT delivery times strongly depend on the individual facility and its characteristics. Proton currents are frequently assumed to be constant or are not specifically addressed. The characteristic of variable currents (dependent on the optimized spot weight distribution) is crucial for some clinics: treatment currents determined by the lowest spot weight per energy layer often result in low currents and prolong treatment times.
With respect to delivery efficiency, certain spot distributions may be more favorable than others for different centers. In this context it should also be mentioned that proton facilities face a lowest spot limit given by the monitor chamber which can lead to undeliverable spot patterns. Post-processing interactions such as manual elimination of the 'responsible spots' can solve the issue but may degrade plan quality (Zhu et al 2010, Howard et al 2014. By incorporating the limitations of the monitor unit into the optimization, deliverable spot distributions are assured while plan quality is maintained (Cao et al 2013).
Here we present an optimization routine for IMPT planning which integrates delivery efficiency into a prioritized optimization scheme in order to reduce treatment times and meanwhile control plan quality. Prioritized treatment planning was first suggested by Wilkens et al (2007) for IMXT plan optimization as a stepwise approach to implementation of optimization goals. It translated the hierarchy of objectives into a mathematical optimization routine. The prioritized efficiency-optimization routine presented here ('PrEfOpt') optimizes the plan quality first and allows treatment time to be reduced via various alternative methods in the final step of the routine.
We demonstrate the feasibility and potential of the efficiency-optimization method by applying it to a clinical astrocytoma case. In order to compare the efficacy of the implemented method and its dependence on proton current properties two generic facility types were implemented, i.e. with constant and variable proton currents.

The prioritized efficiency-optimization routine
The prioritized efficiency-optimization algorithm (PrEfOpt) prescribes a stepwise approach to implementation of optimization objectives. Each step is realized by a corresponding objective function, while the order of steps refers to the clinical relevance of the corresponding goal. PrEfOpt consists of two steps (see figure 1; for mathematical formulations see section 2.2).
Step I comprises two consecutive sub-runs to optimize the dose distribution according to the clinical prescription: first (step Ia) the homogeneous coverage of the PTV is optimized; second (step Ib) the mean dose to selected structures is reduced by considering prior achievements, introduced as constraints. For both sub-runs maximum dose constraints are specified prior to step Ia, and these are maintained throughout the whole optimization routine including the following step II.
The treatment time is optimized in step II by one of four alternative methods. Method A reduces the sum over all spot intensities. Method B maximizes the lowest spot intensity per energy and beam. For both methods prior dosimetric achievements are turned into hard constraints which may be violated within a small specified range given by the so called slipfactor s. Method C and D reduce the spot pattern by eliminating spots and energy layers, respectively. The number of eliminated spots (C) is defined by a variable limiting factor f S (here applied to the median spot weight). Instead of comparing each spot with the median spot weight the criterion may be adapted as one wishes, e.g. to the average PTV dose contribution of each spot. The number of eliminated energy layers per beam (method D) is defined by the parameter l E : the l E lowest weighted spot layers are deleted. For both methods, C and D, a successive re-optimization run of the plan quality using the remaining spots/layers and otherwise identical parameters as in step I is optional.

Mathematical formulation of PrEfOpt
Each step of PrEfOpt comprises an objective function to optimize the vector ω of n spot weights ω j which multiplied by the influence matrix D ij yields the dose ( ) ω D i to the voxel i. OAR and PTV describe the set of voxels of the corresponding structure volumes (total number of voxels: N OAR and N PTV ). Each spot weight ω j belongs to one beam B (number of beams: N B ) and one energy layer.
Let us introduce the following mathematical definitions: The energy of the kth energy layer of the mth beam, with the number of energy slices M m in the mth beam: and the ensemble of spots j in this energy layer: Figure 1. Optimization routine of PrEfOpt: plan quality is obtained in step I according to Wilkens et al (2007); step II offers four alternative methods for reducing the treatment time. The trade-off between plan quality and treatment time, i.e. between steps I and II, is controlled by two types of 'trade-off parameters': the slip factor and the elimination limits.
belongs to th energy layer of the th beam m k , All optimization runs are subject to: and defined hard maximum dose constraints D max to OARs: Step I: Optimization of homogeneous tumor coverage with the prescribed dose D prescr by the standard quadratic objective function (Ia) and reduction of the mean dose D mean to selective OARs (Ib) while the prior optimized PTV coverage is maintained by turning the achievements of step Ia into hard constraints. The introduction of a so-called slip-factor ˜⩾ s 0 allows for small deteriorations to enlarge the solution space. For the discussed case ˜= s 0.1.
Step Ia: Step Ib: subject to: Step II: Optimization of delivery efficiency by one of four alternative methods (A-D).
Objective functions of cases A and B are optimized with respect to previous achievements. A slip factor ⩾ s 0 is introduced to control the trade-off between the plan quality (step I) and treatment time.
Case A: Minimization of the sum over all spot weights: subject to: Note that the slip factor is only applied to the PTV term, not to the achieved D mean of the OAR.
The available number of spot is reduced beforehand by excluding all spots which had zero weight after step Ib: subject to equations (9) and (10) and: Cases C and D reduce the optimization problem by a defined number of spots (limiting factor f S ) (C) or energy layers (limit l E ) (D). Case C: Elimination of spot weights which are lower than the specified threshold given by the median spot weight ω j Ib and the variable factor f S : Case D: Elimination of the l E lowest weighted energy layers of each beam: Case C and D: Optional re-optimization of the resulting optimization problem with the reduced number of spots and the optimization parameter used in step I.

Treatment time implementation
The algorithm for and the potential of each optimization method are demonstrated based on time calculations with two generic facility types. The implemented virtual facility types are based on simplified assumptions and present examples of different current specifications, which primarily serve to demonstrate potential differences in the results derived using the chosen efficiency-optimization method. Treatment times t RT were calculated as the sum over energy switch times ∆t E , the time for the adjustment of each spot position ∆ = t 0.1 ms A (adapted from van de Water et al (2015)) and the actual beam-on time ∆t S . The minimum spot weight was defined as The efficacy of each method was evaluated by relative time savings ∆t RT :

Patient case and treatment planning parameters
The algorithm was implemented in a research planning system for 3D spot scanning protons (Lomax 1999, Nill et al 2004 in the computational environment of CERR (MATLAB) (Deasy et al 2003, Schell andWilkens 2010). All runs were optimized with the commercial toolbox Mosek (www.mosek.com).
PrEfOpt was applied to the case of an astrocytoma patient (figure 2) previously treated with photons at our institute. We selected two opposed fields of 90 and 270 which were simultaneously optimized. The spot spacing within the PTV (plus a margin of 0.5 cm) was ∆ = ∆ = x y 0.4 cm and ∆ = ⋅ z w 1.2 80 (width w at the 80% intensity level of the Bragg peak). Available energies ranged from 50 MeV to 250 MeV in steps of 1 MeV. The beam featured a constant lateral width of = ⋅∆x FWHM 0.5 at the patient surface (FWHM, full width half maximum of a Gaussian profile).
The prescribed dose was 60 Gy ( × 30 2 Gy) to the PTV. Dosimetric quality was analyzed by various dose-volume histogram (DVH) indicators of the PTV and the brainstem, referring to the brainstem excluding the PTV plus 2 mm (figure 2). Minimum and maximum doses were evaluated as ( ) D 1 cm min 3 and ( ) D 1 cm max 3 referring to the dose received by 1 cm 3 of the corresponding volume. Due to the formulation of the algorithm, the maximum dose to OARs was not increased in step II. In order to compare the potential of each method we determined the plan of the shortest treatment time of each facility type which fulfilled the dosimetric criteria of the initial plan up to a maximum change of ∆ = ± D 1% and ±2% of ( . Case-dependent trade-off parameters allowed us to control and steer the compromise between plan quality (in terms of PTV coverage) and treatment time for each optimization method.  Gy. The 100% coverage of the PTV with the 95% isodose was reduced to 99.9% for s = 0.4 and 98.7% for s = 1.2. Due to the defined hard constraints on D max for all runs and on the achieved D mean of step Ib (see equations (4) and (10)), brainstem doses of both efficiency-optimized plans were equally as good or lower than in the initial plan.

Trade-off analysis for different optimization method-facility combinations.
We calculated plans of different trade-offs for each facility-method combination by varying the corresponding parameters s, f S and l E (for the utilized values see section 2.4 and for results see table 1). The correlations between treatment times and plan quality are visualized in trade-off plots (figures 4 and 5). Each plan is represented by three evaluated DVH criteria.

Facility with a constant proton current.
Minimizing the total spot weight sum (figure 4, left) resulted in remarkable time savings, indicated by a gap in treatment times at similar doses between step I and the first plan of step II (s = 0.4) (see also figure 3). The least time was saved via 'maximization of the lowest spot weight' (designed for facility type 2) and spot elimination (C) without re-optimization. The efficacy of the latter was improved by re-optimizing the problem after spot elimination (figure 4, center) which enabled the elimination of a larger number of spots. Dose distributions derived by eliminating energy layers were qualitatively not comparable to the initial plan and required subsequent re-optimization of plan quality     (4) and (10)) with shorter treatment times.
3.1.4. Facility with a variable proton current. By maximizing the lowest spot weight of each energy layer (figure 5, left) times were reduced by, for example, ( ) ∆ = t I 14.8% RT v at a constant coverage of 100% and a decrease in ( ) ∆ < D 1 cm 0.5% min 3 . Similar results were achieved by deleting low-weighted spots (C) ( figure 5, center). A larger number of spots was eliminated with consecutive re-optimization, which led to shorter RT times than without re-optimization. Greater time reductions were obtained by method D with re-optimization (figure 5, right) and by method A. The latter, and the re-optimization of method C, did not continuously reduce RT times with the number of eliminated spots/reduced spot weight sum but indicated that fewer spots or energy layers, as well as smaller sums of spot intensities, do not necessarily result in shorter treatment times.

Varying energy switch times
In order to investigate the impact of energy switch times, delivery times were additionally calculated with ∆ = t 5 E s (see table A3, appendix). Different relative time reductions were derived with a varying impact dependent on the applied optimization method. Compared with the results for ∆ = t 1 E s (see tables 1 and A3 in the appendix) the relative decrease of delivery times was larger for the 'elimination of low-weighted energy layers': delivery times were reduced by ( ) RT c for dose changes of ∆ < D 1% and ∆ < D 2%, respectively. Contrary results were found for all other presented methods (A-C).

Varying geometry: prostate case
PrEfOpt was further applied to the case of a patient with prostate cancer, with the main differences being a larger target volume and higher required energies. Beam angles were fixed at 90 and 270 . For the facility with a constant current, the methods showed similar trends, as suggested by the astrocytoma patient (see table A3  (for plan details see table A2 in the appendix). RT times achieved by minimizing the overall spot weight sum (A) did not strictly monotonically decrease with quality degradation, as observed for the astrocytoma patient (figure 4, left), caused by variations in the number of energy layers. For the facility with a variable current, method A only achieved time reductions for compromises larger than 1% of the minimum PTV dose. The shortest RT times were obtained by the elimination of low-weighted spots without re-optimization (C) with ( ) . A larger number of spots was eliminated with consecutive reoptimization but this led to smaller time savings than without re-optimization (for plan details see table A2 in the appendix).

Discussion
We presented the treatment planning algorithm PrEfOpt, a two-step routine, to optimize IMPT plans dosimetrically and to reduce treatment times while controlling plan quality. PrEfOpt offers different methods for increasing delivery efficiency which were applied to two patients, one with an astrocytoma and one with prostate cancer. For both constant and variable current, reductions in irradiation time were achievable without major compromises in plan quality. A generally valid identification of 'the most effective' optimization method is not possible. Delivery times strongly vary with the implemented currents and energy switch times (see section 2.3). As the initial RT times of facility 1 and facility 2 differed, absolute times are not comparable. The facility simulations served to analyze the potential of the methods for each facility type separately.
To compare the efficacy of the methods we selected the plans with the shortest application times which still fulfilled certain DVH criteria (see table 1). As the number of calculated plans was limited and varied between methods, the proximity of the chosen plan to the exclusion criteria, and thus corresponding qualities, differed slightly. We determined the shortest feasible RT time and the 'first plan' violating the criteria. The 'true' feasible irradiation time lies somewhere in between.
The derived RT times indicated methods A and D to be clearly superior for a constant current for both patients. While these methods also gave the largest time savings for the facility with a variable current for the astrocytoma case, contrary results were found for the prostate cancer case. The minimization of the overall spot weight sum, as well as re-optimizations of quality after spot or energy layer eliminations (methods C and D) do not necessarily lead to shorter RT times for variable currents: (re-)appearing low spot intensities which consequently reduce the current may result in even longer delivery times. Method B was designed to specifically increase the variable current, which, similar to the pure elimination of low-weighted spots, considerably decreased corresponding RT times. The latter is a simple post-processing step that cannot be considered as an actual optimization method. As it allows one to reduce RT times for I v with comparably little effort it is considered to be a valuable alternative to more complex efficiency optimizations.
With respect to methods A and B, the formulations decreased the sum of spot weights and increased the variable currents, respectively, as intended in all cases. A simultaneous influence on distinctive efficiency-determining parameters, such as the number of spots and energy layers or the overall weight sum in case B, is, however, not 'prohibited' and may work Figure 6. Trade-off between quality and treatment times for a prostate cancer case. Presented are the method-facility combinations which achieved the largest time reductions: 'energy layer elimination with re-optimization' for a constant current (left) and 'spot elimination' of low-weighted spots without re-optimization for a variable current (right).
against the intended shorter RT times. The introduction of additional constraints could solve this situation. It should be noted though that a large number of constraints could influence the feasibility of the optimization problem.
Besides the dependence on treatment currents, the selection of the most suitable time reduction method depends on the energy switch times. Obviously the larger the energy switch time, the more impact is given to the 'energy layer elimination' method. The great potential for reducing energy layers to give a decrease in treatment times was published recently (Cao et al 2014, van de Water et al 2015. van de Water et al (2015) presented an iterative optim ization routine to reduce energy layers (and spots) in robust IMPT treatment planning in a multicriterion optimization system. Based on prioritized optimization the algorithm assured the plan quality by introducing dose constraints. The published time savings for energy switches of 1 s and 5 s are consistent with our results (method D with re-optimization; see table A3 in the appendix). Both studies concluded that due to the degenaracy of IMPT plans, times can be reduced without affecting dosimetric plan quality. Our findings underline this conclusion and complement the research by discussing different facility types for which different efficiencyoptimization methods may have greater potential.
We assessed the plan quality by evaluating DVH points of the dose distributions. In IMPT range uncertainties are of great concern; these are reduced by robust planning or incorporated by robust optimization (Unkelbach et al 2007, Pflugfelder et al 2008. We did not perform any robustness analysis here. Each method may have a different impact on how spot weights are (re-)distributed. More pronounced weighting of high proton energies (e.g. by method A) could occur. Since the biological effect of proton radiation increases at the distal edge of the Bragg Peak, these changes in spot patterns may cause undesired, crucial dosimetric consequences in patient treatment, and thus require specific considerations; similarly some methods may derive more homogeneous spot patterns over all energies, which could be beneficial in terms of robustness. To see that weights were not shifted between both beams, such that the majority of the dose is only delivered by one of the beams after efficiency optimization, we evaluated the ratio between the total spot weight sum of both beams; it remained mostly constant with slight variations between cases and methods (see tables A1 and A2 in the appendix). In clinical settings the maintenance of plan robustness has to be assured similar to dosimetric quality.
Further limitation of the treatment planning procedure might be found in the objective functions of step I. The minimization of the mean dose (step Ib) may not lead to the most optimal plan quality for every case.
Step I is merely an example, and can obviously be replaced by different objective functions.
In step II the slip factor was exclusively applied to the results of step Ia, i.e. the PTV coverage. Thus doses to OARs were kept at a constant level or even improved. The introduction of a slip factor on OARs may lead to greater time savings in some cases, and could be considered as an option to modify the PrEfOpt routine.
Provided that a plan has been efficiency optimized (by any method), a method to decide on the best trade-off is required. What degree of dosimetric degradation is still acceptable, and is it acceptable to influence the plan quality at all with the intention of just saving time? These are difficult questions, partly involving ethical issues, which should be considered in this context.
The question of whether anyone should (have to) make a crucial decision on the trade-off between efficiency and quality was discussed earlier by Bortfeld and Webb (2009). We agree completely with their consideration that it would mostly be worth to wait a minute longer for a better plan (Bortfeld and Webb 2009). There may, however, be patients for whom the benefit of shorter treatment times is larger than the impact of minor compromises in the calculated dose distributions. Further, our work and prior publications indicated that large time savings do not necessarily mean a reduction in quality (Coselmon et al 2005, Craft et al 2007, Kang et al 2008, Mittauer et al 2013, Cao et al 2014, van de Water et al 2015. One method for selecting the treatment plan with the best trade-off is by execution of an automatic script which selects the plan with the shortest time that still fulfills specific DVH criteria. This method will always pick the 'worst' plan within the defined interval of acceptable DVH values, and may not always result in the same choice being made as would be selected by a person. Method A for a constant current (figure 4, left) presents a large decrease in RT time between the initial and the efficiency-optimized plan with the best quality, and comparably small time improvements were only achieved with further plan degradation. The 'first plan' would most likely be the 'trade-off choice' made by any person.
Moreover, automated plan selection exclusively takes DVH points into account rather than comparing the actual dose distributions. The best trade-off choice may vary with the derived dose distributions and the patient's background, i.e. the intention of the treatment and the physical condition of the patient.
The creation of a plan database similar to MCO treatment planning could be a useful application for realizing efficiency improving algorithms in clinical practice. PrEfOpt generates plans with different trade-offs which could serve to fill the database. By providing a slider on the planning interface which represents the treatment time along with sliders to control dosimetric objectives, the trade-off between time and quality could be steered interactively by the user. MCO would give control of the the trade-off to the planner or/and responsible physician rather than being hard-coded in a planning system.

Conclusion
Our work has demonstrated the potential of efficiency optimization in IMPT planning via different optimization methods. By prioritizing objectives, PrEfOPt achieved remarkable time reductions while maintaining the previously obtained dosimetric quality. Inevitable trade-offs between quality and delivery time were detected, indicating a limit to the achievement of reasonable time reductions. Even though the application of any type of efficiency optimization does not generally guarantee to decrease RT times, efficiency optimization is still considered to be useful in clinical practice. If a plan cannot be improved with regards to RT time with a similar plan quality the initial plan may still be selected for treatment. Generally the implementation of a routine after plan quality optimization to increase delivery efficiency could be of great potential for clinical treatment planning.
In order to decide on the best trade-off plan, PrEfOpt could serve to calculate plans with different trade-offs to fill a MCO database. The user could then navigate to the preferred tradeoff between plan quality and delivery time.