希腊达夫诺（Dafnero）发现的DFN3-150 Paradolichopithecus aff. arvernensis标本的面部虚拟重建与分析

《Scientific Reports》：Virtual reconstruction and analysis of the face of DFN3-150 Paradolichopithecus aff. arvernensis specimen from Dafnero, Greece

　　**摘要**
欧亚地区晚更新世至早更新世时期的猴类Paradolichopithecus的系统分类一直是灵长类古生物学中长期争论的话题。这一分类单元的系统发育位置对于我们理解疣猴科（Cercopithecidae）的进化具有重要意义：虽然它通常被视为类似狒狒的Macacina亚科成员，但也有学者提出将其归入Papionina亚科的另一种分类方案。这一观点将挑战“狒狒构成非洲特有演化的支系”这一普遍接受的观点。在本研究中，我们利用新颖的虚拟技术及不同的重建方法，对来自希腊下更新世Dafnero遗址的亚成年雌性Paradolichopithecus aff. arvernensis头骨DFN3-150的面部结构进行了重建。通过对这种动物与现存Papionina属物种的形态学特征进行几何形态测量分析，我们发现其面部形态与Papio属物种的相似度高于与猕猴属（Macaca）的相似度。这一结果与近期关于其表型和生态特征的研究结果一致。然而，要直接推断其系统发育关系，还需要更多样化的Papionina属样本数据。

**引言**
疣猴亚科（Cercopithecinae）作为猕猴科（Cercopithecidae）的一个分支，大约在中新世早期至中期（1440万至2060万年前）开始分化1,2,3,4,5,6。该亚科进一步分为两个族：Cercopithecini族和Papionini族，后者在晚中新世（约800万至1100万年前）进一步分化为Macacina亚族（猕猴系）和Papionina亚族（狒狒及其近亲系）2,3,6。它们的系统发育关系因强烈的进化异速生长效应而变得复杂，这种效应影响了头骨形状的变异，从而模糊了基于整体形状的分类指标，可能是因为物种间的体型演化存在高度平行性7,8,9。欧亚地区的Paradolichopithecus（Depéret，1929年命名）及其同时期可能同义的物种Procynocephalus（Schlosser，1924年命名）10,11,12,13，在疣猴科系统分类讨论中占据关键地位。Paradolichopithecus通常被归入Macacina亚科12,14，但也有学者提出将其归入Papionini亚科15,16。根据其系统发育关系的最终确定结果，这一分类归属可能会颠覆“Papio属物种为非洲特有群体”的共识。Paradolichophecus是欧亚地区化石记录中体型最大的疣猴类成员，其生存时间跨度从中新世中期延续至早更新世，目前已确认有四个物种。欧洲发现的物种包括Paradolichopithecus arvernensis（Depéret，1929年命名），以及Paradolichophecus geticus（Necrasov, Samson, and Radulesco，1961年命名，这些物种也可能为同物异名17,18）。Paradolichophecus sushkini（Trophimov，1977年命名）分布于中亚，Paradolichophecus gansuensis（Qui, Deng, and Wang，2004年命名）发现于中国。这些物种分布于从西班牙到中国的整个亚阿尔卑斯山脉地带10,19,20,21,22。
从头部特征来看，Paradolichopithecus与猕猴属相似10,11,17,19,23；而从身体其他部位来看，它们更接近狒狒（Papio属和Theropithecus属），表现出适应陆地行走的特征10,12,24。Paradolichopithecus的一些头骨特征（如轻微凹陷的额鼻部轮廓、宽大的腭部、较小的上颌窦、弱或缺失的上颌窝以及厚实的上颌骨）使其被归入Macacina亚科10。然而，并非所有Paradolichopithecus标本都具有上颌窦25；同时，某些Papio属和Theropithecus标本也具有上颌窦13。另一方面，原始的Papionina属物种也缺乏上颌窝26，这可能是Paradolichophecus的原始特征之一。Paradolichophecus中的一些类似狒狒的头部特征（如泪骨中的泪腔完全被包裹）可能因该物种的原始状态而在化石研究中显得不那么典型，因此系统分类价值有限10,11,13。不过，在Paradolichophecus sushkini中，据报道存在深陷的下颌窝和明显的上颌窝，以及完全位于泪骨中的泪腔，这些特征被认为表明其与狒狒的亲缘关系更近15,16。
DFN3-150是目前已知的为数不多的完整Paradolichopithecus头骨之一10，是研究其系统发育位置的重要标本。先前的研究表明它具有上颌窦这一特征，但这一特征的系统发育意义仍存疑10,13。对其内耳形态的几何形态测量分析27表明，其形态更接近Papionina亚科的基部分支或更接近Papionina亚族的祖先类型。这些发现与其他采用传统和几何形态测量方法研究Paradolichophecus标本的结果相反，后者认为其与猕猴属关系更密切14，尽管Paradolichophecus体型更大且具有适应性强的陆地移动特征28,29。
Papionina亚科成员的体型变异尤为显著7,9,30,31,32,33，其形态受异速生长关系的强烈影响。因此，Paradolichophecus与狒狒的相似性可能是由于相似的体型和头骨大小导致的平行进化结果27,34。进一步复杂化其解释的是，DFN3-150头骨存在破碎和塑性变形，尤其是右侧受影响更为严重10。虽然这种变形现象明显可见，但变形的不对称性（即哪一侧完整）尚未被详细研究和量化35。
DFN3-150头骨是Thessaloniki Aristotle大学地质与古生物学实验室、法国国家科学研究中心（CNRS）的PALEVOPRIM实验室以及Poitiers大学联合团队在希腊西北部Dafnero下更新世化石遗址发现的化石之一10,36,37,38。该标本被归类为亚成年雌性个体（根据第三臼齿的萌出情况判断为4级成熟度39），属于Paradolichophecus aff. arvernensis。结合Lesvos岛Vatera遗址和Achaia地区Karnezeika遗址的发现22,29，这些标本丰富了希腊地区的Paradolichopithecus化石资料，有助于进一步探讨Paradolichophecus/Procynocephalus支系在猕猴科中的系统分类位置10,11,12,13,15,16,19。

**研究方法**
我们应用了两种不同的虚拟重建技术对DFN3-150的面部结构进行了重建，旨在更准确地评估其形态，并分析变形处理方法对其系统分类框架的影响。首先，我们虚拟去除了头骨中的沉积物；接着将头骨的不同碎片或部分重新组合成完整的结构。随后，我们使用两种手动重建方法制作了头骨面部的模板，并对这两个模板及原始头骨应用了几何形态测量分析：(1) Schlager等人（2018年）提出的反变形方法41；(2) Amano等人（2022年）提出的计算机化修复方法42（其中包含和不包含表面半标记）。最终得到了九个虚拟重建模型。随后，我们通过3D几何形态测量分析将这九个模型与猕猴属和Papionina亚族的形态特征进行比较，以评估DFN3-150的面部特征在系统发育上的位置。

**结果**
应用不同的虚拟重建方法得到的结果存在显著差异（见图1）。Amano等人（2022年）的方法重建的模型具有更宽的吻部，无论是否使用表面半标记，都会使吻部在矢状面和冠状面上显得更宽。这种方法还会缩短眼眶半径。这种差异源于重建过程中需要参考对照标本来修复变形的头骨。相比之下，Schlager等人的方法（2018年）则产生了相反的效果。
主成分分析（PCA）显示，在33个解剖标志点的形状空间分布中，体型对标本的整体排序有显著影响（见图2和补充图S6、补充表S2）。前两个主成分（PC1和PC2）解释了66.9%的变异，将标本分为两大类：PC1主要区分了猕猴属和Papio属，后者位于坐标轴的正端。处于该轴正端的标本具有较长的吻部，使得头骨整体更为狭窄且高度更高。这种延长是由于特定标志点（如鼻尖点、鼻棘点、NPM点和prosthion点）的移动导致的（补充表S2）。相比之下，zygo-max、frontomalare、orbitale和temporal等标志点则朝向头骨后部延伸，方向为内侧和前方。
尽管在不同重建方法下吻部的长度相对稳定，但PC2方法（包含表面半标记）下吻部的矩形轮廓凸显了方法选择的影响。各方法在眼眶和鼻孔形状上的差异进一步强调了这一问题。
对于这种特定标本，即使不存在严格意义上的单侧变形，使用多种重建方法仍然是合理的。尽管如此，值得指出的是，所有根据Amano及其同事（2022年）的方案重建的DFN3-150模型都位于两个晚期亚成年Papio样本的附近，其中一个样本Papio anubis也被纳入了重建方案中。尺寸变量的影响：即使在最初的广义Procrustes对齐（GPA）之后，该对齐方法去除了等尺寸效应，非参数的双尾Kendall’s tau（τ）相关性测试仍然显示前两个主成分分数（PC）与质心尺寸之间存在显著关系（PC1：τ = 0.55，PC2：τ = 0.29，p ≤ 0.01；补充表S12），而多变量离散度的同质性测试未发现不同属之间存在显著差异，结果高于接受阈值（p = 0.534）（补充表S13）。这进一步使我们能够解释这两个属之间的异速生长关系，无论是它们共享一个共同的异速生长轨迹，还是每个属遵循自己独特的形状与尺寸关系。对于比较样本的独特异速生长回归分析显示，尺寸与属之间的统计上不显著的交互项仅解释了1.5%的变异（表1；Rsq = 0.015，p = 0.43）。在随后省略了交互项的共同异速生长回归中，尺寸解释了14.7%的变异（Rsq = 0.147，p < 0.001），而在考虑了尺寸因素后，属解释了10%的变异（Rsq = 0.1，p < 0.001）（表1）。这两个主要效应都非常显著。我们用来评估斜率间角度同质性的向量相关性测试进一步支持了共同异速生长模型（Z-score = 0.55，p = 0.29；补充表S15），表明斜率大致平行。这一点尤其重要，因为尺寸与属之间的重叠并不明显。方差膨胀因子（VIF）为1.86，容忍度为0.55，表明尺寸与属之间存在适度的线性关系，验证了它们各自的预测可靠性。即使在尺寸与属之间存在高相关性时也是如此。然而，需要指出的是，样本的性别、年龄和物种状态并未被明确纳入模型，因此可能会进一步影响残差变异。

分组变量的影响：无论共享的异速生长是否反映了真实的生物模式，还是可能受到了较小样本Papio（N = 7）的影响，稀释分析（补充表S16）得出的中位数p值为0.43。在1000次重复实验中，有38次（3.8%）的交互项p值低于α水平0.05。这表明即使样本量平衡，也不存在显著的属与尺寸交互作用。相比之下，在每个稀释集合内随机排列属标签后的空值排列实验中，中位数p值更高（p = 0.82）。在20000次空值排列中，有134次（0.67%）的p值低于平行性零假设的α水平阈值。在所有重复实验中，交互项的Z分数保持较低且为正值，而空值分布的Z分数较低且为负值（补充图S15）。这表明数据具有连贯的结构，当随机化分组分配时，这种结构完全消失。

然而，这些结果与使用系统发育MANCOVA得到的结果不同（补充表S17），在系统发育模型中只有尺寸预测因子对面部形状有显著影响（p < 0.001），而在考虑了共同进化历史后，属的效应不显著（p = 0.99）。鉴于这些变量在系统发育模型中的交互作用最初是显著的（p = 0.018），GLS稀释分析表明这种效应对属间的样本量差异并不稳健；当组大小平衡时，交互项变得不显著（中位数p = 0.88；补充表S18）。在1000次重复实验中，只有15次（1.5%）的p值低于阈值。在20000次空值排列中，有7197次（35.98%）的测试结果表明，当属预测因子被随机排列时，斜率偏差更有可能发生。这一点反映在Z分数分布上（补充图S15），表明面部形状的属级差异在很大程度上是由系统发育分歧造成的。因此，属被验证为反映面部形状进化关系的预测因子，并提出了共同异速生长轨迹的假设。

共同异速生长分析：尺寸和属效应的差异：除了质心尺寸引起的变异之外，不可忽视的还有属特异性差异（Rsq = 0.1，Z = 3.52，p < 0.001）。这表明进化形状变异与整体系统发育有关。然而，如上所述，没有其他潜在的贡献因素被明确纳入模型，因此这些结果应被视为初步的。共同异速生长回归中拟合值的PC1分数（无交互项）解释了74%的变异，当将其与质心尺寸绘制在一起时，显示Paradolichopithecus重建与Papio有显著且一致的接近性（截距差异较小）。当我们根据质心尺寸绘制拟合值的PC2时，情况则相反，PC2解释了20.3%的变异（补充图S16），这表明该轴捕捉到的是仅由尺寸效应引起的形状差异。这一结论主要是基于这样一个事实：在不同属的拟合值PC2轴上存在重叠，同时尺寸对形状的总体影响仍然明显。

图3显示了共同异速生长回归中拟合值的PC1与质心尺寸的关系图。紫色代表Macaca，橙色代表Papio，棕色代表Paradolichopithecus模型。质心尺寸的差异通过符号的大小来图形化显示。为了更清晰的可视化，各种DFN3-150虚拟重建的缩写遵循图1的编号；字母b、c、d对应于虚拟重建协议，罗马数字I、II、III对应于应用这些协议的模板。当检查形状差异及其静态异速生长成分与质心尺寸的关系时（图4），得出了相同的结论。Macaca物种主要分布在回归分数的负值区域，而Papio和Paradolichophecus则分布在正值区域。此外，属之间的质心尺寸重叠程度，从中等大小的成年Macaca样本到较小的Papio（如亚成年个体或雌性个体），支持了继承形状差异的观点。使用Amano及其同事（2022年）的方案重建的Paradolichophecus模型位于最大猕猴和最小到中等大小的Papio范围内，回归分数有轻微的重叠，而使用Schlager及其同事（2018年）的方案重建的模型则没有显示出显著不同的模式。然而，再次强调，使用Amano及其同事（2022年）的方案重建的模型更接近用于恢复DFN3-150的Papio anubis样本。此外，建议在未来的研究中增加对DFN3-150的性别和年龄范围的更密集采样，以便可能澄清潜在的缺陷，并重新校准整体异速生长轨迹。

图4显示了以质心尺寸为唯一预测因子的回归分数（简单异速生长）与质心尺寸的关系图。紫色代表Macaca物种，橙色代表Papio，棕色代表Paradolichophecus模型。紫色和橙色的实线分别代表Macaca和Papio的个体回归线。线周围是95%置信区间。黑色虚线是共同异速生长回归向量。质心尺寸的差异通过符号的大小来图形化显示。为了更清晰的可视化，各种DFN3-150虚拟重建的缩写遵循图1的编号；字母b、c、d对应于虚拟重建协议，罗马数字I、II、III对应于应用这些协议的模板。其他比较材料的缩写已被省略。亚成年Papio样本通过其缩写标识。

当回归分数向正值移动时，形状变化指的是吻部的整体延长和面部的增高，这一点已有文献记录。这种延长是通过吻突、鼻突和鼻后突的向下投影以及齿槽突点的较小向下移动来实现的。同时，颧骨和额骨上的标志点向内侧移动，使得面部显得狭窄（补充图S17）。然而，当我们进一步划分变异并绘制共同异速生长模型的形状分数时，复杂性变得更加明显（补充图S18）。中等范围质心尺寸之间的形状相似性很强，一些Papio个体和所有Paradolichophecus重建与中型到大型猕猴的分数相当。这是属的叠加效应。在这里，Papio anubis以及DFN3-150被归类为雌性亚成年个体，Papio hamadryas样本指的是性别不明的亚成年个体，而另一个Papio样本是成年雌性个体。显然，年龄作为一个因素并未掩盖整体趋势，但同时也突显了属的效应。然而，在考虑属的情况下，一个进化视角出现了，因为不同属在共同异速生长向量上有所区分。

在考虑了共同异速生长后，每个属在固定比较样本平均质心尺寸下的预测平均形状显示出明显的分离（补充图S19）。在这些预测中保持尺寸协变量不变，确保了组间平均形状差异的唯一来源是属效应。对这些预测的最小二乘均值的主成分分析显示，Macaca、Papio和Paradolichophecus重建在前两个主成分上占据不同的形状空间区域。Papio和Paradolichophecus位于PC1的正值区域，这解释了大部分变异（PC1：76.82%，PC2：23.18%）。每个预测均值周围的95%置信椭圆没有重叠，从而表明独立于尺寸的属级形状差异是稳健的。这些结果确认了可以通过属级差异预测形状成分。

与尺寸相关变化中的明显面部延长不同，属级差异的特点是特定标志点的相对位置发生变化。朝向PC1最大值的方向，鼻突、NPM和鼻后突共同使吻部最前端部分向下移动，而其他面部标志点——如颧骨-上颌点和颧骨-下颌点以及颧弓最大点——则向腹部移动，有效收缩了中部面部。与异速生长趋势相反，M3的近端和远端以及M1-M2标志点都向尾部移动，尽管程度不同，M3的远端也略微向腹部移动；值得注意的是，M3的远端位移向量异常长，进一步导致上颌骨结节更加明显。虽然这些结果本身很重要，但所有亚成年个体（即所有DFN3-150重建和两个Papio样本）的M3终端位置的变化程度尚不清楚，需要更多发育多样化的样本来进一步探索发生效应。尽管如此，所有Macaca和大多数Papio样本都是成年个体，因此当前的极性在某种程度上是有效的。另一方面，颅顶上的标志点主要向头部方向移动，尽管整体变化较小（补充图S20）。

在这里，我们展示了多个DFN3-150颅骨的重建模型，以评估其形态亲缘关系，并探讨重建协议对其感知亲缘关系的影响。这些重建模型基于三个手动重建的模板和两种不同的虚拟协议。尽管这些重建尚未通过实验验证其有效性，但它们都在不同程度上纠正了塑性变形（参见图1）。从质量上评估重建协议的性能来看，Schlager及其同事（2018年）提供的算法不应被视为能够可靠地再现DFN3-150的生前形态。我们应用这一协议作为一种无偏的Taxonomic方法来评估双侧变形。然而，由于其设计特点，它不适用于具有单侧变形的样本。在当前研究中，这种误用在形状空间的PCA图（图2）中更为明显。所有三个模板在使用这种方法进行复古变形后，无论是在形状空间还是在回归分析中，都与其他协议42恢复的模型存在不一致。它们可能保留了一定程度的手动重建模板中的变形信号。这在面部的眶部表现尤为明显，而吻部则受影响较小。相反，Amano及其同事（2022年）提出的协议42引入了与参考标本相当大的变异，这使得其应用值得商榷，尤其是在考虑更广泛的Paradolichopithecus讨论10,16,23时。因此，最保守的模型是使用Amano及其同事（2022年）的协议42生成的模型，这些模型没有考虑表面半标记。这与其他成年Paradolichopithecus头骨的总体观点一致，无论雌雄，上颌在冠状面上的轮廓都呈椭圆形13,46,47。这种成熟过程中的结构在考虑性别差异和年龄后仍然存在，因为性别之间的形态变化在成年之前并不显著32。尽管在进化背景下Paradolichopithecus头骨的大小有所增加14，这一事实仍然成立。后一种方法实现了对上颅的精确和可控的恢复，专注于眼眶的形态特征，同时保留了吻部的原始形状。当前研究未明确解决的另一个缺点是关于DFN3-150的发育年龄组的问题。在进化背景下使用成年标本可以解释完整的表型表达33。为了克服这一限制，我们建议从现存的近缘物种中实施一个发育生态量纲向量，可能包括额外的系统发育信息48，同时结合现有的Paradolichopithecus成年形态。通过这种方式，或许可以模拟年轻个体的成年形态，反之亦然。同时，可以使用更严格的系统发育和地质年代学约束，在“正向发育生态投影”49中对DFN3-150进行重建。这种方法还可以专注于头骨的特定区域，反映与参考样本相对应的区域变化。这样我们可以避免对目标物种的成熟形态进行广泛的推断。一种“缩小比例”的区域方法也可能足以区分系统发育信号50，超越同形效应。

为了测试这些差异对DFN3-150重建模型感知的形态亲缘关系的影响，我们在Macaca和Papio的比较框架中探索了整体形态相似性。虽然Paradolichopithecus-Macacina姐妹群关系仍然是Paradolichopithecus系统学中的主导假说10,14,16（及其中的参考文献），但尚未完全排除狒狒在讨论中的位置15,24,27。为了探讨这一点，我们初步检验了Paradolichopithecus群是从欧洲的Macaca祖先种群原地进化的假设14，通过展示尺寸增加的趋势。鉴于质心大小与体重之间的高相关性，至少对于Papionina来说7,30,31,33,51，生态量纲在Papionina的形态变化中起着重要作用。因此，尺寸在很大程度上反映了当前提供的比较物种之间的进化关系。考虑到质心大小的进化重要性和影响Papionini头骨形态的生态量纲关系30,31，当将形态变化建模为尺寸和属的严格函数时，我们的数据显示DFN3-150的重建模型分为两组（图3）：遵循Schlager及其同事（2018年）协议41的重建模型（质心尺寸较大）和根据Amano及其同事（2022年）协议42重建的模型（质心尺寸较小）。尽管如此，DFN3-150模型仍然沿着我们Papio样本的生态量纲轨迹。这表明DFN3-150的重建模型在形态上更接近Papio，而不仅仅是生态量纲的结果。这个假设并不是直接的进化推断，而是一个需要用更具包容性的比较样本（在发育年龄、性别、属和种类的方面，特别是Papionina）进一步检验的假设。为了进一步实现关于生态量纲对Papionina分类亲缘关系影响的研究的可比性，需要一个结合形态测量和非测量数据的综合数据集。此外，我们的研究中没有考虑Papio的发育生态量纲效应33以及整个比较样本中的性别二态性。还需要注意的是，只有两个比较标本代表了接近成年形态的亚成年个体。同样，DFN3-15010也被描述为一个亚成年个体。虽然预期青少年个体的性别二态性不会有显著差异，但我们有意没有按性别合并属，以保持尽可能多的变化。我们的选择基于这样一个事实：归因于性的种内变异是重叠的发育生态量纲的产物52,53,54，而性别之间的差异在后来发育阶段表现为额外的次要轨迹偏差53。因此，在族水平上，属内的明显差异预计不会影响进化生态量纲推断54。此外，Papio的雄性体型范围的扩展预计将进一步增加形状空间的表型变异性，从而允许适当的整体生态量纲轨迹校准。这在比较样本相对有限的情况下尤为重要，正如目前的Papio样本一样。DFN3-150重建模型在多因素回归模型中接近Papio，或者至少它们与Macaca的截距距离较近，这表明它在任何给定尺寸下都与Papio有更密切的形态亲缘关系。特别是考虑到两个属共享的发育轨迹30，这一点尤为重要。后者在图4中得到了体现，我们有限的亚成年Papio样本接近成年猕猴的形态。最大猕猴的质心尺寸，包括代表这一亚族的祖先状态的保守型Macaca sylvanus55,56，与其较小的Papionina成员的质心尺寸有重叠。因此，可以假设DFN3-150可能已经达到了超出成年Macaca的分数范围。特别是当我们考虑到较大的尺寸预计会在Papionini中引入新的形状空间30,31时，这一过程也受到雌性Papio生物学年龄的影响51。因此，亚成年Papio标本中仍在萌出的M3牙齿以及可能的DFN3-150，预计将进一步改变面部形状。这种变化与相应的生态量纲轨迹一致，但保持了M3牙齿相对于颞下颌关节的位置57。无论如何，发育生态量纲可能并不表明同源性，而且之前已经有记录显示Papionini之间存在发育上的差异54。在本研究中，我们确定远端M3牙齿是两个被研究属之间变化最大的标志，显示出不同的方向，与静态的生态量纲效应相反（补充图S20）。

在这项研究中，我们有意排除了其他影响，如性别、年龄以及残余表型表达的更广泛的生态和功能方面31,33,51，以避免生态量纲信号的过度模糊。这可以在共同的生态量纲回归形状分数（补充图S18）中看到，个体在某个范围内变得更加混合，因为变异进一步分离，突出了属之间的形状差异。然而，这种限制是主要问题，使得在更广泛的进化视角下对结果的任何系统评估都是初步的。随着质心大小的增加，生态量纲主要影响的区域是上颌和颧骨前突的相对收缩。这扩展到额骨颧突的根部，较大的标本在前眶区有轻微的凹陷。与这些变化无关的形状信号可能来自其他因素，可能反映了广义上的系统发育关系。无论如何，生态量纲控制的平均形状差异预测突显了DFN3-150与Papio物种的接近性，进而可能与Papionina有关；应进一步评估这一系统学替代方案。额缝（与样本的牙科组特征一致）或其他骨折均未发现；因此，这种结构很可能是由于压缩作用造成的，而非任何解剖学特征所致。额骨在最后右侧的压缩区域被分成了两个部分。最初对所有切片进行分割后得到了36个虚拟段落，经过筛选过程后保留了25个（见图5）。剩余的11个小碎片要么被整合到相邻的较大段落中，要么由于操作困难而被完全省略。面颅由5个部分组成，而脑颅由20个部分组成，其中一部分单独构成了蝶鞍。有些碎片没有对应的另一侧结构，例如鼻骨碎片和突出的头顶骨，它们既不属于额骨部分也不属于顶骨部分。图5 这张图像的替代文本可能是通过AI生成的。全尺寸图像 分割后的颅骨：A) 额骨；B) 尾骨；C) 左侧；D) 右侧；E) 背部；F) 腹侧。未缩放。手动重建 原始标本在背侧和侧面受到了压缩，尤其是右侧受到了明显影响。这种侧向压缩主要表现在两个区域：首先，右侧上颌骨腭突（没有明显的横向腭缝，因此与腭骨合并）重叠在其左侧对应部分上，大致沿着中线腭缝的方向（见图5；补充图S4）；其次，右侧上颌骨的整个额突与前鼻骨碎片重叠（见图5；补充图S4）。这种压缩导致右侧半颅骨相对于左侧半颅骨垂直抬高。我们根据段落间的平滑度标准来验证这些重建的准确性。手动重建方面，我们分别制作了两个模型来校正这些变形：“wideNasal”模型用于校正右侧上颌骨对鼻骨的覆盖（见图1 IIa），“widePalate”模型用于校正右侧到左侧的腭骨覆盖（见图1 IIIa；有关本研究整个流程的概述，请参见补充文件SOM2）。我们将这些模型称为“解压模板”，因为它们看起来比原始标本宽度更大。这些模型实现了颅骨的侧向解压，除了原始化石外，还将作为后续数字重建的模板。对于“widePalate”模板，我们手动平移和旋转了右侧上颌骨段落，使其与左侧对应部分对齐。这一过程通过估计的腭骨水平板远端部分的宽度来控制。这条边的两个相对部分应该被柱状突（Staphylion）地标所分割，该地标位于腭骨水平板后缘切线的中点。当标记腭骨后端的线变得与通过柱状突地标的腭缝垂直时，变换停止。然后，将右侧上颌骨段落的平移和旋转值应用于右侧颅骨的所有其他段落，确保它们整体上平滑地跟随变换，假设只有同一组的均匀向量同时影响右侧各段落。我们对“wideNasal”重建也采用了相同的方法。目标是将右侧上颌骨段落到其正确的解剖位置，重点在于重现其额突的连续轮廓，优先考虑过渡的平滑度而非两部分之间的融合程度。接着，我们将相同的平移和旋转值应用于其他两个面颅段落以及右侧的所有其他碎片，确保它们之间的初始相对位置保持一致。网格清理 在手动重建颅骨的过程中，重建的模板表面网格留下了两个主要缺口（参见补充图S5）。第一个缺口出现在额骨的右侧，大致从颧弓（CZA）地标开始。这个缺口沿着整个额骨延伸，向后到达顶骨之前的位置。这是由于将额骨分成两个独立的部分以抵消上述所说的压缩作用造成的。第二个缺口在调整后的腭骨模型中更为明显，但在“wideNasal”重建中仍然可见。这个缺口是由于右侧上颌骨碎片的移动以及该区域的保存状态造成的。上颌骨缺少朝向鼻部的额突，而这个额突应该与颧骨的额突一起形成一个平滑的圆形轮廓，与眼眶边缘相接。这使得该区域相对于更完整的左侧不完整。我们采用了两种基于地标的方法来填补这两个缺口。对于额骨缺口，我们在边缘沿线采样了两组地标——”wideNasal”模板共97个，“widePalate”模板共69个——这些地标在边缘两侧大致等距离放置（见补充图S5）。基于这些地标，计算出了凸包多边形表面网格。我们选择了凸包拓扑结构而不是Delaunay三角剖分方法，因为它能均匀地填补缺口，避免多余或凹陷，同时确保缺口两侧的高度一致。通过重新网格化（Rvcg::vcgIsotropicRmeshing61）去除了不规则之处，该方法根据定义的结构属性重新计算虚拟表面的三角剖分，从而使凸包网格的分辨率与虚拟标本表面模型的平均值相匹配。然后，将重新网格化的凸包与两种手动重建的表面模型合并。对于鼻骨缺口，目标不仅是填补上颌骨移动留下的缺口，还要重新创建圆形轮廓，重建右侧上颌骨缺失的额突。为此，我们在右侧上颌骨的纵軸上沿缺口采样了三组六个地标。第一组从颧上缘（ZMU）地标开始，沿着上颌骨下边界延伸至鼻前上颌边缘（NPM）地标；第二组覆盖了沿上颌骨近端边界的鼻上颌缝长度；第三组包含八个地标，从鼻骨右侧边缘开始，延伸至鼻骨的前部。然后，将每组地标转换为40个等距点（见补充图S5），沿多项式贝塞尔曲线（bezier::pointsOnBezier62）排列。贝塞尔曲线是根据控制点插值函数定义的，其中第一个和最后一个地标定义了样条的曲率。接下来，我们反转这些新地标，创建了40对垂直于鼻骨缺口的三个地标。使用这些新点，再次定义了贝塞尔曲线，并沿着它们采样了70个等间距的地标。最后，我们将这套新的地标导出并在AVIZO中处理成点云，应用Delaunay三角剖分算法，并进一步处理渲染后的表面，采用与额骨缺口处理相同的步骤。这些协议共同帮助我们重新创建了缺失的轮廓，成功填补了鼻骨和额骨区域的缺口。Delaunay三角剖分算法能够对一组点进行镶嵌，满足Delaunay条件，即没有一个点属于算法生成的三角形的外围超球面。地标采样 对于右侧段落，我们使用了与左侧段落相同的变换值来采样解剖学地标。首先，我们从原始模型中采样了地标（见补充表S2）。然后，我们将上述用于调整右侧段落的变换值应用到每个重建模型的相应地标上，以确保一致性并最小化误差。所有重建模型中的左侧地标保持不变。我们使用相同的坐标系统来保持一致性（无论是“widePalate”模型还是“wideNasal”模型）。虚拟重建 我们采用了Schlager等人2018年描述的逆变形协议41。我们利用了最近提出的表面配准协议65,66,67来获取表面半地标。由于右侧颅骨的变形，使用这种表面配准方法比仅依赖标准方法（如薄板样条（TPS）变形68）更能有效地将表面半地标块从左侧转移到右侧。随后，我们应用了Schlager及其同事（2018年）提出的逆变形协议41（有关该方法的批评，请参见69）。表面配准 这种方法在单层网格上效果最佳，因为它可以避免迭代最近点（ICP）算法与内部或孤立表面的后续不匹配66,67。因此，提取了原始DFN3-150标本的单层三角形网格。左侧面颅部分——从额颧缝到下牙槽缘，以及从颧颞缝到鼻前上颌边缘（NPM）和颧上缘（ZMU）地标（见补充图S6, S7）——在整个表面配准过程中被视为参考，因为在手动重建中这部分保持不变。相比之下，原始模型和两种手动重建的右侧单层网格被视为目标，因为它们的位置不相同。然后我们数字化并等距重新网格化了这些网格，以确保它们具有相同数量的面和分辨率，从而提高了配准算法的性能。在左侧参考网格上，我们使用K均值聚类算法（Mopho::fastKmeans）在网格顶点上采样了120个表面半地标。在这些准备之后，参考网格通过基于先前收集的“锚点”地标和35个等距曲线半地标，使用TPS变形进行镜像和大致对齐。接着，我们应用弹性迭代最近点（ICP）算法来优化参考网格和目标网格之间的配准。这里的“弹性”是指包含高斯平滑位移向量66,67成分。一旦两个网格的顶点被最佳对齐，我们就将左侧参考网格的半地标转移到右侧目标网格上，使用最近匹配方法（Morpho::transferPoint和Morpho::vert2points61）。这个过程产生了240个表面半地标（见补充图S7），覆盖了每个模板的面颅（原始模型和两个手动重建）。逆变形 我们应用了Schlager及其同事开发的“无分类法”协议41。这种工作流程使用Ghosh及其同事设计的“封闭形式解决方案”来对称化双侧变形的物体71，通过对称性来逆变形变形的标本35,72。与原始实现不同，我们在这里引入了半地标，而不仅仅是依赖解剖学定义的地标68。采样半地标可以增加对研究结构的覆盖范围，解决了不完整标本中缺乏解剖学地标的常见问题。为了应用这种方法，我们在标本的八条曲线上采样了34个双侧解剖学地标和88个半地标。我们还使用了240个表面半地标（见补充图S6, S7和补充表S2, S3）。首先，基于地标级别计算对称性（Morpho::retroDeform3d61），然后将其应用于表面模型（Morpho::retroDeformMesh61）。根据Amano等人2018年的方案进行修复42 为了研究不同重建方法对DFN3-150形态的影响，我们还应用了Amano及其同事（2022年）开发的修复方案42。该方案利用地标来描述标本的形状及其变异。其设计旨在利用分类学上的种间形态变异作为参考资料。然后将变形的标本垂直投影到一个表示这种变异的主成分框架超平面上，从而消除标本中的任何变形，使其恢复到与参考颅骨超平面的匹配状态。我们选择了三个雌性参考颅骨，其中两个属于亚成年个体。这些参考标本用于修复DFN3-150的面颅形状。第一步。我们为参考标本数字化了一组具有里程碑意义的62个双边解剖学标志点和124个等距曲线半标志点，以及DFN3-150模板（补充图S6、S8；补充表S1、S2、S3）。对于参考标本Mandrillus sphinx PRICT70，使用TPS算法根据其他两个参考头骨统计重建了P3颊尖上缺失的标志点。同样的方法也被用来计算“original”模板的额突（glabella）和鼻根（nasion），以及“wideNasal”和“widePalate”模板的额突、鼻根和鼻棘（nasion and nasospinale）（Morpho::fixLMtps61）。步骤2：为了最小化由参考头骨结构不对称性引起的形态变异，我们使用反射重标记（Morpho::symmetrize61）校正了偏离对称性的部分，并随后使用TPS对网格进行变形（Morpho::tps3d61），从而为每个标本创建了一个完全对称的形状。由于这种方法的位移值更小，因此比简单地沿着中矢状面镜像更为优越。剩余的修复工作流程是在有和没有表面半标志点的情况下进行的，从而得到了“original”、“wideNasal”和“widePalate”模板的两个不同三元组。通过这种方式，我们考虑了参考头骨对后续步骤可能产生的影响，因为修复过程将遵循这些头骨的变化。

对于表面半标志点，我们通过对Papio anubis PRICT733参考标本的左侧面颅表面网格顶点进行K均值聚类（Morpho::fastKmeans61）来创建了140个半标志点。然后，这些半标志点被使用之前反射和重标记的解剖学标志点建立的中矢状面映射到右侧。为了防止标志点在网格内部错位，我们使用R语言中的一个脚本沿表面法线膨胀标志点，再将它们投影回表面（代码由Costantino Buzi博士和Antonio Profico教授提供）。得到的表面半标志点补丁（280个点）通过TPS变形（Morpho::placePatch61） transferred到了其他比较样本中。总的来说，没有表面半标志点的标志点配置共包含186个解剖学标志点和曲线半标志点，而添加了280个表面半标志点的配置则共有466个点。在这两种情况下，都使用了前三个主成分（PC）得分42，这代表了100%的种间变异。

**误差计算**
我们评估了标志点放置的观察者内误差，这影响了表面注册工作流程的初步步骤（参见材料和方法：I. Schlager等人2018年的协议41）。因此，SK在“wideNasal”模板的表面模型上连续三次数字化所有解剖学和曲线半标志点，总时间为三周。通过在每个标志点位置首先计算到配置中心的欧几里得距离，然后将这些距离从重复实验的平均配置距离中减去来评估误差。后者以毫米（mm）和百分比的形式表示出来（补充表S6）。为了进一步量化标志点放置误差对实际表面注册的影响，我们计算了每个标志点配置试验中所有顶点到生成表面的欧几里得距离，并将其与目标网格表面进行比较（顶点到表面的距离）。我们量化了三个距离阈值的百分比。然后我们计算了欧几里得距离的标准差（Sd）（补充表S7）。接下来，我们计算了黎曼距离（RD）33（也称为Procrustes距离；在预形状空间上两个形状之间的测地距离——一种形状接近性的度量——即部分GPA后的结果空间）73（补充表S7）。该算法成功注册了表面，主要偏差发生在应用区域的边界（补充图S10），正如其他研究66,67所报告的那样。因此，“锚定”标志点的采样并未影响表面注册协议生成的表面半标志点的获取。

**统计分析的测量误差**
我们计算了标志点放置的观察者内误差对统计分析的影响。然而，在这里也考虑了比较样本中数字网格分辨率的影响。SK在3个标本上重复采样了33个解剖学标志点（补充表S1、S2），每个标本代表逐渐降低的分辨率状态，时间间隔分别为：t1 = 0小时，t2 = 1小时后，t3 = 2天后，t4 = 7天后。我们将测量误差表示为每个物种每个标志点的欧几里得距离（ED）与其平均值的标准差（Sd），考虑了四次重复实验。为了突出不同方法在计算此指标（ED）时的不一致性，它是基于原始标志点位置和Boas坐标计算的（补充图S11）。在初步的Procrustes注册之后，我们将Procrustes形状变量乘以它们的质心大小测量值来计算Boas坐标。这逆转了缩放的效果，恢复了配置的原始大小，同时保持了平移和旋转75。这个新空间近似了构型空间76。此外，我们还计算了每个物种试验的成对黎曼距离（RD）（补充表S11）。大多数标志点的欧几里得距离的标准差（Sd）保持在1毫米以下，而在生成Boas坐标过程中的部分注册对这一指标没有重大影响。值得注意的是，在最低分辨率表面上采样的两个zygo-max上部标志点的Sd值超过了2毫米的阈值。其中，t4时期的采样重复性似乎最为变化较大，这从重复实验之间的黎曼距离可以看出。zygo-max上部的Sd值超过1毫米的阈值也适用于最高分辨率的表面模型，仅针对左侧的双侧标志点。

我们使用一种称为双向Procrustes ANOVA（MANOVA）77的程序来评估标志点放置的重复性。在这种设置中，我们将标志点坐标视为因变量，物种和重复次数作为因素。我们使用np-MANOVA的平均平方值计算了相对重复性作为类内相关系数（ICC），范围从零（0）到一（1）74,78。ICC量化了不同组内个体重复测量值的相对重复性。然后计算了欧几里得距离的标准差（Sd）作为测量误差（ME）（补充表S11）。在分析之前，我们对所有标志点配置应用了部分Procrustes注册（Morpho::procSym61），并提取了它们的对称成分79。随后的回归使用了10,000次排列和II型平方和进行残差随机化（RRPP::lm.rrpp.ws80,81,82,83,84）。方差分析显示物种间的形状变异显著，但在采样重复内没有变异（补充表S8）。Rsq值为0.98（p = 0.0002），相对ME为0.35%。我们使用阻塞排列的MANOVA进一步确认了试验效应对标本形状的无显著性。我们在物种因素中限制了10,000次排列，考虑了试验内的非独立性（非交换项），并使用了I型平方和（vegan::adonis285），结果是重复实验的p值为p = 0.97（补充表S9）。对于这两种排列方法，我们都测试了每个物种和试验的标志点配置的方差多变性（vegan::betaspider和vegan::permutest，共10,000次排列85），发现其远高于可接受的α≤0.05水平（补充表S10）。这些额外的误差测量步骤在评估测量误差时至关重要，因为在评估测量误差时没有标准化的阈值。

**数据处理**
鉴于重建的确定性和引入的相关误差，我们决定仅使用在面部骨骼和额骨前部采样的33个解剖学标志点（补充图S6；补充表S2）来研究DFN3-150的面颅亲缘关系。这些标志点属于I型和II型，除了frontomalare temporale属于III型标志点86。由于缺乏骨骼表面，我们使用TPS算法从剩余的比较材料中估计了DFN3-150反变形模型中的鼻棘（nasion）和prosthion标志点，但基于Amano及其同事2022年的恢复协议恢复的两个原始化石模型中的鼻棘除外42。在下面描述的工作流程中，我们不允许这些计算出的标志点发生滑动。一些解剖学标志点（见补充表S2）被允许像表面半标志点一样滑动（见87），以最小化弯曲能量标准88，从而补偿由于网格分辨率低导致的放置不确定性，并确保某些特征（如闭合缝合线）的准确识别。提取调整后的标志点后，我们使用反射重标记程序（Morpho::symmetrize61）对所有比较样本配置进行了对称化。这些标志点配置作为本研究的主要变量。紧接着，我们将所有比较材料的相应网格转换而来，以适应对称化的标志点配置（Morpho::tps3d61）。这些配置和网格，连同最初的DFN3-150模型的标志点配置，将作为本研究的主要形式变量。为了评估标志点数量是否足以捕捉样本间的形态变异，我们应用了LaMBDA::lasec函数89（补充图S12），该函数通过计算每个子集与完整配置之间的Procrustes平方和（PSS）来量化逐渐增加的标志点数量对变异的贡献。对于当前的样本，即使只有30个标志点也能够达到良好的拟合效果。比较材料的标志点配置最初经过了部分Procrustes注册对齐（GPA），这会将标本缩放到单位质心大小并消除位置和方向的差异（讨论见73,90）。标志点配置的质心大小计算为所有标志点到其质心的距离平方和的平方根73,76。之后，我们计算了所有DFN3-150模型的质心大小（Morpho::Csize61），并将其与比较样本的GPA对齐（Morpho::align2procSym61），以确保比较样本不会受到化石标本包含在其平均值计算中的影响，从而得到一个无偏的参考样本。最后，Procrustes坐标被投影到形状切线空间，这是Kendall形状空间的线性近似——两个空间之间的相关性约为0.999（Morpho::regdist61）。

**形状变异和异速生长**
我们将切线空间坐标提交给了主成分分析（PCA），其中化石标本被投影到主成分空间40。在评估了有意义的PC得分的数量（Morpho::getMeaningfulPCs61）之后，我们进行了配对的双尾非参数Kendalls’ tau（τ）44（stats::cor.test91）测试，以研究这些PC轴上观察到的分离是否包含了所有速生效应。通过这个测试，我们检查了每个物种的前三个PC得分与物种质心大小之间的关联。由于前两个PC得分（74.33%的累积比例）和质心大小（MVN::mvn92 bootstrap‘能量距离’多变量正态性测试：E统计量：1.27；p < 0.05）没有多变量正态分布，因此需要进行非参数测试。在分析所有速生效应之前，我们首先使用10,000次排列对每个属组的形状变量进行了多变量方差同质性测试（group dispersion），并使用了空间中位数分析类型（vegan::betaspider和vegan::permutest85）。此外，我们通过进行了Wilcoxon-Mann-Whitney非参数测试93（coin::wilcox_test函数，基于10,000次迭代）来评估尺寸是否独立于属，这指导我们选择了II型或III型平方和进行异速生长回归。我们还计算了它们关联的大小效应大小，使用Wilcoxon测试（rstatix::kruskal_effsize94），得到了一个较大的值。这种非参数测试更受欢迎，因为两个组之间的尺寸方差同质性假设没有得到满足，如Barletts测试（stats::barlett.test91）所示，而通过Shapiro-Wilk测试（stats::shapiro.test91）也表明正态性假设没有违反。为了评估进化速生的统计显著性，我们进行了一系列多变量回归，研究了尺寸和属的回归（RRPP::lm.rrpp.ws80,81,82,83,84，共10,000次排列）。这些额外的误差测量步骤至关重要，因为在评估测量误差时没有标准的阈值，无论试验或网格分辨率如何。最终，我们进行了一项独特的异速生长回归分析，其中包含了交互项（形状 ~ 大小 * 属），并使用第三类平方和（Type III sums of squares）来适当地划分不平衡数据中的方差，以处理交互作用的影响95,96。为了比较异速生长回归中各个预测因子的显著性，我们在系统发育背景下使用MMCOVA（RRPP::lm.rrpp.ws80,81,82,83,84）对比较样本进行了分析，采用系统发育树拓扑结构作为协方差矩阵。作为参考，我们使用了基于本研究中包含的12个物种的两个属的共识时间树的系统发育图（Supplementary Figure S13；https://10ktrees.nunn-lab.org/index.html）。为了检查预测因子之间的共线性（即它们之间的重叠程度），我们使用了方差膨胀因子（VIF）及其倒数（即容忍度）98,99。由于Papio标本数量有限，我们通过对独特异速生长模型进行稀释测试来评估其结果的稳健性，在该测试中，我们通过交互项的统计显著性（F统计量）来评估角度的相似性（同质性）100。这项分析的R代码改编自Ariel Marcy101。我们在普通最小二乘（OLS）和广义最小二乘（GLS）回归模型中应用了这一测试，并使用了系统发育协方差矩阵。实际上，我们是在平衡样本量的情况下检验不同属是否遵循各自的形状-大小关系（不同的斜率）。从这个意义上说，稀释是一种欠采样方法。为此，我们进行了1000次稀释迭代：每次随机选择7个Macaca标本（与Papio的样本量相匹配），将它们与所有Papio标本合并，然后重新计算独特异速生长回归。对于每一组稀释后的样本，我们还通过随机排列这些14个稀释样本的属标签来生成一个零分布，并重新计算20次独特异速生长回归（每次迭代都进行分组）。通过这种方式，我们人为地破坏了将属与大小和形状联系起来的生物学意义。这种方法使我们能够直接评估观察到的属间异速生长斜率差异是否可以仅由抽样效应解释，同时也能检验在属与大小和形状对应关系完全破坏的情况下是否仍然如此。

虽然我们可以使用RRPP::pairwise80,81,82,83函数来评估斜率的相关性，但在我们的分析中只有两个组的情况下，这种选择是多余的。然而，由于方法论的限制82,83，我们无法将相同的测试应用于系统发育回归。对于这两种回归，我们都报告了F统计量的p值和Z分数。Z分数作为效应大小的度量，表示观察到的F统计量与使用RRPP包生成的残差排列分布的标准差之间的差异83,102,103。为了进一步区分大小信号与属相关头骨形状的差异，我们将形状变量建模并预测为属和常数质心大小的函数，该常数质心大小对应于完整比较样本的平均质心大小。这使我们能够在评估属级形状变异时控制异速生长效应。每个属的预测平均形状及其相关的95%置信区间被投影到形状空间的主成分轴上。这种方法能够在保持大小不变的情况下，评估特定属的形状亲缘关系，从而明确那些不能归因于大小相关异速变化的分类学差异。为此，我们使用了stats::predict函数91。

总体而言，本分析中的异速生长关系通过预测线和回归分数进行了可视化104,105。回归分数表示“与回归模型预测的形状变化相关的形状变量，但也包括了该方向上的残差变异”（104；第72页）。为了直观展示每个标志点处的变异以及变化的方向，我们根据landvR::coordinates.difference106、landvR::Procrustes.var.plot106函数、geomorph::shape.predictor80,81,82,83和Morpho::tps3d61等函数生成了弯曲的网格。数据分割、数据收集和Delaunay三角剖分是在Avizo 8.0.0、Avizo 9.2.0 Lite和AVIZO 3D 2023.1.1软件中完成的（?FEI Visualization Sciences Group）。图表是通过rgl（v1.3.31）和ggplo2（v4.0.1）包生成的91。表面半标志点使用fastKmeans函数采样，曲线半标志点则通过Morpho61（v2.13）包的equidistantCurve函数使得等距分布。这些半标志点使用Morpho的relaxLM和slider3d函数进行了滑动调整。网格的清理由Rvcg61包的vcgIsotropicRemeshing函数完成。贝塞尔曲线上的点通过Bezier62（v1.1.2）包的pointOnBezier函数进行采样。在表面配准步骤中，我们使用了Morpho的vert2points、fastKmeans和transferPoints函数，以及mesheR107（v0.4.200213）包的gausMatch函数。DFN3-150模型的反变形使用了Morpho的retroDeform3d和retroDeformMesh函数完成，而在计算机化修复过程中，我们使用了Morpho的fixLMtps、symmetrize、tps3d和placePatch函数。数据处理使用了Morpho的procSym、Csize和alignProcSym函数，以及LaMBDA89（v0.1.10000）包的lasec函数。数据的整体展示使用了landvR106（v0.5.3）包的coordinates.difference和Procrustes.var.plot函数，geomorph80,81,82,83包的shape.predictor函数（v4.0.10）和Morpho的tps3d函数。研究的准备工作包括使用Morpho的getMeaningfulPCs函数、stats91（v4.5.1）包的cor.test、barlett.test、shapiro.test函数，以及MVN92（v6.2）、coin93（v1.4-3）和rstatix94（v0.7.3）包的mvn、wilcox_test和krustal_effesize函数进行统计分析。对于异速生长回归框架，我们使用了RRPP80,81,82,83包（v2.1.2.999）的lm.rrpp.ws和pairwise函数；对于形状预测，使用了stats包的predict函数。统计分析的测量误差通过RRPP的lm.rrpp.ws函数以及vegan85（v2-7.2）包的adonis2、betaspider和permutest函数进行评估。