{'fraction_non_alphanumeric': 0.03128817433108721, 'fraction_numerical': 0.01264924932025062, 'mean_word_length': 4.8758971984255615, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Temporal Registration in Application to In-utero MRI Time Series

We present a robust method to correct for motion in volumetric in-utero MRI time series. Time-course analysis for in-utero volumetric MRI time series often suffers from substantial and unpredictable fetal motion. Registration provides voxel correspondences between images and is commonly employed for motion correction. Current registration methods often fail when aligning images that are substantially different from a template (reference image). To achieve accurate and robust alignment, we make a Markov assumption on the nature of motion and take advantage of the temporal smoothness in the image data. Forward message passing in the corresponding hidden Markov model (HMM) yields an estimation algorithm that only has to account for relatively small motion between consecutive frames. We evaluate the utility of the temporal model in the context of in-utero MRI time series alignment by examining the accuracy of propagated segmentation label maps. Our results suggest that the proposed model captures accurately the temporal dynamics of transformations in in-utero MRI time series.


I. INTRODUCTION
This paper demonstrates a robust temporal image registration framework. Large motion presents significant challenges for time-course analysis in dynamic image series. To correct for motion in medical images, registration that produces voxel correspondences is commonly employed. Unfortunately, registration often fails if the motion is substantial. We take advantage of the temporal smoothness of motion to achieve robust alignment of volumetric MRI time series. This work is motivated by a clinical research application in monitoring placental and fetal health, which we introduce below.

A. Motivating application
In-utero blood oxygenation level dependent (BOLD) magnetic resonance imaging (MRI) is a promising imaging tool for studying functional dynamics of placenta and fetal organs [1]- [4]. Changes in fetal and placental oxygenation levels with maternal hyperoxygenation can be used for charactering and detecting placental dysfunction [5]. Temporal MRI data suffers from serious motion artifacts due to maternal respiration, unpredictable fetal movement, and signal non-uniformities [6], as illustrated in Fig. 1. Different types of motion are presented in these MRI time series. The fetal brains move as rigid objects with a wide range of motion; the placenta moves more locally and deforms non-rigidly. In our study, the image series include about 500 volumes (frames) that are analyzed to study fetal and placental physiology.

B. Related Work
Temporal registration has been previously investigated in time series of cardiac images [7]- [10] and lung images [10]- [13]. Both cardiac motion and breathing patterns are somewhat regular and smooth across time, enabling biomechanical modeling [7], [8], [11]. Electrocardiography (ECG) gating or respiratory gating can be incorporated in the cardiac or lung image motion correction, by either enforcing the transformation models to be periodic or penalizing non-periodic transformations in the cost function. These approaches cannot be generalized to images that contain complex non-periodic motion patterns.
Existing registration techniques for dynamic medical imaging data can be categorized into two distinct groups. The first type of methods performs registration based on the differences between each pair of consecutive frames. Each pair is aligned, and the estimated transformations are concatenated to estimate alignment of all images in the series [13], [14]. In application to long image series like our data, this approach leads to substantial accumulated errors after several concatenation steps. These methods only take the moving image and the reference image into account in each pairwise registration step, and ignore essential temporal information from other images in the series. The second type applies groupwise registration to simultaneously align each individual image to a group template -estimated or selected from one of the input images -aiming to yield acceptable registration results for all images in the series [10], [12], [15]- [21]. Taking this approach in our application, registration fails for a great number of images that are substantially different from the template, necessitating outlier detection and rejection [22]. Alternatively, a global objective function based on pairwise differences between each individual image in the series and an implicitly defined template frame can incorporate the temporal information [10]. This leads to a global optimization process, i.e., the algorithm performs pairwise registration of consecutive frames iteratively 1 2 124 125 Fig. 1. Example twin pregnancy case from the study. The same 2D cross-section from volumetric frames 1, 2, 124, and 125 is shown. Arrows indicate areas of substantial motion of the fetal heads (green and yellow) and highly non-rigid motion of the placenta (red).
until the entire series comes into alignment. This type of global optimization is extremely challenging in large image sets like our time series. The problem of temporal alignment has also been investigated in longitudinal image analysis that aims to characterize temporal variations [20], [23]- [26]. The goal is to examine the behavior of a biological system driven by development, aging, or disease processes [25]. The statistical analysis is performed via temporal shape regression on a Riemannian manifold [20], [23]- [25] to construct a trajectory that characterizes the shape variation across time. Longitudinal data is often sparsely distributed in time, and temporal smoothness constraint is therefore imposed. Typically, the challenges in longitudinal analysis involve growth or degeneration, which requires modeling tools that capture subtle deformation and are distinct from those in our problem of MRI time series alignment in the presence of substantial motion.
Unlike other applications that can require subjects to remain stationary, fetal and placental imaging is an application with irregular motion that cannot be controlled at the time of acquisition.

C. Proposed Temporal Registration Method
To overcome the difficulty of correcting large and irregular motion in a long MRI time series, we assume a Markov structure that induces temporal smoothness and offers an efficient way to estimate the transformations by serially solving pairwise registration problems. We derive filtered estimates of the transformations from the hidden Markov model that represents the MRI series. The resulting sequential procedure estimates the transformations of the template to each frame in the series. We provide a flexible framework for temporal dense alignment of image time series.
Unlike the methods that perform pairwise registration of consecutive frames and concatenate the estimated transformations, we use temporal smoothness explicitly at every estimation step. The sequential estimation framework yields an efficient algorithm that makes alignment of long MRI series feasible. Our approach is related to filtering methods in respiratory motion modeling [11]. We do not model the motion explicitly but rather capture it through transformation of organs of interest to template to adapt to complex motion patterns.
We demonstrate and evaluate the proposed method in the context of the in-utero BOLD MRI study by providing registration-based tracking for organs of interest. We employ the rigid motion model for fetal brains and a non-rigid B-Spline motion model for placentae. The experimental results offer significant improvements over pairwise registration in alignment of fetal brains that tend to move substantially, and robust improvements in alignment of placentae.
This paper presents an extended version of the preliminary algorithm presented at the International Conference on Medical Image Computing and Computer-Assisted Intervention (MICCAI) [27]. In this paper, we derive in detail the instantiation of a hidden Markov model in the context of MRI time series. We model the motion as a latent state in the hidden Markov model, derive a temporal registration algorithm using filtered estimates, and discuss the challenges of developing a smoothing algorithm. Here we employ rigid and B-Spline transformation models on the regions of fetal brains and placentae respectively, instead of one single transformation model for the entire volumetric image used in [27], to adapt to the different motion patterns in the two regions. Finally, extensive experiments investigate the sensitivity of the resulting registration on the model parameters.
This paper is organized as follows. In the next section, we present the generative model for temporal motion and observed images using hidden Markov model, derive the resulting algorithm for filtered estimates, and provide implementation details. In Section III, we introduce the in-utero BOLD MRI data and report the experimental results. In Section IV, we discuss future directions, methodological challenges of developing a backward message passing algorithm for temporal registration, and potential extension for other transformation models. Conclusions follow in Section V.

II. METHODS
We employ hidden Markov model (HMM) to represent the temporal changes in the MRI time series. One image in the MRI time series is chosen as a template, whose transformations represent the rest of the images in the series. We treat the transformation as a hidden state, and each image as an observation that depends on the template and the corresponding hidden state. In the following subsections, we present inference in HMMs using filtered estimates [28], [29] in the context of temporal registration, and describe the implementation details.

A. HMM and Filtered Estimates
The temporal dynamics of latent (hidden) states in a HMM exhibit the Markov property, i.e., the current state depends on the history only through the previous state. Thus, inference based on the entire series of past observations is equivalent to integrating the information from the previous state and from the corresponding observation, leading to an efficient sequential estimation algorithm.
In our application, one image is chosen as the template I to be treated as a global parameter shared by all observed images. We assume that the template I deforms at each time point, and the transformation ϕ n of the template is the latent state at time n ∈ {1, ..., N }, where N is the total number of images in the series. The observed image J n at time n is generated by applying the anatomical transformation ϕ n to the template I independently of all other time points. Fig. 2 illustrates this model.
Our aim is to estimate the latent state variable ϕ n from the observations J 1:n , where we use J k:m and ϕ k:m to denote sub-series {J k , J k+1 , ..., J m } of the MRI time series {J 1 , J 2 , ..., J N } and the associated transformations {ϕ k , ϕ k+1 , ..., ϕ m }. Formally, the estimation of the latent variable ϕ n is reduced to maximization of its posterior distribution p(ϕ n |J 1:n ; I) parameterized by the template I. This posterior distribution, also known as the filtered estimate of the state, can be efficiently computed using forward message passing [28], [29], also known as sequential estimation. The posterior distribution p(ϕ n |J 1:n ; I) is determined by integrating the previous posterior distribution p(ϕ n−1 |J 1:n−1 ; I) with the temporal dynamics p(ϕ n |ϕ n−1 ) and the likelihood p(J n |ϕ n ; I) of the current observation: where Eq. (2) follows from the Markov property: Applied recursively, the forward message passing produces the posterior distribution p(ϕ n |J 1:n ; I) for each time point n in the number of steps that is linear with n. Similarly, backward message passing enables inference of the posterior distribution p(ϕ n |J 1:N ; I), which is based on all data, often referred to as smoothing. In this paper, we incorporate temporal smoothness in temporal registration by filtering. We discuss the challenges of developing a smoothing algorithm for temporal registration in Section V.
The estimate ϕ * n is then used to compute the optimal transformation ϕ * n+1 at time n + 1, and so on until we reach the end of the time series.

C. Model Setup
We model the likelihood term p(J n |ϕ n ; I) in Eq. (5), also known as emission probability, as the exponential likelihood on image similarity: where Dist(·, ·) is a measure of dissimilarity (distance) between images. This probability measures how different the transformed template and the currently observed image are given the template I and the transformation ϕ n . The state transition probability p(ϕ n |ϕ n−1 ) in Eq. (5) encourages temporal smoothness of transformations and spatial smoothness of ϕ n : where · is an appropriate norm that encourages transformation ϕ n to be similar to ϕ n−1 , and Reg(·) is the regularization term that encourages spatial smoothness and other desired properties. The regularization term Reg(ϕ n ) is often used in image registration to restrict the estimated motion to a specific transformation group. λ 1 and λ 2 are regularization parameters that control the amount of temporal and spatial smoothness.
We manipulate Eq. (5) to obtain ϕ * n ≈ arg max ϕn p(J n |ϕ n ; I) p(ϕ n |ϕ * n−1 ) and observe that this optimization problem is reduced to serial pairwise image registration of the template I and the observed image J n . The algorithm proceeds as follows. Given the transformation estimate ϕ * n−1 of the template I to the image J n−1 , we apply the registration algorithm to I and J n by minimizing the cost function in Eq. (8), while also using ϕ * n−1 as an initialization, resulting in the estimate ϕ * n . This HMM-based temporal transformation model can be easily generalized to different types of similarity metrics between moving images and fixed images, temporal regularization, and transformation models. The model can be readily augmented to account for signal non-uniformities and a latent template that is estimated jointly with the transformations, similar to prior work in groupwise registration [10], [12], [20], [21]. In the following subsection, we describe our implementation for aligning in-utero BOLD MRI time series.

D. Implementation Details
We implemented our method using the Insight Segmentation and Registration Toolkit (ITK) [30]. Negative crosscorrelation [31] is employed as the metric of image dissimilarity Dist(·, ·) as cross-correlation can handle images with signal non-uniformities that are present in the in-utero BOLD MRI data. Cross-correlation is computed over a 5×5×5 voxel volumetric window. We implemented the rigid transformation model to account for the fetal brain motion and the B-Spline transformation model to characterize the placental transformation. We apply rigid transformation model for fetal brain alignment by defining a region-of-interest (ROI) mask in the template that is used to compute the registration cost function. For placentae, the B-Spline transformation model corrects the non-rigid deformation in addition to motion. For the B-Spline model, the grid spacing was set to be 10×10×10 voxels, in the volume of 150×150×120 voxels of our MRI data.
Our choice of the transformation model offers a natural implementation of spatial smoothness. For the rigid model, the third term in Eq. (8) disappears, and for B-Spline model, it becomes a sum-of-squares of the B-Spline coefficients. The temporal regularization is implemented as the L 2 norm of the difference between the current transformation and the previously estimated transformation parameters. Specifically, the six parameters in 3-D rigid transformation model are split into two groups: three parameters representing rotation angles and three parameters representing translation. The weight of the second term in Eq. (8) is then represented by λ R 1 and λ T 1 respectively for the two groups of parameters. λ R 1 and λ T 1 are set to 1 and 0 in the rigid registration. λ 1 is set to 0.005 in the B-Spline registration.
Note that the template can be any image in the series, and the algorithm can proceed in either ascending or descending order of the series. In our experiments, we chose the image in the middle of each series as the template I = J 250 . In each series, we register images J 1:249 to the template J 250 in a descending order, and register images J 251:N to the template J 250 in an ascending order.

III. EXPERIMENTAL RESULTS
We evaluated our method on in-utero BOLD MRI time series acquired as part of a placental function study [32].

A. Data
Ten pregnant women were consented and scanned on a 3T Skyra Siemens scanner (single-shot GRE-EPI, 3 × 3mm 2 in-plane resolution, 3mm slice thickness, interleaved slice acquisition, TR= 5.8 − 8s, TE= 32 − 36ms, FA= 90 o ) using 18-channel body and 12-channel spine receive arrays. Each series contains around 300 volumetric images. This study included three singleton pregnancies, six twin pregnancies, and one triplet pregnancy, between 28 and 37 weeks of gestational age.
To eliminate the effects of slice interleaving, we resampled odd and even slices of each volumetric image onto a common isotropic 3mm 3 image grid. We split each interleaved volumetric image into even and odd slices, and linearly interpolated the acquired slices to complete the volume. The number of images in each MRI series was thus doubled to around 600. To enable quantitative evaluation, we manually delineated the fetal brains (total of 18) and placentae (total of 10) in six frames (1, 50, 100, 150, 200, 250) in each series.

B. Evaluation
To evaluate the temporal model, we compare it to a variant of our method that does not assume the temporal structure and instead aligns each image in the series to the reference frame using the same registration algorithm used by our method. This approach is commonly used for fetal MRI series alignment [22], [33]. Algorithmically, this corresponds to setting λ 1 = 0 in Eq. (8) and initializing the registration step with an identity transformation instead of the previously estimated transformation ϕ * n−1 . The other parameter settings in the temporal model are the same with those in the pairwise model.
To quantify the accuracy of the alignment, we transform the manual segmentation labels (fetal brains and placentae) in the template I = J 250 to the five segmented frames (J 1 , J 50 , J 100 , J 150 , J 200 ) in each series using the estimated transformations. We employ Dice coefficient [34] to quantify volume overlap between the transferred and the manual segmentations. In our application, the goal is to study the average temporal signal for each ROI, and therefore delineation of an ROI provides an appropriate evaluation target. In addition, we also investigate how the alignment accuracy changes with the weight λ 1 of the temporal smoothness term in Eq. (8). For B-Spline model, we additionally evaluate the log determinants of Jacobian of the estimated transformation fields. The top row shows statistics of 18 fetuses for our method (red), pairwise registration to the template frame (green), and no motion correction (blue). The 18 fetuses in the study are reported in the increasing order of volume overlap of original segmentation labels (no motion correction). The two plots in the bottom row show the Dice coefficients of the 90 images from temporal registration and pairwise registration and Dice coefficients from temporal registration and no motion correction. The volume overlap for temporal registration is 9.0%(±16.4%) on average higher than that from pairwise registration. In subjects with the pairwise volume overlap below 0.75 (34 out of 90), the average volume overlap improvement is 31.2%(±20.9%). Subfigure (b): Volume overlap as a function of the weights of temporal regularization. The plot shows the volume overlap in the region of the fetal brain with rigid transformation model as a function of λ R 1 , when λ T 1 = 0 (red), λ T 1 = 0.1 (green), and λ T 1 = 1 (blue).

C. Results
in the regions of fetal brains and placentae to represent the first frame J 1 in the series that is substantially different from the template. The delineations achieved by transferring manual segmentation labels from the reference frame to the coordinate system of the selected frame are in good alignment with the manual segmentation outlines for that selected frame. Fig. 4 reports volume overlap statistics (Dice coefficient) for the fetal brains in this study. We observe that temporal alignment improves volume overlap in important ROIs over the baseline pairwise registration. We note that temporal alignment offers particularly substantial gains in cases with a lot of motion, i.e., low original volume overlap, in fetal brains. Fig. 4 also shows how alignment accuracy changes with the weights of temporal smoothness regularization in Eq. (8). For the rigid model, the Dice coefficients are more sensitive to λ T 1 , and the registration achieves the best accuracy for λ T 1 = 0. Note that the temporal smoothness of translation is still somewhat maintained by initializing the current registration with the previous transformation estimation, even if λ T 1 = 0. In our tests, the temporal rigid registration yields the best result when λ T 1 = 0 and λ R 1 = 1 for fetal brains. Fig. 5 reports volume overlap statistics for the placentae in this study. Our temporal alignment offers consistent improvement for placentae for all cases over pairwise registration to the reference frame. Fig. 5 also reports the sensitivity of the alignment accuracy with respect to the weight of temporal smoothness regularization in Eq. (8). For B-Spline registration on placenta, the algorithm yields the best results when λ 1 = 0.005. Finally, Fig. 5 shows the mean and the standard deviation of the log determinants in each case (in total, 10 cases). The mean values of the log determinants in each subject are all between −0.15 and 0.15. A log determinant value above 0 indicates volume expansion of the transformation, and a log determinant value below 0 indicates volume compression. All the determinants of Jacabian of our estimated B-Spline transformations for placentae are positive.

IV. DISCUSSION
Our temporal registration approach improves alignment of fetal brains and placentae in the in-utero MRI times series. The temporal model offers substantial improvement for fetal brain alignment in the presence of large motion, which presents significant challenges for the pairwise registration method. The incorporation of temporal smoothness improves the alignment of the placenta for every subject. In the study, the log determinants of Jacobian of the estimated B-Spline transformations are close to 0, which indicates no over-expansion or overcompression in the resulting placenta transformation estimates.
One important downstream use of this temporal registration method is modeling of ROI-specific intensity changes in the in-utero MRI time series [22], [32]. Currently the signal dynamics in fetal brains, which tend to have a wide range of motion, include substantial fluctuations, due to spatial signal non-uniformities that affect the estimated signal if the fetus moves substantially in the volumetric image. The crosscorrelation measure of image similarity used in our temporal registration captures local image differences, but tracking the intensity changes in corresponding regions across images is still affected by the signal non-uniformities. Future work will obtain robust estimates of the in-utero MRI signal time-courses by augmenting the temporal method with a model of ROIspecific intensity changes.
In addition, the registration of each frame to the template frame can be based on all the images in the time series by including a backward message pass. One of the challenges is to find a template that is close to both the first observation (image) and the last observation (image) in the HMM. Otherwise, estimating the first or the last hidden state (transformation) might fail due to the substantial difference between the moving image and the template image, since there is no temporal information passed to the very first step of the sequential estimation .
Our temporal registration method naturally accepts different types of transformation models besides the rigid model and the B-Spline model. For example, it can be extended to diffeomorphic transformation. The diffeomorphic demons algorithm [35] adapts convolutional kernels on the incremental vector field and the updated displacement field at each iteration acts as a regularization in the cost function, which can naturally be implemented with our temporal registration method to incorporate temporal smoothness.

V. CONCLUSIONS
We present a HMM-based registration method to align images in MRI time series. Forward message passing incorporates the temporal model of motion into the estimation procedure. The filtered estimates are therefore based on not only the present image frame and the template, but also on the previous frames in the image series. Our method can be easily extended on other forms of transformation models and image similarity metrics. The experiments on in-utero BOLD MRI time series demonstrate the promise of our approach in this novel and challenging application.