解剖几何结构与纤维分布对跟腱力学特性的影响：基于有限元建模的方法

《Journal of Biomechanics》：Influence of anatomical geometry and fiber dispersion on Achilles tendon Mechanics: A finite element modelling approach

　　yeswanth S Pydi | Maria Pierantoni | Anna Gustafsson | Hanna Isaksson
瑞典隆德大学生物医学工程系

摘要
肌腱会根据机械负荷调整其结构和机械性能，其中机械刺激（如应变）被认为起着关键作用。最近使用有限元建模的研究提供了关于不同负荷条件下肌腱结构中应变分布的详细见解。然而，大多数研究人员采用了高度简化的几何模型，这些模型无法捕捉肌腱中的固有应变异质性，从而强调了开发解剖学上准确的模型的必要性。
本研究通过使用大鼠跟腱的解剖学详细有限元模型，探讨了几何形状和纤维方向如何影响肌腱的应变异质性。为此，利用相位对比增强同步辐射X射线断层扫描图像来获取肌腱的几何形状和纤维方向。分析了基于解剖学分割的两种不同类型的几何模型以及具有直线纤维分布的简化圆柱模型，以了解几何形状在循环负荷下的影响。还对使用平均全局和局部映射纤维分布的解剖模型进行了相同的参数分析，以了解纤维方向和分布的影响。
结果表明，几何形状和纤维方向显著影响肌腱的机械性能。解剖学模型显示出更大的应变异质性，最大主应变分别为（0.070 ± 0.020）与圆柱模型（0.058 ± 0.015）相比，而包括全局和局部分散的纤维方向后，在负荷过程中从肌腱中挤出的流体速度增加。这可能对肌腱功能的力学生物学预测具有重要意义。

1. 引言
跟腱是人体最大的肌腱，在将小腿肌肉的收缩力传递到跟骨方面起着重要作用（Wang, 2006）。肌腱具有机械敏感性，过高的或过低的机械刺激都会损害其生物力学特性，导致刚度下降和肌腱病变（Ottani et al., 2002, Screen, 2008）。通过特定的康复锻炼施加的机械刺激旨在引起肌腱的应变，从而促进其适应和恢复（Bohm et al., 2015, Escriche-Escuder et al., 2022）。因此，理解肌腱复杂的力学特性及其所经历的应变区域变化对于理解其适应性至关重要。最近关于人类肌腱的实验研究表明，在初始拉伸负荷期间，肌腱近端区域会出现高横向应变（Nagelli et al., 2022）。然而，在峰值负荷下，肌腱中部显示出最高的纵向应变。这说明了肌腱在变形过程中复杂的非均匀应变分布。

表1. 优化过程中保持不变的基质本构参数（Notermans et al., 2019）

有限元模型（FE）可以提供有关肌腱在各种负荷条件下的机械性能的关键见解（Funaro et al., 2025, Kinugasa et al., 2016）。像大鼠这样的小型动物模型对于研究力学生物学原理至关重要（Fang and Lake, 2016, Ganji et al., 2023, Wu et al., 2017）。现有文献探讨了大鼠跟腱在循环拉伸负荷下的各向异性力学性能（Khayyeri et al., 2016），以及强调了应变场的速率依赖性粘弹性响应（Notermans et al., 2019）。此外，研究还强调了异质应变场直接影响细胞反应和胶原蛋白生成的重要性（Notermans et al., 2021）。最近，我们的团队研究了跟腱中的胶原纤维如何适应体内机械负荷（Janssen et al., 2025）。利用同步辐射X射线断层扫描数据构建了FE模型，这些模型采用了全局（所有元素中纤维方向相同）和局部（纤维方向空间变化）的方法（Janssen et al., 2025, Pierantoni et al., 2023）。纤维方向的变化可以解释体内未负荷肌腱的机械响应受损现象。
然而，大多数这些研究使用的跟腱解剖模型过于简化，通常假设为圆柱形。因此，先前的研究未能捕捉到解剖学上的复杂性和可能的应变异质性，从而忽略了其真实的动态特性。因此，为了为通常由空间变化的应变刺激驱动的计算力学生物学研究铺平道路，有必要开发能够考虑整体几何变化的解剖学特定肌腱FE模型。
本研究的目的是评估特定受试者的解剖学几何形状和纤维方向对大鼠肌腱机械性能的影响。我们团队之前的研究已经对大鼠肌腱进行了实验表征并分析了其纤维方向（Pierantoni et al., 2023）。在本研究中，使用了（Pierantoni et al., 2023）中的相位对比同步辐射X射线断层扫描数据来开发解剖学FE模型。此外，还在FE模型中实现了特定样本的纤维方向和分布，既采用全局平均分布，也基于（Janssen et al., 2025）的元素特定分布的局部映射，以评估它们在应变分布和流体速度中的作用。

2. 材料与方法
本研究使用来自三只大鼠的相位对比同步辐射X射线断层扫描图像，开发了简化的圆柱形和更复杂的解剖学几何模型，并将其平均机械响应与实验数据进行比较。解剖学模型和圆柱模型都假设胶原纤维方向是直的。随后，根据实验数据分析了具有直线和分散胶原纤维方向的解剖学模型。

2.1. 大鼠肌腱的体外测试
在之前的一项研究中，使用瑞士保罗谢勒研究所（Villigen）的TOMCAT光束线（SLS）的相位对比增强同步辐射断层扫描技术，对10-12周龄雌性Sprague-Dawley大鼠的跟腱进行了成像（Pierantoni et al., 2023）。扫描的视场为4.2 × 3.5 mm，各向同性体素大小为1.63 μm。在本研究中，我们使用来自三只完全负荷大鼠肌腱的图像来获取特定样本的几何形状和空间胶原纤维分布。之前的研究还使用配备Interface SMT1-250 N载荷传感器的Instron 8511加载框架对另一组九只大鼠的肌腱进行了拉伸机械负荷（Pierantoni et al., 2023）。每个肌腱 preloading 到1 N，然后施加6%的循环拉伸应变，加载率为0.1 mm/s，持续10个周期。在本研究中，我们使用了该实验数据集的前3个周期。

2.2. 图像分析、分割和肌腱FE模型
通过分割高分辨率同步辐射图像（Pierantoni et al., 2023）为三个样本开发了FE模型，如图1a所示。分割工作使用3D Slicer 5.2.2（美国波士顿）半自动完成（见补充图1）。得到的二值化体积用三角形元素网格化，并导出为STL文件。然后使用ANSA v2024（希腊塞萨洛尼基BETA CAE）生成具有4节点四面体元素的实体网格。为了创建圆柱模型，在解剖学FE肌腱的顶部和底部测量了横截面积，并构建了与解剖学模型长度和横截面积匹配的圆柱模型（见图1b）。三个样本的长度从5.7到6.7 mm不等。然而，为了保持一致性，所有三个样本的远端和近端区域都被修剪，以便在所有三种情况下模拟相似长度的肌腱（约5.5 mm）。每个解剖学FE模型包含大约32,000个元素（元素大小约为0.1 mm），而相应的圆柱FE模型包含约45,000个元素（元素大小为0.1 mm）。支持这一选择的网格收敛性分析见补充图2。

2.3. 本构材料模型
肌腱被建模为一个横向各向同性的纤维增强多孔粘弹性材料，基于（Khayyeri et al., 2016）提出的模型，该模型基于（Wilson et al., 2004）的原始材料模型。该模型由充满间隙液体的多孔非纤维基质和代表胶原纤维的纤维成分组成，胶原纤维是组织的主要结构和承重元素。材料模型的更多细节见补充材料。

2.4. 优化
使用来自实验数据（Pierantoni et al., 2023）的边界条件，对具有直线纤维分布的圆柱模型进行了材料属性优化（见图1）。由于胶原是承重最多的成分，承担了肌腱总应力的95%（Khayyeri et al., 2016），因此优化了描述胶原纤维的标准线性固体模型中的五个参数（E1, E2, k1, k2, η），而其他参数保持不变（见Notermans et al., 2019）。使用移动渐近线（MMA）优化算法（Svanberg, 2007）来最小化实验和模拟力-时间响应之间的RMS误差（公式（8）：
(8) RMSerror = 1/n ∑i=1n Fiexp - Fisim
其中Fiexp和Fisim分别对应实验和模拟的负荷值，n表示数据点数量。

2.5. 纤维方向
为了研究纤维方向的影响，使用了三种纤维配置：直线纤维、具有全局和局部映射纤维方向的分散纤维。最初，在解剖学和圆柱模型上都实现了直线纤维配置，以了解几何形状的影响。随后，通过分析分段后的同步辐射X射线断层扫描图像（Pierantoni et al., 2023）（图2a）确定了肌腱的分散纤维方向。初始图像体积（像素大小1.625 μm）被缩放到各向同性像素大小3.25 μm。然后使用Otsu阈值法（Otsu, 1979）对缩放后的图像进行二值化，以区分肌腱组织与其周围环境。纤维方向分析按照（Pierantoni et al., 2023）中的方法进行，使用标准差为1.5的两个高斯滤波器生成颜色图图像，以量化每个点的3D方向（图2b）。

2.6. 下载：下载高分辨率图像（143KB）
2.7. 下载全尺寸图像

图1. (a) 来自高分辨率同步辐射图像的肌腱分割。(b) 解剖学和圆柱FE模型。(c) 应用于FE模型的位移边界条件。
解剖学模型预先加载0.1 mm，持续5秒，以确保它们完全笔直。然后将变形后的配置的几何形状进行循环加载（见图1c），该加载与实验协议一致（6%的应变，持续6.6秒，应变速率为0.91%/s）。圆柱模型在没有预加载的情况下进行了相同的循环加载，因为它们没有轴向偏移且是笔直的。

2.8. 优化
使用来自实验数据的边界条件（Pierantoni et al., 2023）优化了具有直线纤维分布的圆柱模型（见图1）。由于胶原是承重最多的成分，承受了肌腱总应力的95%（Khayyeri et al., 2016），因此优化了描述胶原纤维的标准线性固体模型中的五个参数（E1, E2, k1, k2, η），而其他参数保持不变（见Notermans et al., 2019）。这些纤维定向直方图是根据全局定向方法为每个元素生成的。最后，在保持相同纤维密度的基础上，将特定元素的纤维定向映射回有限元（FE）几何结构中，从而得到了一个具有空间变化（局部）纤维定向的模型（见图3b）。下载：下载高分辨率图像（198KB）下载：下载全尺寸图像

图3. 负载-时间响应的比较。(a) 平均实验值和圆柱形FE模型。(b) 具有直纤维的圆柱形和解剖学FE模型。(c) 具有分散纤维的解剖学FE模型（全局和局部）。

3. 结果
3.1. 优化后的材料参数
使用具有直纤维配置的圆柱形FE模型，根据实验数据优化了纤维材料参数。优化后的材料参数见表2。模拟数据和实验数据的负载-时间响应吻合良好（见补充图3），整体均方根误差为0.35 N。

表2. 优化的肌腱胶原材料参数。
材料参数 E1(MPa) E2(MPa) k1 k2 η(MPa·s)
值 24.9 20.0 0.49 74.3 89.3 3.0

3.2. 几何形状和纤维定向的影响
3.2.1. 负载-时间变化
对于圆柱形模型，力-时间响应与实验结果非常吻合。圆柱形模型的最大峰值载荷为4.5 N，峰值载荷的总体均方根误差为0.3 N（见图3a）。所有三种具有直纤维和分散纤维的解剖学模型的平均力-时间响应分别显示在图3b和图3c中。具有直纤维的解剖学模型的最大峰值载荷为2.7 N，而具有全局纤维的模型为2.5 N，具有局部纤维的模型为1.2 N。在所有模型中，粘弹性响应导致第一周期和第三周期之间的峰值载荷减少了：圆柱形模型减少了20.2%，具有直纤维的解剖学模型减少了18.8%，具有全局纤维的模型减少了18.0%，具有局部纤维的模型减少了13.3%。

3.2.2. 应变异质性
请参见图4，该图展示了在最大峰值载荷（t = 6.6 s）时，圆柱形模型（图4a）、具有直纤维的解剖学模型（图4b）、具有全局纤维的解剖学模型（图4c）以及具有局部纤维的解剖学模型（图4d）的最大主应变（εmax）的空间分布。圆柱形模型的εmax最大值为0.11，具有直纤维的解剖学模型的εmax为0.16，具有全局纤维的模型的εmax为0.11，具有局部纤维的模型的εmax为0.17（见图4d）。
下载：下载高分辨率图像（306KB）下载：下载全尺寸图像

图4. 在最大峰值载荷（t = 6.6 s）时，(a) 具有直纤维的圆柱形解剖学FE模型，(b) 具有直纤维的模型，(c) 具有全局纤维的模型，以及(d) 具有局部纤维的模型的最大主应变的空间变化。(e) 所有上述FE模型在第一个峰值载荷（t = 6.6 s）时的最大主应变比较。

从t = 6.6 s时的全宽半高（FWHM）分析来看，圆柱形模型的平均εmax为0.069 ± 0.015，而具有直纤维的解剖学模型的εmax为0.073 ± 0.020，具有全局纤维的模型的εmax为0.054 ± 0.018，具有局部纤维的模型的εmax为0.063 ± 0.023（见图4e）。从空间分布和FWHM分析来看，解剖学模型的应变异质性高于圆柱形模型。

3.2.3. 流体速度
流体速度主要在第一周期的峰值载荷（6.6 s，见图5a）和第一周期完全卸载后（13.2 s，卸载见图5b）进行分析。在加载过程中，施加在顶面的拉伸位移产生了孔隙压力梯度，将流体从肌腱中推出，导致流体速度矢量向外，圆柱形模型的最大值为约0.0036 mm/s（见图5c），具有直纤维的解剖学模型的最大值为0.0029 mm/s（见图5d）。在具有全局纤维（图5e）和局部纤维（图5f）的解剖学模型中，这种向外流动的流体进一步放大，显著增加到0.0092 mm/s。流体运动的空间分布也因几何形状而异：圆柱形模型中的最高流体流速发生在靠近施加位移的顶部边界附近（见图5c），而在解剖学模型中则集中在肌腱的中部（见图5d和图5e）。此外，全局和局部纤维模型产生的峰值流体速度幅度几乎相同，但由于局部模型的纤维定向，其空间变化略大（见图5f）。尽管存在这些区域差异，所有模型中的流体速度矢量的加载-卸载方向性保持一致，分散纤维增强了流体速度的幅度。在解剖学模型2和3中也观察到了类似的流体速度和最大孔隙压力现象，这些内容已添加到补充图8中。

下载：下载高分辨率图像（529KB）下载：下载全尺寸图像

图5. (a) 圆柱形和(b) 解剖学模型中肌腱横截面的流体速度矢量。在6.6 s时的(c) 圆柱形、(d) 解剖学模型、(e) 直纤维和(f) 局部纤维的轴向流体速度的空间分布。

4. 讨论
本研究全面了解了解剖学FE模型在模拟大鼠跟腱生物力学响应方面与传统圆柱形模型的差异，并强调了纤维定向对载荷分布、应变异质性和流体速度的重要影响。

4.1. 负载-时间响应的影响
结果清楚地表明，肌腱的几何形状对反应力有影响。由于在优化过程中使用了圆柱形模型，实验数据和模拟的负载-时间响应非常吻合，均方根误差为0.32 N。然而，在解剖学模型的情况下，峰值载荷减少了36%。这种减少可以归因于解剖学模型中较小的横截面积（见补充图4），尤其是在肌腱的中部区域（Pierantoni等人，2023年；Reeves和Cooper，2017年）。进一步分散的纤维，无论是全局的还是局部的，在决定机械响应方面也起到了关键作用。与直纤维配置相比，具有分散纤维的模型使得峰值载荷减少了15.8%（全局纤维配置）和52%（局部纤维配置）。与解剖学模型中的直纤维不同，分散纤维并未与加载轴完全对齐，从而降低了整体刚度和峰值载荷。在局部纤维配置中，这种效应更为明显，因为元素级别的纤维直方图引入了更大的空间异质性和错位，进一步降低了整体刚度（Janssen等人，2025年）。还包括纤维分散也引起了显著的平面内剪切应力。对于直纤维的圆柱形模型，峰值（S13）达到了0.010 MPa，而由于其不规则的几何形状，解剖学模型的峰值增加了到0.013 MPa。引入分散纤维（全局和局部）进一步增强了剪切耦合，最大S13约为0.06 MPa（见补充图5）（Federico和Gasser，2010年；Gasser等人，2006年）。此外，解剖学几何形状在加载边界和固定边界表面之间引入了空间偏移，这进一步促进了剪切应力的产生（见补充图6）。

4.2. 应变异质性
最大主应变显示，与圆柱形模型相比，解剖学模型表现出更明显的应变异质性（最大εmax=0.11，平均εmax=0.058±0.015）（见补充图7）。这可以归因于解剖学肌腱的复杂几何形状，它引入了局部化的变形和应变分布。值得注意的是，在器官层面，大鼠肌腱的横截面积在其长度方向上动态变化，变化幅度为6%至16%（Pierantoni等人，2023年）。这种变形是不对称的，包括横向收缩和伸长，导致肌腱中的应变场不均匀。对于具有直纤维的解剖学模型，最大εmax达到了0.16，平均εmax为0.073±0.020。而在具有全局纤维的模型中，最大εmax降低到0.11，平均εmax为0.054±0.018，表明峰值应变明显减少。相比之下，局部纤维配置产生的最大εmax为0.17，平均εmax为0.063±0.023，显示出比直纤维和全局纤维模型更不均匀的应变场。这些差异源于纤维与肌腱几何形状的相互作用。直纤维与加载轴完全对齐，在横截面变化的区域产生突然的不连续性，导致局部应变集中。具有35个φ-θ定向的全局纤维模型引入了各向异性，使应变过渡更加平滑，并更均匀地分布载荷。然而，在局部纤维配置中，由于纤维的元素级定向，存在更大的空间变化（见图5f）。尽管存在这些区域差异，流体速度矢量的加载-卸载方向性在所有模型中保持一致，分散纤维增强了流体速度的幅度。在样本2和3的解剖学模型中也观察到了类似的流体速度和最大孔隙压力现象，这些内容已添加到补充图8中。

4.3. 流体速度
横观各向同性的肌腱模型由多孔的非纤维基质组成，其中充满间质液和纤维状胶原网络，后者是主要的承重组成部分。在轴向张力作用下，固相由于泊松效应而发生径向收缩，这减少了基质中的孔隙空间并产生了孔隙压力增加。这种压力升高根据达西定律驱动流体传输，其中渗透率取决于变化的孔隙比。由于在肌腱表面施加了零孔隙压力边界条件，因此最高压力梯度出现在边缘，导致流体速度向外部边缘增加。因此，空间变化的应变分布，结合肌腱组织的 porelastic 性质，导致肌腱中心的流体速度最小，而边缘的最大流体速度最大。从圆柱形模型到解剖学模型，流体速度增加了19.44%。在圆柱形模型中，峰值流体速度主要出现在加载边界附近，这是由于施加的位移边界条件造成的。相比之下，解剖学模型的最大流体速度出现在中部区域，这是由于其几何形状引入的应变异质性改变了流体流动模式。当结合了分散纤维定向（全局和局部）时，流体速度的变化更大，与直纤维配置相比增加了5.34倍。在直纤维定向中，纤维与加载轴对齐，使它们能够承受最大应力。因此，周围的非纤维基质受到的应力最小，导致孔隙压力降低，从而流体速度降低。然而，分散纤维定向通过φ和θ的分布引入了各向异性。这降低了纤维承受最大载荷的能力，从而增加了非纤维基质中的应力。基质应力的增加导致了孔隙压力的增加，进而增加了流体速度（见补充图8）。尽管全局和局部纤维配置都显著增加了流体速度，但总体流动幅度在两者之间相似。然而，局部纤维配置引入了空间变化的纤维定向，导致应力传递和孔隙压力梯度的局部差异。因此，局部模型中的流体速度场变得更加空间异质（见图5）。当前研究中的最大流体速度范围为3至18 μm/s。一般来说，文献报道跟腱的流体速度范围为5至50 μm/s，这可以触发细胞反应（Archambault等人，2002年）。然而，据我们所知，目前文献中还没有关于跟腱中间质液速度的直接体内测量数据。

这种对流体速度的理解在肌腱的生物物理刺激中起着关键作用。先前的研究（Khayyeri等人，2017年；Khayyeri等人，2015年）已经考虑了纤维的重定向现象；而本研究的结果表明，分散的纤维取向会增强流体的运动，从而能够更好地预测通过肌腱传递的机械力。这种增强的流体运动可以激活机械敏感的信号通路，如COX-2和TGF-β（Lavagnino和Arnoczky，2005年），并有助于细胞在适应和重塑过程中的存活（Butler等人，2008年）。关于应变异质性和载荷传递的这些机制发现是适用于所有物种的。无论在大鼠还是人类中，水分约占肌腱总重量的65-75%，但流体的分布以及组织对压力的响应会随体型和机械负荷的不同而变化。尽管如此，我们的研究结果表明，由结构异质性驱动的应变和流体速度在未来的肌腱生物力学研究中可能具有关键作用。然而，这些定量结果仅适用于大鼠，不应直接用于推断人类的肌腱情况。

5. 局限性
由于样本数量有限，本研究没有进行应变和流体速度的统计比较。此处报告的定量数据应用于理解肌腱的生物力学特性，而非估计尺寸效应。虽然具有分散纤维的解剖模型可能最能代表实验样本，但每个解剖模型由于其独特的几何结构和纤维取向，都会产生一套不同的材料参数。这使得无法得出通用的材料参数，同时也增加了整体计算成本。为了避免这一问题，本研究采用了简化的圆柱模型来确定优化的材料参数。此外，有限元模型可以纳入肌腱下结构的细节，从而更准确地模拟肌肉骨骼系统的负荷。未来研究将探讨肌腱下结构（包括肌腱下结构之间的接触和滑动）及其对不同步态负荷模式下局部应变和流体速度的影响。

6. 结论
总之，本研究将几何形状和纤维取向纳入计算框架中，以更好地理解肌腱的力学特性以及应变异质性的作用。研究结果表明，解剖几何结构是导致肌腱模型中应变异质性的主要因素。此外，加入分散的纤维取向会提高肌腱模型中的流体速度。这些发现表明，由结构异质性驱动的应变和流体速度在未来的肌腱生物力学研究中可能具有关键作用。总体而言，这项工作提出了一种全面的建模方法，强调了解剖细节和纤维结构在准确模拟肌腱生物力学中的重要性。

作者贡献声明：
Yeswanth S Pydi：撰写初稿、验证、方法论制定、数据分析、概念构思。
Maria Pierantoni：撰写修订稿、监督、方法论制定、数据管理。
Anna Gustafsson：撰写修订稿、监督、方法论制定、概念构思。
Hanna Isaksson：撰写修订稿、监督、资源管理、项目协调、资金筹措、概念构思。