Introduction to constrained registration trough diffeomorphic metric mapping.
Anton FRANÇOIS
Plan
Introduction
Shapes: Many Modalities
mesh
MRImages
spatial transcriptomic
Diffusion Tensor Imaging (DTI)
Diffeomorphic image registration:
Finding a smooth one-to-one (non-linear) deformation to be biologically plausible (no holes, shearing, tearing...)
Image registration
J. Sassen, 2024
Recent examples
Recent examples
J. Sassen, 2024
Recent examples
J. Sassen, 2024
Recent examples
Shape Spaces
Point clouds
- \(\Omega\) : space and \(x \in \Omega\)
- \( \{c_k\}_{k \in \llbracket 1, N\rrbracket } \) : control points
- \( \{\alpha_k\}_{k \in \llbracket 1, N\rrbracket } \) : momenta or vector weights, \(\forall k, \alpha_k \in \mathbb R ^d\)
- \(K \) : Radial basis function
- \(\Phi(x) = x + v(x)\) : Deformation
\[v(x) = \sum_{k=1}^N K(x, c_k)\alpha_k \]
\(c_1\)
\(\alpha_1\)
Deformation parametrized trough vector fields
Norm on vectors fields
Given two vectors fields
\[v_1 = \sum_i K(\quad , c_i)\alpha_i;\qquad v_2 = \sum_j K(\quad, c_j') \beta_j \]
we define the scalar product
\[\langle v_1, v_2 \rangle_V = \sum_i \sum_j \alpha_i^\top K(c_i, c'_j) \beta_j \]
and the norm \(V\):
\[\|v\|_V^2 = \langle v, v \rangle_V ^2 \]
Let \(V\) be the Reproducing Kernel Hilbert Space (RKHS) of vector fields whose kernel is \(K\). We denote \(L: V \rightarrow V^*\) the Riesz operator such that \[\|v\|^2_V = ( Lv,v).\]
V is an admissible RKHS of vector fields.
Ensuring the diffeomorphism
To this end, we make the control points \(c_k\) and weights \(\alpha_k\) to depend on a “time” \(t \in [0,1]\) that plays the role of a variable of integration.
\[ v(t,x) = \sum_{k=1}^{N} K(x, c_k(t)) \alpha_k(t) \]
a particule \(x\) follows the integral curve according to the PDEs
\[\dot x(t) = \frac{dx_t}{dt} = v(t,x(t)), \qquad x(0) = x_0\]
\[\dot c(t) = K(c(t), c(t)) \alpha (t), \qquad c(0) = c_0 \]
the transformation is uniquely caracterized by the set of initials condition \(\{x_0, c_0, \alpha_0\}\).
Varifold Norm
\[\langle S, S' \rangle_{W^*} = \sum_i \sum_j K(c_i, c_j') \frac{(n_i ^\top n_j')^2}{|n_i| |n'_j|}\]
\[d_W(S, S') = \|S - S'\|_{W^*}^2 = \langle S, S \rangle_{W^*} + \langle S', S' \rangle_{W^*} - 2 \langle S, S' \rangle_{W^*}\]
Varifold Norm
\[\langle S, S' \rangle_{W^*} = \sum_i \sum_j K(c_i, c_j') \frac{(n_i ^\top n_j')^2}{|n_i| |n'_j|}\]
\[d_W(S, S') = \|S - S'\|_{W^*}^2 = \langle S, S \rangle_{W^*} + \langle S', S' \rangle_{W^*} - 2 \langle S, S' \rangle_{W^*}\]
\(c_i\)
\(n_i\)
\(n_i\)
\(c_i\)
\[\langle S, S' \rangle_{W^*} = \sum_i \sum_j K(c_i, c_j') \frac{(n_i ^\top n_j')^2}{|n_i| |n'_j|}\]
\[d_W(S, S') = \|S - S'\|_{W^*}^2 = \langle S, S \rangle_{W^*} + \langle S', S' \rangle_{W^*} - 2 \langle S, S' \rangle_{W^*}\]
Varifold Norm
\[\langle S, S' \rangle_{W^*} = \sum_i \sum_j K(c_i, c_j') \frac{(n_i ^\top n_j')^2}{|n_i| |n'_j|}\]
\[d_W(S, S') = \|S - S'\|_{W^*}^2 = \langle S, S \rangle_{W^*} + \langle S', S' \rangle_{W^*} - 2 \langle S, S' \rangle_{W^*}\]
Varifold Norm
Text
Time Series - LDDMM
Time Series - LDDMM
Graph : \(G(s) := \{(t,s(t)) : t \in I\}\)
Solved the problem by choosing a well crafted kernel:
\[ K_G((t,x),(t',x')) = \left(\begin{matrix} c_0 K_{\mathrm{time}} & 0 \\ 0 & c_1K_{\mathrm{space}} \end{matrix}\right)\]
Time Series - LDDMM
Images
Anton François - 23/05/2023
Lagrangian Formulation
$$x'(t) = v(t,x)$$
Lagrangian scheme
Not suitable for images
Eulerian scheme
\(I_{t+\delta_t} = I_t - \delta_t(v_t \times\nabla I_t)\)
Eulerian Formulation
$$\partial_t I(t,x) = - v(t,x) \cdot \nabla I(t,x)$$
Integration of the geodesics equations
Anton François - 23/05/2023
\(G\) is a group of diffeomorphisms
The deformation defined from the ordinary differential equation:
\[\partial_t \varphi_t = v_t \circ \varphi_t; \quad \varphi_0 = \mathrm{Id}\]
with \(v_t \in V,\forall t\in [0,1]\), is a diffeomorphism. We note \(\varphi^v \in G\) such a diffeomorphism.
\(G\) := Space of Deformations
Deformation as a flow of vector field
Anton François - 23/05/2023
\[\partial_t \varphi_t = v_t \circ \varphi_t; \quad \varphi_0 = \mathrm{Id}\]
Image transport (advection): \[\partial_t I_t = v_t \cdot \nabla I_t \]
Deformed image
\(I_t = I_0 \circ (\varphi^v)^{-1}\)
The deformation acts on images.
Geodesic := Shortest path for the exact matching:
\[\mathrm{inf}_v \int_0^1 \|v_t\|_V^2 dt\]
Anton François - 23/05/2023
LDDMM can reach only images with same topology
Metamorphosis deforms images and adds intensity.
Making the registration of images B and C possible.
\[\partial_t I_t = v_t \cdot \nabla I_t\]
\[+ \mu z_t\]
Topology and appearances variations
Semi-Lagrangian scheme principle
\(\partial_t I_t = v_t \cdot \nabla I_t + \mu z_t\)
Semi-Lagrangian scheme
\(I_{t+\delta_t} = \mathrm{Interp}(I_t,\varphi^{v_t}) + \delta_t \mu z_t\)
Metamorphosis
$$\left\{\begin{array}{rl}v_t &= -\frac{\rho}{\mu} K_\sigma \star (z_t \nabla I_t)\\ \dot z_t &= -\quad \nabla \cdot (z_t v_t) \\ \dot I_t &= -\langle v_t , \nabla I_t\rangle + \mu z_t\end{array}\right.$$
Advection equation with sourceContinuity equation
By computing the Euler-Lagrange equation of the exact matching cost :
and doing the variation with respect to \(I\) and \(v\) we obtain this set of geodesic equations:
$$ E_\mathrm M (I,v) = \frac12\int_0^1 \|v_t\|_V^2 + \rho \|z_t\|_{L_2}^2 dt \quad s.t.\ I_0 = S, I_1 = T$$
Metamorphosis
\[\left\{\begin{array}{rl} v &= - \sqrt{ \rho } K_{V} (p \nabla I)\\ \dot{p} &= -\sqrt{ \rho } \nabla \cdot (pv) \\ z &= \sqrt{ 1 - \rho } p \\ \dot{I} &= - \sqrt{ \rho } v_{t} \cdot\nabla I_{t} + \sqrt{ 1-\rho } z.\end{array}\right. \]
Advection equation with sourceContinuity equation
Given the image evolution model
\[\dot I_t = - \sqrt{\rho}v_t \cdot \nabla I_t + \sqrt{1 - \rho} z_t\]
and the Hamiltonian :
$$H(I,p,v,z) = - (p |\dot{ I}) - \frac{1}{2} \|v\|^2_{V} - \frac{1}{2}\|z\|^2_{2}. $$
we get the optimal trajectory as a geodesic of the form:
Constrained Metamorphosis
Anton François - 23/05/2023
Registration with topological changes
Healthy brain
Brain with a glioblastoma
Anton François - 23/05/2023
On real data: BraTS 2021
- Big data-set
- 4 modalities
- Glioblastoma ground truth segmentation
Registration validation using manually segmented ventricles (40 subjects).
Registered ventricles overlap measured with DICE score.
Constrained Metamorphosis
Let \((M_t)_{t\in [0,1]}\) and \((P_t)_{t\in [0,1]}\) be two continuous temporal masks and \(w_t \in V\), an admissible vector field.
By computing the Euler-Lagrange equation of the exact matching cost:
$$ E_{\mathrm{CM}}(I,v) = \int_0^1 \|v_t\|_V^2 + \rho \langle z_t, M_t z_t \rangle_{L^2} +\gamma\| P_t(v_t - w_t) \|_V^2 dt $$
and doing the variation with respect to \(I\) and \(v\) we obtain this set of geodesic equations:
Oriented & Weighted
\(s.t., I_0 = S, I_1 = T\)
Theorem:
$$\left\{\begin{array}{rl} v_t &= - \frac{\rho}{\mu (1 + \gamma P_t)} K \star (z_t\nabla I_t) + \frac{\gamma P_t }{1 + \gamma P_t} w_t\\ \dot z_t &= -\quad \nabla \cdot (z_t v_t) \\ \dot I_t &= - v_t , \cdot\nabla I_t + \mu M_t z_t\end{array}\right.$$
\(\frac{\rho}{\mu(1+\gamma P_t)}\)
\(\frac{\gamma P_t}{1+\gamma P_t}\)
Constrained Metamorphosis
Results on ToyExamples
LDDMM
CM
WM
\(I_t\)
\(I_t\) & \(\varphi^{v_t}\)
\(z_t\)
Constrained Metamorphosis
Results on ToyExamples
Constrained Metamorphosis
Weighted Metamorphosis
Qualitative results
Anton François - 23/05/2023
Quantitative Results:
Lower the better
Higher the better
Modulated Metamorphosis
The standard \(n\)-Simplex is the subset of \(\mathbb R^{n+1}\) given by
$$\Delta^n = \left\{p \in \mathbb R^{n+1} : \forall i \in \{1,\cdots,n+1\}, p_i > 0; \sum_{i=1}^{n+1} p_i = 1\right\}.$$
A path on the 2-Simplex fig from [Jasminder, 2020]
Data from :
segmented MRI
Data Nature
Rigid alignment
- Cancerous Growth (red - orange)
- White Matter Intensity reallocations (orange)
- Ventricule expansion (gray)
Registration difficulties
Control allocation matrix
Table with the values of \(\rho\) to transition from one class to another.
Up to now, the evolution of the image was described by:
\[\dot{q}_t = -\sqrt{1- \rho} \, v_t \cdot \nabla q_t + \sqrt{\rho} \, z_t., \rho \in \mathbb R\]
\(R =\)
\((\)
\()\)
Up to now, the evolution of the image was described by:
\[\dot{q}_t = -\sqrt{\rho} \, v_t \cdot \nabla q_t + \sqrt{1 - \rho} \, z_t., \rho \in \mathbb R\]
Let \( R \) be the control allocation matrix (CAM), such that each entry \( R(i,j) \) indicates the value of \( \rho \) for the transition from class \( i \) to class \( j \).
For \( f, g \in \Delta^d \), we define:
\[\rho(x, f, g) = \sum_{i \in \mathcal{C}} \sum_{j \in \mathcal{C}} f^i g^j R(i,j) = \left\langle f, gR \right\rangle.\]
We then introduce a new dynamics:
\[\dot{q}_t = -\sqrt{\tilde{\rho}_t} \, v_t \cdot \nabla q_t + \sqrt{1 - \tilde{\rho}_t} \, z_t,\]
where \( \tilde{\rho}_t = \rho(x, q_t, T) \).
Adapative re-allocation control
Conclusion
deck
By Anton FRANCOIS
