Skip to main content
UKPMC Funders Author Manuscripts logoLink to UKPMC Funders Author Manuscripts
. Author manuscript; available in PMC: 2020 Sep 1.
Published in final edited form as: IEEE Trans Med Imaging. 2020 Feb 18;39(9):2750–2759. doi: 10.1109/TMI.2020.2974844

Deformable Slice-to-Volume Registration for Motion Correction of Fetal Body and Placenta MRI

Alena Uus 1, Tong Zhang 1, Laurence H Jackson 1, Thomas A Roberts 1, Mary A Rutherford 1, Joseph V Hajnal 1, Maria Deprez 1
PMCID: PMC7116020  EMSID: EMS88512  PMID: 32086200

Abstract

In in-utero MRI, motion correction for fetal body and placenta poses a particular challenge due to the presence of local non-rigid transformations of organs caused by bending and stretching. The existing slice-to-volume registration (SVR) reconstruction methods are widely employed for motion correction of fetal brain that undergoes only rigid transformation. However, for reconstruction of fetal body and placenta, rigid registration cannot resolve the issue of misregistrations due to deformable motion, resulting in degradation of features in the reconstructed volume. We propose a Deformable SVR (DSVR), a novel approach for non-rigid motion correction of fetal MRI based on a hierarchical deformable SVR scheme to allow high resolution reconstruction of the fetal body and placenta. Additionally, a robust scheme for structure-based rejection of outliers minimises the impact of registration errors. The improved performance of DSVR in comparison to SVR and patch-to-volume registration (PVR) methods is quantitatively demonstrated in simulated experiments and 20 fetal MRI datasets from 28-31 weeks gestational age (GA) range with varying degree of motion corruption. In addition, we present qualitative evaluation of 100 fetal body cases from 20-34 weeks GA range.

Index Terms: MRI, Motion correction, Fetal motion, Slice-to-volume registration, Deformable registration

I. Introduction

OVER the past two decades development of fast acquisition sequences along with advanced motion compensation techniques [1] has gradually allowed incorporation of MRI into clinical practice for imaging of fetal pathologies [2].

Single shot fast spin echo (ssFSE) sequences allow acquisition of each slice in less than one second, which minimises the impact of fetal motion artefacts on image quality. However, in 3D stacks, inter-slice motion still exists leading to either minor misalignments (Fig. 1.a) or complete loss of volumetric information (Fig. 1.b).

Fig. 1.

Fig. 1

Fetal MRI: examples of stacks with minor (a) and severe (b) motion corruption acquired during the same fetal exam and under the same orientation at different time points.

Slice-to-volume registration in combination with super-resolution (SR) reconstruction is considered to be an efficient motion correction approach since it resolves out-of-plane motion [3]–[5]. The fact that the volumetric region of interest (ROI) is oversampled at different stack orientations ensures consistency of reconstructed volumes. Recent validation of SVR for fetal brain reconstruction showed strong correlation between 2D and 3D biometry values [6].

Since the SVR is based on rigid registration, its application has primarily focused on the brain as it undergoes only rigid motion. In general, fetal body and placenta are affected by local non-rigid deformations [7], for example due to bending of fetal abdomen (Fig. 2.a) or maternal breathing, and the use of rigid registration leads to deteriorated reconstruction quality. Slices that are inconsistent with the most prevalent body position are either rejected as outliers [5] or contribute as an error to the reconstructed volume, resulting in blurring of local features (e.g., spine) and loss of texture information (Fig. 2.b).

Fig. 2.

Fig. 2

Example of a sequence of slices (a) acquired during the change of fetal trunk position and the corresponding SVR reconstructed volume (b). Note that bending deforms the shape of the spine and internal organs in both in-plane and through-plane directions.

A. Related work

The original concept of SVR for reconstruction of fetal brain from motion-corrupted MRI stacks was to interleave slice-to-volume registration with scattered data interpolation based on weighed sum of Gaussian kernels [8] or multilevel B-splines [9]. The SVR reconstruction framework was gradually extended with SR reconstruction [3], [4], edge-preserving regularisation [4], outlier rejection [3], [5], intensity matching [5], total variation regularisation [10], sinc PSF model and GPU-parallelisation [11]. More recent works proposed to employ convolutional neural network (CNN) for fetal brain segmentation from MRI slices [12], [13] and prediction of transformations of slices to canonical atlas space [14], [15].

Though SVR has been developed primarily to reconstruct the fetal brain MRI, it has also been applied to the reconstruction of fetal thorax [16], [17] under the assumption of approximately rigid motion, due to the embedding in the rib cage. Given the known cardiac phases of each of the slices, SVR can also be employed for 4D fetal cardiac reconstruction from dynamic MRI [18]. Rigid SVR has also been applied to reconstruct MRI of placenta [19], which undergoes deformation due to the maternal breathing. Here the quality of the reconstructed image relies on the robust statistics to exclude the slices which have been deformed compared to the most common shape of the placenta present in the data. In [20], rigid SVR was further extended to patch-to-volume registration (PVR) which addresses the reconstruction of deformable organs by performing motion compensation in piecewise rigid fashion, thus allowing the reconstruction of the entire uterus.

With respect to application of deformable SVR for motion correction, the existing solutions primarily focus on registration of intra-operative slices with a pre-operative planning volume [21], [22] or multimodal registration (e.g., histology to MRI) [23], [24]. The majority of monomodal methods are based on rigid SVR for global alignment followed by Free Form Deformation (FFD) registration [25] for correction of non-rigid shape changes. Recently, [26] formalised deformable graph-based SVR approach validated on a 3D heart MRI dataset. However, the existing implementation is limited to in-plane deformations only. Model-based SVR methods integrating biomechanical models for physics-based regularisation were proposed in works of [27], [28].

B. Contributions

In this paper we present a novel approach for non-rigid motion correction in 3D volumes based on an extension of the rigid SVR reconstruction method [5] with hierarchical deformable registration scheme and structure-based outlier rejection. The proposed DSVR method allows for both in-plane and out-of-plane correction of local non-rigid deformations of fetal body and placenta resulting in high quality volumetric reconstruction. The structure-based outlier rejection ensures that the residual registration errors are not propagated to the reconstructed volume.

The performance of DSVR is evaluated in simulated experiments, and using 20 fetal MRI datasets from 28-31 weeks GA range with varying degrees of motion corruption, by comparison to the rigid SVR and PVR. Furthermore, we present qualitative evaluation of DSVR reconstruction quality on 100 iFIND1 cases from 20-34 weeks GA range.

II. Background

A typical fetal MRI dataset consists of multiple motion corrupted stacks S = {Sl}l=1,…,L covering a specific ROI and acquired under different orientations. In the context of motion correction, SVR is used for iterative recovery of a high resolution ’motion-free’ isotropic volume X from an array of K low resolution motion-corrupted slices Y = {Yk}k=1,…,K. It is formalised as follows [5]:

Yk*=MkX,yjk*=skebjkyjk, (1)

where k is a slice index and j is a pixel index within each slice. The matrix Mk describes spatial relationship between the acquired slices Yk = {yjk}j=1,…,Nk and the reconstructed volume X = {xi}i=1,…,N, where each row {mijk}i=1,,N represents the acquisition PSF for voxel yjk, transformed according to the estimated motion parameters and sampled on the grid of the high resolution volume X. Intensity corrected slices are denoted by Yk*={yjk*}j=1,,Nk, while Bk = {bjk}j=i,…,Nk and sk are the slice-dependent bias fields and scaling factors, correspondingly.

The algorithm proceeds by interleaving slice-to-volume registration and super-resolution reconstructions for a fixed number of iterations (q = 1,…,Q). At each SVR iteration q, the current estimation of X is registered to slices {Yk} and the transformed PSFs {mijk(q)}q=1,,Q are updated according to the resulting spatial transformations. Next, as initialisation of the SR reconstruction loop, the weighted Gaussian interpolation [8] is performed for estimation of X (0,q). Then, we perform super-resolution reconstruction [3] of the volume X by iterative gradient descent optimisation based on minimisation of the sum of squared errors jkejk2+λR(X), where λR(X) is an edge preserving regularisation term [5] and ejk are the errors between the original {yjk*} and simulated Ȳk = {ȳjk}j=1,…,Nk slices

ejk=yjk*y¯jk,y¯jk=imijk(q)xi (2)

Each SR reconstruction iteration (n=1,,NSR(q)) includes expectation-maximization (EM) robust statistics scheme for rejection of outliers and estimation of {bjk} and sk

xi(n+1,q)=xi(n,q)+αkjpjkpkslicemijk(q)ejk+αλxiR(X), (3)

where pjk and pkslice are posterior probabilities of a voxel or a slice not being an outlier.

III. Method

A. Overview of the algorithm

The proposed DSVR method is an extension of the rigid SVR SR reconstruction framework [5] described in Sec. II and it is summarised in Fig. 3 with the novel modules highlighted by red outline.

Fig. 3.

Fig. 3

Proposed DSVR reconstruction algorithm. The novel elements are highlighted by red outline.

Following global rigid and FFD registration of all stacks, the template stack is used for initialisation of the registration target for the DSVR loop. At each DSVR iteration (Sec. III-C), the current estimation of X (n,q) is registered to each of the slices {Yk}. This is followed by the Weighted Gaussian reconstruction and iterative SR reconstruction (Sec. III-B) with additional structure-based outlier rejection (Sec. III-D) step. The pipeline is performed for a predefined number of DSVR (Q) and SR (NSR) iterations. The rest of this section describes the novel DSVR modules in detail.

B. Incorporating deformations into SR reconstruction

In super-resolution SVR (Sec. II), the forward problem is modelled by applying the transformed PSFs {mijk}i=1,,N to high-resolution volume X to simulate (intensity-corrected) voxel ȳjk of the acquired slice k (2). The underlying PSF in the space of acquired data is a continuous function fjk (u) = f(uujk) where u is a location in the space of the acquired stack, ujk is the position of the voxel yjk and f is the shape of the acquisition PSF, which we model by a 3D Gaussian with zero mean and principal axis aligned with the axis of coordinate system of the acquired image. The transformation Tk between locations u in space of acquired stack and the anatomical locations v is estimated by registration of the acquired slice Yk and the volume X, defining the transformed PSFs by mijk=fjk(Tk1(vi)), where vi is the location of the voxel xi in the anatomical space. In case of rigid SVR, the transformed PSFs are re-oriented Gaussians. To correct for deformation of the fetal body and placenta the transformations Tk(u) need to be deformable. Therefore in DSVR the PSFs are deformed using non-rigid transformations Tk(u) and their shapes become non-Gaussian. The high resolution volume X is estimated iteratively using (3) for both rigid and deformable cases.

C. Hierarchical motion correction

SVR of the fetal brain can be well constrained by acquiring several stacks {Sl}l=1,…,L in different orientations and assuming rigid motion of the region of interest, to recover the ’true’ shape of the fetal brain. The fetal body and placenta, on the other hand, undergoes continuous deformation in time and DSVR is under-constrained in comparison to SVR. In order to overcome this limitation, a hierarchical scheme for gradual refinement of transformation during slice-to-volume registration is proposed.

1. Volumetric registration

At first, one of the stacks is selected as a template (Stemplate) for initialisation of the global registration target. In order to eliminate the impact of global rotations and translations between stacks, the initial step includes 3D-to-3D rigid registration of the input stacks {Sl} to a masked ROI in the template stack. The resulting transformations GlR are then used for initialisation of 3D-to-3D global deformable registration of stacks to Stemplate. The template stack cropped with a bounding box of the mask acts as the preliminary initialisation of the reconstructed volume Xinit.

For the purpose of formalisation of DSVR registration steps, we define a deformable registration operator as 𝔇(Itarget, Isource, Tinit, d), where Itarget and Isource are the source and target images, d represents the resolution of the deformable transformation and Tinit is the input transformation. Then the global deformable stack registration is expressed as:

GlD=𝔇(Sl,Xinit,GlR,dinit) (4)

In order to avoid over-fitting to the motion corrupted features of stacks, global transformation with low resolution dinit is chosen and all stacks are smoothed using Gaussian blurring. Therefore, the output transformations GlD provide estimation of only large range (global) deformations between the trunk positions in the stacks and the template. The trunk mask is then transformed to all stacks and they are cropped to large bounding box ROIs.

2. Slice-to-volume registration

We refine the motion parameters by iteratively aligning the reconstructed volume to each individual slice. The volume is registered to the slices (i.e., 3D-to-2D) rather than the slices to the volume to ensure that both in-plane and out-of-plane deformable motion is resolved. The first iteration of DSVR is performed by deformable registration of the smoothed template stack Xinit to each of the slices with low resolution transformations (5).

Tk(0)=𝔇(Yk,Xinit,Gk,lD,d(0)) (5)

The B-spline control point spacing d (q) is decreased at every DSVR iteration q, thus refining the resolution of the transformation Tk(q).

Tk(q+1)=𝔇(Yk,X(n,q),Tk(q),d(q+1)) (6)

This is coupled with decreasing SR regularisation parameter λ (3) to prevent overfitting to the residual motion in the early stages and progressively allowing more localised deformations as the features in X (n,q) become better defined. Fig. 4 illustrates an example of the refinement of 3D transformation (in our case implemented by a B-spline control point grid) of Tk(q) with respect to DSVR iteration q and reconstructed volume X (n,q) used as a template. It was identified experimentally that rigid registration should not be used during slice-to-volume registration steps since it disrupts the local deformation fields obtained during hierarchical refinement. The spatial relationship coefficients Mk(q) are computed during each iteration q after slice-to-volume registrations were completed, by transforming 3D Gaussian PSFs fjk using Tk(q).

Fig. 4.

Fig. 4

An example of refinement of B-spline control point grid (d (q)) at each DSVR iteration (sagittal plane with respect to uterus).

D. Structure-based outlier rejection

1. Global structure-based outlier rejection

During each iteration q, after the registration and prior to the Gaussian reconstruction step (see Fig. 3), misregistered slices are removed to ensure that global registration errors are not propagated into the initial estimation and further SR reconstruction loop. The quality of registration is assessed as global normalised cross correlation (NCC) between the original slice Yk* and the current estimation of the output volume X (n,q) transformed with Tk(q) within the masked slice ROI. The slices with transformations resulting in correlation lower than a threshold value TNCC are excluded. The corresponding slice outlier criteria are computed as:

wkG={1,if1Nkj=1Nk(yjk*μy*)(xjkTμxn)σy*σx>TNCC0,otherwise, (7)

where {xjkT=X(n,q)(Tk(q)(ujk))}j=1,,Nk is the transformed high resolution volume resampled on the grid of the slice Yk. The values μ y*, μx, σ y* and σx are the corresponding intensity means and standard deviation of {xjkT} and Yk*.

2. Local structure-based outlier rejection

DSVR of deformable objects can be prone to regional misregistrations and overfitting. Therefore, an additional step for regional outlier rejection is introduced. At each SR iteration n (see Fig. 3), the regions of simulated slices with low structural similarity are excluded from contribution to the reconstructed volume. It is based on local structural similarity (SSIM) maps {smjk} between the simulated and original slices:

smjk=(2μr*μr¯+c1)(2σr*r¯+c2)(μr*2+μr¯2+c1)(σr*2+σr¯2+c2), (8)

where r* and r¯ are the regions in the original Yk* and simulated Ȳk (n) slices centered around voxel j with circular window and μr*,μr¯,σr*2 and σr¯2 are the corresponding average and variance intensity values and σr*r¯ is the covariance of r* and r¯. As defined in [29], values c 1 and c 2 used in order to balance the division with weak denominator are computed as (k 1 L)2 and (k 2 L)2, where L is the dynamic range of intensities in r* and k 1 and k 2 are equal to 0.001 and 0.003, correspondingly. The outlier criteria of slice voxels with similarity < TSSIM are set to zero:

wjkL={1,ifsmjk>TSSIM0,otherwise (9)

3. Incorporating structure-based outlier rejection into SR reconstruction

Outlier rejection is incorporated into SR reconstruction (3) by scaling the contribution of the reconstruction error ejk by a weight wjk:

xi(n+1,q)=xi(n,q)+αkjwjkmijk(q)ejk+αλxiR(X) (10)

In our previous work [5] this weight was set based on EM robust statistics on voxel and slice intensities to pjkpkslice. The proposed total structural outlier criteria for a voxel defining its contribution to the reconstructed volume is computed as:

wjkS=wkGwjkL (11)

In Section V-D we will show that the combination of the two schemes by setting wjk=pjkpkslicewjkS outperforms the individual outlier rejection schemes.

IV. Implementation

A. Input data requirements

The input dataset includes stacks of different orientations, an approximate mask covering the ROI (i.e., fetal body or placenta) and a selected template, which can be either one of the stacks or a scout scan. The global structure of the target object needs to be preserved in the template stacks, as it is used for initialisation of the registration target. Minor to average degree of motion corruption of the template is acceptable and is resolved by Gaussian blurring. The severity of motion is visually assessed with respect to the degree of loss of volumetric structural information (see Fig. 1).

Rigid SVR requires masking of the ROI in the template stack in order to eliminate the impact of the independent motion of the mother and the fetus. On the other hand, FFD registration does not require precise masking and in our experience preforms better for larger ROIs.

Variable orientations of input stacks help to prevent over-fitting to a particular motion-corrupted stack. The minimal requirement for SR reconstruction of an isotropic volume from multiple stacks is two sufficiently different orientations. It was identified experimentally that using 5 to 8 stacks is sufficient for good quality reconstruction depending on the amount of motion corruption, resolution and SNR level of the original volumes. Similarly to the capture range limitation of SVR [11], gradient-descent based deformable registration methods are not capable of resolving motion involving large degree rotations or excessive bending, which should be taken into account with respect to selection of input stacks.

B. Deformable registration

The B-spline FFD registration [25] with NMI similarity measure was chosen for both deformable SVR and global registration steps due to the lower computational requirements compared to the diffeomorphic registration such as FFD parameterised by stationary velocity (SV) fields [30]. Although SV FFD ensures invertibility of transformations, it did not lead to an indicative improvement of reconstruction results while significantly increasing processing time, which made it not feasible for our application.

Typically, the reconstruction pipeline requires 3 SVR iterations (Q = 3) each of which is followed by 10 to 30 SR (NSR(q)) iterations with gradually refined regularisation parameters. The resolution of B-spline FFD transformation is controlled by changing the resolution of the B-spline control point (CP) spacing d. We choose resolution scheme with B-spline control point spacings d (0) = dinit, d (1) = 2/3 · dinit, d (2) = 1/3 · dinit. As we show in Sec. V-C, it was identified experimentally that 15mm → 10mm → 5mm CP refinement (dinit = 15mm) produces the optimal reconstruction quality for our fetal test cases and 0.85mm output resolution. The corresponding optimal regularisation λ (q) values were chosen as 0.1 → 0.05 → 0.02 similarly to the SVR settings in [5].

C. Structure based outlier rejection parameters

Analysis of the choice of the structural similarity thresholds showed that the optimal values corresponding to adequate registration quality are TNCC = 0.75 for global and TSSIM = 0.6 for local regions. Using lower values might lead to inclusion of regions that were erroneously overfitted.

The 20 mm diameter for SSIM kernel was experimentally identified as optimal for the feature sizes in 28-31 weeks GA range subjects, e.g., the transverse diameter of fetal kidneys for this GA range varies within 15 – 25mm [31].

D. Software packages and hardware requirements

DSVR framework was implemented based on MIRTK2 library with multi-CPU parallelisation of registration and reconstruction steps. The structure and functionality of the core reconstruction steps follow the IRTK-based3 implementation of the original SVR reconstruction method [5]. The code is available online as a part of SVRTK4 package.

The major advantage of MIRTK registration library is the use of conjugate gradient descent optimisation [32] that significantly increases computational efficiency of FFD registration that constitutes the most time-consuming part of DSVR pipeline. Depending on the ROI size (related to GA of the subjects), number of stacks, output resolution and the system configuration, the reconstruction time can typically vary between 15 to 60 minutes.

V. Experiments and Results

We evaluate DSVR based on the comparison to the rigid SVR method [5] that was recently reported to produce the best results for placenta reconstruction [19] and its GPU version [11] was employed for 3D fetal cardiac reconstruction in [17]. Furthermore, DSVR is compared to the recently introduced PVR method designed for motion correction in large FoV regions [20]. The reconstruction quality is evaluated using intensity and structural similarity metrics.

A. Fetal MRI data

The fetal MRI data used for evaluation contains 20 iFIND T2-weighted datasets of fetuses from 28-31 weeks GA range. This particular GA range was selected due to the lower amplitude of movement and higher prevalence of bending and stretching [7], since this study focuses on correction of local non-rigid deformations of organs rather than global change of fetal body position.

The iFIND acquisitions were performed on a 1.5T MRI using ssFSE sequence with TR = 15000 ms, TE = 80 ms, voxel size = 1.25 x 1.25 x 2.5 mm, slice thickness 2.5 mm and slice spacing 1.25 mm. The stacks were acquired under different orientations, with 100-160 slices per stack, depending on GA and orientation. Each of the datasets contains 6 stacks with minimum 3 different orientations without major SNR loss.

The datasets were divided into 2 groups with respect to severity of motion. The minor motion group contains 10 cases that include stacks only with minimal loss of structural information. In 10 cases from the severe motion group, the majority of stacks have severe misalignment of slices. The severity of motion corruption was visually estimated by an operator with respect to the consistency of volumetric information in all three planes (Fig. 1) as well as the changes of the global fetal body position between the stacks. Template stack selection was performed manually based on the degree of motion corruption and the position of the fetal body similar to the other stacks.

The amount of motion in each stack was also assessed as the average NCC between sequential slices for a masked ROI. Tab. I demonstrates that the severe motion datasets have lower average NCC range than the minor motion datasets.

Table I. Stack Motion Corruption Assessment: Sequential Slice NCC.

Minor motion datasets Severe motion datasets
0.688 ± 0.071 0.475 ± 0.130

B. Simulated experiment

In order to assess the general capability of DSVR to recover consistent volumetric information and local anatomy features, we perform a phantom experiment with simulated non-rigid motion. At first, a high quality volume reconstructed from a minimal motion dataset is selected as a reference. Next, five sets of slice transformations (incorporating both local non-rigid and global rigid and non-rigid components) extracted from other existing reconstruction cases are used to generate motion-corrupted stacks from the reference volume.

In this experiment, each of the five simulated datasets contains six generated stacks that have different orientations and a mask defining trunk ROI in the template stack. The default DSVR reconstruction pipeline is executed for all datasets. In addition, rigid SVR [5] and PVR [20] reconstructions are performed for comparison to the state-of-the-art methods.

Fig. 5 shows coronal plane view of the original reference volume (X), one of the stacks with simulated deformable motion (S′), DSVR reconstructed volume (XDSV R) along with its difference with the reference (XDSV RX) and SVR (XSV R) and PVR (XPV R) results. Prior to the analysis of the results, in order to avoid possible impact of the global change of the body position, the reconstructed volumes were aligned to the reference using FFD registration with large CP spacing (15 mm). All three methods successfully reconstructed the major anatomy structures including topology of kidneys and spine. However, in SVR and PVR outputs, misregistrations due to non-rigid deformations led to blurring of texture of local features and higher errors.

Fig. 5.

Fig. 5

Simulated experiment: original reference volume (X), one of the generated motion-corrupted stacks (S′), SVR (XSV R), PVR (XPV R) and DSVR (XDSV R) reconstruction results and their difference with the reference (coronal plane).

The corresponding quantitative comparison of the motion-free reference volume and five DSVR, SVR and PVR reconstructions in terms of normalised root mean square error (NRMSE), peak signal-to-noise ratio (PSNR) and NCC computed for the same masked ROI covering fetal trunk is given in Tab. II. All results are statistically significant with p < 0.001. There is a strong correlation between the original and DSVR volumes. The worse results for the SVR outputs indicate that the impact of non-rigid deformations on texture cannot be resolved by rigid registration even with rejection of outliers. The reconstructed PVR volumes also have lower similarity to the reference volumes.

Table II. Simulated Experiment: SVR, PVR and DSVR vs. Reference.

Method NRMSE PSNR NCC
SVR 0.119 ± 0.002 28.916 ± 0.479 0.938 ± 0.003
PVR 0.182 ± 0.009 25.602 ± 0.383 0.882 ± 0.009
DSVR 0.078 ± 0.008 32.560 ± 0.859 0.973 ± 0.006

We also calculated target registration error (TRE) to evaluate average displacement error of the estimated non-rigid transformations compared to the the original ones used for generation of motion-corrupted stacks. Comparison of TRE for DSVR and SVR reconstruction is presented in Tab. III. DSVR TRE values are significantly lower in comparison to SVR and statistically significant with p < 0.001.

Table III. Simulated Experiment: SVR vs. DSVR TRE [mm].

SVR TRE [mm] DSVR TRE [mm]
2.279 ± 0.509 0.797 ± 0.158

C. Reconstruction of fetal body

In fetal body reconstruction, due to the absence of the ground truth as well as the constantly changing shape of the trunk organs, assessment of the quality of reconstructed volumes is challenging. In [5], leave-one-out analysis was proposed for evaluation of SVR results. It is based on the comparison of the original {Yk*} to simulated {Ȳk} slices for a stack that was registered in SVR step but excluded from SR reconstruction thus not contributing to the output volume.

Tab. IV presents quantitative comparison of SVR, PVR, DSVR and DSVR with structural outlier rejection (DSVR+S) results for the same masked ROI of the excluded stack for 20 datasets. The values for PVR results are given primarily for a reference, since, due to the differences in implementation, comparison is performed for simulated and original patches with overlapping regions rather than for slices. Furthermore, PVR employs different SR reconstruction pipeline and does not provide an option for stack exclusion.

Table IV. Fetal Body Reconstruction. Leave-one-out Analysis: SVR, PVR, DSVR and DSVR+S.

Method NRMSE PSNR NCC
Minor motion group (10 datasets):
SVR 0.229 ± 0.023 24.328 ± 0.920 0.780 ± 0.084
PVR* 0.264 ± 0.050 22.438 ± 1.371 0.712 ± 0.120
DSVR 0.177 ± 0.021 26.764 ± 1.139 0.863 ± 0.038
DSVR+S 0.174 ± 0.025 26.764 ± 1.249 0.867 ± 0.043
Severe motion group (10 datasets):
SVR 0.275 ± 0.025 23.219 ± 0.861 0.646 ± 0.068
PVR* 0.280 ± 0.046 21.873 ± 1.436 0.633 ± 0.107
DSVR 0.222 ± 0.032 25.422 ± 0.619 0.831 ± 0.043
DSVR+S 0.214 ± 0.033 25.746 ± 0.597 0.844 ± 0.040
*

PVR comparison was performed on the patch level.

The results for both minor and severe motion datasets show that DSVR surpasses SVR and PVR for both intensity and structural characteristics. Additional structural outlier rejection (DSVR+S) produces a significant improvement only for the severe motion datasets. This is expected since minor motion assumes high NCC values of registration output and DSVR+S should produce only minimal impact. All results apart for comparison of DSVR and DSVR+S for minor motion cases are statistically significant with p < 0.005.

For one of the minor motion cases shown in Fig. 6 (sagittal plane view) when the trunk positions in all stacks are approximately aligned and there are no severe non-rigid deformations, SVR successfully reconstructs the global trunk topology. However, due to bending motion, there is a noticeable loss of structure in the spine region as well as the general degradation of texture. PVR allows reconstruction of the large ROI and partially resolves these artefacts improving definition of the spine. However, the introduced smoothing lowers image quality in terms of interpretation and resolution of small features. On the other hand, DSVR results are characterised by high definition of the local anatomy structures. It also has to be noted that there is a noticeable change in the position of the trunk between different reconstruction methods caused by the different approaches for initialisation of the registration target (Xinit). SVR and PVR use the average of all stacks after global rigid stack registration and DSVR uses the selected template stack. Therefore, the SVR and PVR solutions converge to an intermediate averaged state, whereas DSVR converges to the trunk shape in the template stack.

Fig. 6.

Fig. 6

Example of motion correction for a minor motion dataset: motion corrupted stack, SVR, PVR and DSVR reconstructions (sagittal plane).

A typical example of failed rigid SVR due to non-rigid motion is given in Fig. 7 where one of the original slices Yk is compared to the corresponding simulated slices Ȳk from SVR and DSVR reconstructions. The kidney and bladder regions are segmented in order to assess the registration accuracy. In this case, SVR could not correct the impact of spine bending thus converging to an average position with displaced kidney and resulting in large errors {ejk}. On the other hand, FFD registration improves the mapping (Tk) between Yk and X (n,q). The deformation of ROI boundaries in DSVR output indicates the high degree of non-rigid deformation.

Fig. 7.

Fig. 7

Comparison of SVR and DSVR in presence of non-rigid motion: original acquired slice (Yk) vs. slices simulated from SVR and DSVR (Ȳk). The yellow isolines delineate the structure in the original slice, and show misalignment with SVR reconstruction due to limitation of rigid motion correction. The problem is resolved by non-rigid motion correction in DSVR.

For one of the severe motion dataset results shown in Fig. 8 (coronal plane view), large slice misregistration errors lead to a severe degradation of local features in SVR. PVR resolves this producing a clear trunk structure, however, similarly to the previous example, the smoothed texture lowers the quality of definition of abdominal organs. Although there is an improvement in DSVR vs. SVR output, a significant amount of artefacts due to misregistrations still remains. As mentioned in Sec. III, the employed gradient-descent FFD method is not capable of resolving large bending and rotations therefore leading to misregistrations.

Fig. 8.

Fig. 8

Example of motion correction for a severe motion dataset: motion corrupted stack, SVR, PVR and DSVR reconstructions (coronal plane).

Structural outlier removal (DSVR+S) improves the output by minimising the contribution of registration errors to reconstruction. It also increases the proportion of rejected slices, which, for severe motion datasets, can vary between 20 - 50 % of the total slice number. In comparison, the original EM robust statistics [5] results in only 10 - 15 % slices being rejected, which seems to be insufficient in major motion cases.

D. Parametric study

Regularisation parameters that control the smoothness of the transformations in DSVR have significant impact on the quality of reconstruction. Transformations with large CP spacings allow for correction of the global body shape, but are not efficient for recovery of local features. On the other had, too small CP spacing will lead to overfitting and unrealistic deformations. We experimentally determined that the optimal multi-resolution scheme for fetal body dimensions is 15mm → 10mm → 5mm for 3 iterations.

We evaluated the regularisation schemes with the CP spacings d (0), d (1) = 2/3 · d (0) and d (2) = 1/3 · d (0) for the subsequent DSVR iterations. The d (0) value was varied from 30mm to 5mm for five severe motion datasets. We calculated NCC between simulated {Ȳk} and original {Yk*} slices for the masked fetal body ROI in all stacks. Fig. 9 shows the average NCC values over the five subjects for different initial transformation resolutions d (0). The average SVR output value is provided for the reference. We can observe that d (0)=15mm results in optimal performance and further refinement does not improve the results due to overfitting to the motion artifacts.

Fig. 9.

Fig. 9

CP spacing analysis: NCC between the original (Yk*) vs. simulated (Ȳk) slices.

In addition, we performed quantitative assessment of the impact of the structural outlier rejection steps on the reconstruction quality and number of excluded slices. Tab. V presents the results for 5 severe motion datasets with the default processing pipeline and different combinations of outlier rejection methods: EM-based method proposed in [5], global structural slice-level rejection based on global NCC (G-STR), and local region rejection method based on SSIM maps (L-STR). We report NCC between original Yk and simulated slices Ȳk of the excluded stack within the masked fetal body ROI, and the average proportion of the excluded slices. It can be observed that combination of all three outlier rejection steps results in best performance in severe motion cases.

Table V. Outlier Rejection Scheme Assessment: NCC Between (Yk) and (Ȳk) and Proportion of Excluded Slices.

Method NCC % of excluded slices
EM 0.828 ± 0.028 13.41 ± 2.88 %
EM + G-STR 0.839 ± 0.037 38.51 ± 15.28 %
EM + L-STR 0.834 ± 0.031 18.85 ± 5.17 %
EM + L-STR + G-STR 0.841 ± 0.036 35.59 ± 12.29 %
L-STR + G-STR 0.837 ± 0.039 30.19 ± 11.06 %

E. Reconstruction of placenta

Furthermore, the performance of SVR, PVR and DSVR for reconstruction of placenta was compared for 10 iFIND fetal MRI cases. Each of the datasets contains 5 orthogonal stacks covering the entire uterus and reconstructions are performed for the masked uterus ROI similarly to [19]. The corresponding results of the leave-one-out analysis are presented in Tab. VI in terms of NCC between the original Yk and simulated Ȳk slices of the same masked ROI of an excluded stack. The proposed DSVR method outperforms both other methods. All results are statistically significant with p < 0.005. Fig. 10 shows an example of motion-corrected outputs of SVR, PVR and DSVR of placenta. Placenta is primarily affected by respiratory motion, and is subject to stretching and bending. As a result, DSVR provides better reconstruction quality compared to the alternative methods.

Table VI. Placenta Reconstruction. Leave-one-out Analysis: SVR, PVR and DSVR (NCC).

SVR PVR* DSVR+S
0.639 ± 0.085 0.726 ± 0.045 0.792 ± 0.072
*

PVR comparison was performed on the patch level.

Fig. 10.

Fig. 10

Example of motion correction for placenta: motion-corrupted stack, SVR, PVR and DSVR+S reconstructed volumes (coronal plane).

F. Qualitative analysis

To assess the applicability of DSVR across gestational ages, we performed qualitative evaluation for 100 iFIND fetal body cases randomly selected from 20-34 weeks GA range. The reconstructed volumes were graded by trained clinicians with respect to the image quality in [0; 4] range (4 corresponding to high quality). Volumes with grades ≥ 2 were considered to have sufficient quality for further analysis and interpretation. The distribution of grades for each week of GA is presented in Fig. 11. Average grades per week of GA varied within 2.5 - 3.5 range (3.09 ± 0.78). The primary causes of lower grades were motion for younger subjects and low SNR for older subjects. Only 6% of all cases failed (grades < 2) due to severe motion which could not be resolved by gradient descent based FFD registration. Further details of the analysis are presented in Supplementary material.

Fig. 11.

Fig. 11

Quality of DSVR reconstructions of fetal body region for 100 reconstructed iFIND cases vs. GA. Squares represent the average grade per week of GA, and bars represent the range of values.

VI. Discussion and Conclusions

We proposed and implemented a novel DSVR method for compensation of non-rigid motion in fetal MRI. It allows reconstruction of high resolution 3D volumes from multiple stacks of slices affected by non-rigid motion. Therefore, DSVR extends application of slice-to-volume reconstruction to fetal body and placenta.

Unlike the conventional rigid SVR methods, DSVR is capable of correction of local deformations of organs caused by bending and stretching. Correction of both in- and out-of-plane non-rigid motion is ensured by registration of the volume to slices rather than slices to volume. The challenge of the absence of a ’stable’ shape is addressed by hierarchical FFD SVR scheme initialised by one of the stacks that gradually converges to a stable state. The fact that the stacks have different orientations helps to prevent overfitting to motion artefacts. In addition, structure-based outlier rejection step is introduced to minimise the impact of misregistration errors on the reconstructed volume.

The method was quantitatively evaluated on 20 fetal MRI datasets from 28-31 weeks GA range. This age range was chosen due to high incidence of stretching and bending, whereas large rotation and translation motion is less prevalent compared to younger subjects [7]. Comparison to the state-of-the-art solutions showed that DSVR surpasses both SVR and PVR methods for minor and severe motion datasets. This was further confirmed by an additional experiment with simulated non-rigid motion. DSVR reconstructed 3D fetal body volumes are characterised by well defined features and texture of the spine, heart and abdominal organs.

The current implementation of DSVR is not designed for correction of large amplitude motion, including, large rotations and bending, due to the limited capture range of the employed gradient descent based optimisation methods, and therefore relies on outlier rejection in such cases. The qualitative study of DSVR reconstruction of fetal body for 100 iFIND cases from 20-34 weeks GA range showed that large amplitude motion primarily affects datasets of subjects under 25 weeks GA. In future, this limitation can be addressed by application of CNN-based methods, as already proposed for rigid SVR fetal brain reconstruction [14], [15].

We further demonstrated that DSVR outperforms both SVR and PVR for 3D placenta reconstruction. In future, introducing additional decoupling of maternal motion would potentially improve correction of large amplitude motion.

Although DSVR reconstructed volumes can be used for qualitative analysis, the question of volume-preservation for DSVR reconstruction still remains open. Quantitative measurements performed on DSVR reconstructed volumes can be influenced by various factors such as the position and shape of the fetal body in the template, number of stacks or CP spacing. This limitation should be addressed in future by introducing further model-based constrains on deformation fields along with automation of template selection and masking steps.

Supplementary Material

Supplementary material presents the details of qualitative analysis of DSVR reconstruction quality for 100 iFIND cases, the full versions of Fig. 5, 6, 8, DSVR reconstruction of the whole uterus for twin pregnancy case, and visual representation of the DSVR algorithm.

Supplementary material

Acknowledgment

Thank you to Matthew Fox, Joanna Allsop, Ana Gomes and Elaine Green for their oversight during the scanning of volunteers and patients. Thank you to Jacqueline Matthew, Milou van Poppel and Johannes Steinweg for grading DSVR reconstruction quality for the investigated iFIND cases.

The iFIND project data used in this research were collected subject to the informed consent of the participants. The MIRTK library was used under the Apache License, Version 2.0. The original IRTK reconstruction code was used under the creative commons public license from IXICO Ltd. The views expressed are those of the authors and not necessarily those of the NHS, the NIHR or the Department of Health.

This work was supported by the NIH Human Placenta Project grant [1U01HD087202-01], the Wellcome EPSRC Centre for Medical Engineering at King’s College London (WT 203148/Z/16/Z), the Wellcome Trust and EPSRC IEH award [102431] for the iFIND project and by the National Institute for Health Research (NIHR) Biomedical Research Centre based at Guy’s and St Thomas’ NHS Foundation Trust and King’s College London.

Footnotes

References

  • [1].Malamateniou C, et al. Motion-compensation techniques in neonatal and fetal MR imaging. American Journal of Neuroradiology. 2013;34(6):1124–1136. doi: 10.3174/ajnr.A3128. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • [2].Story L, Rutherford M. Advances and applications in fetal magnetic resonance imaging. The Obstetrician & Gynaecologist. 2015;17(3):189–199. [Google Scholar]
  • [3].Gholipour A, Estroff JA, Warfield SK. Robust super-resolution volume reconstruction from slice acquisitions: Application to fetal brain MRI. IEEE Transactions on Medical Imaging. 2010;29(10):1739–1758. doi: 10.1109/TMI.2010.2051680. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • [4].Rousseau F, Kim K, Studholme C, Koob M, Dietemann J-L. On Super-Resolution for Fetal Brain MRI. Med Image Comput Comput Assist Interv. 2010;13(Pt 2):355–362. doi: 10.1007/978-3-642-15745-5_44. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • [5].Kuklisova-Murgasova M, Quaghebeur G, Rutherford MA, Hajnal JV, Schnabel JA. Reconstruction of fetal brain MRI with intensity matching and complete outlier removal. Medical Image Analysis. 2012;16(8):1550–1564. doi: 10.1016/j.media.2012.07.004. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • [6].Kyriakopoulou V, et al. Normative biometry of the fetal brain using magnetic resonance imaging. Brain Structure and Function. 2017;222(5):2295–2307. doi: 10.1007/s00429-016-1342-6. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • [7].Nowlan NC. Biomechanics of foetal movement. European Cells and Materials. 2015;29:1–21. doi: 10.22203/ecm.v029a01. [DOI] [PubMed] [Google Scholar]
  • [8].Rousseau F, et al. Registration-Based Approach for Reconstruction of High-Resolution In Utero Fetal MR Brain Images. Academic Radiology. 2006;19(9):1072–1081. doi: 10.1016/j.acra.2006.05.003. [DOI] [PubMed] [Google Scholar]
  • [9].Jiang S, Xue H, Glover A, Rutherford M, Rueckert D, Hajnal JV. MRI of moving subjects using multislice Snapshot images with Volume Reconstruction (SVR): Application to fetal, neonatal, and adult brain studies. IEEE Transactions on Medical Imaging. 2007;26(7):967–980. doi: 10.1109/TMI.2007.895456. [DOI] [PubMed] [Google Scholar]
  • [10].Tourbier S, Bresson X, Hagmann P, Thiran JP, Meuli R, Cuadra MB. An efficient total variation algorithm for super-resolution in fetal brain MRI with adaptive regularization. NeuroImage. 2015;118:584–597. doi: 10.1016/j.neuroimage.2015.06.018. [DOI] [PubMed] [Google Scholar]
  • [11].Kainz B, et al. Fast Volume Reconstruction from Motion Corrupted Stacks of 2D Slices. IEEE Transactions on Medical Imaging. 2015;34(9):1901–1913. doi: 10.1109/TMI.2015.2415453. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • [12].Salehi S, et al. Real-Time Automatic Fetal Brain Extraction in Fetal MRI by Deep Learning. ISBI 2018; 2018. pp. 720–724. [Google Scholar]
  • [13].Ebner M, et al. An Automated Localization, Segmentation and Reconstruction Framework for Fetal Brain MRI. MICCAI 2018; 2018. pp. 313–320. [Google Scholar]
  • [14].Hou B, et al. 3D Reconstruction in Canonical Co-ordinate Space from Arbitrarily Oriented 2D Images. IEEE Transactions on Medical Imaging. 2018;37(8):1737–1750. doi: 10.1109/TMI.2018.2798801. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • [15].Salehi SS, Khan S, Erdogmus D, Gholipour A. Real-Time Deep Pose Estimation With Geodesic Loss for Image-to-Template Rigid Registration. IEEE Transactions on Medical Imaging. 2019;38(2):470–481. doi: 10.1109/TMI.2018.2866442. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • [16].Kainz B, et al. Motion corrected 3D reconstruction of the fetal thorax from prenatal MRI. MICCAI 2014; 2014. pp. 284–291. [DOI] [PubMed] [Google Scholar]
  • [17].Lloyd DFA, et al. Three-dimensional visualisation of the fetal heart using prenatal MRI with motion corrected slice-volume registration. The Lancet. 2019;393(10181):1619–1627. doi: 10.1016/S0140-6736(18)32490-5. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • [18].van Amerom JF, et al. Fetal cardiac cine imaging using highly accelerated dynamic MRI with retrospective motion correction and outlier rejection. Magnetic Resonance in Medicine. 2018;79(1):327–338. doi: 10.1002/mrm.26686. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • [19].Torrents-Barrena J, et al. Fully automatic 3D reconstruction of the placenta and its peripheral vasculature in intrauterine fetal MRI. Medical Image Analysis. 2019;54:263–279. doi: 10.1016/j.media.2019.03.008. [DOI] [PubMed] [Google Scholar]
  • [20].Alansary A, et al. PVR: Patch-to-Volume Reconstruction for Large Area Motion Correction of Fetal MRI. IEEE Transactions on Medical Imaging. 2017;36(10):2031–2044. doi: 10.1109/TMI.2017.2737081. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • [21].Tadayyon H, Lasso A, Kaushal A, Guion P, Fichtinger G. Target motion tracking in MRI-guided transrectal robotic prostate biopsy. IEEE Transactions on Biomedical Engineering. 2011;58(11):3135–3142. doi: 10.1109/TBME.2011.2163633. [DOI] [PubMed] [Google Scholar]
  • [22].Xu H, Lasso A, Fedorov A, Tuncali K, Tempany C, Fichtinger G. Multi-slice-to-volume registration for MRI-guided transperineal prostate biopsy. Int J Comput Assist Radiol Surg. 2015;10(5):563–572. doi: 10.1007/s11548-014-1108-7. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • [23].Osechinskiy S, Kruggel F. Slice-to-Volume Nonrigid Registration of Histological Sections to MR Images of the Human Brain. Anatomy Research International. 2010;2011:1–17. doi: 10.1155/2011/287860. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • [24].Rivaz H, Karimaghaloo Z, Fonov VS, Collins DL. Nonrigid registration of ultrasound and MRI using contextual conditioned mutual information. IEEE Transactions on Medical Imaging. 2014;33(3):708–725. doi: 10.1109/TMI.2013.2294630. [DOI] [PubMed] [Google Scholar]
  • [25].Rueckert D, Sonoda LI, Hayes C, Hill DLG, Leach MO, Hawkes DJ. Nonrigid Registration Using Free-Form Deformations: Application to Breast MR Images. IEEE Transactions on Medical Imaging. 1999;18(8):712–721. doi: 10.1109/42.796284. [DOI] [PubMed] [Google Scholar]
  • [26].Ferrante E, Paragios N. Graph-based slice-to-volume deformable registration. Int J Comput Vision. 2018 Jan;126(1):36–58. [Google Scholar]
  • [27].Nir G, et al. Model-based registration of ex vivo and in vivo MRI of the prostate using elastography. IEEE Transactions on Medical Imaging. 2013;32(7):1349–1361. doi: 10.1109/TMI.2013.2269174. [DOI] [PubMed] [Google Scholar]
  • [28].Kuklisova-Murgasova M, et al. Distortion correction in fetal EPI using non-rigid registration with a Laplacian constraint. IEEE Trans Med Imaging. 2018;37(1):12–19. doi: 10.1109/TMI.2017.2667227. [DOI] [PubMed] [Google Scholar]
  • [29].Wang Z, et al. Image Quality Assessment: From Error Visibility to Structural Similarity. IEEE Transactions on Image Processing. 2004;13(4):600–612. doi: 10.1109/tip.2003.819861. [DOI] [PubMed] [Google Scholar]
  • [30].Schuh A, et al. Construction of a 4D brain atlas and growth model using diffeomorphic registration. 1st ed. Vol. 8682. Springer International Publishing Switzerland; 2015. pp. 27–37. [Google Scholar]
  • [31].Van Vuuren SH, et al. Size and volume charts of fetal kidney, renal pelvis and adrenal gland. Ultrasound in Obstetrics and Gynecology. 2012;40(6):659–664. doi: 10.1002/uog.11169. [DOI] [PubMed] [Google Scholar]
  • [32].Modat M, et al. A parallel-friendly normalized mutual information gradient for free-form registration. SPIE Medical Imaging; 2009. 72 590L–1–8. [Google Scholar]

Associated Data

This section collects any data citations, data availability statements, or supplementary materials included in this article.

Supplementary Materials

Supplementary material

RESOURCES