通过分层Wasserstein核快速且可解释地量化生物形状的异质性

《PLOS Computational Biology》：Fast and interpretable quantification of biological shape heterogeneity via stratified Wasserstein kernel

　　摘要
现代成像技术产生了大量的细胞和亚细胞结构数据，这需要原理性的方法来实现在个体和群体之间的形状比较。我们引入了分层Wasserstein框架，该框架将每个形状视为一个非结构化的点云，并通过排名的局部距离分布将其嵌入到欧几里得空间中。这种嵌入产生了一个等距不变的欧几里得距离和一个正定核，用于群体分析，同时提供了一个基于样本的估计器，能够在接近 quadratic 的时间内处理大型数据集。通过利用核方法，该框架能够执行统计上严谨的任务，如非参数假设检验，提供理论保证和可解释性。我们展示了该框架在大型生物数据集上的适用性。通过分析2D癌细胞轮廓，我们量化了群体层面的差异，并识别了对观察到的差异贡献最大的代表性细胞。利用3D细胞膜和细胞核体积，我们揭示了跨细胞群体的形态变化模式。这些结果为群体层面的生物形状分析提供了一个简单而原理性的工具，可能在计算成像和数据科学的不同领域产生重要影响。

作者总结
生物结构有多种形态，从整个组织到单个细胞，甚至蛋白质构象。现代成像技术现在产生了大量的这些形状数据，使我们有机会研究结构在不同条件下的变化或随时间的演变。然而，快速且可解释地比较大量复杂形状仍然很困难。许多现有方法依赖于手工选择的特征或标志点，而其他方法则太慢，无法应用于生物学中常见的大型数据集。在这项工作中，我们介绍了一种广泛适用的方法，该方法直接从原始几何形状进行比较。我们的方法使用其内在距离来总结每个形状，而不依赖于标志点或对齐，并且足够高效，可以处理大型数据集。由于该方法基于坚实的统计原理，它使我们能够确定两组形状在哪些方面存在有意义的差异，并识别出对这些差异贡献最大的例子。我们还展示了该方法可以揭示平滑的形态趋势，使研究人员能够追踪形状沿生物轨迹的演变。总体而言，我们的框架为分析不同生物系统中的形状变化提供了一个通用、可扩展且基于统计的工具。

引言
高通量成像的最新进展产生了从整个组织到单个细胞、亚细胞结构和蛋白质构象的大量生物形状数据集。这些数据集使得在条件变化下进行群体层面的形态分析成为可能，但同时也暴露了一个方法论缺口：我们缺乏既可解释又在计算上高效的工具，用于比较单个对象和对象分布层面的形状。例如，关于细胞周期中形态如何演变[1]或不同物种之间蛋白质结构如何不同[2]等问题，需要尊重内在几何的比较方法，能够扩展到大型数据集，并自然地与统计分析和机器学习联系起来。传统的形状比较方法面临众所周知的限制。许多方法依赖于将数据简化为特定的形状特征，如体积、高度、拓扑数据分析中的持久条形码[3]、球谐函数[4]，其中感兴趣的特征是预先确定的。基于标志点的分析流程通过特定标志点的对应关系来量化成对差异，这些标志点是手动放置的，可能劳动密集、主观且难以复制[5]，适用于大型数据集。其他经典参数方法针对特定领域进行了定制[6,7]，并且/或者需要大量的预处理，包括下采样、对齐或插值[8,9]。在包含数千到数百万个形状的当代成像流程中，这些限制成为了瓶颈。尽管大多数现有工作也集中在单个层面的成对距离或形状对齐上，但群体层面的形状集合之间的度量仍然很少被探索。从历史上看，这种差距可以归因于形状数据集的相对较小规模，分析仅限于几十到几百个样本。随着现代成像流程现在产生数千到数百万个形状，对可扩展且基于统计原理的群体层面距离的需求变得紧迫。已经提出了一些群体层面的方法。例如，在特定背景下考虑了Fréchet均值和到所有形状的距离[10]、线性子空间嵌入和Kullback–Leibler散度[11]，但它们都受到特定假设的约束，没有一个提供通用的、基于理论的、计算高效的方法来比较任意形状的群体。在这个背景下，最优运输理论提供了一种原理性的方法，用于比较作为非结构化点云或概率直方图提供的形状对象。特别是，Gromov–Wasserstein距离通过它们的内部成对距离来比较形状，避免了显式对齐[12]。相关的工作通过距离分布来总结形状，无论是在全局还是围绕每个点[12–15]。然而，众所周知，这样的计算可能是难以处理的[16]且成本高昂[17]，因此很难应用于大规模数据集。此外，很难直接从这些维度大于一的距离构建正定核，这限制了下游统计工具，如核假设检验[18–20]和表示学习[21,22]。现有的工作[23,24]通过切片度量并将运输问题简化为一维来解决这个问题，但这需要切片的采样角度，并且可能无法充分利用整个形状的信息。我们通过引入分层Wasserstein距离来应对这些挑战，该框架简单而高效，将每个非结构化点云嵌入到欧几里得空间中。我们展示了这种距离是厄米的，并产生一个正定的核，在温和的规则条件下，该核是特征的。这种构造在设计上是等距不变的，并且在等距意义上是单射的。在计算上，该方法每个形状的点数几乎是二次方的，并且在通过熵正则化加速后，其性能优于Gromov–Wasserstein（后者是三次方的）。我们在2D和3D形状数据集上展示了该框架的实用性，重点关注个体和群体层面的形状分析。

结果
我们首先介绍了分层Wasserstein框架的概述，描述了其构造和关键属性。这为本文的其余部分提供了概念和理论基础。然后，我们在二维和三维生物成像数据集上展示了其性能，突出了其处理多种形状类型和支持群体层面统计分析（如假设检验和依赖性检测）的能力。

分层Wasserstein框架概述
我们提出了分层Wasserstein框架，该框架将每个形状（表示为非结构化点云）嵌入到欧几里得空间中，并便于在该空间中使用核方法进行群体层面的量化任务。与现有的形状距离（如Gromov–Wasserstein及其下界[12]（定义在方程M3、M4和M5中（方法部分））相比，所提出的框架在区分能力方面相似，同时在计算上更高效，复杂度几乎是点数的二次方，而实验运行时间通常低于Gromov–Wasserstein方法所需的1%。诱导的核是特征的，因此通过核方法得到的群体层面统计保持了标准的统计保证，包括一致性和对替代方案的检验能力[18,19,25]。图1提供了我们程序的示意图，该程序使用局部距离分布来产生形状嵌入并推导群体层面统计。这里我们描述了主要的构造并总结了其关键属性；详细的陈述和证明被推迟到方法部分和补充信息中。

下载：PNG（较大图像）、TIFF（原始图像）
图1. 我们程序的概述。
给定一个表示各种类型形状的2D或3D非结构化点云，我们首先计算每个形状的内在距离，并将其局部化。使用这些局部距离，可以对大量形状执行各种下游任务，如降维、聚类、假设检验和特征选择。详细的方法论在后续部分中提出，并进行了严格的理论基础讨论。

首先，我们用一个点云和内在距离来表示每个形状，例如表面上的测地距离或环境空间中的欧几里得距离。对于形状上的每个点，我们计算到所有其他点的距离分布，这捕捉了其局部几何结构。然后根据这些局部距离分布的连续函数对点进行排序，从而对形状进行一维分层。在每个层次内，我们通过记录距离分布的一组分位数来总结局部几何结构。这定义了形状到由分层变量和分位数级别索引的二维函数空间的嵌入。两个形状之间的分层Wasserstein距离仅仅是它们嵌入之间的距离，在这个二维域上计算得出。在这项工作中，我们展示了分层Wasserstein距离继承了许多Gromov–Wasserstein的理想属性，同时显著提高了计算效率。由于局部距离分布的不变性属性，它在等距变换下自然是对称的。在对形状进行分解的排名功能上满足温和的规则条件时，距离在等距意义上是单射的：两个形状之间的距离为零当且仅当它们是等距的。因为距离是由嵌入空间中的范数诱导的，所以基于它的标准核（如高斯核和拉普拉斯核）是正定的。当距离是单射时，这些核也是特征的，意味着它们可以区分任何两个形状的分布。基于有限点云的距离估计器在统计上是一致的：随着每个形状的点数增加，以及相应地细化离散化，估计的距离趋近于群体距离。计算成本几乎与每个形状的点数成二次方，这比Gromov–Wasserstein方法高效得多。最后，因为核是特征的，从中构建的标准群体层面统计具有强大的理论保证。例如，最大平均差异（MMD）[18]测量两个群体之间的任何形式的差异，而Hilbert–Schmidt独立性准则（HSIC）[19]测量变量之间的依赖强度。这两种统计方法都是一致的，并且对所有固定替代方案都有检验能力，这使得可以进行严格的双样本非参数检验和形状分布与外部协变量之间的依赖性检测。我们提供了表1，将我们提出的分层距离与其他相对替代方案进行了比较。与其他距离相比，分层Wasserstein是希尔伯特的，可以从有限样本中一致估计，并且实现了最低的渐近计算复杂度。如果选择排序和分箱的预言函数得当，它是条件单射的，并产生一个特征核。在下一部分中讨论的数字示例中进一步展示了区分能力和计算复杂度之间的平衡，这些示例涉及真实的生物形状数据集。

在简单合成形状上的这些距离的经验比较可以在S3和S4图中找到。

下载：PNG（较大图像）、TIFF（原始图像）
表1. 以2D/3D中的点云表示的形状之间的距离比较。对于Gromov–Wasserstein和Wasserstein之间的局部距离，分别显示了精确计算（左）和熵近似（右）的计算成本。N = 每个形状的点数。*精确的GW被广泛认为是NP难的[16]；?熵GW的复杂性取决于（正则化）和（容忍度）[26]，两者都是小的正值，对于准确的近似；§关于计算成本的系统研究可以在[27–29]中找到。§我们的框架在假设（H1）和（H2）下具有条件单射和特征核，详细信息在补充信息中。

乳腺癌细胞轮廓形状来自荧光显微镜
我们将我们的框架应用于从荧光显微镜获得的2D癌细胞形状[32,33]。细胞图像被二值化，其边界被提取出来形成由细胞轮廓的2D坐标给出的离散曲线。数据集包括来自三种不同乳腺癌细胞系的癌细胞形状：（1）MCF10A（228个细胞）：非肿瘤性的人类乳腺上皮细胞系，是正常细胞的经典模型；（2）MCF7（225个细胞）：转移潜力相对较低的乳腺癌细胞系；（3）MDA-MB-231（缩写为MDA，224个细胞）：高度侵袭性的三阴性乳腺癌细胞系，常用于作为转移进展的模型。给定这些细胞轮廓的群体，目标是测试不同细胞系定义的群体之间的形状分布是否存在显著差异。理解这些形状差异可能有助于揭示细胞形态是否以及如何在不同的癌症类型中编码功能行为。图2A和2B展示了原始图像示例以及每个细胞系中的随机细胞样本，更多图像可以在S1图中找到。实施细节可以在“方法”部分和表2中找到。下载：PNG（较大图像）TIFF（原始图像）。表2. 2D和3D形状分析的实验参数和实施细节总结。https://doi.org/10.1371/journal.pcbi.1014254.t002 下载：PNG（较大图像）TIFF（原始图像）。图2. 通过轮廓表示的2D癌细胞形状群体的案例研究。A：[32]整理的细胞系的荧光显微镜图像示例。B：分割后的细胞轮廓，每个群体中随机选取了8个形状进行可视化。C：通过Fréchet平均形状和在SRV度量下每个群体内到平均形状的距离检查，使用Geomstats [10]计算，表明MCF10A和MDA群体之间的相似性。D：为了进一步探讨MCF10A和MDA之间的差异，我们使用我们提出的分层核计算了基于MMD的统计量及其p值，利用了所有可用的细胞形状。E：每个群体中得分最高的8个细胞突出显示了与其群体大多数形状差异最大的细胞。F：在小型（15）到中型（60）样本量下，通过五种独立试验平均的II型错误的经验概率和计算时间（误差条：±1 SD）。https://doi.org/10.1371/journal.pcbi.1014254.g002 相对MMD测试揭示了三阴性癌细胞和非癌细胞之间的相似性：在测试三个群体之间的相对MMD时，我们观察到：尽管MDA是癌细胞系，但它与MCF10A（正常细胞）的相似性大于与MCF7（低转移性癌症）的相似性。平方MMD的经验估计显示MMD2(MCF10A,MDA)=0.0189，远小于MMD2(MCF7,MDA)=0.1409。为了测试这种差异是否具有统计学意义，我们对以下情况进行了相对MMD测试：p值是使用排列测试计算的，在零假设下，两个群体与MDA的距离相等。具体来说，我们排列MCF7/MCF10A的标签，并重新计算每次排列后的MMD与MDA的差异。经过1000次排列，该测试得出的p值为0.001（S1C图），表明MDA（高转移性癌症）在形态上更接近MCF10A（正常细胞系）而不是MCF7（低转移性癌症）。这一观察结果与基因表达研究[34]的发现一致，该研究报告称MCF10A和MDA都属于基底型亚型，而MCF7属于管腔型。我们的基于形状的结果表明，细胞形态与分子亚型（基底型与管腔型）的关系更紧密，而不是与癌症状态的关系。这一结果也与之前的方法[35]一致，该方法使用geomstats [10]包在平方根速度（SRV）度量下计算Fréchet平均形状。所有3个群体的处理后的轮廓的平均形状显示在S1B图中，其中MCF10A和MDA的平均形状都是细长的，并且几乎相同。绝对MMD测试表明群体在尾部存在差异：鉴于图2C中显示的平均形状和到平均距离的分布几乎相同，看起来MDA和MCF10A的两个形状群体在平均值和方差上是一致的。一个自然的后续问题是，这两个形状群体的生成分布是否真的相同。为此，我们放大这个意外的配对，并对这些群体进行了MMD测试，检测到分布的任意阶差异，超出了二阶：经过1000次排列，结果显示p值为0.002（图2D），表明这两个群体之间存在显著差异。为了定位形态差异，我们检查了MMD证人函数识别的极端形状，该函数突出了形状分布差异最大的区域。图2E中显示了每个群体中MMD证人函数值最高的8个细胞。MDA细胞的极端子集表现出高度不规则和突出的形态，这些形态是区分MCF10A的关键因素。我们假设这种丰富的形状可能对应于MDA已知的侵袭性行为[36]。值得注意的是，MMD的统计差异不是由平均形状或到平均形状的距离（图2C）驱动的，而是由MDA中的这个尾部群体驱动的，这只能使用核MMD这样的非参数测试来捕捉。这些测试在计算上高效且在统计上具有足够的功效：为了评估绝对和相对测试的功效，我们进行了随机实验以量化经验I型和II型错误率。我们将我们的核算法与Zhang等人[37]建议的框架进行了比较，该框架使用了平方根速度（SRV）度量及其变体（Elastic）进行MMD。为了确保与需要更平滑、参数化曲线的方法进行公平比较，将涉及SRV和弹性度量的基准方法应用于预处理数据，其中每个形状都已被插值到2000个点，并事先对齐到相同的参考形状。相比之下，在我们的方法中，我们直接在原始数据上进行测试，无需任何处理。我们通过MMD对样本量在15到60之间变化的样本进行了测试。为了系统地量化错误，我们对每个样本量重复实验1000次，通过排列测试计算p值，并在图2F中报告犯错误的相对频率。对于II型错误，我们从每个组中抽取一个样本，进行测试，并计算无法拒绝形状分布相等的相对频率。在我们提出的距离中，Gromov–Wasserstein一致地实现了最低的II型错误率，其次是分层Wasserstein距离和Gromov–Wasserstein的第二下界（SLB）。所有Wasserstein类型的距离在样本量小于40时都优于基于SRV的度量，这与Wasserstein类型距离作用于内在距离分布并且对曲线参数化和离散化效应不太敏感的事实一致。对于I型和II型错误率，也进行了类似的实验和观察，并在S1图中显示。我们还注意到错误率对超参数的选择具有鲁棒性，敏感性结果在S2图中报告。为了评估计算效率，我们记录了计算每个子样本中所有成对距离所需的墙时间，平均了五次独立运行。这些结果显示在图2中。虽然GW提供了最佳的统计功效，但其计算成本相当高。我们的分层Wasserstein实现了有利的平衡：虽然在II型错误控制方面略逊于Gromov–Wasserstein，但它比Gromov–Wasserstein及其第二下界更高效，并且仅次于简单聚合和排序所有成对距离的汇总Wasserstein方法。总的来说，我们的非参数、分布级方法直接在原始轮廓形状上操作，检测尾部效应，并定位差异的驱动因素，同时在计算上仍然实用。艾伦研究所的3D细胞和细胞核形状 在这个例子中，我们使用Viana等人的[1]的3D细胞成像数据展示了距离矩阵的分位数编码了有意义的形态信息。我们分析了原始作者用于训练分类器的5764个细胞的子集，包括反映六个有丝分裂阶段（M0、M1M2、M3、M4M5、M6M7早期和M6M7中期）以及三种异常类型（斑点、死亡、错误）的细胞类型注释。我们方法的实施细节可以在“方法”部分和表2中找到。现有的降维方法要么依赖于手工制作的特征（例如，细胞和细胞核的体积和高度）[1]，要么依赖于来自深度模型（如变分自编码器[38]）的学习到的潜在表示，然后使用UMAP。然而，这些方法要么针对特定应用高度专业化，要么需要大量的训练工作。形状的降维揭示了相关的形态变化。我们通过侵蚀每个二值图像来提取表面点，然后进行下采样以保留每个细胞200个点。对于每个细胞，我们将其分为100个区间，计算每个区间内细胞形状和细胞核形状的100个分位数，将它们连接起来，并对得到的20,000维向量应用UPAM。图3A显示了使用分层Wasserstein距离获得的2D UMAP嵌入。为了基准测试我们的结果，我们使用了4种额外的方法进行了降维，并在图3B中展示了它们，这些方法基于（1）原始二值图像中的欧几里得距离，（2）通过细胞和细胞核体积的主导特征，（3）使用预训练的PointNet模型的中间潜在层[39]，以及（4）由于计算成本过高而将每个细胞下采样到100点的Gromov-Wasserstein距离。为了公平比较，我们还将提出的分层Wasserstein嵌入应用于下采样到100点的形状。如图S4所示，得到的嵌入在质量上与使用200点得到的嵌入相似。我们的降维最好地保留了形状动态的循环和连续性，这与基于特征和PointNet的嵌入一致，但轨迹更平滑，反映了Viana等人[1]对手动注释模糊性的评论。下载：PNG（较大图像）TIFF（原始图像）。图3. 3D细胞和细胞核形状群体的案例研究。A：使用分层局部距离分布，通过Viana等人[1]的注释有丝分裂数据集的UMAP嵌入，以及周期中的代表性细胞（蓝色）和细胞核（红色）形状。B：标准方法在细胞级别上的嵌入无法揭示动态。原始二值图像之间的欧几里得距离和降采样点云上的Gromov–Wasserstein距离没有显示出有意义的模式，而基于特征的嵌入（体积或神经网络）捕获了异质性但没有捕获有丝分裂的循环连续性。C：使用UMAP显示的细胞特征允许解释每个阶段内的形状变化。D：使用所提出的MMD进行的群体级别嵌入揭示了与A一致的循环过程。E：个别细胞的HSIC证人分数，显示了在有丝分裂后期附近具有混合正负值的区域。我们选择了3个具有高正证人分数的细胞（黄色十字）和3个具有高负证人分数的细胞（洋红色星形），并在旁边显示了细胞核形状。https://doi.org/10.1371/journal.pcbi.1014254.g003 这种嵌入提供了有关有丝分裂周期中形态变化的生物学洞察。在图3C中，我们根据形态特征对细胞进行着色，以显示它们如何随周期演变。从M1M2阶段开始，细胞和细胞核的体积显著增加，随后在M3和M4M5阶段细胞核高度显著升高。这些趋势与已知的有丝分裂阶段一致，其中体积增加反映了DNA和细胞器的复制，而结构转变（例如，染色质凝聚和核包膜破裂）是后期阶段的标志[40]。异常细胞的体积和细胞核形状都特别小，可能反映了不完全的有丝分裂、异常状态或分割伪影。群体的嵌入尊重了有丝分裂周期的进展。为了更好地总结对应于不同有丝分裂阶段的形状群体之间的差异，我们使用MMD作为群体之间的距离度量对细胞群体进行了嵌入。图3D显示了得到的距离矩阵，该矩阵呈现出带状结构：相邻阶段的群体（例如，M0和M1M2）之间的距离较低，而异常细胞类型则与有丝分裂群体相距较远。由于MMD定义了有效的群体级别距离，我们进一步应用了多维缩放（MDS）将群体嵌入到2D空间中，揭示了有丝分裂阶段的明显进展。我们注意到，两种替代的群体度量，能量距离（和Wasserstein距离），并没有像基于MMD的嵌入那样表现良好，这可能是由于核方法通过平滑对噪声更加抵抗，这对于真实成像数据至关重要。与其他群体级别距离的候选者相比，可以在S3图中找到。HSIC证人分数像在基于MMD的测试中使用的证人函数一样，提供了一种解释个别样本如何贡献于变量之间整体依赖性的方法。我们排除了异常细胞类型，并使用从0到1的等间隔值来更好地捕捉周期性，其中使用了von Mises核。将这些评分应用到细胞周期中的细胞和核形状上，可以突出显示出形状对有丝分裂阶段具有特别信息价值或无信息价值的特定区域（图3E）。为了更好地理解这些区域，我们可视化了具有最强正负见证评分的代表性细胞。有趣的是，尽管所有选定的细胞都属于有丝分裂的后期阶段，但它们的核形态却有很大的不同。具有正见证评分的细胞（黄色十字）倾向于显示双叶或分段的核形状，这与活跃的分裂一致；而具有负评分的细胞（洋红色星星）则表现出平滑、未分裂的核。这表明HSIC见证评分能够捕捉到细微的形态差异，这些差异反映了形状与有丝分裂进程的匹配程度，有助于识别具有不同发育状态的阶段意识原型。总体而言，我们的框架在揭示已知的有丝分裂轨迹的同时，还能够找到罕见的、特定阶段的形态，展示了揭示数据内部未知异质性的潜力。

讨论：为了提供在个体和群体层面上量化形状数据异质性的工具，我们提出了一个分层Wasserstein框架，将形状数据嵌入到欧几里得空间中，并利用其中的核方法。这种构建在等距变换下是自然不变的，并且可以通过最大均值差异（MMD）为形状-形状距离和群体间距离提供一致的样本估计。在2D乳腺细胞轮廓和3D有丝分裂细胞及核的分割上，该框架支持降维、聚类和非参数假设检验，其性能与最先进的方法相当或更优，同时大大减少了运行时间。虽然我们关注的是平滑的平面曲线和表面，但该公式可以扩展到其他度量数据类型，例如，当配备了测地距离和适度规则性时，可以用于建模为度量图的神经元树[14,41]。由于该方法直接作用于大小可变的非结构化点云，并且不需要地标或全局对齐，因此它为跨领域的形状分析提供了一种统一且可扩展的途径。更广泛地说，我们的框架在涉及形状/图形数据量化的领域中具有潜在应用，如脑成像、蛋白质密度图和社交网络。与基于变形的形状模型不同，分层Wasserstein框架并不是为了提供个体形状的可逆参数化或定义形状之间的测地变形，而是为了支持基于统计原理且计算效率高的形状和形状群体比较。接下来，我们概述了几项限制和改进的潜在途径。首先，在实践中，我们根据点的局部距离分布的平均值对它们进行排序；在近似对称或距离简化的形状上（S3文本中的图B），这种统计量可能具有较低的区分能力。尽管其性能与Gromov-Wasserstein和随机多边形的第二下界相当（S3文本中的图C），但在具有对称性的T形形状上，分层框架下的形状距离被显著高估了（S3文本中的图D）。一个原则性的补救措施是使用局部距离律的前k个矩（或分位数）进行字典排序；我们展示了这样的排序在有限样本上是可逆的，但选择一个小而数据适应的k值以及相应的排序函数仍然是一个未解决的问题。其次，箱数和分位数的选择是用户自定义的，没有具体的指导原则。为分层估计器推导出有限样本误差界限将有助于数据驱动的默认值，以平衡离散化偏差和箱间方差。第三，目前核带宽是通过中位距离启发式方法或手动调整设置的；在当前成像中常见的低信噪比情况下，更原则性的规则（例如，最大化估计的MMD测试功率或核对齐，或噪声意识插件）可以提高性能。解决这些问题将进一步增强其鲁棒性和易用性。所提出的框架有潜力以多种方式促进新的计算工具和生物学发现。通过将形状空间与核联系起来，该框架解锁了一个广泛的推理和表示工具箱。(i) 结合多模态测量（例如，形状和转录组学），基于核的条件独立性（KCI）[20]可以产生p值，以识别解释形态的驱动基因，同时控制混杂因素，有助于推断形态变化的基因调控。(ii) 结合核表示学习[22]和功能数据分析（例如，嵌入上的核PCA或高斯过程），可以学习低维替代物，捕捉形状变化最快的离散表型或连续轨迹。(iii) 由于其计算效率，该框架可以集成到体内扰动筛选中，以优先考虑引起最大形状分布变化的条件。总体而言，这些分析使该框架成为假设检验、因果推理、扰动发现等实用工具。通过解决上述挑战并进一步扩展框架到其他任务，分层Wasserstein有可能成为一个更强大的工具，从而更全面地洞察多种生物环境中高度异质的形状空间。

方法：形状表示和内在距离 我们将每个形状表示为一个非结构化的点云。每个形状都配备了一个内在距离和一个统一的度量。根据应用的不同，可以选择不同的距离，包括欧几里得距离、k最近邻图上的测地距离[41]或扩散距离[42]。内在距离矩阵为...。本工作中一个核心的形状统计量是局部距离分布，定义如下：对于每个点，局部距离分布是前向度量，其中表示上的均匀概率度量，#表示前向运算符。在点云的离散设置中，这简化为经验度量（M1）。为了比较，人们还可以定义形状的全局距离分布为前向（M2），它捕捉形状上成对距离的总体分布。在离散设置中，这对应于所有成对距离的经验分布。众所周知，与全局距离分布相比，局部距离分布具有更好的区分能力[12,14]，因为它们能够保留点级的几何信息，并且无法仅从...中恢复。为了尺度不变性，形状被重新缩放，使得中间成对距离等于一。Gromov-Wasserstein及相关距离 对于两个度量度量空间，i=1,2，p-Gromov-Wasserstein（GW）距离[12,43]定义为（M3）。在离散设置中，这对应于最小化点之间的耦合矩阵。人们普遍认为精确计算是NP难的[16]，通常使用熵正则化来获得近似解，其复杂度为[26]。在所有数值实验中，我们使用了Python Optimal Transport[44,45]提供的实现。全局Wasserstein距离。一个更简单的下界比较了成对距离的全局分布[13,15,46]。设为方程（M2）中定义的经验全局距离的分位数函数。这些分布之间的p-Wasserstein距离为（M4）。这个距离可以在...中高效计算，但它不是单射：不同的形状可能会产生相同的全局距离分布。S3文本中的图A给出了一个说明性的例子[12]。第二下界（SLB）。一个更强的度量比较局部距离分布[12,14]：（M5），其中是点x处的局部距离分布，是点之间的耦合。在规则性假设下，SLB严格强于全局Wasserstein距离，并且对于平滑封闭形状是单射的[14]，但它仍然计算要求高且不是Hermitian的。分层Wasserstein距离 所提出的分层Wasserstein距离结合了GW类型距离的几何区分能力和基于分位数的嵌入的效率。对于每个点，我们计算其局部距离分布，并进一步通过一个标量排序函数对其进行总结，该函数以反映点在形状中的几何上下文。在实践中，我们使用局部距离分布的平均值，这倾向于将靠近形状重心的点与靠近外围或突出区域的点分开。为了使这个描述符在形状之间可比并且对单调变换不变，我们将转换为归一化排名，其中表示的累积分布函数。这种转换将形状上的点映射到单位区间[0,1]，从而对形状进行一维分层。分层保留了按几何上下文排序的点序，同时丢弃了绝对尺度信息，使得不同形状之间的局部几何能够一致聚合。然后沿...分解形状度量，其中表示层次u处的条件律。我们通过评估在级别处的分位数函数来定义函数Q：其中表示的（左连续的）广义逆。对于每个层次和分位数级别，我们计算...总的来说，我们定义了一个嵌入。两个形状之间的分层Wasserstein距离是它们嵌入之间的距离：在排序函数的适度规则性条件下，该距离在等距变换下是单射的（S1文本中的定理2）。由于它是由L2范数诱导的，因此从D2构建的高斯和拉普拉斯核是正定的，如果嵌入是单射的，这些核是特征的。基于分箱和分位数的经验估计器在每个形状的采样点数增加且离散化细化时，在联合情况下具有统计一致性。在排名函数的规则性下，这确保了排名的稳定性，并且在每个层次有足够样本的情况下，估计的距离趋于群体距离（见S2文本）。其计算复杂度接近于点的数量的一次方，相比GW和SLB节省了大量成本。我们框架中的详细数值实验设置 用于比较的基准 在2D平面曲线之间的弹性度量。在2D癌细胞轮廓形状的任务中，为了比较，我们使用了具有两个参数选择的弹性度量[47-49]。设为平滑的参数化曲线，和为沿的两个切向量场。具有参数a, b>0的弹性度量定义为...，其中s表示的弧长参数，是关于弧长的h的导数。平方根速度的a=1/2且b=1，而Geomstats[10]实现的弹性度量的默认值为a=1和b=1/2。对于2D细胞系数据集，我们通过对齐后的曲线进行弹性度和平方根速度度量进行基准测试，每条曲线插值后包含2,000个点，可以从https://github.com/wxli0/dyn/tree/main%4092c7a58/dyn/datasets/breast_cancer/aligned/projection_rescale_rotation_reparameterization访问。3D图像的数据预处理。通过使用scipy.ndimage.binary_erosion进行单次二值侵蚀，将二值体转换为表面点云，并使用默认设置，然后提取边界体素。对于分层Wasserstein框架，结果表面点均匀随机下采样到每个细胞200个点，对于Gromov-Wasserstein基准，则下采样到每个细胞100个点。3D图像之间的Pointnet距离。作为一个现代基于深度学习的基准，我们使用了在ModelNet40分类任务上预训练的pointnet++模型的特征编码层[39,50]，可以通过https://guochengqian.github.io/PointNeXt/modelzoo/访问。每个细胞和核的图像首先下采样到2,048个点，并嵌入到维度为1,024的特征空间中。我们使用这种嵌入进行降维，并在其中使用欧几里得距离进行MMD，作为我们3D细胞图像案例研究的基准。3D细胞形状的特征。我们计算了3D细胞形状的特征，用于基准测试和可视化目的，具体来说，包括体积、高度和主要倾斜度。每张图像都以质心为中心，并应用主成分分析来确定其主导变化轴。体积定义为核心轴上的垂直范围（最大值和最小值之间的差异）。体积是为点的凸包计算的。主要倾斜度是通过主导主轴与第三垂直轴之间的角度计算的。形状群体的核方法 从分层Wasserstein距离构建的正定核使得可以使用核方法进行群体级别的统计分析。我们关注最大均值差异（MMD）和Hilbert-Schmidt独立性准则（HSIC），这些是广泛使用的具有强理论保证的非参数统计量。最大均值差异（MMD）[18]是在再生核希尔伯特空间（RKHS）中两个分布的核均值嵌入之间的平方距离，与核k相关：给定样本和，（有偏的）经验估计器是...通过组标签的排列测试来评估显著性。当k是特征性的（如分层Wasserstein度量上的高斯或拉普拉斯核）时，如果且仅当，MMD等于零。该检验对于所有固定替代假设都是一致的，这意味着随着样本量的增加，其检验功效会收敛到1 [18,51]。MMD见证函数：MMD见证函数[18]可以识别出两个分布差异最大的形状空间区域。通过评估样本形状来突出显示对观察到的差异有贡献的代表性形状。希尔伯特-施密特独立性准则（HSIC）：HSIC[19]使用每个域上的核函数来量化随机变量之间的依赖性，并可以解释为联合分布与其边缘分布乘积之间的最大均值差异（MMD）。对于成对数据，如果形状的特征核为k，协变量的特征核为?，则HSIC定义为相关核函数空间（RKHSs）中交叉协方差算子的平方希尔伯特-施密特范数。当Gram矩阵和中心化矩阵满足条件时，只要两个核都是特征核，（有偏的）估计量的HSIC值为零，这表明变量Y和X是独立的。HSIC对于所有替代假设都是一致的，并且能够检测到任意非线性依赖性[19,21]。HSIC见证函数类似于MMD，它也允许一个见证函数来定位各个形状-协变量组合对整体依赖性的贡献。给定在RKHS中固定的核k和?，经验HSIC见证函数为：其中 和 分别是中心化的核函数。正值表明一对形状-协变量支持观察到的依赖性，而负值则表示相反。这个函数实际上是将MMD见证函数应用于联合分布及其边缘分布的乘积得到的。这个函数可以用来解释特定形状或协变量值的贡献。

统计分析：所有群体层面的比较都是使用基于分层Wasserstein距离和其他基准距离的核方法进行的。形状群体之间的双样本检验使用了最大均值差异（MMD）统计量[18]。对于每个检验，我们使用1,000次排列来估计零假设分布，并根据超过观察值的排列统计量的比例来报告p值。MMD见证函数用于识别对组间显著差异有贡献的代表性形状。为了评估形状分布与外部协变量（例如发育阶段）之间的依赖性，我们使用了希尔伯特-施密特独立性准则（HSIC）[19]，在形状和协变量域上都使用了Laplace核，带宽通过中位数距离来选择。为了清晰起见，本节中的所有检验都是使用Laplace核在形状之间进行的，带宽的选择基于中位数启发式方法：其中 表示形状空间中的距离，如Gromov-Wasserstein及其替代方法。对于3D有丝分裂细胞形状数据，有丝分裂阶段变量包括六个有序类别（M0, M1M2, M3, M4M5, M6M7_early, M6M7_half），这些类别被编码为整数。有丝分裂阶段标签被编码为一个角度变量，并映射到单位圆上的圆坐标系。对这些圆形嵌入应用高斯核，相当于（在一个常数因子下）von Mises核。带宽参数 是通过中位数启发式方法选择的，即设置为单位圆上阶段嵌入点之间平方距离的中位数。这种编码方法保留了有丝分裂过程的循环拓扑结构，同时允许阶段变量与形状之间的欧几里得距离自然交互。通过对协变量标签进行1,000次排列试验来评估显著性。HSIC见证函数用于可视化对依赖性贡献最大的形状-协变量对。

计算细节：所有实验都是用Python（v3.13.7）实现的。数值计算使用了NumPy（v2.3.4）。最优传输计算使用了POT库（v0.9.6）。基于核的统计检验依赖于scikit-learn（v1.7.2），降维使用了UMAP（v0.5.9），而基准测试中使用的弹性形状分析则采用了Geomstats（v2.8.0）。计时实验是在配备Apple M4 Pro芯片和48 GB统一内存的MacBook Pro（14英寸，2024年款）上进行的，操作系统为macOS Sequoia（v15.7.3）。所有方法都是在相同的硬件和软件条件下评估的；因此报告的运行时间用于相对比较而非绝对基准测试。

支持信息：参考文献