Cardiac motion estimation from medical images: a regularisation framework applied on pairwise image registration displacement fields

Accurate cardiac motion estimation from medical images such as ultrasound is important for clinical evaluation. We present a novel regularisation layer for cardiac motion estimation that will be applied after image registration and demonstrate its effectiveness. The regularisation utilises a spatio-temporal model of motion, b-splines of Fourier, to fit to displacement fields from pairwise image registration. In the process, it enforces spatial and temporal smoothness and consistency, cyclic nature of cardiac motion, and better adherence to the stroke volume of the heart. Flexibility is further given for inclusion of any set of registration displacement fields. The approach gave high accuracy. When applied to human adult Ultrasound data from a Cardiac Motion Analysis Challenge (CMAC), the proposed method is found to have 10% lower tracking error over CMAC participants. Satisfactory cardiac motion estimation is also demonstrated on other data sets, including human fetal echocardiography, chick embryonic heart ultrasound images, and zebrafish embryonic microscope images, with the average Dice coefficient between estimation motion and manual segmentation at 0.82–0.87. The approach of performing regularisation as an add-on layer after the completion of image registration is thus a viable option for cardiac motion estimation that can still have good accuracy. Since motion estimation algorithms are complex, dividing up regularisation and registration can simplify the process and provide flexibility. Further, owing to a large variety of existing registration algorithms, such an approach that is usable on any algorithm may be useful.

Exclusion table for cardiac motion challenge data is demonstrated in table 2. Since Volunteer 7,8 and 11 (V7,8,11) have more than 50% of the landmarks excluded, therefore the entire volunteers' data would be excluded. Observer 1 data of V13 would be excluded for the same reason, while observer 2 of V13 data is still acceptable.

Consistency Correction's Gradient Descent
The initial value of F * is obtained using C init values, and updated based on the following formulation: Where Where F n includes F at all frequency modes, rearranged as a vector. d n is the descent direction at iteration n, H is the pseudo-Hessian matrix, λ is the damping parameter, initialised with λ = 0.001 and updated according to 1 with an update factor of 5. (1-β ) is the weight of the momentum term that further dampens the change in d based on d n−1 .
Convergence is achieved when the change in F n over an iteration is small, where ε = 10 −5 : (2)

Lagrangian Tracking Gradient
Descent X re f will descent from its known coordinate at time t. This value is updated based on formulation in equation (1), but with F replaced by X re f . λ is kept to 10 −5 < λ < 1 with an initial value of 0.01 and an update factor of 10, β set to be 0.8, and d n−1 is set to zeros when n = 1.
During the descent, a measure of relative descent speed is used to detect local minima, plateau and oscillations, which we found to be quite common. The relative descent speed (Cost n rel ) is defined as the absolute change in exponential moving average cost function (Cost n−1 ave ) over the iteration, normalised against the last iteration's averaged cost: Cost n ave = Cost n−1 ave +Cost n 2 At the start of the iteration, Cost 0 ave can be set to any large number, such as 100. When a point is stuck at a local minima or plateau Cost n rel will be low due to diminishing gradient. Cost n rel is also designed to detect oscillations since Cost n ave remained relatively constant under oscillation. We detect these situations with the criteria that Cost n rel < 0.4, and if found, a new starting point is randomised within ±2 b-spline grid spacing from the old starting point and the descent is repeated. This imparts a slight stochastic property into the algorithm. Convergence is considered achieved when Cost( X re f ) < 10 −6 , which took less than 50 iterations in most cases. -3.6 × 10 −2 2.9 × 10 −3 2.1 × 10 −4 * values as mean ± standard deviation * * p-values are obtained from one-sided paired t-test of various methods against CC

1.6Êu of Symmetric Log Domain Demons (SLDD) and Free Form Deformation (FFD) Registration
We calculated average Euclidean distance errors of the simple Eulerian and Lagrangian marching without our regularization, and errors after our regularization (the BSF Consistency Correction). Table 1 below shows the results, demonstrating that our proposed regularization framework could reduce errors (Êu) from both registration methods, regardless of their initial accuracies.
Further, optimisation of the demon's result were performed by varying Gaussian standard deviation were tested, with (SLDD) or without (LDD) symmetric forces and cost function. Whichever method, Consistency Correction showed reduction in Euclidean error.

Jacobian of the Transformation
From the CMAC data, the segmentation of myocardium were already given. We used the segmentation as a mask on the image. The Jacobian were obtained at points that were sampled at regular interval of 1/4 the width of the voxel in all 3 dimensions within the myocardium. This were repeated for 40 regularly spaced time points per cardiac cycle.