1 Introduction
Image registration is a common task in MR image analysis and has been an active research topic in many areas, which aligns the image space into a common anatomical space. The deformable registration strategy calculates the dense correspondence between two images and generates the affine transformations for global alignment, which are usually much slower and have higher degrees of freedom. There are two types of methods in deformable registration, including classical registration methods (nonlearningbased) and learningbased methods. Classical registration methods solve the optimization problems on spatial of deformation (the space of displacement vector field), including elastictype models [1], bsplines [7], statistical parameter mapping [17], Demons [9] or discrete methods [5]. Diffeomorphic transformations have undergone extensive advancement, resulting in stateoftheart tools such as LDDMM, DARTEL and SyN. But these tools demand amounts of time and computational resources for a given image pair on a CPU.
The emergence of deep learning enables neural networks to penetrate into medical image registration. Recently, some methods are proposed to train neural networks that map a pair of input images to output deformation. Balakrishnan et al. [11] performed a rigorous analysis of their method and presented the results on a complete MR volume. They studied MR brain image registration based on 3D slices and the computational efficiency of this method is also higher than the existing general methods. They proposed an unsupervised learning network called VoxelMorph[15], using a pair of 160x192x224 3d images as input. They tested their method on eight publicly available MRI brain image datasets. They also present a probabilistic generative model and derive an unsupervised learningbased inference algorithm [3][4] to provide diffeomorphic guarantees and uncertainty estimates.
In this paper, we propose a new deformable medical image registration method based on differential geometric information which trained on original network in [4] and get the stateoftheart performance among else methods before to the best of our knowledge. The constructions of this paper are: 1. we propose a new deformable medical image registration method based on differential geometric information including the Jacobian determinant(JD) and the curl vector(CV) of diffeomorphic registration field on the VoxelMorph CNN architecture in the paper[4].
2.we also build a registration template for 7 MRI image samples according to average of differential geometric transformation.
2 Method
2.1 The Deformation Method
In deformable registration problem, it is essential to find diffeomorphic registration field . Here, we introduce the deformation method that can generate JD and CV by the grid generation (namely registration field ), which will be used in our later methods. Diffeomorphism is an active research topic in differential geometry [13]. JD and CV play an important role in determining a diffeomorphism. Since JD has the same direct physical meaning in grid generation which represents the size of grid cell and CV represents the grid cell rotations, the deformation method was applied successfully to grid generation and registration problems of MRI brain image[2].
We use MATLAB to generate the grid images (diffeomorphic registration field ) using the deformation method on the gradient and gray value of original image. And we extract the images formed by JD and CV information from the grid images (Fig.1(c)(d)). The following Fig.1(b)(c)(d) show the grid image, JD and CV image of MRI brain respectively, which present geometric deformation features of the image, especially highlight the change of morphological features of CSF, GM, WM. According to manifold learning [10], CNN can extract the manifold features from images. So, we extract JD and CV images from generated known and use them as input of CNN for better extracting manifold features.
2.2 Variational Method for Image Registration
Here, we use variational method to construct diffeomorphisms , through controlling directly the Jacobian determinant and the curl vector of a transformation, which shows the important role of the JD and the CV in determining a diffeomorphism. It is called the Variational Method[13]. Let , , ,, we define the similarity measure in 3D as:
(1) 
Here, is the prescribed Jacobian determinant monitor function, and is the prescribed curl vector monitor function. So, the problem we study now is to construct numerically a transformation , namely, .
For two different images and , image registration process is to fnd a transformation such that . Here, a new registration algorithm based on variational method is developed. Given the moving image T and the fixed image R in 3D, we determine a registration transformation by iteratively minimizing a similarity measure:
(2) 
3 Proposed Framework
3.1 Training Strategy Based on Average JD and CV
The following Fig.2 shows the proposed framework and the training strategy of our method. Given a same atlas or reference volume, N volumes of our training set
were used as the input of the network, and output the approximate posterior probability parameters including the velocity field mean
and variance
. The velocity field z is sampled and transformed to a diffeomorphic deformation field using squaring and scaling layers. Finally, we use a spatial transform function warps to . At registration stage, we generate the image of deformation filed by evaluating during registration between original each two volumes (x,). As for known generated registration filed , we extract JD and CV images (note CV images have 3 volumes) in the first training stage and use them as channels of original networks for second train. So, we can generate N pairs of , for each moving image subject of train set, and combine them with corresponding JD and CV generated from as 5 channels, and for fixed image subject, we combine it with the average image of , namely as 5 channels, finally combined moving and fixed images are input into Voxelmorph to train again and then we test registration between pairs of unseen subjects from testing set. The loss can also be the original composed of two terms by estimating the variational lower bound using SGD methods in MIT’s paper. To make comparison in following experiment, we also get JD and CV from diffeomorphisms which generated directly using deformation method on the gradient and gray value of input moving volumes described in section 2.1, but not from the registration filed in the first training stage, then perform the same process as above.The following Fig.3 shows the deformation filed of initial training stage and generated JD and CV image from . Fig.3 (a) shows the RGB image of and Fig.3 (b)(c)(d) shows each channel of Fig.3(a). Fig.3 (e) shows the generated JD of . It is important to note that the curl of the image has three components(channels) coming from different axes in 3D. Fig.3(f) represents the RGB image of generated CV and each channel of it is shown in Fig.3(g)(h)(i). In this work, we use generated JD and CV as channels of the training model to train again based on the VoxelMorph CNN architecture and we get better performance than train only once in MIT’s paper [21] which we describe in the following section.
3.2 Average Transformations
We use average a set of transformations[16] by averaging their JD and CV. Given a series of transformations for , we construct their average in the following steps:
(1)Let and as in Eq.(1).
(2)Use Eq. (1) to calculate a transformation such that , and .
(3)We define this as the average of for .
The average deformation has certain geometric meaning: the local size change ratios modeled by
and the local rotations modeled by are averaged to average transformation. So, we apply average transformations to medical image registration. This is also the reason why we combine average of JD and CV with Fixed image as 5 channels in the above framework to predict better registration filed . Because it synthesizes all the information of which are obtained in the first training stage and it can improve the accuracy of medical image registration.3.3 Construction of Atlas(Template)
Here, we want to construct an unbiased template using average of transformations according to Fig.4.
Step1: for our dataset with N subjects, we take one subject out as the initial template, then register to all N images for to get all deformation filed with Voxelmorph CNN.
Step2: find the average transformation = of all the deformation filed by the Average of Transformations method described in above. This step will get N average transformations for .
Step3: Warp on the average transformation to get the temporary templates = for . This step will get N temporary templates .
Step 4: repeat Step1 to Step3 on biased temporary templates to reduce their bias. So we get a new group of average transformation and new temporary templates .And we will do above steps until the average transformation is close to the identify map Id (unit orthogonal grid which corresponds to unit orthogonal map).
4 Experiments
4.1 Datasets, Preprocessing and Evaluation Criteria
4.1.1 Datasets and Preprocessing
We validate our method on two datasets including ADNI[8] and MRBrainS18. ADNI We use one dataset from ADNI including 199 T1weighted brain MRI scans. And we split the dataset into 159,20 and 20(8:1:1) volumes for train, validation, and test sets respectively.
MRBrainS18 We apply the average transformation method on this dataset. Seven brain MRI scans (size:240x240x48) of MRBrainS18 with manual segmentations are provided for training and no data are provided publicly for testing. The data can be downloaded from http://mrbrains18.isi.uu.nl/.
Before the CNN training stage, standard preprocessing steps for structural brain MRI include the following steps, including skull stripping, resampling and affine spatial normalization using FreeSurfer[6]. Segmentation maps including 29 anatomical structures obtained using FreeSurfer for evaluation.
4.1.2 Evaluation Criteria
In this paper, Dice Score(DS) is the most used criteria for volume overlap of anatomical segmentations. We expect the regions in and according to the same anatomical structure to overlap well. Let , be the voxels for and F respectively.DS is computed by:
(3) 
Dice score is between 0 and 1, and Dice score is closer to 1 indicates that the structures are more identical. And vice, a score of 0 indicates that there is no overlap. We also evaluate the regularity of the registration field . The Jacobian matrix captures the local properties of .We compute the number of all nonbackground voxels for which , where is not diffeomorphic.
4.1.3 Implementation
Our implementation uses Keras[18] with a Tesorflow backend[12] and the Adam optimizer[14] with a learning rate of
. We used MATLAB to generate the images formed by JD and CV information based on brain MRI, and saved it as nii image format. We set the epochs as 1500, batchsize as 1, steps of per epoch as 100 using one GeForce GTX 1080 Ti GPU. Each training batch consists of one pair of volumes to reduce the memory usage. We choose the model that optimizes Dice on validation set and get results on test set. Algorithm runtime were computed for a GeForce GTX 1080 Ti GPU and an Intel Core i76800k CPU.
5 Results
5.1 Experiments on Test Set
We present the following experiments demonstrating the results of our proposed method compared with stateofthearts methods. Specifically, we register each scan to an atlas computed using external data [6] which is also used on MIT’s paper. Table 1 gives a summary of the results on test sets. According to Fig.2, we get JD and CV from two different , one is directly generated using deformation method on the gradient and gray value of input moving volume described in section 2.1 and another is the registration filed in the first training stage, which we refer to two methods as VoxelMorphdeformation and VoxelMorphregistration. We compare them with stateofthearts methods the original VoxelMorph and VoxelMorphdiff from MIT’s work [15][4]. The first method VoxelMorphdeformation get the best average Dice with 0.764. The latter three method have comparable runtimes but all faster and yield uncertainty estimates than original VoxelMorph. However, our second method, VoxelMorphregistration produces better diffeomorphic registration fields(having lower numbers of nonnegative Jacobian locations) than other methods.
Method  Avg. Dice  GPU sec  CPU sec  Percentage() of  Uncertainty 

VoxelMorph(CC)[15] 
0.752(0.142)  0.53(0.01)  85.1(1.21)  0.562(0.112)  No 
VoxelMorphdiff[4] 
0.753(0.138)  0.45(0.01)  56(0.13)  0.511(0.182)  Yes 
VoxelMorphdeformation 
0.764(0.125)  0.46(0.02)  57(0.12)  0.504(0.186)  Yes 
VoxelMorphregistration  0.758(0.134)  0.48(0.14)  59(0.24)  0.492(0.235)  Yes 
VoxelMorphatlas 
0.756(0.141)  0.41(0.23)  58(0.22)  0.487(0.194)  Yes 

Fig.5 shows representative results. We present the example MR slices (Two anatomical planes) of input moving image, atlas, and moved images for different experiments, with overlaid boundaries of ventricles (yellow). From the red box, according to atlas, our latter two methods get clearer details of moved images and better registration performance than VoxelMorphdiff. Fig.6 shows the registration fild for different experiments (RGB image and warped grid), the method VoxelMorphregistration get smoother diffeomorphic registration field than else two method.
Based on FreeSurfer segmentation, we used three methods already trained to test the anatomical scans from the unseen anatomical structures of ADNI. This dataset contains expert manual delineations of anatomical structures, which is used for our evaluation. The following experiments shows the results of anatomical segmentation for different experiments. We can see VoxelMorphdeformation method performs better segmentation results than other methods from the different color regions.
5.2 Analysis of Velocity Sampling and Uncertainty
In this section, we evaluate generated velocity sampling and uncertainty for different experiments. From the bottom row, we can see the VoxelMorphregistration method learns smaller values for approximation and yields smoother velocity sample than other two methods, namely the diagonal covariance has the less potential to add noise to
, so, the loss function coupled with the squaring and scaling layers lead to smoother and more accurate registration field
as Fig.6 shows. And this also testify the result that it can generate better diffeomorphic registration fields which shown in Table 1.5.3 Experiments on Atlas
We test the average transformation method using MRBrainS18 dataset to construct a template(atlas) for 7 samples and use it as fixed images for medical registration based on above original VoxelMorphdiff network. We do average transformation method for two iterations, the two average registration field are shown in Fig.9 (a)(b) respectively. Because the average transformation is close to the identify map Id, so we stop the iteration step and use one of the Template as fixed image for original VoxelMorphdiff network. From the Table 1, we refer this experiment as VoxelMorphatlas, compared with VoxelMorphdiff, it has comparable average Dice scores, runtimes and yield uncertainty estimates but produces better diffeomorphic registration fields(having lower numbers of nonnegative Jacobian locations, for VoxelMorphatlas and for VoxelMorphdiff). Moving, fixed image(atlas), fixed image and generated registration field are shown in Fig.9(c)(d)(e)(f) respectively.
6 Conclusions
In conclusion, we propose a new deformable medical image registration method based on differential geometry and VoxelMorph CNN architecture. We compute the differential geometric information including the Jacobian determinant (JD) and the curl vector(CV) of diffeomorphic registration field and use them as channels of the model to train again. We get JD and CV from two sources, generated diffeomorphisms directly by original MR images and the registration filed in the first training stage, we test our method on two datasets including ADNI dataset and the results show our method VoxelMorphdeformation has better Dice and VoxelMorphregistration has better diffeomorphic registration fields than MIT’s result. However, our method needs amounts of iterations and mathematical computation when generating JD and CV.
We also build a registration template for 7 testing samples according to average transformation method based on MRBrainS18 Challenge dataset, and obtain excellent improvement on medical image registration with nonnegative Jacobian locations compared with MIT’s VoxelMorphdiff. However, it also needs some iteration and computing resource which depends on the volume numbers of dataset. It may be useful to construct an enormous atlas dataset for medical registration. We believe the proposed method can advance the performance in medical image registration and clinical diagnosis.
References
 [1] Bajcsy R. & Kovacic, S. Multiresolution elastic matching. Computer Vision, Graphics, and Image Processing.46,1–21(1989).
 [2] Zhu Y.P., Zhou Z.C., Liao G.J., Yang Q.X., Yuan K.H. Effects of Differential Geometry Parameters on Grid Generation and Segmentation of MRI Brain Image. IEEE Access.7(1),68529–68539(2019)
 [3] Dalca A.V., Balakrishnan G., Guttag J., Sabuncu M.R. Unsupervised Learning for Fast Probabilistic Diffeomorphic Registration. arXiv preprint arXiv: 1805.04605, 2018.
 [4] Dalca A.V., Balakrishnan G., Guttag J., Sabuncu M.R. Unsupervised Learning of Probabilistic Di?eomorphic Registration for Images and Surfaces. arXiv preprint arXiv: 1903.03545, 2019.
 [5] Adrian V.D., Andreea B., Natalia S.R., Polina G. Patchbased discrete registration of clinical brain images.In MICCAIPATCHMI Patchbased Techniques in Medical Imaging.Springer(2016).
 [6] Fischl B. Freesurfer. Neuroimage.62(2),774–781(2012)
 [7] Daniel R., Luke I.S., Carmel H., Derek L.G.H., Martin O.L., David J. H. Nonrigid registration using freeform deformation: Application to breast mr images.IEEE Transactions on Medical Imaging.18(8),712–721(1999).
 [8] Susanne G.M., Michael W.W., Leon J.T., Ronald C.P., Cli?ord R.J., William J., John Q.T., Arthur W.T., Laurel B. Ways toward an early diagnosis in Alzheimers disease: the Alzheimers Disease Neuroimaging Initiative (ADNI).Alzheimer?s & Dementia.1(1),55–66(2005).
 [9] Thirion J.P.Image matching as a di?usion process: an analogy with maxwell?s demons.Medical Image Analysis.2(3),243?260(1998).
 [10] Lei N., Luo Z.X., Yau S.T., Gu X.F. Geometric Understanding of Deep Learning. arXiv preprint arXiv: 1805.10451,2018.

[11]
Balakrishnan G., Zhao A., Sabuncu M. R., Guttag J., Dalca A. V.
An unsupervised learning model for deformable medical image registration.In
Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition
. pp. 9252–9260(2018)  [12] Martin A., Ashish A., Paul B., Eugene Brevdo, Zhifeng Chen, Craig Citro, Greg S Corrado, Andy Davis, Je?rey Dean, Matthieu Devin, et al. Tensor?ow: Largescale machine learning on heterogeneous distributed systems. arXiv preprint arXiv:1603.04467, 2016.
 [13] Chen, X., Liao, G. New variational method of grid generation with prescribed jacobian determinant and prescribed curl. Computer Science.2–6(2015)
 [14] Diederik P.K., Jimmy B. ADAM: A method for stochastic optimization. arXiv preprint arXiv:1412.6980, 2014.
 [15] Balakrishnan G., Zhao A., Sabuncu M.R., Guttag J., Dalca A.V. VoxelMorph: A Learning Framework for Deformable Medical Image Registration. IEEE Transactions on Medical Imaging.99,1–1(2019)
 [16] Zhou Z.C., Hildebrandt B., Chen X., Liao GJ.Computational Technologies for Brain Morphometry.arXiv preprint arXiv: 1810.04833, 2018.
 [17] Ashburner J., Friston K.Voxelbased morphometrythe methods.Neuroimage.11,805–821(2000).
 [18] Francois C.Keras. https://github.com/fchollet/keras, 2015
Comments
There are no comments yet.