Recalculation of dose for each fraction of treatment on TomoTherapy

Objective: The VoxTox study, linking delivered dose to toxicity requires recalculation of typically 20–37 fractions per patient, for nearly 2000 patients. This requires a non-interactive interface permitting batch calculation with multiple computers. Methods: Data are extracted from the TomoTherapy® archive and processed using the computational task-management system GANGA. Doses are calculated for each fraction of radiotherapy using the daily megavoltage (MV) CT images. The calculated dose cube is saved as a digital imaging and communications in medicine RTDOSE object, which can then be read by utilities that calculate dose–volume histograms or dose surface maps. The rectum is delineated on daily MV images using an implementation of the Chan–Vese algorithm. Results: On a cluster of up to 117 central processing units, dose cubes for all fractions of 151 patients took 12 days to calculate. Outlining the rectum on all slices and fractions on 151 patients took 7 h. We also present results of the Hounsfield unit (HU) calibration of TomoTherapy MV images, measured over an 8-year period, showing that the HU calibration has become less variable over time, with no large changes observed after 2011. Conclusion: We have developed a system for automatic dose recalculation of TomoTherapy dose distributions. This does not tie up the clinically needed planning system but can be run on a cluster of independent machines, enabling recalculation of delivered dose without user intervention. Advances in knowledge: The use of a task management system for automation of dose calculation and outlining enables work to be scaled up to the level required for large studies.


INTRODUCTION
Patients having curative treatment with radiotherapy receive radiation dose over a number of separate attendances known as "fractions" of radiotherapy. The radiotherapy treatment is planned using a kilovoltage (kV) CT image set, onto which a clinician outlines target volumes and organs at risk. The internal anatomy of the patient varies from day to day, therefore the dose delivered each day will vary from the dose determined at planning. Where a patient is being treated with daily image-guided radiotherapy on a TomoTherapy® system (Accuray Inc., Madison, WI), a megavoltage (MV) image is taken before treatment to determine the position of the patient; this is matched against the planning kV CT images, and the patient's position is adjusted by means of couch movements and gantry rotations (shifts in x/y/z and roll) to minimize geometric uncertainties.
For prostate radiotherapy, a number of authors have found variation in daily rectal position and differences between planned and delivered doses to the rectum. [1][2][3][4][5] These studies support the idea that accumulated dose over the course of treatment may differ from the planned dose. 6 When determining the relationship between effects such as normal tissue toxicity and dose, it is essential to have a means of estimating the dose actually delivered.
This work forms part of the VoxTox Programme, 7 which is funded by Cancer Research UK and is part of the National Institute for Health Research Clinical Research Network portfolio. The ultimate aim of VoxTox is to establish the accumulated dose and its relationship with toxicity, in approximately 2000 participants treated for prostate cancer, head and neck cancer or a central nervous system tumour. All patients have been treated using helical tomotherapy using daily MV image guidance. 8 The study was granted ethics approval by the Essex Research Ethics Committee in February 2013. All participants had given informed consent.
The TomoTherapy Planned Adaptive® module 9 allows dose to be recalculated based on the MV CT. However, this application has a graphical user interface and ties up the use of a planning terminal that is in daily use for routine clinical work. For the purpose of this study, we require recalculation of dose for typically 20-37 fractions per patient, for approximately 2000 patients, which in total therefore requires in excess of 60,000 fractions to be recalculated. For our preliminary study on ten patients, 5 it was possible to perform this by manually reoutlining and using an interactive dose calculation system. Scaling up to hundreds of patients requires a non-interactive interface that permits batch calculation with multiple computers. We have developed a system for automatic recalculation of daily delivered dose, which can be integrated within a batch calculation framework derived from the one used in particle physics. 10 We describe results of calculation times applying this system to the first 151 patients with prostate cancer in this study.

METHODS AND MATERIALS
Patients were treated on two TomoTherapy HiArt® machines. The dose calculation algorithm used is based on that described by Thomas et al. 11 Data are extracted from the TomoTherapy archive using methods developed in-house. The files extracted include the digital imaging and communications in medicine (DICOM) RTPLANs, the kV CT images, the MV images acquired before treatment on each fraction and the shifts (x, y, z and roll) applied for each fraction. All data are converted into the DICOM RT format. The shifts that are applied for each fraction of radiotherapy are determined on the basis of the matching, performed by the radiographers, between the MV CT images and the planning kV CT images to ensure that the correct target is covered by the radiation distribution. The system produces an initial automatic registration, which is manually adjusted by the radiographers to ensure soft-tissue matching. For recalculation of delivered dose, it is necessary to use the final moves applied, by applying the shift values from the coordinates of the in-plane isocentre values and applying the superior-inferior (SI) shift value to the couch position of the MV images. The roll is corrected for in-treatment by offsetting the gantry angle of each projection by the roll. This has the effect of rotating the dose distribution relative to the MV image.
The length of the MV scan was chosen to cover the length of the prostate, to enable daily image-guidance whilst minimizing imaging dose and time. As a result, the MV image set was shorter than the planning kV set. Imaging was usually carried out using the "coarse" settings which gives slices at 6-mm spacing.
The scan circle of the MV CT is 38.6 cm in diameter, which is frequently too small to cover the edge of the patient, as shown in Figure 1. It is therefore necessary to utilize kV CT in areas outside the MV CT scan circle, masked and transformed in a manner that allows for the shifts and roll that are applied.
All Hounsfield unit (HU) values within this circle are set to 21000 (corresponding to zero density) in the kV image. Figure 1 shows the process required to combine the two images. Where image matching required a roll to be applied, we offset the gantry angle of all control points in the dose calculation by this roll, exactly mimicking the situation in treatment delivery. Since the purpose of image matching is to ensure that dose is not rotated relative to the planning kV CT image, it is necessary to rotate the kV CT to compensate for this roll.
The ray-tracing algorithm described by Thomas et al 11 is carried out twice, once for each image. The radiological path length to each point is the sum of the radiological path lengths from the MV CT and from the masked kV CT images. The two images have very different pixel sizes, and require different HU to electron density conversions. Producing a composite image, by interpolation and replacement of HU values, would require more computational resources (calculation time and storage space) than is required to ray trace twice through separate images. For MV images, the HU to electron density conversion varies over time. The calibration of MV CT images on Tomo-Therapy has not been stable over time. [12][13][14] Prior to 2011, the TomoTherapy HiArt used a rotating target. The target required frequent replacement, which caused step changes in HU calibration with each target change, in addition to random fluctuations. TG-148 15 recommends that monthly HU calibration tests should verify that water-equivalent materials vary by ,30 HU from the calibrated image value to density table to maintain a dosimetric uncertainty of 2% or less. Since 2011, two changes introduced by Accuray Inc. have improved this: the use of a nonrotating target with a much longer lifetime and the introduction of a weekly user-run MV CT recalibration procedure. Measurements were performed using the density inserts (Gammex RMI Middleton, WI) inside the cylindrical phantom from TomoTherapy Inc. commonly referred to as the "cheese phantom". 15 We have repeatedly measured the phantom with density inserts over an 8-year period, at intervals varying between monthly and quarterly, and recorded the mean HU within each of a set of regions corresponding to the density inserts. The inserts had electron density relative to water of 0.30, 0.49, 0.942, 0.979, 1.018, 1.053, 1.080, 1.104, 1.112, 1.274, 1.472 and 1.694. Two additional inserts contained air and distilled water.

Implementation
The files extracted from the archive were tokenized to remove patient identifiers before exporting from the hospital to the university department of physics. The files were stored in a hierarchical structure to facilitate their processing using the computational task-management system GANGA, 10 which was originally developed for use on the ATLAS and LHCb experiments at CERN. 16 GANGA is an extensible system for defining and managing arbitrary computing tasks (called "jobs"). Applications consist of a set of algorithms run in sequence; an application may run the same algorithm multiple times, with different configurations. A job definition consists of the algorithm to be run, the data set specifying patient data for processing, the backend specifying where to perform processing, plus handling options for splitting jobs into subjobs and for merging outputs of subjobs.
The files extracted include the treatment plan data, the kV CT scan (downsampled by TomoTherapy from 512 3 512 to 256 3 256, and with the CT couch replaced by the TomoTherapy couch), the MV images acquired before treatment on each fraction, the positional corrections (x, y and z shifts and roll) applied for each fraction and the planned dose data. All data were saved in the DICOM RT format following the Tomo-Therapy DICOM conformance statement. 17 The shift and roll information were added to the MV image headers using private DICOM tags. Where an MV image does not immediately precede a treatment (e.g. in the event of a machine breakdown or where excessive rectal filling necessitated the patient going to the lavatory and being rescanned), a tag is added to the header indicating that this image should not be used for recalculation. Also stored are the RTSTRUCT objects (DICOM RT structure sets) of outlines drawn in the planning process, plus any RTSTRUCT objects associated with each MV image set, produced either by manual outlining or automatic segmentation. 5 The automatic segmentation tool was written in MATLAB® (MathWorks®, Natick, MA), using MATLAB'S two-dimensional implementation of the Chan-Vese algorithm. 18 The HU values were rescaled to highlight the critical regions, with regions of air identified separately and rescaled to match the HU of rectal material. The starting contour for each slice was based on the equivalent contour in the planning scan. At the inferior end of the rectum, the contrast between the rectum and muscles is not sufficient to distinguish these consistently, even if viewed by an experienced oncologist. This lack of contrast leads the algorithm to overcontour on these slices, by including muscles near the rectum. This is detected by identifying the most superior slice with an unexpectedly large contour and replacing contours from the inferior end to this slice with the planning scan contours. In a small proportion of cases, incorrect outlining occurs on single slices; to reduce the incidence of this, contours were automatically compared with a smoothed three-dimensional surface formed from all contours in a scan; where a large difference was found, the contours were replaced by the corresponding smoothed contour. The use of the smoothing would typically improve conformity indices (measured against a test set of manually drawn contours) from 0.7 to 0.8. Key parameters in the autocontouring algorithm were tuned using a training set of data, using manual scanning as the gold standard for comparison and a testing set of data used to verify the accuracy of the autocontouring. Further details of the algorithm, including results of agreement with manually outlined contours, are given by Sutcliffe et al. 19 Figure 2 gives an overview of the data flows in VoxTox. The dose calculation program was written as a command-line application using MATLAB. Interactive user input was replaced by automatically produced text files to deploy the program using GANGA. The calculation is performed over a set of points

RESULTS AND DISCUSSION
Megavoltage image calibration It can be seen that up to 2011, there was considerable variation between measurements, although still within the 630 HU recommended by TG148. 15 Following the introduction of the fixed In each case, results are from processing data for a single fraction for each of the 33 patients. Processing was performed on machines that have 24-core Intel® Xeon X5650 CPUs (Intel, Santa Clara, WI), with a clock speed of 2.67 GHz. A downsampling factor of 2 gives an in-plane spatial resolution of approximately 1.5 mm. target in 2011 and the weekly HU calibration procedure in 2012, this variation considerably decreased, such that it is now ,610 HU. Hence for results from 2008 to 2011, variations in HU calibration could be introducing uncertainties in dose of up to 2% whilst since 2011 this has reduced to ,1%. If recent calibration data were to be used to recalculate doses from early in the life of the machine, they would erroneously show a systematic underdosing, since they would be systematically overestimating densities.
The change in weekly HU calibration procedure introduced in July 2012, as well as making the results more stable, caused a systematic reduction in the slope of the line of electron density vs HU as shown in Table 1. It should also be noted that there is considerable variation in slope between two machines. This is in contrast to the results observed for diagnostic kV CT scanners, 25 where for low-Z materials, an offset of zero and a slope of 1000 appeared to fit a wide range of CT scanners. Table 2 shows the central processing unit time to calculate a dose cube for a single fraction of radiotherapy, memory requirements and size of output file for downsampling factors from 4 to 1. A downsampling factor of 2 has been chosen for subsequent calculations. Figure 3 shows close agreement between the doses calculated using this software and those calculated using the TomoTherapy planned adaptive module. For a typical patient with prostate cancer, the mean dose difference was 0.36%, and .98% of points agreed between the two calculations with a 3%/3 mm box index, 26 excluding points ,5 mm deep in the patient (Figure 4).

Dose calculation
To test the ability of the software to cope with large changes in patient density, a prostate plan was delivered to the "cheese" phantom with and without the insert plugs. Ionization chamber measurements were made at the two points shown in Figure 5.
The plan was delivered three times, once with unit density plugs in all the holes, once with the inner ring of plugs removed and once with all the plugs removed. Table 3 shows the agreement between calculation and measurement. Even for such an extreme case, agreement was within 3.1%.

Time saving
Daily dose calculation For the ten patients in the study by Scaife et al, 21  Equivalent uniform dose to rectum Figure 6 shows the results for the calculations of the EUD to the surface of the rectum. The delivered EUD is shown relative to the planned EUD for the 151 patients calculated. Although dose tracking at the voxel level remains an important goal, in practice there are many barriers in identifying structures with sufficient accuracy to achieve this. Summation of dose surface maps achieves many of the objectives of dose tracking, with fewer uncertainties. By splitting the surface at a point posterior of the centroid of each image, we ensure that the central column of the DSM corresponds to the middle of the anterior surface of the Figure 5. Measurement of the impact of changes in patient densities. A clinical prostate plan was recalculated on the cheese phantom, in the section of it with removable plugs. The plan was delivered three times, once with unit density plugs in all holes, once with the inner ring removed and once (as shown) with all the plugs removed. The figure shows a dose distribution superimposed on a megavoltage CT image of the phantom. rectum. This is in contrast to the splitting at the most posterior point of each outline as described by Murray et al, 22 which we found led to large variations from slice to slice. The anterior surface of the rectum is stable relative to the prostate and hence relative to the high dose volume in situations where daily imageguided radiotherapy is used to position the high dose volume to the prostate. Work is in progress identifying features in the delivered DSMs and linking them to toxicity. Figure 6 shows that delivered EUD appears to be on average 3% lower than the planned EUD, in addition to large differences between patients. Further work is needed to determine whether this systematic difference results from clinical factors such as changes in rectal shape during the course of treatment or whether it results from differences between the methods of outlining on the planning CT scans and the on-treatment MV CT scans.

CONCLUSION
We have developed a system for automatic dose recalculation of TomoTherapy dose distributions. The system does not tie up the clinically needed planning system but can be run on a cluster of independent machines. The use of the computational task-management system GANGA 10 enables automation of the process of recalculating delivered doses, dose-volume histograms and DSMs without user intervention.
The automation, using a framework adapted from particle physics, has enabled us to scale up the dose calculation and image segmentation and to process data for hundreds of patients, rather than the tens of patients that can be achieved by manual means. Use of automatically generated input files guaranteed consistency, reproducibility and freedom from user input errors for the result data set. In conjunction with the ongoing collection of patient toxicities, this method of automatic recalculation will enable the goals of the VoxTox programme 7 to be achieved, enabling us to compare planned dose and delivered dose with observed toxicity.