非刚性3D形状到图像匹配的组合解决方案

τ 0i, τ 0jτi, τj10000非刚性3D形状到图像匹配的组合解决方案0Florian Bernard 1,2,3 Frank R. Schmidt 3 Johan Thunberg 1 Daniel Cremers 301卢森堡系统生物医学中心，卢森堡大学，卢森堡2卢森堡中央医院，卢森堡3慕尼黑工业大学（TUM），德国0（a）模板（3D表面网格）0（d）匹配的模板0（b）数据项（每个三角形）（c）平滑项（相邻三角形）0每个三角形的刚体变换03D图像0（低成本）τ0i∈SE(3)0Eij(τi,τj) Ei(τi)0三角形邻域0低成本0高成本0三角形Fiτi∈SE(3)0τi(Fi) τ0i(Fi)0E(τ) =0i = 1 Ei(τi) + 0（i,j）∈EEij(τi,τj)0图1. 3D形状到图像的匹配。 (a)给定一个三角网格形状，我们为每个三角形Fi关联一个刚性变换τi∈SE(3)，以使数据项和平滑项的总和最小化。(b)数据项Ei(τi)衡量变换后的三角形τi(Fi)与体积图像的匹配程度。(c) 平滑项Eij(τi,τj)惩罚变换后的三角形τi(Fi)和τj(Fj)之间的差异。(d)优化E给出了形状到图像的匹配。0摘要0我们提出了一个组合解决方案来解决非刚性匹配3D形状与3D图像数据的问题。为此，我们将形状建模为一个三角网格，并允许该网格的每个三角形进行刚性变换，以实现与图像的合适匹配。通过对相邻三角形之间的距离和相对旋转进行惩罚，我们的匹配在图像和形状信息之间进行权衡。在本文中，我们解决了两个主要挑战：首先，我们用适当的图论方法解决了由此产生的大规模和NP困难的组合问题。其次，我们提出了对无界的6维Lie群SE(3)进行高效离散化的方法。据我们所知，这是非刚性3D形状到图像匹配的第一个组合公式。与现有的局部（梯度下降）优化方法相比，我们得到的解决方案不需要良好的初始化，并且在最优解的范围内。我们在非刚性3D形状到形状和非刚性3D形状到图像配准的两个问题上评估了所提出的方法，并证明它提供了有希望的结果。01. 引言0在计算机视觉和图像分析中，将形状模板与图像匹配是一个经过深入研究的问题。它产生了许多应用，包括图像分割和物体检测。图像中线条和参数化曲线的检测的早期方法是基于投票的Hough变换[14]，后来被推广到任意形状的检测[1]。虽然Hough变换考虑了刚性形状，但在图像分割任务中利用形状信息也已经在非刚性情况下得到了解决，包括基于主动形状模型[10]、水平集[12]、凸形状空间[13]、多相图割[58]或统计形状模型[19,63]的方法。对于从单个2D图像重建物体形状，基于模板的形状方法旨在通过3D到2D投影将给定的3D模板与图像匹配[50, 38,46]。一些作者已经考虑了非刚性形状到图像匹配问题的组合形式。1.2. Main ContributionsThe main contribution of this paper is to present for thefirst time a combinatorial formulation of the non-rigid 3Dshape to 3D image matching problem. Whilst our problemis a natural extension to the afore-mentioned “dimensionone” matching approaches [11, 16, 53, 30], a generalisationto (intrinsic) dimension two problems is more intricate. Ourmain contributions are:• By using a surface mesh transformation model thatmakes use of per-triangle rigid transformations, we10010某些形状类别的匹配问题已经得到了解决。对于匹配轮廓[11,53]或2D弦图多边形[16]，可以全局解决优化问题。然而，将这些方法推广到3D形状是非平凡的，我们不知道有关此问题的先前工作。本文的目的是通过提出一个组合公式来填补这个空白，用于将3D形状模板非刚性匹配到3D图像。为此，我们将形状建模为一个三角网格，并允许该网格的每个三角形通过刚性变换τi∈SE(3)进行独立变换。使用对无界的6维Lie群SE(3)进行离散化，我们将匹配任务形式化为一个流形值多标记问题，可以将其作为能量最小化问题来解决。0E(τ) =0i = 1 Ei(τi) +0(i,j) ∈ E Eij(τi,τj). (1)0在这里，数据项Ei(τi)考虑了图像信息，而平滑项Eij(τi,τj)衡量了观察到的形状与建模形状先验之间的差异。通过惩罚相邻三角形之间的距离和相对旋转，我们的匹配在图像和形状信息之间做出了妥协。一般来说，最小化（1）中形式的函数是NP难的[26]。01.1. 相关工作0据我们所知，本文是第一篇考虑非刚性三维形状到三维图像匹配问题的组合形式的论文。接下来，我们将总结与我们的工作最相关的方法。连续优化：在许多情况下，我们自然地假设图像或形状的变形是空间连续和平滑的。通常，这类问题是在微分同胚空间上的优化问题来表述的[15, 40, 2,41]。通常，使用梯度下降等方法来获得（局部）问题的（通常是非凸的）最优解。然而，这些方法的一个主要缺点是良好的初始估计至关重要，而且通常对解的最优性没有界限。为了处理在以最优传输为基础的二维形状到图像匹配问题中的非凸性，[52]中的作者提出使用分支定界方案。最短路径和动态规划：与连续局部优化方法相反，许多视觉问题可以以离散的方式表述，从而可以使用基于图算法和动态规划（DP）的解决方案[17]。由于曲线本质上是一维的，各种曲线匹配形式也可以简化为在特定图中寻找最短路径。此外，基于递归公式的方法可以将曲线匹配问题简化为图中的最短路径问题。0通过使用易于解决的子问题进行匹配，具有树结构模板的匹配问题通常可以通过DP来解决。对于将开放轮廓变形匹配到二维图像，[11]提出了一种基于DP的全局解决方案。在[16]中，作者基于DP提出了一种解决弦图多边形变形匹配到二维图像问题的方法。在[53]中，作者提出了一种基于轮廓和图像的乘积图中的循环的全局最优匹配方法。最近在[30]中引入了一种基于乘积图的相关形式，用于可变轮廓到三维形状匹配。图割：众所周知（参见[4]），图的任何割都可以解释为在环境空间中找到一个共维数为1的闭流形（例如，二维中的闭曲线，三维中的闭曲面等）。一个例子是从一组稀疏的三维点重建三维形状，其中后者在离散的三维网格上表示[35]。标签问题：标签问题在计算机视觉中无处不在，既出现在连续设置中，也出现在离散设置中[62]。流行的马尔可夫随机场（MRF）框架提供了贝叶斯处理[37]。此外，还研究了MRF的线性规划松弛[59]。多标签的连续方法包括各种凸松弛[47, 32, 55,18]，具有在流形上取值的函数的总变差正则化的多标签问题[33]，以及子标签精确的凸松弛[42,31]。在离散多标签方法中，前面提到的图割可以用于找到某些二元标签问题的全局解，包括具有次模二元成本的问题[25]。对于一类多标签问题，也可以找到全局解[25]。这个子类包括在完全有序标签方面是凸的成本[22]。此外，还提出了一些高效的算法来找到一般多标签问题的局部最优解[6,27]，甚至具有理论最优性保证。关于可以使用图割进行优化的能量函数的更详细描述，请参见[26, 24, 25]。formulate the 3D shape to 3D image matching problemin terms of a manifold-valued multi-labelling problem.• We introduce a pairwise term that defines a metricon the label space SE(3), which itself is a high-dimensional Lie group. With that, our energy func-tion is amenable to be minimised by the α-expansionalgorithm [6], which has been shown to work well inpractice, is efficient even for very large label spaces,and has theoretical optimality guarantees.• In contrast to continuous optimisation methods thatuse gradient descent-like algorithms, our combinato-rial method does not require a good initialisation.• In order to deal with the computationally challengingdiscretisation of SE(3), we propose to use a coarse-to-fine discretisation of the Lie group.2. Non-Rigid 3D Shape-to-Image MatchingIn this section we first specify our objective, followed bya description of the data term and smoothness term. Afterintroducing the combinatorial problem, we describe the dis-cretisation of the label space and we discuss the algorithmicsolution of the problem.2.1. Objective FunctionIn the following, we assume that a 3D shape S ⊂ R3is given as a triangular mesh. This means we have n ∈ Ntriangles F1, . . . , Fn ⊂ R3 such thatS =n�i=1Fi.(2)We use the set E ⊂ {1, . . . , n}2 to define the neighbour-hood between pairs of (different) triangles. We assume thatfor all (i, j) ∈ E the neighbouring triangles Fi and Fj arenon-disjoint and that the intersection Fi ∩ Fj results eitherin a common edge or a common vertex. Also, w.l.o.g. weassume that for each (i, j) ∈ E it holds that i < j, i.e.(i, j) ∈ E ⇒ (j, i) /∈ E.Our objective is it now to match the 3D shape S ontoa volumetric image I : Ω → Rc, where Ω ⊂ R3 denotesthe compact image domain and c ∈ N describes the amountof image channels. While we are interested in a non-rigidshape-to-image matching, we like to favour matchings thatare as-rigid-as-possible, similar to the approach in [54] thatapplies (locally regularised) rigid transformations to eachvertex. However, in our case this is done by applying toeach triangle Fi a rigid transformationτi = (˜τi,⃗τi) ∈ SE(3) = SO(3) ⋉ R3,(3)where ˜τi ∈ SO(3) ⊂ R3×3 represents the rotational partand ⃗τi ∈ R3 represents the translational part of τi. Thetask of finding the best matching τ = (τ1, . . . , τn) can beformulated as minimising the energyE(τ) =n�i=1Ei(τi) +�(i,j)∈EEij(τi, τj).(4)In Section 2.2 we define the data term Ei(τi) that evalu-ates how well the transformed triangle τi(Fi) fits to the im-age data. In Section 2.3 we define the smoothness termEij(τi, τj) that measures the geometric dissimilarity be-tween the shape model S and the transformed shapeτ(S) :=n�i=1τi(Fi).(5)Using the proposed piecewise rigid transformation modelwe may end up with a model τ(S) having (small) gaps orintersections between neighbouring triangles. We will lateraddress this issue and present a simple yet effective way ofdealing with this irregularity.J[F] :=�FJ(x)dx.(6)Ei(τi) = −J[τi(Fi)].(7)dSO(3)(˜τi, ˜τj) =�12log(˜τ Ti ˜τj) F(8)100202.2. 数据项0数据项  E i ( τ i )  ∈  R  衡量变形三角形  τ i ( F i )与图像数据  I  的拟合程度。为此，我们引入得分图像  J  :Ω  →  [0 ,  1] ，该图像源自图像  I（例如梯度幅度图像，或基于神经网络的更高级的预测器）。对于三角形  F  �  Ω ，我们定义0这样，值  J [ F ]  表示三角形  F与图像数据的拟合程度，其中得分图像的高值表示良好的拟合。数据项可以表示为0在离散设置中，数据项  E i  通过对三角形上的函数值  − J (x )  进行加权求和来计算。权重考虑了变形三角形 τ i ( F i )在图像中的光栅化。02.3. 光滑项0配对项  E ij ( τ i , τ j )  ∈  R + 0 惩罚相邻三角形  F i  和  Fj  在经过  τ i  和  τ j变换后的不一致性。为了定义配对项，我们首先引入适当的距离。表达式0是旋转  ˜ τ i  和  ˜ τ j 在 SO(3) 上的测地距离 [ 21 ]，其中矩阵对数  log( ∙ ) 。对于  q idSE(3),X(τi, τj) = maxx∈X ∥τi(x)−τj(x)∥2 ,(9)Eij(τi, τj) =(10)λBdSO(3)(˜τi, ˜τj) + λSdSE(3),Fi∩Fj(τi, τj).minτ∈SE(3)nE(τ).(11)(i) dSE(3),X(τi, τj) ≥ 0, dSE(3),X(τi, τi) = 0,10030其中  q i  和  q j  是  ˜ τ i  和  ˜ τ j的四元数表示，可以有效地计算距离，即  d SO(3) = 2 cos− 1 ( |� q i , q j �| ) ，其中  �∙ ,  ∙�  是四元数内积 [ 21]。为了定义相邻三角形之间的距离，我们利用了群作用的概念。更具体地说，我们定义0其中群  SE(3)  作用于非空紧致集合 X  �  R 3。在我们的情况下，我们使用  X  =  F i  ∩  F j ，因此 dSE(3) ,F i ∩ F j ( τ i , τ j )  可以看作是变形三角形  τ i ( F i )和  τ j ( F j )  之间的距离。在这种情况下，  d SE(3) ,F i ∩ Fj ( τ i , τ j )  的最大值在  F i  和  F j的公共顶点处达到，从计算的角度来看是有吸引力的。利用引入的距离，我们将配对项定义为其加权和，即0弯曲项的目的是通过权重  λ B  >  0确保相邻三角形的旋转相似。拉伸项的权重为  λ S  >  0，确保相邻三角形保持紧密在一起。02.4. 组合公式0将形状  S  与图像  I  的匹配  τ  给出了优化问题的解0由于可行集  SE(3) n  的非凸性，问题（ 11）是非凸的。这种非凸性使得直接在无界连续空间  SE(3) n上求解问题变得困难。我们的方法是在搜索空间的离散化上进行优化。通过这样做，我们得到了一个多标签问题，可以使用高效且有效的算法来解决。为了离散化  SE(3)，我们利用了它是  SO(3)  和  R 3的乘积空间的事实。因此，我们定义 L �  SO(3) �R 3 = SE(3)为包含  ℓ  =  |L|  个  SE(3)元素的（有限）流形值标签空间。平移：由于平移部分由  R3  编码，因此李群  SE(3)  是非紧致的。然而，由于图像域Ω是紧致的，图像大小为离散化平移提供了自然的界限。设  nx , n y , n z  是图像  I  的体素数，ℓ x , ℓ y , ℓ z  是  x ,  y ,  z平移的标签数。为方便起见，我们假设我们的模板相对于图像域的中心定义，即模板的重心与0图像的中心。此外，不失一般性，我们假设我们正在寻找这样一个匹配，使得（变换后的）模板的大部分位于图像内部1。让我们定义Zm(n) = {−n02 }为包含m个均匀间隔元素的集合，中心为0，其中m是一个奇数正整数。直径n定义为最大元素和最小元素之间的差。平移的离散化由集合�L = Zℓx(nx) × Zℓy(ny) ×Zℓz(nz)给出，其中|�L| = ℓx ∙ ℓy ∙ℓz。旋转：之前已经提出了与SO(3)离散化相关的各种工作。这些工作包括刚体路径规划的采样策略[29]，基于向量距离的SO(3)邻域的近似[36]，或者对3D旋转的各种度量进行分析[21]。我们对SO(3)的离散化基于Hopf纤维化，它用圆S1和2-球S2来描述SO(3)。这种方法的直觉是将S1和S2的离散化转移到旋转空间。对于详细的描述，我们将感兴趣的读者引用到[61]。让˜L表示包含˜ℓ =|˜L|个均匀采样的SO(3)的集合。通过在标签空间Ln上优化E，我们现在得到组合优化问题如下：0min τ ∈Ln E(τ). (12)03. 算法0为了解决问题（12），我们使用α-扩展[6，026 ,  5]，它贪婪地一次只更新一个标签。我们注意到，α-扩展也有潜在的替代方法（例如，对于非度量的两两项[27,55]或融合移动[34]）。虽然α-扩展要求两两项是度量，但从实际和理论角度来看，它具有吸引力。更具体地说，它是高效的，对初始化具有鲁棒性，获得的局部最优解保证在全局最优解的因子范围内，而且支持平滑项在线计算的高效实现也是可用的[6, 26,5]，这对于我们解决的问题的规模是至关重要的。我们现在证明我们的两两项是度量，因此α-扩展是适用的。0引理1：设X � R3是一个非空紧集，τi，τj，τk ∈SE(3)。拉伸项0d SE(3) ,X ( τ i , τ j ) = max x ∈ X ∥ τ i (x) - τ j (x)∥ 2 ,0是一个拟度量，即满足0（ii）对称性：dSE(3)，X(τi，τj) = dSE(3)，X(τj，τi)，以及01 如果不是这种情况，可以相应地增加图像大小。̸τ (s)i= (˜τ (s)i,⃗τ (s)i) ∈ ˜L (s)i× ⃗L (s)i⊂ SO(3) ⋉ R3, (13)⃗L (s+1)i:=(14)��,˜L[r]p := ˜L[r] ∩ BSO(3)ǫ(Id),(15)˜L (s+1)i:= {τ ˜τ (s)i: τ ∈ ˜L[s+1]p⊆ SO(3)},(16)10040（iii）三角不等式：dSE(3)，X(τi，τk) ≤ dSE(3)，X(τi，τj)+ dSE(3)，X(τj，τk)。0证明：(i)和(ii)直接根据定义得出。三角不等式成立，因为0dSE(3)，X(τi，τk) = max x ∈ X ∥ τi(x) - τk(x) ∥ 20= max x ∈ X ∥ τ i (x) - τ k (x) + τ j (x) - τ j (x) ∥ 20≤ max x ∈ X ∥ τi(x) - τj(x) ∥ 2 + ∥ τj(x) -0≤ max x ∈ X ∥ τ i (x) - τ j (x) ∥ 2 + max x ∈ X ∥ τ j (x) - τ k (x) ∥ 20= dSE(3)，X(τi，τj) + dSE(3)，X(τj，τk)。■0命题1：对于λS，λB > 0和(i, j) ∈E，定义在(10)中的两两项Eij(∙, ∙)是一个度量。0证明：根据第2.1节的假设，对于(i, j) ∈ E，有X := Fi ∩ Fj≠�是紧集。因此，dSE(3)，X(∙，∙)是一个拟度量（引理1）。虽然dSO(3)(∙，∙)在SO(3)上被认为是度量，但在SE(3)上只是一个拟度量。由于Eij(∙，∙)是两个拟度量的正线性组合，Eij(∙，∙)也是一个拟度量。为了证明Eij(∙，∙)是一个度量，我们证明Eij(τi，τj)=0意味着τi=τj。对于Eij(τi，τj)=0，有dSO(3)(˜τi，˜τj)=0，这意味着˜τi=˜τj。此外，由于Eij(τi，τj)=0和˜τi=˜τj，有dSE(3)，X(τi，τj)= max x ∈ X  ∥ (˜τi(x)+�τi) - (˜τj(x)+ �τj) ∥ 2 = ∥ �τi - �τj ∥ 2 =0，这意味着�τi=�τj。因此τi=τj。■03.1. 粗到细的处理0实际上，对于一个相当大的标签数ℓ=˜ℓ∙ℓx∙ℓy∙ℓz，直接求解问题（12）是不可行的。为了解决这个问题，我们提出使用粗到细的策略，在标签空间的不同级别s（近似）解决问题（12）。设s=0表示最粗（初始）级别，s=smax≥0表示最细（最终）级别。一旦在级别s得到解τ(s)，为了在级别s+1运行算法，标签初始化为τ(s)，并相应地更新标签空间。为了计算效率，在粗到细的方法中，每个三角形都有自己的可行标签空间。设L(s)i表示三角形i在级别s的可行标签空间。最初，在基本级别s=0时，每个三角形的标签空间相同，即L(0)i=L(0)。获得L(s+1)i的一般思想是在级别s上考虑对变换τ(s)i的邻域进行（均匀）离散化，其中邻域的半径随着级别的减小而减小。让我们在SO(3)上引入一个邻域：0定义1（SO(3)上的�-球）：以半径为�、中心为˜τ的SO(3)球定义为B SO(3)�(˜ τ)={τ∈SO(3):d SO(3)(τ,˜τ)<�}。0接下来，我们根据标签空间的乘积空间性质描述粗到细的标签空间结构。由于每个标签都可以写成0我们可以独立考虑平移和旋转。注意，˜ L(s)i×�L(s)i�SE(3)不再是一个群。通过强制要求τ(s)i∈L(s+1)i，级别s的解τ(s)也包含在新的标签空间中。因此，能量从级别s增加到s+1。平移：定义� L(0):=�L（参见第2.4节）。为了获得三角形i在级别s+1的平移集合，级别s+1的新平移网格以�τ(s)i为中心。此外，从级别s开始的直径缩小了一半，导致0�τ(s)i+0Zℓx(n(s)x2)×Zℓy(n(s)y2)×Zℓz(n(s)z2)0其中向量集的加法是逐元素的。最初，n(0)x=nx，n(0)y=ny和n(0)z=nz。旋转：设˜L[r]是分辨率为r的（整个）SO(3)的离散化，包含˜ℓ[r]个元素，其中˜ℓ[r]随着r的增加而增加。不要将˜L[r]与˜L(s)i混淆，后者是级别s上第i个三角形的（旋转）标签空间，将在下面定义。根据[61]中的构造，我们获得了r=0,...,4的SO(3)离散化˜L[r]的5个分辨率。在将恒等式包含在˜L[r]中后，元素的数量范围从˜ℓ[0]=577到˜ℓ[4]≈2∙106。让我们定义一个包含p个元素的集合，这些元素来自于˜L[r]，它们“最接近”恒等式，即0其中对于每个分辨率r，使用满足|˜L[r]p|≥p的最小�作为半径。现在，我们定义˜L(0):=˜L[0]，并且0它是由˜τ(s)i与˜L[s+1]p中的所有旋转组合而成的集合。对于所有级别s，我们使用p=˜ℓ[0]=577。对于预定义的r=0,...,4的SO(3)网格化，总是存在一个�，使得上述不等式是紧的，即|˜L[r]p|=˜ℓ[0]=577。通过将恒等式包含在˜L[s+1]p中，我们确保˜τ(s)i∈˜L(s+1)i。由于0∈Z∙(∙)，因此�τ(s)i∈�L(s+1)i。因此，我们有τ(s)i∈L(s+1)i，确保能量在从级别s移动到s+1时不会增加。100503.2. 实际考虑因素0在本节中，我们描述了将所提出的方法应用于实践的一些方面。03.2.1 网格连通性0将模板网格的每个三角形应用单独的刚体变换后，通常得到的网格可能存在相邻三角形之间的间隙或交叉。然而，由于引入了正则化项，这些间隙或交叉可以预期较小。为了恢复原始的网格拓扑结构，我们将模板网格中具有相同位置的顶点子集替换为它们的重心。03.2.2 内存需求0为了运行我们的算法，我们预先计算数据项，需要 O ( n ∙ ℓ)的内存（在线计算具有恒定的内存需求，但会导致显著增加的运行时间）。由于预先计算完整的成对项需要 O ( n 2 ∙ ℓ2 ) 的内存，我们只预先计算弯曲项  d SO(3) ，需要 O ( ˜ ℓ2 ) 的内存。拉伸项  d SE(3) ,X  在线计算。04. 结果0为了评估我们的方法，我们重点展示了我们方法的普适性。通过将点集配准[3, 48, 8, 36, 44, 43, 23, 20, 49,39]和3D形状的相关对应问题[56, 60, 28,7]重新构建为形状到图像匹配问题，我们的方法可以解决这些问题。因此，在下面的评估中，除了3D图像分割，我们还考虑了可变形网格配准的情况。04.1. 可变形网格配准0在第一组实验中，我们展示了我们的方法可以用于进行可变形网格配准。为了强调我们的方法对初始化不敏感，我们将其与一种基于期望最大化的广泛使用的点云配准方法Coherent Point Drift (CPD)算法[44,43]进行了比较。模板和目标：为了评估，我们使用低分辨率的斯坦福兔子网格作为模板（n =498），如图2（左上角）所示。共使用20个具有随机姿态的兔子网格的变形版本作为配准目标。为此，我们基于一个立方体网格上的8个控制点上定义的随机3D位移向量，合成了高分辨率兔子网格的变形版本（图2，右上角）。然后，将这些位移向量传递给网格，使用基于样条的方法。0为了实现平滑和非线性变形，我们使用插值方法。图2显示了其中一个变形后的网格（左下角）。最后，对变形的形状应用随机姿态变换，如图2（右下角）所示。0图2.兔子网格。左上：模板。右上：高分辨率目标。左下：变形目标（与原始目标重叠）。右下：具有随机姿态的变形目标。0图3.兔子模板与变形目标进行配准的定性结果（参见图2右下角）。左上：CPD结果（形状破坏）。右上：CPD误差。左下：我们的结果。右下：我们的误差。0评分图像：为了将我们的方法用于网格配准，我们为每个目标网格创建一个评分图像，然后将模板网格拟合到评分图像上。对于 d : Ω → R + 0，我们将位置 x到目标网格边界的距离表示为 d ( x )。现在，我们将评分图像定义为 J ( x ) = exp( − d ( x)0β ) ，其中我们使用了 β =2 。050100150200020406080100Jo(x) = wo exp�−−d2o(x)2σ2o�.(17)Yo = {x ∈ R3 : Iuneto(x) = 1}(18)10060误差0%顶点0我们的CPD0图4. 垂直轴上具有小于或等于水平轴上数值的顶点的百分比。0max({nx,ny,nz})和λB = 27 ∙ g02|E|考虑了问题规模。正数cmax是单个三角形数据项可能的最大绝对值的上界（参见公式（6）），我们通过将分数图像的最大值乘以最大三角形的面积来计算。标签空间的大小为|L(0)| = 93 ∙ 577 =420,633，分数图像的尺寸范围从208^3到262^3，导致每次在MacBookPro（2.5GHz，16GB）上的注册的平均处理时间约为92分钟。结果：在CPD的情况下，我们首先解决刚性配准问题，然后解决非刚性配准问题（直接进行非刚性配准效果更差）。由于CPD对初始化非常敏感，在20个评估案例中有17个失败，其中一个代表性的失败案例如图3（顶部行）所示。这种极端的配准失效数量强调了对初始化鲁棒性的必要性。相比之下，在所有20个案例中，我们的方法都能够实现良好的配准，如图3（底部行）所示为一个代表性结果。在图4中，我们呈现了定量评估。讨论：虽然我们不声称对网格配准方法进行了全面评估，但我们以概念验证的方式证明了我们的方法对初始化的不敏感性。我们方法的一个优点是，我们既不需要兼容的网格拓扑，也不需要兼容的网格离散化，因为目标是以分数图像的形式表示的。由于分数图像是目标形状表面的离散表示，我们的方法相当于基于表面的配准，而不是像CPD那样基于点的配准，后者倾向于将点对齐得尽可能接近。此外，分数图像还提供了更多的灵活性，因为可以集成其他信息（例如不确定性、网格纹理、形状特征等）。04.2. 分割0在第二组实验中，我们将我们的方法应用于16个多模态3T磁共振图像中的四个脑结构（单个对象的黑质和丘脑核以及红核，均为双侧）。即使对于人类来说，划分丘脑核也是一个具有挑战性的任务[51]。0主要困难包括图像对比度较弱和脑结构的尺寸较小（图5中显示的结构包含在约6×3×2.5cm3的边界框中，MRI图像覆盖约20^3cm3的体积）。模板：为了捕捉脑结构之间的相互关系，我们使用一个多对象模板（n =379），如图5所示。该模板通过（退化的）三角形连接相邻的脑结构，称为“幻影三角形”，这些三角形仅用于平滑项，并且在数据项方面是“自由”的。0图5. 脑结构模板。0max({nx,ny,nz})和λB = 90 ∙ g0π ，其中 g 的定义与之前相同。标签空间的大小为 |L(0)| =113 ∙ 577 =767,987，所有分数图像的维度为364×436×364，导致每次拟合的平均处理时间约为58分钟。分数图像：为了使用我们的方法进行图像分割，我们使用基于最近提出的3DU-Net CNN[9]的数据项。对于所有16个图像，我们以留一法的方式训练网络，用于预测体积分割。图7的中间列显示了三个预测的体积分割示例。对于这个具有挑战性的分割任务，U-Net能够识别出脑结构的（粗略）位置，但在许多情况下不能产生类似脑结构形状的输出（图7的前两行）。因此，我们使用我们的方法补充U-Net分割的几何信息。对于每个脑结构o∈{1,2,3,4}，我们使用一个独立的分数图像Jo。给定脑结构o的二进制U-Net分割Iuneto: Ω →{0,1}，我们首先使用形态学运算提取（预测的）边界。对于d o: Ω → R+0，我们用do(x)表示位置x到提取的边界的距离。然后，我们使用高斯核将分数图像定义为0权重  w o  和带宽  σ o用于将对U-Net分割的置信度纳入考虑。为此，对于表示脑结构  o  的U-Net分割  I unet o ，我们使用以下公式：MAD(Yo) =0.20.40.600.20.4DSCour − DSCunet10070对于表示分割体素的一级集合  I unet o，我们使用（每个坐标的）中位数绝对偏差（MAD）的平均值作为鲁棒分散度度量，计算如下：03 ∥ MAD ( Y o ) ∥ 1 ，(19)0对于  MAD ( Y o ) = median( |Y o  −  median( Y o ) | )  ∈  R 3 。0中位数是以每个坐标为基础理解的，集合向量差异和绝对值在元素上理解。现在，对于给定模板的结构  o  的分割体素的点云  Y ′ o ，平均MAD差异的绝对值由  h o  =  | MAD (Y o )  −  MAD ( Y ′ o ) | 给出。因此，我们定义  σ o  =  ρ ( h o +1) ，其中 ρ =3缩放。因此，如果U-Net分割和模板的平均MAD相等，则带宽对应于 ρ，而分散度差异较大导致较大的带宽，考虑到U-Net分割中的更多不确定性。此外，我们定义  w o  = 10h o +1，这样脑结构  o的U-Net分割中的不确定性增加会导致数据项的权重减小。结果：为了评估分割，我们使用Dice相似系数（DSC）作为体积重叠度量，定义为 2 |Y∩Y ′ |0|Y| + |Y ′ | ，其中  Y  和  Y ′分别是分割体素的点云。我们计算每个单独脑结构的DSC，然后报告四个DSC值的平均值。图6中的定量结果比较了普通U-Net分割和我们得到的分割结果。左侧的箱线图显示，总体上我们的方法在16个案例中实现了更高的体积重叠。此外，右侧的排序DSC差异图表明，我们的方法在大多数情况下改善了DSC（值大于零），只有少数情况下略微降低。在图7中，展示了三个分割案例的定性结果，分别对应图6中的最好、中位数和最差案例（右侧）。0our unet0DSC0最好 中位数 最差0图6.左：我们方法与U-Net分割的DSC箱线图。右：16个案例的排序DSC差异（值大于零表示对U-Net的改进）。0定性结果显示了三个分割案例，分别对应图6（右侧）中的最好、中位数和最差案例。在图7的第一行中可以看到，在某些情况下，我们的方法甚至能够在较差的U-Net分割基础上实现合理的分割。这是通过更强调形状信息相对于数据项来实现的，这也会对结果产生偏差。0方法朝着形状信息的方向发展。关于模型对数据拟合方法在运动结构中的偏倚性的讨论可以在[45]中找到。groundtruth U-Net our0最好0中位数0最差0图7.脑结构分割实验的定性结果。每一行显示了不同结果的实例，从上到下依次为最好、中位数和最差的DSC差异（参见图6右侧）。05. 结论0我们引入了第一个用于非刚性匹配3D形状和3D图像的组合方法。关键思想是将3D形状表示为三角网格，并在三角形集合上解决一个流形值多标签问题。我们通过最小化成本函数确定与每个三角形相关联的刚体变换的分配，其中一元项编码图像中的局部匹配成本，而成对项惩罚变形中的