Continuous real time 3D motion reproduction using dynamic MRI and precomputed 4DCT deformation fields

Abstract Radiotherapy of mobile tumors requires specific imaging tools and models to reduce the impact of motion on the treatment. Online continuous nonionizing imaging has become possible with the recent development of magnetic resonance imaging devices combined with linear accelerators. This opens the way to new guided treatment methods based on the real‐time tracking of anatomical motion. In such devices, 2D fast MR‐images are well‐suited to capture and predict the real‐time motion of the tumor. To be used effectively in an adaptive radiotherapy, these MR images have to be combined with X‐ray images such as CT, which are necessary to compute the irradiation dose deposition. We therefore developed a method combining both image modalities to track the motion on MR images and reproduce the tracked motion on a sequence of 3DCT images in real‐time. It uses manually placed navigators to track organ interfaces in the image, making it possible to select anatomical object borders that are visible on both MRI and CT modalities and giving the operator precise control of the motion tracking quality. Precomputed deformation fields extracted from the 4DCT acquired in the planning phase are then used to deform existing 3DCT images to match the tracked object position, creating a new set of 3DCT images encompassing irregularities in the breathing pattern for the complete duration of the MRI acquisition. The final continuous reconstructed 4DCT image sequence reproduces the motion captured by the MRI sequence with high precision (difference below 2 mm).

in the same position and with the same constraints as during the treatment, several studies have shown that during treatment, the motion induced by the breathing can differ significantly from the motion captured by the 4DCT in the treatment planning step, 6 resulting in potential errors or suboptimal treatment.
Different methods have been developed to address the motionrelated issues and increase confidence in tumor localization. Gating based on an external surrogate can be used but relies on the basic assumption that the correlation between internal and external motion has been the same at planning and in the treatment room. 7 Greater accuracy can be achieved by tracking fiducials implanted in the tumor, 8 but this last solution requires heavy and risky invasive preintervention. Other methods have been developed to reduce the impact of breathing-related motion during treatment such as abdominal compression, 9 audio coaching, 10 or mechanically assisted ventilation. 11 These options yield lower motion amplitude or a more regular breathing pattern that can result in a first reduction of the motion margins, or a better gating or tracking precision.
The ultimate reduction in the motion margins would entail adapting the treatment in real-time, based on precise tracking of the 3D anatomical structures' volumes. This would make it possible to transform a complex 4D treatment into a sequence of more precise 3D treatments synchronized with the anatomical motion. To achieve this, the real time positions of the target and of the surrounding organs must be known throughout treatment delivery.
With the recent developments of hybrid Linac-MRI solutions in standard photon based radiotherapy, 12 continuous real time MR imaging during treatment is now possible and more research is focusing on online adaptive treatment [13][14][15] and tumor tracking. 16 In proton therapy, no hardware solution is as yet available, but is under serious consideration. 17 MRI is the ideal imaging modality for this application. It can give good soft tissue contrast and it is not irradiant, which is crucial for an unconstrained use during treatment.
Unfortunately, MRI devices do not allow real-time full 3D acquisition for now and therefore only 2D slices can be acquired in real-time during treatment to follow the motion (this is usually called cine or dynamic MRI). Full 3D motion must then be derived from this limited amount of geometric information using motion models 8 Different methods have been developed to use 2D slices efficiently to drive 3D patient-specific motion models built on vector fields coming from nonrigid deformation algorithms. The diaphragm position can be used as a navigator, 19 a PCA-based similarity metric can be derived from the vector fields, 20,21 or more recently ROI have been used to drive a motion model. 22 All these options are used to transfer the 2D motion information from dynamic 2D MRI slices to complete 4D motion retrospectively. The treatment delivery quality can then be controlled between fractions and the accumulated dose for both the tumor and organs-at-risk can be computed, taking the patient's real breathing motion into account. However, these methods are not fast enough to be used in real-time. Tumor tracking is sometimes possible using cine MRI 2D slices, 23 but even with MRI soft tissue contrast the target is not always visible inside soft tissue such as in the liver.
While MRI is the ideal solution for continuous imaging during treatment, tissue density related images such as CT's are still necessary for dosimetric quantification. Recent results have shown that for standard radiotherapy and proton therapy alike, MRI-only workflows can be precise enough to generate virtual 4DCT images derived from 4DMRI acquisitions. This method can replace 4DCT acquisition in the planning phase and be used for dose calculation and image guidance. 24,25 Another important issue is controllability. To be accepted as a real-time treatment guiding tool, such a method has to give real-time feedback to allow a treatment to be controlled and stopped in realtime if necessary. To address these issues, we propose a method to transfer, in real time and continuously, the breathing-induced anatomic motion tracked on 2D dynamic MRI to a virtual 4DCT sequence. Like the previously mentioned motion models, the proposed workflow is based on precomputed deformation fields, but the driving mechanism is a simple multimodal MRI-CT interface tracking method. It is designed in such a way that is reliable and controllable through observation of the results in real-time, allowing the practitioner in the treatment room to take immediate actions if necessary. It is also scalable in precision and can benefit from any future improvement in fast MR imaging. The method could be used as a treatment verification or guiding tool and as part of a real-time dose accumulation observation method for photons or protons (see Discussion section).

| METHOD
The initial step of our method relies on the acquisition or creation of a 4DCT (referenced as 4DCTo for original 4DCT) before treatment delivery. Prior to the real-time application, deformation fields between the 4DCTo phases are computed using the diffeomorphic morphons algorithm 26 from the open source platform OpenReggui. 27 This algorithm computes 3D to 3D deformation fields that are consistent with the anatomy of the organs (diffeomorphism allows elimination of unrealistic artifacts in the deformation field) and can be summed and scaled easily while preserving rotations.
Using the deformation from the morphons, we also compute the midposition (MidP) image 28 offline prior to the real-time application. Two deformation fields are saved for each phase N of the 4DCTo: • The deformation field between the MidP phase and the phase N.
• The deformation field between the phase N and the phase N + 1 (the last phase is registered on the first one to create a cycle).
The continuous motion is captured through 2D dynamic MRI sequences. The acquisition frequency of the dynamic MRI is chosen high enough to be able to follow the breathing motion of the patient, creating a 2D video of the breathing pattern (See material section). The motion will be tracked on those image series using interfaces between moving anatomical structures. DASNOY-SUMELL ET AL.

| 237
While the 2D acquisition is running, the last step before launching the real-time tracking is to find the position of the dynamic MRI plane inside the 4DCTo. This way, the corresponding 2DCT slices are extracted from the 4DCTo to form a short 2DCT video of 10 frames in the same plane position as the 2D MRI Fig. 1. This is done using the first few MRI frames of the sequence, based on the matching of nonmoving bony structures such as vertebrae. This can be done using multimodal rigid registration, but for the 15 patients in this study we did it manually.
In a real treatment situation, it would be done right before the start of treatment.
The general workflow of our approach is summarized in Fig. 2.
Using the two data sets, namely the continuous 2D MRI and the 4DCTo with the associated deformation fields, our method transfers the positions of tracked anatomical structures from the 2D MRI frames to a new set of virtual 3DCT phases in real-time. A 3DCT image is generated for each MRI frame by inter-or extrapolating the deformation fields to constitute the continuous virtual 4DCT (4DCTc). Note that by using precomputed deformation fields to deform 3DCT images, our method relies on the hypothesis that breathing-unrelated anatomical changes such as stomach fullness, bowel gas, tumor growth/shrinkage and patient weight loss stay small between the 4D image set acquisition and the application of our method, at least in the treatment path. Under this condition, the entire 3D motion description based on a small set of 2D slices is reliable. If it is not the case, the quality of our results might drop but the validation step which is also part of the real-time workflow will detect the difference and let the user to take action if necessary (see discussion). The different steps in the process are detailed in the following subsections.

2.A | Interface tracking
The breathing-induced organ motion and positions are captured by tracking the interfaces between two tissues. To track an interface, the user can place a navigator, a small line in the direction of the motion and crossing the interface to track. The tracking works with a simple image processing pipeline, Fig. 3 occurs in real-time (~2 ms per navigator) and is robust to noise. The name "navigator" is borrowed from the MRI navigator echoes, a 1D acquisition MRI scanners can acquire and use to sort images in post processing or trigger acquisitions. 29,30 In this case, the 1D vector is a vector of pixels manually selected on the frame, but the outcome is similar. We chose to position them manually for reasons linked to reliability and controllability by the practitioner and are developed in the discussion section.
When an interface is chosen on one image series, a navigator pair is automatically created to track the same interface on both image series using the registration of the MRI position on the 4DCTo. The tracked position is given as the distance between the navigator end and the tracked interface in mm. The interface motion measured on the MRI frames is compared to the one measured on the 10 frames of the 2DCT video Fig. 4.
As the aim is to reproduce with fidelity the motion of the tumor, of surrounding tissues and of anatomical structures in the path of a photon or proton beam, these regions are our main concern and were chosen in priority on all the sequences to place the navigators. A set of two to five navigator pairs was used for each MRI slice position, depending on the quantity of interesting tissue interfaces that were visible on both the MRI and CT images. In most cases, tumor borders were not visible enough to be tracked on either the MRI or CT image series. We chose interfaces with which the tumor motion seemed highly correlated, but we did not compute the correlation coefficient for this work.

2.B | Phase selection method
Once the interface position is measured on the 10 phases of the 4DCTo and at the same position for the new MRI slice, both motion signals are combined to determine which phase has to be created with the following steps: 1. Check the current breathing state.

2.D | 3DCT phase creation
The final step is to generate the new 3DCT phase. The velocity field between two phases, N and N + 1 in interpolation cases or MidP and N in extrapolation cases, is multiplied by the interpolation or extrapolation ratio. Then, in order to use it to deform an image, the velocity field has to be converted into a deformation field using field exponentiation. Finally, the phase N is deformed using this interpo-

2.E | Validation
The evaluation of the quality of the results is also part of the realtime process. To see if the motion of the MRI is well-matched, pairs of navigators are applied again, this time on the current MRI frame and on the 3DCT we have created for this frame. As the MRI contains the ground truth anatomy that we aim to replicate with the 4DCTc, this difference should be as small as possible. This comparison is done for two navigator sets:    When measured on the navigators of the UN set, the absolute value of the motion signals difference between the MRI and the 4DCTc was usually under 2 mm in average during the 2 min of acquisition (Example in Fig 9). The average and maximum difference for all patients are reported in Table 1. The same results for the NN navigator set are reported in Table 2  The same results are reported in Table 1  This is caused by the fact that extrapolation also multiplies the errors by the extrapolation phase ratio. In some extreme cases, it can create a big difference between the MRI frame amplitude and the reconstructed 4DCTc phase Fig 11. Of all the 74 MRI positions used, the worst results shown in Table 1 and Table 2  As our next step, we are working on a way to apply this method even in the case of interfraction changes between the anatomy of the day and the 4DCTo. In any event, it is important to focus on motion in the beam's path to create the 4DCTc to get the best results for a treatment verification/guiding tool.

4.B | Manual tracker positioning
We chose to use a manual tracker positioning for different reasons. They have the advantage of being fast to acquire and process.
For a real-time use, several applications have to be considered: • In order to be used as an offline tool to recompute the fraction dose afterwards, only the MRI acquisition needs to be fast to capture any irregularities in the breathing pattern. The method presented here can be used offline.
• Using this tool as a real-time anatomy-based guiding tool with preplanned treatment means that tissue interface tracking and breathing state evaluation have to be done in real time, but these steps take only a few milliseconds. The main issue is the MRI acquisition and reconstruction time. method with such speed is difficult to achieve, but fast methods already exists 36,37 and may become even faster in the future.

| CONCLUSION
The proposed method can generate in real time a set of 3DCT images that imitates well the motion observed on dynamic MRI. Our method relies on a simple real-time interface tracking method and on precomputed 3D deformation fields. It could be used to recompute the delivered dose between treatment fractions more precisely, and eventually trigger a plan adaptation. It could also be used to compute the delivered dose in real time, either as a real-time verification tool or in the future as a real-time guiding tool for an on-thefly adapted treatment method. Positioning the trackers manually requires a human intervention, but it also makes the method more robust to imaging artifacts or noise and gives the physician final control. The resulting 4DCTc quality is verified in real-time using the same tracking tool. Our next step will be to use it to compare different treatment strategies to reduce the motion margins. The different treatment strategies outcomes can be compared by recomputing the delivered dose on the 4DCTc with consideration of the real motion tracked on the MRI sequences.

CONF LICT OF I NTEREST
Kevin Souris is now involved in a research collaboration with IBA s.a. and is funded by the Walloon Region (Belgium).