This paper was converted on www.awesomepapers.org from LaTeX by an anonymous user.
Want to know more? Visit the Converter page.

11institutetext: School of Data Science, Fudan University, Shanghai, China 22institutetext: Department of Electrical and Computer Engineering, University of Washington, Seattle, United States

BInGo: Bayesian Intrinsic Groupwise Registration via Explicit Hierarchical Disentanglement ††thanks: X. Wang and X. Luo contributed equally to this work.

Xin Wang 1122    Xinzhe Luo 11    Xiahai Zhuang 11
Abstract

Multimodal groupwise registration aligns internal structures in a group of medical images. Current approaches to this problem involve developing similarity measures over the joint intensity profile of all images, which may be computationally prohibitive for large image groups and unstable under various conditions. To tackle these issues, we propose BInGo, a general unsupervised hierarchical Bayesian framework based on deep learning, to learn intrinsic structural representations to measure the similarity of multimodal images. Particularly, a variational auto-encoder with a novel posterior is proposed, which facilitates the disentanglement learning of structural representations and spatial transformations, and characterizes the imaging process from the common structure with shape transition and appearance variation. Notably, BInGo is scalable to learn from small groups, whereas being tested for large-scale groupwise registration, thus significantly reducing computational costs. We compared BInGo with five iterative or deep-learning methods on three public intrasubject and intersubject datasets, i.e. BraTS, MS-CMR of the heart, and Learn2Reg abdomen MR-CT, and demonstrated its superior accuracy and computational efficiency, even for very large group sizes (e.g., over 1300 2D images from MS-CMR in each group).

1 Introduction

Multimodal groupwise registration aims to align multimodal images into a common structural space. Unlike conventional pairwise registration which aligns moving images separately to a fixed image, groupwise registration can ameliorate the bias from designating a reference image, by estimating a common space to which all the images are co-registered. Therefore, it has become an essential task in multivariate image analysis, including longitudinal research, atlas construction, motion estimation, and population studies [14, 10, 18, 6].

Conventional iterative methods usually estimate the desired spatial transformations by optimizing intensity-based similarity measures [13, 20, 25, 22, 16]. Recently, unsupervised deep-learning methods attempt to realize groupwise registration in an end-to-end fashion [4, 7]. For instance in [4], the authors devised a network to optimize the conditional template entropy (CTE) introduced in [22]. These learning-based models estimate conventional similarity measures developed for iterative registration using stochastic gradient methods, and hence may inherit the same suboptimal performance due to insufficient modeling of the image similarity. Besides, these measures rely on correctly characterizing the joint intensity profile over the entire image group, which may be computationally prohibitive and applicability-limited due to large group sizes.

In this work, however, we shift attention from devising similarity metrics to learning the underlying generative imaging process. Specifically, we establish a probabilistic generative model for the observed images, in which transformations and the common structure are disentangled as latent variables. In this way, groupwise registration can be achieved by unsupervised variational auto-encoding: the encoder extracts structural representations of images and then infers the spatial deformations; the decoder imitates the imaging process by reconstructing original images from estimated common structure and the inverse deformations. The contributions of this work can be summarized as follows:

  1. (1)

    We propose BInGo, a theoretically-grounded unsupervised Bayesian framework for groupwise registration, which learns intrinsic structural similarity of multimodal images by a principled variational distribution, and is capable of disentangling structural representations from image appearance.

  2. (2)

    BInGo is scalable to be trained with small image groups while being applied to large-scale and variable-size test groups, significantly improving applicability and computational efficiency.

  3. (3)

    We demonstrate the superiority of BInGo over similarity-based methods on three multimodal intrasubject and intersubject datasets.

To the best of our knowledge, this is the first work that realizes scalable multimodal groupwise registration by disentangled representation learning.

2 Methodology

Let the original image group be 𝑿=(Xm)m=1M\bm{X}=(X_{m})_{m=1}^{M}, with MM the number of modalities and Xm:ℝd⊃Ωm→ℝX_{m}:\mathbb{R}^{d}\supset\mathrm{\Omega}_{m}\rightarrow\mathbb{R}. Groupwise registration aims to find spatial transformations ϕ=(ϕm:ℝd⊃Ω→Ωm)m=1M\bm{\phi}=(\phi_{m}:\mathbb{R}^{d}\supset\mathrm{\Omega}\rightarrow\mathrm{\Omega}_{m})_{m=1}^{M}, such that the aligned images (Xm​ ​ϕm)m=1M(X_{m}\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m})_{m=1}^{M} share a common structural representation 𝒁\bm{Z}.

We assume XmX_{m} is generated from 𝒁\bm{Z} and ϕm−1\phi_{m}^{-1} by a transformation-equivariant imaging process fmf_{m} that contains modality-specific appearance information, i.e.

Xm=fm​(𝒁​ ​ϕm−1)=fm​(𝒁)​ ​ϕm−1.X_{m}=f_{m}(\bm{Z}\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m}^{-1})=f_{m}(\bm{Z})\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m}^{-1}. (1)

Therefore, a 3-step unsupervised auto-encoding scheme can be formed: I) By modeling fm−1f_{m}^{-1} we could encode individual structural representation for XmX_{m} as 𝒁​ ​ϕm−1=fm−1​(Xm)\bm{Z}\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m}^{-1}=f_{m}^{-1}(X_{m}), from which ϕ\bm{\phi} could be estimated more easily. II) 𝒁\bm{Z} could be then inferred from 𝒁=fm−1​(Xm​ ​ϕm)\bm{Z}=f_{m}^{-1}(X_{m}\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m}) following the equivariance assumption. III) By modeling fmf_{m} we could reconstruct XmX_{m} from the estimated 𝒁\bm{Z} and ϕm−1\phi_{m}^{-1}.

The following establish BInGo to realize this scheme by variational and disentangled auto-encoding, thus learning ϕ\bm{\phi} in a unified and end-to-end fashion.

2.1 Hierarchical Bayesian Inference

Refer to caption
Figure 1: The proposed Bayesian framework. (a) Graphical model of the imaging process. (b) Inference model for the velocity fields 𝒗\bm{v} (orange), the common structural representation 𝒛\bm{z} (green), and the reconstruction (blue). Random variables are in circles, deterministic variables are in double circles, and observed variables are shaded.

We first formalize the estimation of spatial deformations ϕ=(ϕm)m=1M\bm{\phi}=(\phi_{m})_{m=1}^{M} through Bayesian inference, as illustrated in Fig. 1. Specifically, let 𝒙=(xm)m=1M\bm{x}=(x_{m})_{m=1}^{M} be a sample of 𝑿\bm{X}. We decompose the latent variables generating 𝒙\bm{x} into two independent subgroups: 1) the common structural representation 𝒛\bm{z}, and 2) the stationary velocity fields 𝒗=(𝒗m)m=1M\bm{v}=(\bm{v}_{m})_{m=1}^{M} that parameterize ϕ\bm{\phi} [1]. Thus, following the variational Bayes framework, the objective function to maximize is the evidence lower bound (ELBO) of the log-likelihood, which takes the form

ℒ​(𝒙)≜\displaystyle\mathcal{L}(\bm{x})\triangleq 𝔼q​(𝒛,𝒗|𝒙)​[log⁡p​(𝒙|𝒛,𝒗)]−DKL​[q​(𝒛,𝒗|𝒙)∥p​(𝒛)​p​(𝒗)]\displaystyle\ \mathbb{E}_{q(\bm{z},\bm{v}\,|\,\bm{x})}\big{[}\log p(\bm{x}\,|\,\bm{z},\bm{v})\big{]}-D_{\text{KL}}\big{[}q(\bm{z},\bm{v}\,|\,\bm{x})\;\|\;p(\bm{z})\,p(\bm{v})\big{]} (2)

where q​(𝒛,𝒗|𝒙)=q​(𝒗|𝒙)​q​(𝒛|𝒗,𝒙)q(\bm{z},\bm{v}\,|\,\bm{x})=q(\bm{v}\,|\,\bm{x})\,q(\bm{z}\,|\,\bm{v},\bm{x}) defines the variational posterior distribution, and DKLD_{\text{KL}} is the Kullback-Leibler (KL) divergence. After optimizing the ELBO, the desired velocity fields 𝒗^\widehat{\bm{v}} can be inferred via maximum a posteriori.

Furthermore, we express the latent variables with LL hierarchical levels [24], i.e. 𝒛=(𝒛l)l=1L\bm{z}=(\bm{z}^{l})_{l=1}^{L} and 𝒗m=(𝒗ml)l=1L\bm{v}_{m}=(\bm{v}_{m}^{l})_{l=1}^{L}, and a higher level indicates a finer-scale (larger) resolution. Thus, the total deformation ϕm=ϕm1​ ​⋯​ ​ϕmL\phi_{m}=\phi_{m}^{1}\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\cdots\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m}^{L} can be deterministically computed by exponentiation of velocity fields ϕml=id+exp⁡(𝒗ml){\phi}_{m}^{l}=\mathrm{id}+\exp({\bm{v}}_{m}^{l}), where ∂∂t​ϕml​(𝝎,t)=𝒗ml​(ϕml​(𝝎,t))\frac{\partial}{\partial t}\phi_{m}^{l}(\bm{\omega},t)=\bm{v}_{m}^{l}(\phi_{m}^{l}(\bm{\omega},t)), ∀𝝎∈Ω\forall\,\bm{\omega}\in\mathrm{\Omega}. The hierarchical strategy allows to model a complex deformation by several simpler and easier-to-learn ones.

To simplify the KL, we introduce additional independence assumptions: 1) both the prior and posterior of the velocities factorize among different images, i.e. p​(𝒗l|𝒗<l)=p​(𝒗l)=∏m=1Mp​(𝒗ml)p(\bm{v}^{l}\,|\,\bm{v}^{<l})=p(\bm{v}^{l})=\prod_{m=1}^{M}p(\bm{v}_{m}^{l}) and q​(𝒗l|𝒙,𝒗<l)=∏m=1Mq​(𝒗ml|𝒙,𝒗<l)q(\bm{v}^{l}\,|\,\bm{x},\bm{v}^{<l})=\prod_{m=1}^{M}q(\bm{v}_{m}^{l}\,|\,\bm{x},\bm{v}^{<l}), where <l<l denotes the group of latent variables in levels less than ll, and 2) the common structure 𝒛l\bm{z}^{l} can be inferred directly from the observed images and the estimated velocity fields, i.e. q​(𝒛l|𝒙,𝒗,𝒛<l)=q​(𝒛l|𝒙,𝒗)q(\bm{z}^{l}\,|\,\bm{x},\bm{v},\bm{z}^{<l})=q(\bm{z}^{l}\,|\,\bm{x},\bm{v}), since {𝒙,𝒗}\{\bm{x},\bm{v}\} can determine registered images, thus containing all information about 𝒛l\bm{z}^{l} for any ll.

Taking into account these assumptions, the KL divergence can be written as

DKL​[q​(𝒗|𝒙)​q​(𝒛|𝒗,𝒙)∥p​(𝒗)​p​(𝒛)]\displaystyle D_{\text{KL}}\big{[}q(\bm{v}\,|\,\bm{x})\,q(\bm{z}\,|\,\bm{v},\bm{x})\;\|\;p(\bm{v})\,p(\bm{z})\big{]} (3)
=\displaystyle= ∑m=1M[∑l=1L𝔼q​(𝒗m<l|𝒙)​{DKL​[q​(𝒗ml|𝒙,𝒗<l)∥p​(𝒗ml)]}](i)\displaystyle\sum_{m=1}^{M}\left[\sum_{l=1}^{L}\mathbb{E}_{q(\bm{v}_{m}^{<l}\,|\,\bm{x})}\Big{\{}D_{\text{KL}}\big{[}q(\bm{v}_{m}^{l}\,|\,\bm{x},\bm{v}^{<l})\;\|\;p(\bm{v}_{m}^{l})\big{]}\Big{\}}\right]\qquad\qquad\text{(i)}\
+𝔼q​(𝒗|𝒙)[∑l=1L𝔼q​(𝒛<l|𝒙,𝒗){DKL[q(𝒛l|𝒙,𝒗)∥p(𝒛l|𝒛<l)]}](ii),\displaystyle+\mathbb{E}_{q(\bm{v}\,|\,\bm{x})}\left[\sum_{l=1}^{L}\mathbb{E}_{q(\bm{z}^{<l}\,|\,\bm{x},\bm{v})}\Big{\{}D_{\text{KL}}\big{[}q(\bm{z}^{l}\,|\,\bm{x},\bm{v})\;\|\;p(\bm{z}^{l}\,|\,\bm{z}^{<l})\big{]}\Big{\}}\right]\ \text{(ii)},

where we prescribe p​(𝒛1|𝒛<1)≜p​(𝒛1)p(\bm{z}^{1}\,|\,\bm{z}^{<1})\triangleq p(\bm{z}^{1}) and q​(𝒗m<1|𝒙)=q​(𝒛<1|𝒙,𝒗)≜1q(\bm{v}_{m}^{<1}\,|\,\bm{x})=q(\bm{z}^{<1}\,|\,\bm{x},\bm{v})\triangleq 1 for simplicity. The key idea here is: the overall KL divergence is decomposed w.r.t. (i) the velocity fields 𝒗\bm{v}, and (ii) the common structure 𝒛\bm{z}. The former serves as regularization to ensure diffeomorphism, fulfilled by the constraint introduced in [5]; the latter is intended to estimate the structural similarity among the images, as is detailed in the next subsection.

2.2 Intrinsic Similarity over Structural Representations

For registration, it is crucial to efficiently measure the structural similarity, which is achieved in this work by learning instead of a pre-defined metric. To this end, we propose to learn multilevel “expert” distributions qm​(𝒛l|xm,𝒗m)q_{m}(\bm{z}^{l}\,|\,x_{m},\bm{v}_{m}), which serve as fm−1f_{m}^{-1} (i.e., the inverse of imaging processes) to extract the structural representations of warped images. We further assume q​(𝒛l|𝒙,𝒗)q(\bm{z}^{l}\,|\,\bm{x},\bm{v}) and p​(𝒛l|𝒛<l)p(\bm{z}^{l}\,|\,\bm{z}^{<l}) to be the geometric mean [15] and arithmetic mean [23] of the experts, respectively:

q​(𝒛l|𝒙,𝒗)∝[∏m=1Mqm​(𝒛l|xm,𝒗m)]1M,p​(𝒛l|𝒛<l)≜1M​∑m=1Mqm​(𝒛l|xm,𝒗m).q(\bm{z}^{l}|\,\bm{x},\bm{v})\propto\left[\prod_{m=1}^{M}q_{m}(\bm{z}^{l}|\,x_{m},\bm{v}_{m})\right]^{\frac{1}{M}},\ p(\bm{z}^{l}|\bm{z}^{<l})\triangleq\frac{1}{M}\sum_{m=1}^{M}q_{m}(\bm{z}^{l}|\,x_{m},\bm{v}_{m}). (4)

Therefore, the KL in Eq. 3(ii) essentially measures the intrinsic (dis)similarity, i.e., the dissimilarity of the experts (intrinsic structural representations), whose minimization, as a part of the maximization of the ELBO, encourages the experts to be identical, thus forcing the multilevel posteriors q​(𝒛l|𝒙,𝒗)q(\bm{z}^{l}\,|\,\bm{x},\bm{v}) to represent the common structure. Meanwhile, the velocity could be learned in tandem.

We model the experts qm​(𝒛l|xm,𝒗m)q_{m}(\bm{z}^{l}\,|\,x_{m},\bm{v}_{m}) as Gaussians 𝒩​(𝝁ml,𝚺ml)\mathcal{N}(\bm{\mu}_{m}^{l},\bm{\Sigma}_{m}^{l}), and thus so is the joint posterior, i.e. q​(𝒛l|𝒙,𝒗)=𝒩​(𝝁l,𝚺l)q(\bm{z}^{l}\,|\,\bm{x},\bm{v})=\mathcal{N}(\bm{\mu}^{l},\bm{\Sigma}^{l}) with

𝚺l=M⋅[∑m=1M(𝚺ml)−1]−1,𝝁l=𝚺lM⋅∑m=1M𝝁ml​(𝚺ml)−1.\bm{\Sigma}^{l}=M\cdot\left[\sum_{m=1}^{M}\left(\bm{\Sigma}_{m}^{l}\right)^{-1}\right]^{-1},\quad\bm{\mu}^{l}=\frac{\bm{\Sigma}^{l}}{M}\cdot\sum_{m=1}^{M}\bm{\mu}_{m}^{l}\left(\bm{\Sigma}_{m}^{l}\right)^{-1}. (5)

In light of the computational intractability of the KL divergence involving Gaussian mixture distributions, we further exploit its convexity to obtain

DKL[q(𝒛l|𝒙,𝒗)∥p(𝒛l|𝒛<l)]⩽1M∑m=1MDKL[q(𝒛l|𝒙,𝒗)∥qm(𝒛l|xm,𝒗m)].D_{\text{KL}}\big{[}q(\bm{z}^{l}\,|\,\bm{x},\bm{v})\;\|\;p(\bm{z}^{l}\,|\,\bm{z}^{<l})\big{]}\leqslant\frac{1}{M}\sum_{m=1}^{M}D_{\text{KL}}\big{[}q(\bm{z}^{l}|\,\bm{x},\bm{v})\;\|\;q_{m}(\bm{z}^{l}|\,x_{m},\bm{v}_{m})\big{]}. (6)

Hence, the minimization of the KL in Eq. 3(ii) is actually approximated by minimizing the right-hand side of Eq. 6, which has a closed-form expression.

2.3 Explicit Disentanglement with Neural Networks

The ELBO has been decomposed into three terms: reconstruction of the original images (the first term in Eq. 2), the KL divergence for velocity fields to ensure diffeomorphism (Eq. 3(i)), and the KL divergence for the intrinsic (dis)similarity to optimize the registration (Eq. 3(ii)), with weights (hyperparameters) (λi)i=13(\lambda_{i})_{i=1}^{3}. To estimate the ELBO, we propose a dedicated hierarchical variational auto-encoder (VAE) as the inference model in Fig. 1(b), whose architecture with the number of hierarchy L=3L=3 is depicted in Fig. 2. The network learns p​(𝒙|𝒛,𝒗),q​(𝒗ml|𝒙,𝒗<l)p(\bm{x}\,|\,\bm{z},\bm{v}),q\left(\bm{v}_{m}^{l}\,|\,\bm{x},\bm{v}^{<l}\right) and qm​(𝒛l|xm,𝒗m)q_{m}\left(\bm{z}^{l}\,|\,x_{m},\bm{v}_{m}\right), and the expectations in the ELBO are estimated by Monte Carlo sampling, similar to traditional VAEs.

Refer to caption
Figure 2: Architecture (L=3L=3) of BInGo for explicit disentanglement. The channel numbers are indicated around particular feature maps (arrows).

The motif of BInGo is to explicitly disentangle the common structure, spatial transformations and appearance information from multimodal images, thus realizing Eq. 1. To this end, the network comprises three types of cooperative modules: 1) MM encoders that extract multi-level transformation-equivariant structural representations 𝒄ml\bm{c}_{m}^{l}, corresponding to the inverse imaging function fm−1f_{m}^{-1}, 2) multi-level registration (Reg) modules that infer spatial transformations from individual structural representations, corresponding to estimating ϕ\bm{\phi} from {𝒛​ ​ϕm−1}m=1M\{\bm{z}\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m}^{-1}\}_{m=1}^{M}, and 3) MM decoders, with modality-specific appearance information embedded as parameters, which reconstruct registered images from learned common structure, corresponding to the imaging function fm​(𝒛)f_{m}(\bm{z}).

Particularly, the Reg module at level ll assumes the coarser-scale deformations ϕm<l≜ϕm1​ ​⋯​ ​ϕml−1\phi_{m}^{<l}\triangleq\phi_{m}^{1}\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\cdots\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m}^{l-1} have been resolved, and thus takes the individual representations 𝒄ml​(xm)​ ​ϕm<l\bm{c}_{m}^{l}(x_{m})\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m}^{<l} as inputs, which are like “partially distorted” version of the expert distributions, since the experts are representations of fully registered images xm​ ​ϕmx_{m}\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m}. Thus, in the same spirit of Eq. 4, the Reg module computes their geometric mean as the “partially common” structure, and then concatenates the mean and 𝒄ml​(xm)​ ​ϕm<l\bm{c}_{m}^{l}(x_{m})\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m}^{<l} to infer q​(𝒗ml|𝒙,𝒗<l)q(\bm{v}_{m}^{l}\,|\,\bm{x},\bm{v}^{<l}) for each mm separately.

Therefore, BInGo can realize the aforementioned 3-step unsupervised scheme induced by Eq. 1, as follows:

I. Bottom-up Inference of Velocity Fields. The encoders first produce multilevel feature maps 𝒄ml​(xm)\bm{c}_{m}^{l}(x_{m}) as the individual structural representations of original images. Then, up from the bottom (ll increasing from 1 to LL), given that ϕm<l\phi_{m}^{<l} have been inferred (ϕm<1≜id\phi_{m}^{<1}\triangleq\mathrm{id}), the feature maps are warped to obtain 𝒄ml​(xm)​ ​ϕm<l\bm{c}_{m}^{l}(x_{m})\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m}^{<l}. Then, the Reg module at level ll can infer q​(𝒗ml|𝒙,𝒗<l)q(\bm{v}_{m}^{l}\,|\,\bm{x},\bm{v}^{<l}), from which 𝒗ml\bm{v}_{m}^{l} is sampled and ϕml\phi_{m}^{l} is finally computed through integration.

II. Top-down Inference of Common Structures. As the total deformations ϕm=ϕm1​ ​⋯​ ​ϕmL\phi_{m}=\phi_{m}^{1}\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\cdots\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m}^{L} have been inferred, the encoders are fed again with warped images 𝒙​ ​ϕ≜(xm​ ​ϕm)m=1M\bm{x}\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\bm{\phi}\triangleq(x_{m}\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m})_{m=1}^{M} to produce the experts qm​(𝒛l|xm,𝒗m)=𝒄ml​(xm∘ϕm)q_{m}(\bm{z}^{l}\,|\,x_{m},\bm{v}_{m})=\bm{c}_{m}^{l}(x_{m}\circ\phi_{m}) in a top-down (ll decreasing from LL to 11) manner. The posterior q​(𝒛l|𝒙,𝒗)q\left(\bm{z}^{l}\,|\,\bm{x},\bm{v}\right) is then computed by Eq. 4, from which the modality-invariant common structural representation 𝒛=(𝒛l)l=1L\bm{z}=(\bm{z}^{l})_{l=1}^{L} is sampled. Note that we can avoid using encoders twice by warping 𝒄ml​(xm)\bm{c}_{m}^{l}(x_{m}) to compute the experts as 𝒄ml​(xm)∘ϕm\bm{c}_{m}^{l}(x_{m})\circ\phi_{m} (thanks to equivariance), but this requires more computational cost.

III. Disentangled Auto-Encoding. Based on the common structure 𝒛\bm{z} inferred by encoders, the decoders reconstruct the registered images 𝒙^=(x^m)m=1M\widehat{\bm{x}}=(\widehat{x}_{m})_{m=1}^{M}. Then reconstructed original images 𝒙^​ ​ϕ−1≜(x^m​ ​ϕm−1)m=1M\widehat{\bm{x}}\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\bm{\phi}^{-1}\triangleq(\widehat{x}_{m}\,\raisebox{1.0pt}{ \leavevmode\hbox to3pt{\vbox to3pt{\pgfpicture\makeatletter\hbox{\hskip 1.5pt\lower-1.5pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{ }\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{ }\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{ }{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{ }\pgfsys@setlinewidth{0.6pt}\pgfsys@invoke{ }{}\pgfsys@moveto{1.2pt}{0.0pt}\pgfsys@curveto{1.2pt}{0.66273pt}{0.66273pt}{1.2pt}{0.0pt}{1.2pt}\pgfsys@curveto{-0.66273pt}{1.2pt}{-1.2pt}{0.66273pt}{-1.2pt}{0.0pt}\pgfsys@curveto{-1.2pt}{-0.66273pt}{-0.66273pt}{-1.2pt}{0.0pt}{-1.2pt}\pgfsys@curveto{0.66273pt}{-1.2pt}{1.2pt}{-0.66273pt}{1.2pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@stroke\pgfsys@invoke{ } \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}}\,\phi_{m}^{-1})_{m=1}^{M} are obtained by inverse deformations. We utilize (negative) L1 loss between the original images and their reconstructions to estimate the first (reconstruction) term in Eq. 2. In addition, to better disentangle structure and appearance, the convolutional layers of the network are shared across modalities, while domain-specific batch normalization layers (BNs) encode appearance information for each modality [3].

2.4 Towards Large-Scale and Variable-Size Groupwise Registration

Most conventional methods optimize similarity measures defined over the entire image group, making the computation for large groups formidable. Besides, learning-based models that rely on these measures require groups to have the same size. These drawbacks significantly reduce their applicability.

However, BInGo is scalable to handle multimodal image groups with various sizes for training and test. As shown in Fig. 3, it can be trained with either bimodal image pairs or image groups of complete modalities, referred to as partial or complete learning, respectively. Trained in either way, BInGo can be flexibly applied to larger unseen test groups with arbitrary numbers of images, such that computational efficiency of training is substantially boosted.

Refer to caption
Figure 3: Partial and complete learning schemes of our model for large-scale and variable-size groupwise registration.

The scalability of our model is conducive in two scenarios:

Intrasubject Images with Missing Modalities. The goal is to co-register multimodal scans for each patient. In practice, there could be many patients with different missing modalities. BInGo can be trained with these heterogeneous datasets via partial learning, while performing groupwise registration for test subjects with complete modalities.

Intersubject Populations. Intersubject groupwise registration, which is crucial in atlas construction and population analysis, could involve plenty of images for every modality. Conventional methods to this problem rely on iterative optimization and suffer from computational complexity scaling with image group size, while BInGo can greatly relieve training burden and be subsequently applied to test groups of complete populations in one shot.

3 Experiments and Results

3.1 Datasets and Preprocessing

BraTS-2021. The dataset provides 3D pre-operative T1, T1Gd, T2 and T2-FLAIR MR scans of patients with glioblastoma [17, 2]. We randomly selected 300/50/150 patient sets for training/validation/test. The volumes were downsampled into 2×2×2​mm32\times 2\times 2\text{mm}^{3} with ROI of size 80×96×8080\times 96\times 80. As the images of each patient are pre-registered, we use synthetic free-form deformations (FFDs) with different control point spacings to simulate misalignments.

MS-CMR. The MS-CMRSeg challenge [26] provides cardiac MR sequences LGE, bSSFP, and T2 from 45 patients, which exhibit complementary information of the cardiac structure and pathology. The images were preprocessed by affine co-registration, ROI cropping, and slice selection, producing 39/15/44 slices for training/validation/test. We simulated worse misalignments by applying additional FFDs on original images to better demonstrate model efficacy.

Learn2Reg Abdomen MR-CT. This dataset [11, 8] collects 3D MR and CT volumes. The images were resampled into 3×3×3​mm33\times 3\times 3\text{mm}^{3} with size of 128×107×128128\times 107\times 128.

3.2 Experimental Setups

Compared Methods. Three types of unsupervised methods for multimodal groupwise registration were compared on the intrasubject BraTS and MS-CMR datasets: 1) similarity-based iterative methods using information-theoretic metrics CTE, APE or 𝒳\mathcal{X}-CoReg [22, 25, 16], 2) similarity-based deep-learning models that optimize CTE or APE using an attention residual U-Net (AttResUNet) [9, 19] as the backbone, 3) BInGo, the proposed model learning intrinsic similarity through hierarchical disentanglement. For intersubject population groupwise registration on Learn2Reg, the models in 2) are not applicable since they can only be trained with complete image groups, whereas BInGo can be trained partially. In addition, for BraTS and MS-CMR, all baseline methods used single-level velocity fields as the transformation model; for Learn2Reg, iterative baselines performed rigid, affine and FFD registration successively for optimal accuracy, due to the severe initial misalignments.

Implementation Details. Training was done through the Adam optimizer [12] (learning rate: 10−310^{-3}; batch size: 20/120/1 for MS-CMR/other datasets). The experiments were conducted with PyTorch [21] on an NVIDIA RTXTM\text{RTX}^{\textnormal{TM}} 3090 GPU.

Evaluation Metrics. We reported the Dice similarity coefficient (DSC) averaged over all pairwise combinations of the segmentation masks for registered images. For the BraTS dataset, as the ground-truth misalignments were available, we also reported the groupwise warping index (gWI) [16], which would reduce to zero if the misaligned images were perfectly co-registered.

3.3 Results

Table 1: Results on intrasubject image groups from BraTS and MS-CMR. The mean values and standard deviations are presented for the gWI (in voxels) and DSC, with top 2 bolded. The number of parameters (in millions) of each model for each dataset are also reported. Statistically significant improvement (p<0.05p<0.05 for one-sided paired tt-tests) of BInGo with complete† or partial∗ learning was marked with daggers or asterisks, respectively.
Method BraTS MS-CMR
DSC ↑\uparrow gWI ↓\downarrow #Params. DSC ↑\uparrow #Params.
None 0.610±0.150†∗0.610\pm 0.150^{\dagger*} 1.430±0.644†∗1.430\pm 0.644^{\dagger*} — 0.722±0.101†∗0.722\pm 0.101^{\dagger*} —
\hdashline APE [25] 0.726±0.078\bm{0.726\pm 0.078} 0.596±0.149\bm{0.596\pm 0.149} 7.373 0.811±0.072†∗0.811\pm 0.072^{\dagger*} 0.154
CTE [22] 0.561±0.148†∗0.561\pm 0.148^{\dagger*} 1.087±0.411†∗1.087\pm 0.411^{\dagger*} 7.373 0.816±0.077†∗0.816\pm 0.077^{\dagger*} 0.154
𝒳\mathcal{X}-CoReg [16] 0.707±0.089†0.707\pm 0.089^{\dagger} 0.697±0.212†0.697\pm 0.212^{\dagger} 7.373 0.840±0.077†∗0.840\pm 0.077^{\dagger*} 0.154
\hdashline APE+AttResUNet     0.693±0.078†0.693\pm 0.078^{\dagger} 0.757±0.153†∗0.757\pm 0.153^{\dagger*} 22.955 0.846±0.048†∗0.846\pm 0.048^{\dagger*} 8.036
CTE+AttResUNet 0.659±0.096†∗0.659\pm 0.096^{\dagger*} 0.916±0.210†∗0.916\pm 0.210^{\dagger*} 22.955 0.874±0.043†∗0.874\pm 0.043^{\dagger*} 8.036
\hdashline BInGo (complete†) 0.717±0.068\bm{0.717\pm 0.068} 0.596±0.132\bm{0.596\pm 0.132} 13.429 0.887±0.033\bm{0.887\pm 0.033} 4.516
BInGo (partial∗) 0.693±0.0750.693\pm 0.075 0.709±0.1720.709\pm 0.172 13.429 0.877±0.042\bm{0.877\pm 0.042} 4.516
Table 2: Results on intersubject population groups from Learn2Reg. The means and standard deviations of DSCs versus the numbers of images in each group are reported.
Method #Images 2 4 8 16
None 0.396±0.1680.396\pm 0.168 0.386±0.0700.386\pm 0.070 0.319±0.0070.319\pm 0.007 0.306±0.0000.306\pm 0.000
\hdashline APE [25] 0.586±0.3760.586\pm 0.376 0.574±0.2610.574\pm 0.261 N/A N/A
CTE [22] 0.609±0.3060.609\pm 0.306 0.088±0.0360.088\pm 0.036 N/A N/A
𝒳\mathcal{X}-CoReg [16] 0.675±0.3290.675\pm 0.329 0.567±0.2240.567\pm 0.224 N/A N/A
\hdashline BInGo (partial∗) 0.781±0.108\bm{0.781\pm 0.108} 0.715±0.122\bm{0.715\pm 0.122} 0.677±0.059\bm{0.677\pm 0.059} 0.645±0.000\bm{0.645\pm 0.000}
Refer to caption
Figure 4: Results (mean values with one standard deviation bands) before (blue) and after (orange) registration by trained BInGo for images groups with different sizes.
Refer to caption
Figure 5: Example image groups (with segmentation contours overlaid) before or after registration using the best baselines or BInGo. Deformations and example features (extracted from registered images) output by BInGo are also presented. Green and orange arrows indicate some areas of significant improvement (compared to the best baselines) and relatively poor performance, respectively.

Intrasubject Image Groups. Section 3.3 presents registration accuracy on the BraTS and MS-CMR datasets. BInGo (with complete† or partial∗ learning) worked consistently better than similarity-based iterative or learning approaches, except for the APE-based iterative method on BraTS. Notably, BInGo could achieve superior performance to deep-learning baselines with only half of training parameters, even though their backbone AttResUNet is more advanced than the vanilla UNet-like structure used in BInGo. Furthermore, the proposed partial learning strategy could perform better than most baselines with significantly lower computational demands, which potentiates large-scale end-to-end groupwise registration.

Intersubject Population Image Groups. For Learn2Reg, we merged every KK test MR-CT pairs to form larger groups (each containing M′=2​KM^{\prime}=2K images to co-register) as test sets of intersubject populations. Section 3.3 presents the results versus different M′M^{\prime}. All iterative approaches failed when the groups became large due to excessive GPU cost, while our partial learning strategy worked consistently better and maintained a good performance for large groups. Note that the deep-learning baselines are not capable of handling test groups with varying sizes, and thus not presented.

Additional Scalability Tests. For MS-CMR and BraTS, we merged test groups in a similar way to evaluate BInGo (trained via complete learning) on image groups with different sizes. As shown in Fig. 4, with M′M^{\prime} increasing, co-registration becomes significantly more difficult (indicated by the worse initial DSC/gWI), whereas BInGo maintained decent performance for ultra-large groups (e.g., over 1300 images from MS-CMR), and achieved even better accuracy for larger 3D image groups from BraTS, showing remarkable robustness on large-scale groupwise registration.

Qualitative Results. We visualized the results from the best baselines (APE/ CTE+AttResUNet/𝒳\mathcal{X}-CoReg for BraTS/MS-CMR/Learn2Reg) and BInGo in Fig. 5. BInGo could achieve better alignment for both large-scale anatomy and local fine structures in most cases, and the predicted deformations reached great smoothness and diffeomorphism. The feature maps from BInGo were nearly modality-invariant and shared similar structures, illustrating successful disentanglement of the common structure, spatial transformations and appearance information. Besides, the relatively poor performance in certain local areas may be due to the too severe initial misalignment, which made different tissues with similar intensities appear in the same location. Still, BInGo performed no worse than the baselines in such regions, and achieved satisfactory accuracy for the entire foreground.

4 Conclusion and Discussion

In this work we have presented BInGo, a generative Bayesian framework for unsupervised multimodal groupwise registration. This new formulation of image registration has achieved comparable performance with similarity-based iterative and unsupervised methods. In particular, we demonstrated that equipped with unique scalability, BInGo could reduce the computational burden of groupwise registration without compromising accuracy. This opens up the possibility to realize learning-based multimodal groupwise registration on a large scale and with various group sizes. A potential limitation of our work is that there may exist performance drop on specific datasets, e.g. the highly challenging abdominal images. Future work includes investigation into the negative factors that may inhibit BInGo from generalizing to large image groups.

References

  • [1] Ashburner, J.: A fast diffeomorphic image registration algorithm. Neuroimage 38(1), 95–113 (2007)
  • [2] Baid, U., Ghodasara, S., Mohan, S., Bilello, M., Calabrese, E., Colak, E., Farahani, K., Kalpathy-Cramer, J., Kitamura, F.C., Pati, S., et al.: The rsna-asnr-miccai brats 2021 benchmark on brain tumor segmentation and radiogenomic classification. arXiv preprint arXiv:2107.02314 (2021)
  • [3] Cang, W.G., You, T., Seo, S., Kwak, S., Han, B.: Domain-specific batch normalization for unsupervised domain adaptation. In: Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition. pp. 7354–7362 (2019)
  • [4] Che, T., Zheng, Y., Sui, X., Jiang, Y., Cong, J., Jiao, W., Zhao, B.: Dgr-net: Deep groupwise registration of multispectral images. In: International Conference on Information Processing in Medical Imaging. pp. 706–717. Springer (2019)
  • [5] Dalca, A.V., Balakrishnan, G., Guttag, J., Sabuncu, M.R.: Unsupervised learning of probabilistic diffeomorphic registration for images and surfaces. Medical image analysis 57, 226–236 (2019)
  • [6] Geng, X., Christensen, G.E., Gu, H., Ross, T.J., Yang, Y.: Implicit reference-based group-wise image registration and its application to structural and functional MRI. Neuroimage 47(4), 1341–1351 (Oct 2009)
  • [7] He, Z., Chung, A.C.: Unsupervised end-to-end groupwise registration framework without generating templates. In: 2020 IEEE International Conference on Image Processing (ICIP). pp. 375–379. IEEE (2020)
  • [8] Hering, A., Hansen, L., Mok, T.C.W., Chung, A.C.S., Siebert, H., Häger, S., Lange, A., Kuckertz, S., Heldmann, S., Shao, W., et al.: Learn2reg: comprehensive multi-task medical image registration challenge, dataset and evaluation in the era of deep learning. IEEE Transactions on Medical Imaging (2022)
  • [9] Hu, Y., Modat, M., Gibson, E., Li, W., Ghavami, N., Bonmati, E., Wang, G., Bandula, S., Moore, C.M., Emberton, M., et al.: Weakly-supervised convolutional neural networks for multimodal image registration. Medical image analysis 49, 1–13 (2018)
  • [10] Joshi, S.C., Davis, B.C., Jomier, M., Gerig, G.: Unbiased diffeomorphic atlas construction for computational anatomy. NeuroImage 23, S151–S160 (2004)
  • [11] Kavur, A.E., Selver, M.A., Dicle, O., Barış, M., Gezer, N.S.: CHAOS - Combined (CT-MR) Healthy Abdominal Organ Segmentation Challenge Data (Apr 2019). https://doi.org/10.5281/zenodo.3362844
  • [12] Kingma, D.P., Ba, J.: Adam: A method for stochastic optimization. In: 3rd International Conference on Learning Representations (2015)
  • [13] Learned-Miller, E.G.: Data driven image models through continuous joint alignment. IEEE Transactions on Pattern Analysis and Machine Intelligence 28(2), 236–250 (2005)
  • [14] Liao, S., Jia, H., Wu, G., Shen, D.: A novel framework for longitudinal atlas construction with groupwise registration of subject image sequences. Neuroimage 59(2), 1275–1289 (Jan 2012)
  • [15] Lorenzen, P., Prastawa, M., Davis, B., Gerig, G., Bullitt, E., Joshi, S.: Multi-modal image set registration and atlas formation. Medical image analysis 10(3), 440–451 (2006)
  • [16] Luo, X., Zhuang, X.: X-metric: An n-dimensional information-theoretic framework for groupwise registration and deep combined computing. IEEE Transactions on Pattern Analysis and Machine Intelligence (2022)
  • [17] Menze, B.H., Jakab, A., Bauer, S., Kalpathy-Cramer, J., Farahani, K., Kirby, J., Burren, Y., Porz, N., Slotboom, J., Wiest, R., et al.: The multimodal brain tumor image segmentation benchmark (brats). IEEE transactions on medical imaging 34(10), 1993–2024 (2014)
  • [18] Metz, C.T., Klein, S., Schaap, M., van Walsum, T., Niessen, W.J.: Nonrigid registration of dynamic medical imaging data using nd+ t b-splines and a groupwise optimization approach. Medical image analysis 15(2), 238–249 (2011)
  • [19] Oktay, O., Schlemper, J., Folgoc, L.L., Lee, M., Heinrich, M., Misawa, K., Mori, K., McDonagh, S., Hammerla, N.Y., Kainz, B., Glocker, B., Rueckert, D.: Attention u-net: Learning where to look for the pancreas. In: Medical Imaging with Deep Learning (2018)
  • [20] Orchard, J., Mann, R.: Registering a multisensor ensemble of images. IEEE Transactions on Image Processing 19(5), 1236–1247 (2009)
  • [21] Paszke, A., Gross, S., Chintala, S., Chanan, G., Yang, E., DeVito, Z., Lin, Z., Desmaison, A., Antiga, L., Lerer, A.: Automatic differentiation in pytorch (2017)
  • [22] Polfliet, M., Klein, S., Huizinga, W., Paulides, M.M., Niessen, W.J., Vandemeulebroucke, J.: Intrasubject multimodal groupwise registration with the conditional template entropy. Medical image analysis 46, 15–25 (2018)
  • [23] Shi, Y., Paige, B., Torr, P., et al.: Variational mixture-of-experts autoencoders for multi-modal deep generative models. Advances in Neural Information Processing Systems 32, 15718–15729 (2019)
  • [24] Vahdat, A., Kautz, J.: Nvae: A deep hierarchical variational autoencoder. arXiv preprint arXiv:2007.03898 (2020)
  • [25] Wachinger, C., Navab, N.: Simultaneous registration of multiple images: Similarity metrics and efficient optimization. IEEE transactions on pattern analysis and machine intelligence 35(5), 1221–1233 (2012)
  • [26] Zhuang, X., Xu, J., Luo, X., Chen, C., Ouyang, C., Rueckert, D., Campello, V.M., Lekadir, K., Vesal, S., RaviKumar, N., et al.: Cardiac segmentation on late gadolinium enhancement mri: a benchmark study from multi-sequence cardiac mr segmentation challenge. Medical Image Analysis 81, 102528 (2022)