Numerical Frameworks for Elastic Shape Analysis
using Second Order Sobolev Metrics
Emmanuel Hartman^1, Yashil Sukurdeep^2, Emery Pierson^3,
Eric Klassen^1, Mohamed Daoudi^{4,5}, Martin Bauer^1, Nicolas Charon^4
^1Department of Mathematics, Florida State University
^2Center of Imaging Sciences, Johns Hopkins University
^3LIX, Ecole Polytechnique
^4Univ. Lille, CNRS, Centrale Lille, Institut Mines-Télécom, CRIStAL
^5IMT Nord Europe, Institut Mines-Télécom, Univ. Lille, Centre for Digital Systems
^6Department of Mathematics, University of Houston
AMS Joint Mathematics Meeting
6 January 2024
Elastic Shape Analysis of Surfaces
The goal of elastic shape analysis is to define a Riemannian metric on the space of unparameterized surfaces.
Let M be a smooth oriented, Riemannian 2-manifold.
Parameterized Shapes:
\qquad \qquad \qquad \Imm =\operatorname{Imm}(M,\R^3)/\text{translations}
Reparameterization Group:
\qquad \qquad \qquad \Diff = \operatorname{Diff}_+(M)
Unparameterized Shapes:
\qquad \qquad \qquad \Shape = \Imm/\Diff\qquad [q]=\{q\circ\phi|\phi\in\Diff\}
We equip \Imm with a metric G that is invariant under the action of \Diff and equip \Shape with the quotient metric.
To solve geodesic boundary value problems between [q_0],[q_1]\in\Shape, we find a path between q_0 and q_1\circ\phi which minimizes the path energy functional.
Relaxed Matching Frameworks
We utilize a relaxed matching framework which considers an optimization problem of the form:
\inf\limits_{\{q\in C^\infty([0,1],\Imm)\,|\,q(0)=q_0\}} \displaystyle\int_0^1 G_{q(t)}(\partial_t q(t),\partial_t q(t)) dt + \lambda \Gamma([q(1)], [q_1])
where \Gamma is a data loss function on \Shape with \Gamma([p],[q])=0 if p\in[q].
An optimizer, q, minimizes the path energy functional and has q(1)\in [q_1].
Thus, q is a geodesic and orthogonal to the orbits of \Diff.
This approach was implemented in [1] for the SRNF metric.
q_0:
q(0):q(1):
q_1:
^1. Bauer, Charon, Harms, and Hsieh, “A numerical framework for elastic surface matching, comparison, and interpolation”.
Second Order Sobolev Metrics
Let q\in \Imm and h,k\in T_q\Imm. We let g_q be the pullback metric of the Euclidean metric on \R^3. A second order Sobolev metric is then given by
\begin{equation}
\label{eq:H1metric}
G_q(h,k)=\int_M \langle h,k \rangle +g_q^{-1}(dh,dk)+ \langle\Delta_q h,\Delta_q k\rangle\vol_q.
\end{equation}
where
+ dh and dk are viewed as a vector valued one forms,
+ and \Delta_q is the surface Laplacian.
Note
For fixed coordinate we can view g_q, dh, and dk as matrix fields and g_q^{-1}(dh,dk)=\operatorname{tr}(dh\cdot g_q^{-1}\cdot dk^T)
Following the construction of Su et. al.^1 we will further decompose the first order term, in four different terms which each have a geometric interpretation. Therefore, we write \begin{equation} dh= dh_m+dh_++dh_\perp+dh_0. \end{equation}
^1. Su, Bauer, Preston, Laga, and Klassen, “Shape analysis of surfaces using general elastic metrics”.
Splitting Second Order Sobolev Metrics
A straight-forward calculation shows that these terms are orthogonal with respect to the inner product. Thus, we can decompose the first order term of our second order metric producing the following family of second order Sobolev metrics \begin{multline*} G_q(h,k)=\int_M\bigg( a_0 \langle h,k \rangle + a_1 g_q^{-1}(dh_m,dk_m) +b_1g_q^{-1}(dh_+,dk_+)+\\ c_1g_q^{-1}(dh_\bot,dk_\bot)+ d_1 g_q^{-1}(dh_0,dk_0) +a_2 \langle\Delta_q h,\Delta_q k\rangle\bigg)\vol_q. \end{multline*}
Theorem
Assume a_0,a_1,b_1,c_1,d_1,a_2>0. Then, \Imm equipped with the metric G has a non-vanishing geodesic distance and is invariant under the action of \Diff.
Varifold Representations^{1,2,3,4}
For [q]\in \Shape we can associate a varifold \mu_q\in\mathcal{M}(\mathbb R^3\times S^{2}).
In particular, \mu_q:= (q,n_q)_*\operatorname{vol}_q where n_q is the unit normal map of q
We take a norm \|\cdot\|_{V} on \mathcal{M}(\mathbb R^3\times S^{2}) where for [q],[p]\in\Shape the scalar product between the associated varifolds \mu_{q} and \mu_{p} can be written as: \begin{equation*}\label{equ:norm_var} \langle \mu_{q},\mu_{p}\rangle_{V}=\iint_{M \times M}e^{-\alpha||q(x)-p(y)||^2}\langle n_q(x), n_p(y)\rangle^2\vol_{q}(x) \vol_{p}(y). \end{equation*}
Thus we define our choice of relaxation term as,\Gamma([q],[p]):=||\mu_p-\mu_q||^2_{V}=\langle \mu_{p},\mu_{p}\rangle_{V}+\langle \mu_{q},\mu_{q}\rangle_{V}-2\langle \mu_{p},\mu_{q}\rangle_{V}.
^1Charon & Trouvé. "The varifold representation of nonoriented shapes for diffeomorphic registration."
^2Kaltenmark, et al. "A general framework for curve and surface comparison and registration with oriented varifolds."
^3Feydy et al. "Optimal transport for diffeomorphic registration"
^4 Roussillon & Glaunès. "Representation of surfaces with normal cycles and application to surface registration."
Solutions to Geodesic Boundary Value Problem
Thus, our relaxed matching energy between q_0 and q_1 becomes
\inf\limits_{q\in C^\infty([0,1],\Imm)|q(0)=q_0} \displaystyle\int_0^1 G_{q(t)}(\partial_t q(t),\partial_t q(t)) dt
\qquad\qquad\qquad\qquad\qquad+ \lambda ||\mu_{q_1}-\mu_{q(1)}||^2_{V}
with choices of the following hyper-parameters: \alpha,\lambda,a_0,a_1,b_1,c_1,d_1, and a_2.
We discretize this energy for piecewise linear paths via finite differences and triangular meshes via discrete differential geometry.
Surface Matching Results
The Varifold Norm and Partial Matching
What to do when there are parts of q_0 that are not matched to q_1?
 |  | 
|
Consider \omega:M\to [0,1] interpreted as the probability that q_0(x) is matched to q_1 and use a \omega-weighted varifold representation of [q(1)] denoted \mu_{q(1),\omega}.
Thus, our relaxed partial matching problem becomes
\inf\limits_{\{\omega:M\to [0,1]\}}\inf\limits_{\{q\in C^\infty([0,1],\Imm)\,|\,q(0)=q_0\}} \displaystyle\int_0^1 G_{q(t)}(\partial_t q(t),\partial_t q(t)) dt + \lambda ||\mu_{q_1}-\mu_{q(1),\omega}||^2_{V}
Partial Matching Results
Statistical Analysis Pipeline
In addition to the relaxed matching framework, with respect to our family of second order sobolev metrics we implement methods for:
- solving geodesic initial value problems.
- computing Frechet means.
- performing tangent space principle component analysis.
- performing parallel transport via Schild's (or Pole) ladder schemes.
Basis Restricted Elastic Shape Analysis (BaRe-ESA)
One drawback of the previous framework is that geodesics with respect to the chosen metric may not correspond to biologically relavent or physically feasible deformations.
To overcome this we restrict our solutions to the affine subspace of \Imm which can be reached by linearly deforming a template surface, \mathcal{T}, by a basis of allowable deformations, \{h_i\}_{i=1}^n.
Specifically we consider surfaces in the image of the map \begin{align*}F:\mathbb{R}^n \to \Imm \text{ defined via }\\ (w_1,...,w_n)\mapsto \mathcal{T}+\sum_{i=1}^n w_ih_i .\end{align*}
Basis Restricted Elastic Shape Analysis
We perform our optimizations in the lower dimensional space of coeffecients equipped with the pullback of G via F denoted by \overline{G}.
To find the coeffecients mapped via F as close as possible to the orbit of an arbitrary mesh q we optimize
\inf\limits_{\alpha\in C_0^\infty([0,1],\R^n)\,|\, \alpha(0)=0} \displaystyle\int_0^1 \overline{G}_{\alpha(t)}(\partial_t \alpha(t),\partial_t \alpha(t)) dt + \lambda ||\mu_{q}-\mu_{F(\alpha(1))}||^2_{V}.
To solve a geodesic boundary value between [q_0] and [q_1] in our restricted space we relax both enpoints of our path and optimize
\inf\limits_{\alpha\in C^\infty([0,1],\R^n)} \displaystyle\int_0^1 \overline{G}_{\alpha(t)}(\partial_t \alpha(t),\partial_t \alpha(t)) dt
\qquad\qquad\qquad+ \lambda_0 ||\mu_{q_0}-\mu_{F(\alpha(0))}||^2_{V} + \lambda_1 ||\mu_{q_1}-\mu_{F(\alpha(1))}||^2_{V}.
Conclusion
We propose a relaxed matching framework for elastic shape analysis of 3D surfaces equipped with a family of elastic second-order Sobolev metrics. Furthermore we discussed:
+ extending this framework to handle to partial matching tasks.
+ statistical analysis methods under the proposed family of metrics.
+ restricted subspaces of surfaces which can be expressed with respect to a basis of phisically feasible or biologically relevant deformations.
This talk is based on:
"Elastic Shape Analysis of Surfaces with Second-Order Sobolev Metrics", IJCV. ,
"BaRe-ESA: A Riemannian framework for unregistered human body shapes", ICCV. ,
"Basis restricted elastic shape analysis on the space of unregistered surfaces".
The code for these projects is available at:
github.com/emmanuel-hartman/H2_SurfaceMatch and
github.com/emmanuel-hartman/BaRe-ESA
These slides are available at:
www.math.fsu.edu/~ehartman/Slides/JMM2024_Slides/talk.html
Thank you for your attention!