Nonrigid motion registration, required for 3D liver perfusion imaging, assume spatially smooth motion parameters to avoid physically implausible motion and local minima. However, this can yield inaccurate estimation of discontinuous motion fields, particularly for the sliding motion of the liver against the abdominal wall. Various methods^1-4 have been proposed to adjust the smoothness constraints, e.g. Total Variation regularization¹ and Markov Random Field (MRF) based regularization on an image-derived Minimum Spanning Tree (MST).² However, performance remains suboptimal without prior segmentations^3-4, which require manual intervention. In this work, we investigated the use of truncated L1 regularization and automated motion based clustering. At the boundaries of these clusters, no spatial smoothness of the motion parameters is assumed, allowing discontinuous motion.

The motion field was modeled using multi-level 3D B-spline meshes. At each level, the registration problem was formulated as a MRF labelling problem: $${\bf u}=\substack{{\tt\large arg\,min}\\{{\bf u}_p,\,p\in\mathcal P}}\left(\sum_{p\in\mathcal P}\|I\left({\bf x}_p\right)-J\left({\bf x}_p+{\bf u}_p\right)\|_1+\alpha\sum_{\left(p,q\right)\in\mathcal N}\min\left(\frac{\|{\bf u}_p-{\bf u}_q\|_1}{\|{\bf x}_p-{\bf x}_q\|_2},\tau\right)\right),$$ where $$$p,\,q$$$ were two neighboring control nodes in the mesh, $$${\bf x}_p,\,{\bf x}_q$$$ were their coordinates, $$${\bf u}_p,\,{\bf u}_q$$$ were their motion vectors (MVs), $$$I\left({\bf x}_p\right)$$$ was the patch partitioned to node $$$p$$$ in the target image $$$I$$$, $$$J\left({\bf x}_p+{\bf u}_p\right)$$$ was the corresponding patch in the source image $$$J$$$, $$$\tau$$$ was a threshold, and $$$\alpha$$$ was a weighting parameter. In this MRF model, each vertex represented a control node, and its label was a quantized MV. The source image was upsampled to support sub-voxel MV search. The first term in this cost function represented the matching error between the target and the transformed source, and the second pair-wise term enforced smoothness of the MVs of neighboring control nodes. To allow for motion discontinuity, we used a truncated L1 regularization. Each vertex was initially connected with its 6 immediate neighbors. The optimization was initialized using the L2-norm of differences between neighboring patches in the target image as edge costs and a single MST was constructed using Prim’s⁵ algorithm. After that, we used the L2-norm of finite MV differences between neighboring vertices as edge costs to construct a MST and used the Maximum Standard Deviation Reduction (MSDR)⁶ algorithm to partition the vertices into clusters. Next, a new MRF cost function was constructed using only edges completely within these clusters and solved by fast message passing⁷. A multi-resolution approach was used to speed up convergence. Cluster guided interpolation was used to generate the initial motion field for the next resolution level and the label space was successively reduced using incremental MV search. At each level, motion clustering and regularization were alternated for several iterations.

Data were acquired in one volunteer using a 1.5T scanner, 8-channel cardiac coil and LAVA-Flex sequence with IRB approval. Imaging parameters were: voxel size 0.7813x0.7813x2 mm and matrix size 512x512x88. Five breath-holds M0,…,M4 at different exhalation levels were obtained, approximately evenly distributed between maximum exhalation at residual volume (M0) to maximum inhalation at total lung capacity (M4). Registration of M1,…,M4 to M0 was performed using the proposed method as well as two previously published nonrigid registration algorithms Elastix⁸ and deedsMST² for comparison.

Fig. 1 shows M0 in example axial, coronal and sagittal sections, overlaid with the computed motion fields and segmentation contours generated by the proposed algorithm. The contours show that the sliding organs were correctly separated from the surrounding anatomy. The Mean Square Error (MSE) and Normalized Mutual Information (NMI) between M0 and transformed M1~M4 are compared in Fig. 2. In all cases, the proposed method and deedsMST performed better than Elastix. Compared to deedsMST, the proposed method achieved significantly better registration for larger sliding motion.The preliminary data presented here show the feasibility of achieving high quality nonrigid registration in the presence of sliding motion, such as that of the liver against the abdominal wall. The method iterates between determining a segmentation of the liver based on a prior motion estimation and determining motion based on a prior segmentation. Compared to two other registration algorithms, it performed significantly better, especially for larger respiratory motion. The segmentation obtained by the algorithm correctly identified the liver and surrounding organs as undergoing motion patterns that are different from the surrounding anatomy, including the aorta and abdominal and thoracic wall.Motion clustering for automated segmentation to improve nonrigid motion estimation in the presence of sliding motion has been shown to be feasible and superior to prior nonrigid registration algorithms. This will positively impact and ease the processing of 3D liver perfusion analysis by removing the need for manual segmentation, previously required for optimal motion registration.