Four-dimensional deformable image registration using trajectory modeling

Abstract

A four-dimensional deformable image registration (4D DIR) algorithm, referred to as 4D local trajectory modeling (4DLTM), is presented and applied to thoracic 4D computed tomography (4DCT) image sets. The theoretical framework on which this algorithm is built exploits the incremental continuity present in 4DCT component images to calculate a dense set of parameterized voxel trajectories through space as functions of time. The spatial accuracy of the 4DLTM algorithm is compared with an alternative registration approach in which component phase to phase (CPP) DIR is utilized to determine the full displacement between maximum inhale and exhale images. A publically available DIR reference database (http://www.dir-lab.com) is utilized for the spatial accuracy assessment. The database consists of ten 4DCT image sets and corresponding manually identified landmark points between the maximum phases. A subset of points are propagated through the expiratory 4DCT component images. Cubic polynomials were found to provide sufficient flexibility and spatial accuracy for describing the point trajectories through the expiratory phases. The resulting average spatial error between the maximum phases was 1.25 mm for the 4DLTM and 1.44 mm for the CPP. The 4DLTM method captures the long-range motion between 4DCT extremes with high spatial accuracy.

Figure 1. Displacement versus Lagrangian path

The Lagrangian coordinate of the particle located at ξ for t = 0, is given by the function x(ξ, t), which represents the trajectory of the particle originally located at ξ, through Ω as a function of time. The displacement vector, d, is the vector difference x(ξ, t_final) − x(ξ, 0).

Figure 2. 4D Landmark point trajectories

4D landmark point sets were utilized to test the adequacy of the trajectory model and spatial accuracy of the DIR algorithms studied. a) The 4DCT image sets used in this study consisted of the 6 images spanning the expiratory phases from maximum inhalation (T00) to maximum exhalation (T50). 75 landmark point sets were identified on 10 cases as shown for the example point. Each 4D landmark point was identified (yellow arrow) for phases T00 through T50 as shown. b) A sample 4D trajectory of the landmark point depicted is plotted. Note that the T30 and T40 points overlay each other.

Figure 3. Trajectory modeling with polynomials

4D landmark point sets were utilized to test the adequacy of the trajectory modeling across the 6 images spanning the expiratory phases (from T00 to T50). a) A linear model of the z-displacement is shown versus a sample landmark point displacement. b) A cubic model of the z-displacement is shown versus a sample point. c) A quintic (5^th order polynomial) model of the z-displacement is shown versus a sample landmark point displacement.

Figure 4. Landmark points and DIR errors

Manually identified landmark point sets were utilized to compare the spatial accuracy of the DIR algorithms studied. a) 419 manually determined displacement vectors are shown in anterior (top row) and lateral (bottom row) projection for a sample case. The lung silhouette is in gray and the gross tumor volume is shown in red. Residual error vectors are also shown for b) 4DLTM and c) CPP DIR algorithms. Each error vector points from the manually delineated feature location in the target image to that determined from the respective DIR transformation.

Figure 5. Spatial error versus displacement magnitude

The absolute distances between the reference landmark displacement vectors and the a) 4DLTM and the b) CPP algorithms were tallied for the set of 8832 landmarks versus size of the displacement in 4 mm increments. Though the complete set of landmark measurements was used to generate the box plots shown, outlier data points have been removed from the figure for clarity.

Figure 6. Phase-step registration errors

A box plot is shown illustrating the range of 4DLTM cubic magnitude registration errors at each phase increment. The 75 sampled 4D reference trajectories for each case were combined to pool the measured errors, resulting in 750 error measurements for each phase bin.

Figure 7. Calculated temporal motion sequences

For each reference case, an example T00 sagittal view is shown with a random sample of the corresponding in-plane trajectories calculated using the cubic 4DLTM algorithm. The plotted trajectories are color-coded to indicate their temporal sequence. The initial T00 →T10 displacements are shown in blue, while each subsequent displacement gradually changes shade towards dark green. The calculated motion sequences are seen to vary widely across the 10 reference cases, both in time and space. The corresponding quantitative error assessment for each case is shown in Table 4.