在不同研究中纤维的材料特性、体积比例以及髓核的横截面积比例差异很大，而这些因素对椎间盘力学行为的影响尚不明确。通过建立不同参数的椎间盘有限元模型研究椎间盘在压力、压力与角位移或压力与扭矩等载荷下的压力、高度变化、转动角度、应力和应变。结果发现，纤维材料特性对椎间盘力学行为有很大影响，特别是对转动角度的影响尤为显著。当纤维体积比例较小时，其变化对椎间盘的旋转角影响很大。髓核的面积比例对椎间盘力学行为的影响相对较小。在合理范围内通过适当调整纤维的弹性模量或体积比例，可以得到能够模拟正常椎间盘力学行为的模型。在椎间盘有限元模型中髓核的面积比例在25%～50%内是合理的。本文为椎间盘有限元模型创建和椎间盘退化机制研究提供了指导和参考。
引用本文： 董瑞春, 刘忠, 郭云强, 安钰坤, 石舟, 时明. 纤维的体积与材料特性以及髓核横截面积对椎间盘力学行为的影响研究. 生物医学工程学杂志, 2024, 41(1): 144-151. doi: 10.7507/1001-5515.202305001 复制
0 引言
椎间盘由髓核和纤维环组成，其中纤维环又包含基质和纤维。椎间盘的老化和退行性变通常会伴随着上述组织成分和几何结构的改变[1-2]，组织成分的改变将导致材料属性的改变，而材料属性和几何结构都会影响椎间盘的生物力学特性和运动范围。然而，在不同的研究中，纤维的材料特性有很大的差别。例如，有研究发现，由基体和纤维束组成的纤维环试样的等效模量在49～95 MPa之间[3]；在部分研究中，从单个片层获得的实验样品也是纤维环基质和纤维束的组合，发现这些样品的弹性模量在28～136 MPa之间[1-2, 4]；另有研究表明，纤维的弹性模量与其在纤维环中的位置有关，在纤维环外部大而在内部小，在纤维环前部大而在后部小[1-2]；许多研究中纤维的弹性模量取值为500 MPa[5-7]，但也有研究取值为175 MPa[8]。另一方面，在不同文献中，纤维相对于周围基质的体积含量也不同。例如，一项研究认为纤维在纤维环基质中的体积含量在外层的23%到内层的5%之间变化[9]；而另一项研究中，纤维的总体积假定为纤维环基质体积的19%[10]。此外，不同文献中髓核横截面积与椎间盘横截面积的比例也不同，多在25%～50%之间变化[5, 11-16]。上述研究中纤维弹性模量和体积以及髓核横截面积的差异会影响椎间盘有限元模型的准确性。然而，这些因素的变化对椎间盘生物力学的影响规律尚不明确。准确描述椎间盘材料特性和结构的模型对于理解健康、损伤和退化椎间盘的生物力学特性非常重要[17]，因此本研究的目的是探讨纤维材料性质和体积以及髓核横截面积对椎间盘力学行为的影响机制，为有限元建模和椎间盘退变机制的研究提供指导和参考。
1 方法
1.1 椎间盘有限元模型
根据50岁健康女性（身高160 cm，体重60 kg）的腰椎间盘L3-L4几何数据建立椎间盘有限元模型。采用六面体单元对髓核、纤维环基质和终板进行模拟。纤维采用只承受拉力的桁架单元（truss单元）模拟。纤维取向θ从纤维环内侧的±43°下降到纤维环外侧的±28°（相对于横切面或水平面）[18]。不同文献中髓核横截面积占椎间盘总横截面积的比例分别为：0.28[11]，0.36[15]，0.38[12]，0.40[13]，0.42[19]，0.45[16]，0.47[14]和0.50[5]。本研究中创建的五个椎间盘模型，其髓核面积比例分别为0.20、0.29、0.40、0.46和0.53，分别标记为P1、P2、P3、P4和P5，如图1所示。
图1
不同髓核横截面积比的椎间盘L3-L4有限元模型
AO：前外侧；AI：前内侧；LO：左侧外侧；LI：左侧内侧；RO：右侧外侧；RI：右侧内侧；PO：后外侧；PI：后内侧。
AO: anterior-outer; AI: anterior-inner; LO: left lateral-outer; LI: left lateral-inner; RO: right lateral-outer; RI: right lateral-inner; PO: posterior-outer; PI: posterior-inner.
如图1所示，纤维环前部、侧部和后部纵切面的面积（SAF）大小不同，因此，纤维环被分为前部、两个侧部和后部等四个部分，每个部分内部的纵切面面积认为近似不变。各部分的纤维体积比Vp可通过式(1)近似估算：
 |
其中，S是单根纤维的横截面积；n是纤维环纵切面内的纤维根数。假设Vp和n在四个部分中是一致的。纤维环中采取了四种体积比例：5%、19%[10]、33%和100%，分别标记为S1或Vp1、S2或Vp2、S3或Vp3和S4或Vp4，如表1所列，则不同位置单根纤维的横截面积可以通过式(2)获得:
表1
单根纤维横截面积（mm2）
Table1.
The cross-sectional area of single fiber (mm2)
 |
1.2 材料特性
在有限元模型中，将纤维和终板模拟为线性弹性材料。纤维和纤维环基体分别赋予不同的材料属性[3]。把文献中纤维的材料归属为两大类：第一类材料（标记为M1，见表2），纤维的弹性模量从最外层的58 MPa下降到最内层的14 MPa[1, 11]；第二类材料（标记为M2，见表2），其弹性模量从最外层的550 MPa下降到最内侧的357 MPa[6, 10, 20-24]。终板模量为500 MPa，泊松比为0.25[25]。
表2
不同区域纤维的弹性模量（MPa）
Table2.
Young modulus of the fibers in different regions (MPa)
髓核和纤维环基质采用可压缩的Mooney-Rivlin超弹性材料模拟[11, 26]。Mooney-Rivlin材料的应变函数由公式(3)描述：
 |
其中I1和I2是Cauchy-Green右变形张量的第一和第二不变量，J是变形梯度张量的行列式。根据文献[11, 27]，对于髓核，c1、c2和k的值分别为0.05、0.01和50 MPa；对于纤维环基质，c1、c2和k的值分别为0.2、0.01和6 MPa。
1.3 模拟方案
在软件ABAQUS 6.14（法国达索模拟系统公司）中进行模拟计算。由于通过单独创建椎间盘模型来研究椎间盘力学特性的文献数据较少，为了与文献结果对比验证模型，本文基于文献研究选择载荷方案：载荷1，第一步先在椎间盘的上终板施加520 N或0.48 MPa的压力，第二步再对椎间盘上终板进行绕Z轴±5°的旋转[11]，以模拟获得椎间盘最上端高度变化和纤维应变；载荷2，第一步先施加720 N或0.64 MPa的压力，第二步再施加绕Z轴5.5 Nm的扭矩[28]，以模拟获得髓核平均压力、椎间盘高度和转动角度以及纤维环基质应力。所有载荷在每个分析步骤中线性增加或减少。椎间盘下终板在各个方向的自由度都是固定的。对计算求得的椎间盘压力、高度变化、转动角度、应力和应变，分别与文献实验和模拟结果进行对比以验证模型。
本文对不同组合模型进行模拟研究：第一类，研究两种纤维材料下不同纤维体积比例（S1、S2、S3和S4）对椎间盘力学行为的影响，选用的组合模型为M1S1P4、M1S2P4、M1S3P4、M1S4P4、M2S1P4、M2S2P4、M2S3P4和M2S4P4；第二类，研究不同髓核面积比例（P1、P2、P3、P4和P5）对椎间盘力学行为的影响，因纤维材料与体积比例为M1S4和M2S2组合时模拟结果与实验结果相近，所以选用的组合模型为M1S4P1、M1S4P2、M1S4P3、M1S4P4、M1S4P5、M2S2P1、M2S2P2、M2S2P3、M2S2P4和M2S2P5。
2 结果
2.1 髓核压力
本部分结果（见图2）为所有组合模型在载荷2（720 N压力与5.5 Nm扭矩）下的髓核压力。所有模型的髓核平均压力均在体外测量结果范围内[27]，说明本研究的模拟压力结果是正确的。此外，发现纤维材料为M2时髓核压力高于纤维材料为M1时的髓核压力。在荷载2的第2个分析步下，髓核平均压力随轴向扭矩增加很小。髓核平均压力随纤维体积比的增加而显著增加，而随髓核横截面积的增加先增加后减小，并且变化较小。
2.2 椎间盘高度变化
本部分结果（见图3）为所有组合模型分别在载荷1（520 N压力与±5°转动角）和载荷2（720 N压力与5.5 Nm扭矩）下的高度变化。不同载荷下纤维材料为M2的椎间盘高度变化量几乎都大于纤维材料为M1的椎间盘高度变化量。当轴向扭矩在0～5.5 Nm范围内变化时（在载荷2分析步2下），椎间盘高度略微增加。椎间盘轴向应变随纤维体积比的增加而降低，随髓核横截面积的增加而略有降低。部分模型的轴向应变（椎间盘高度变化量除以初始高度）在文献的结果范围内（0.08 ± 0.01）[11, 28]，说明部分参数的组合模型可以模拟正常椎间盘的力学行为。
2.3 椎间盘转动角度
本部分结果（见图4）为所有组合模型在载荷2（720 N压力与5.5 Nm扭矩）下的转动角度。椎间盘轴向转动受到材料和纤维体积比例的影响很大。然而，髓核的横截面积对椎间盘轴向转动的影响相对较小。部分参数组合模型的椎间盘旋转在体内测量结果范围内[27, 29]，M2S2P4、M2S3P4和M2S4P4模型的椎间盘旋转在体外测量结果范围内[27]，这进一步说明部分参数的组合模型可以模拟正常椎间盘的力学行为。
通过拟合结果可以看出，当纤维体积比例较小时，其变化对椎间盘的旋转角影响很大；纤维材料为M1时，当体积比例超过60%后，其变化对椎间盘的旋转角影响较小；纤维材料为M2时，当体积比例超过20%后，其变化对椎间盘的旋转角影响很小。此外，可以观察到髓核面积的比例对椎间盘旋转角度的影响相对较小。
2.4 椎间盘应力、应变与突出变形
本部分结果（见图5）为模型M2S2P4在载荷1（520 N压力与±5°转动角）和载荷2（720 N压力与5.5 Nm扭矩）下的最大主应变和应力。在520 N压力下纤维最大主应变发生在纤维环内侧靠近终板处，数值为0.07，在520 N压力与5°转动角下，纤维最大主应变在纤维环后外侧与端面连接区域，数值为0.16，结果与文献[11]中的数据一致。纤维环基质在单独720 N压力以及720 N压力与5.5 Nm扭矩的组合下应力分布是相同的，最大应力都在纤维环内部和靠近端板处，与文献研究结果也相同[30]，进一步说明模型可以模拟正常椎间盘的力学行为。
图5
模型M2S2P4中纤维最大主应变、纤维环基质最大应力分布与椎间盘中间横截面最外侧变形
Figure5.
Maximum principal strain of fibers, maximum stress distribution of fiber ring matrix and the bulge of the outermost region of the IVD in the mid-transverse plane of model M2S2P4
在本文建立的椎间盘有限元模型中，纤维环后部较薄，而且材料的弹性模量较小，因此，在轴向压力作用下，纤维环后部的凸起最大。所以，建立患者特定结构特性的有限元模型，对于准确研究椎间盘的生物力学特性具有重要意义[17]。
3 讨论与结论
通过上述模拟结果分析与文献对比，可以认为本研究创建的腰椎间盘有限元模型有效性得到了验证，可以用于进一步研究。在本研究中，我们发现纤维的材料性质和体积比例对椎间盘的生物力学特性有很大的影响，特别是对椎间盘的旋转角度（见图4）。例如，当纤维采用材料M1、纤维体积比为S1时，计算的椎间盘旋转角是实验结果[27, 29]的10倍。但是，它们对椎间盘的压力和高度的影响相对较小（见图2～3）。因此，对椎间盘有限元模型进行验证时，应观察椎间盘的旋转角度是否在实验范围内。此外，我们发现，髓核的面积比例对椎间盘的生物力学影响不大（见图2～4），因此，在建立椎间盘的有限元模型时，将髓核的面积比例控制在25%～50%[3, 31]是合理的。
虽然文献中使用的纤维材料性能有很大的不同[1, 9, 11, 32]，但通过调整纤维的体积比例（或横截面积），模拟结果也可以在实验范围内（见图4）。一方面这说明不同文献中纤维的体积可能不同，因为只有很少的文献[3, 7]说明有限元模型中纤维的体积比例或横截面积，如果模型中的纤维体积较小，则使用的弹性模量较大。这可以解释不同文献中纤维弹性模量的差异，特别是实验研究与有限元研究中的差异。因为在纤维力学实验中，样品是单一的片层，是纤维环基体和纤维束的组合，而不是单一的纤维[1-2, 4]。因此，有限元模型中纤维的体积可能小于实验体积，则模型中纤维的弹性模量将大于实验弹性模量。另一方面实验中纤维弹性模量的差异可能来源于椎间盘不同程度的退化，椎间盘的退化可能导致其弹性模量增大，但这需要进一步实验验证。
在轴向压缩和轴向扭矩（载荷2）共同作用下，所有模型计算的平均髓核压力均在体外测量结果[27]范围内（见图2），但是只有部分模型计算的椎间盘高度变化和旋转角度在体外模拟[11, 27]和实验测量结果范围内[28]（见图3～4）。这说明两点：第一，验证了本文模型的有效性；第二，可以通过调整纤维的弹性模量或体积比，使模型更加真实合理。为了给模型标定提供指导或参考，本文总结了两种材料特性下纤维体积比与椎间盘旋转角度之间的关系，如图4所示。结果表明，当纤维体积比例较小时，其变化对椎间盘的旋转角影响很大；而当纤维体积比例增大时（M1纤维材料体积比例超过60%，M2纤维材料体积比例超过20%），其变化则对椎间盘的旋转角影响很小。此外，可以观察到髓核面积的比例对椎间盘旋转角度的影响相对较小（见图4）。
不同方向载荷对椎间盘应力或应变分布的影响不同。例如在轴向压力或在轴向压力与较小扭矩的组合作用下，纤维环后部鼓起最大，说明纤维环后部在轴向压缩载荷下易受损伤，在轴向转动载荷下纤维环后外侧的位移或变形最大，说明后外侧纤维环在旋转载荷下易受损伤（见图5）。因此，可以通过患者椎间盘损伤位置判断损伤原因。需要注意的是在相同载荷下纤维和纤维环基质的应力和应变分布是相同的，所以在研究椎间盘运动范围和应力分布时可以用基质代替纤维（适当增大基质弹性模量）简化模型。但是由于纤维和纤维环基质材料特性不同，在研究椎间盘断裂时需要创建详细的纤维模型。另外，加载速率可影响椎间盘的生物力学[3, 33-34]，但因为加载速率通常不被视为模型参数[35]，因此本研究未考虑加载速率。如果进行更广泛的详细研究，则应考虑加载速率[35]。
重要声明
利益冲突声明：本文全体作者均声明不存在利益冲突。
作者贡献声明：董瑞春负责方案设计、模拟计算、数据整理分析和正文撰写；刘忠和时明负责方案设计、结果分析与讨论；郭云强负责模拟计算、数据整理分析和撰写草稿；安钰坤和石舟负责模拟计算和结果分析。
0 引言
椎间盘由髓核和纤维环组成，其中纤维环又包含基质和纤维。椎间盘的老化和退行性变通常会伴随着上述组织成分和几何结构的改变[1-2]，组织成分的改变将导致材料属性的改变，而材料属性和几何结构都会影响椎间盘的生物力学特性和运动范围。然而，在不同的研究中，纤维的材料特性有很大的差别。例如，有研究发现，由基体和纤维束组成的纤维环试样的等效模量在49～95 MPa之间[3]；在部分研究中，从单个片层获得的实验样品也是纤维环基质和纤维束的组合，发现这些样品的弹性模量在28～136 MPa之间[1-2, 4]；另有研究表明，纤维的弹性模量与其在纤维环中的位置有关，在纤维环外部大而在内部小，在纤维环前部大而在后部小[1-2]；许多研究中纤维的弹性模量取值为500 MPa[5-7]，但也有研究取值为175 MPa[8]。另一方面，在不同文献中，纤维相对于周围基质的体积含量也不同。例如，一项研究认为纤维在纤维环基质中的体积含量在外层的23%到内层的5%之间变化[9]；而另一项研究中，纤维的总体积假定为纤维环基质体积的19%[10]。此外，不同文献中髓核横截面积与椎间盘横截面积的比例也不同，多在25%～50%之间变化[5, 11-16]。上述研究中纤维弹性模量和体积以及髓核横截面积的差异会影响椎间盘有限元模型的准确性。然而，这些因素的变化对椎间盘生物力学的影响规律尚不明确。准确描述椎间盘材料特性和结构的模型对于理解健康、损伤和退化椎间盘的生物力学特性非常重要[17]，因此本研究的目的是探讨纤维材料性质和体积以及髓核横截面积对椎间盘力学行为的影响机制，为有限元建模和椎间盘退变机制的研究提供指导和参考。
1 方法
1.1 椎间盘有限元模型
根据50岁健康女性（身高160 cm，体重60 kg）的腰椎间盘L3-L4几何数据建立椎间盘有限元模型。采用六面体单元对髓核、纤维环基质和终板进行模拟。纤维采用只承受拉力的桁架单元（truss单元）模拟。纤维取向θ从纤维环内侧的±43°下降到纤维环外侧的±28°（相对于横切面或水平面）[18]。不同文献中髓核横截面积占椎间盘总横截面积的比例分别为：0.28[11]，0.36[15]，0.38[12]，0.40[13]，0.42[19]，0.45[16]，0.47[14]和0.50[5]。本研究中创建的五个椎间盘模型，其髓核面积比例分别为0.20、0.29、0.40、0.46和0.53，分别标记为P1、P2、P3、P4和P5，如图1所示。
图1
不同髓核横截面积比的椎间盘L3-L4有限元模型
AO：前外侧；AI：前内侧；LO：左侧外侧；LI：左侧内侧；RO：右侧外侧；RI：右侧内侧；PO：后外侧；PI：后内侧。
AO: anterior-outer; AI: anterior-inner; LO: left lateral-outer; LI: left lateral-inner; RO: right lateral-outer; RI: right lateral-inner; PO: posterior-outer; PI: posterior-inner.
如图1所示，纤维环前部、侧部和后部纵切面的面积（SAF）大小不同，因此，纤维环被分为前部、两个侧部和后部等四个部分，每个部分内部的纵切面面积认为近似不变。各部分的纤维体积比Vp可通过式(1)近似估算：
|
其中，S是单根纤维的横截面积；n是纤维环纵切面内的纤维根数。假设Vp和n在四个部分中是一致的。纤维环中采取了四种体积比例：5%、19%[10]、33%和100%，分别标记为S1或Vp1、S2或Vp2、S3或Vp3和S4或Vp4，如表1所列，则不同位置单根纤维的横截面积可以通过式(2)获得:
表1
单根纤维横截面积（mm2）
Table1.
The cross-sectional area of single fiber (mm2)
|
1.2 材料特性
在有限元模型中，将纤维和终板模拟为线性弹性材料。纤维和纤维环基体分别赋予不同的材料属性[3]。把文献中纤维的材料归属为两大类：第一类材料（标记为M1，见表2），纤维的弹性模量从最外层的58 MPa下降到最内层的14 MPa[1, 11]；第二类材料（标记为M2，见表2），其弹性模量从最外层的550 MPa下降到最内侧的357 MPa[6, 10, 20-24]。终板模量为500 MPa，泊松比为0.25[25]。
表2
不同区域纤维的弹性模量（MPa）
Table2.
Young modulus of the fibers in different regions (MPa)
髓核和纤维环基质采用可压缩的Mooney-Rivlin超弹性材料模拟[11, 26]。Mooney-Rivlin材料的应变函数由公式(3)描述：
|
其中I1和I2是Cauchy-Green右变形张量的第一和第二不变量，J是变形梯度张量的行列式。根据文献[11, 27]，对于髓核，c1、c2和k的值分别为0.05、0.01和50 MPa；对于纤维环基质，c1、c2和k的值分别为0.2、0.01和6 MPa。
1.3 模拟方案
在软件ABAQUS 6.14（法国达索模拟系统公司）中进行模拟计算。由于通过单独创建椎间盘模型来研究椎间盘力学特性的文献数据较少，为了与文献结果对比验证模型，本文基于文献研究选择载荷方案：载荷1，第一步先在椎间盘的上终板施加520 N或0.48 MPa的压力，第二步再对椎间盘上终板进行绕Z轴±5°的旋转[11]，以模拟获得椎间盘最上端高度变化和纤维应变；载荷2，第一步先施加720 N或0.64 MPa的压力，第二步再施加绕Z轴5.5 Nm的扭矩[28]，以模拟获得髓核平均压力、椎间盘高度和转动角度以及纤维环基质应力。所有载荷在每个分析步骤中线性增加或减少。椎间盘下终板在各个方向的自由度都是固定的。对计算求得的椎间盘压力、高度变化、转动角度、应力和应变，分别与文献实验和模拟结果进行对比以验证模型。
本文对不同组合模型进行模拟研究：第一类，研究两种纤维材料下不同纤维体积比例（S1、S2、S3和S4）对椎间盘力学行为的影响，选用的组合模型为M1S1P4、M1S2P4、M1S3P4、M1S4P4、M2S1P4、M2S2P4、M2S3P4和M2S4P4；第二类，研究不同髓核面积比例（P1、P2、P3、P4和P5）对椎间盘力学行为的影响，因纤维材料与体积比例为M1S4和M2S2组合时模拟结果与实验结果相近，所以选用的组合模型为M1S4P1、M1S4P2、M1S4P3、M1S4P4、M1S4P5、M2S2P1、M2S2P2、M2S2P3、M2S2P4和M2S2P5。
2 结果
2.1 髓核压力
本部分结果（见图2）为所有组合模型在载荷2（720 N压力与5.5 Nm扭矩）下的髓核压力。所有模型的髓核平均压力均在体外测量结果范围内[27]，说明本研究的模拟压力结果是正确的。此外，发现纤维材料为M2时髓核压力高于纤维材料为M1时的髓核压力。在荷载2的第2个分析步下，髓核平均压力随轴向扭矩增加很小。髓核平均压力随纤维体积比的增加而显著增加，而随髓核横截面积的增加先增加后减小，并且变化较小。
2.2 椎间盘高度变化
本部分结果（见图3）为所有组合模型分别在载荷1（520 N压力与±5°转动角）和载荷2（720 N压力与5.5 Nm扭矩）下的高度变化。不同载荷下纤维材料为M2的椎间盘高度变化量几乎都大于纤维材料为M1的椎间盘高度变化量。当轴向扭矩在0～5.5 Nm范围内变化时（在载荷2分析步2下），椎间盘高度略微增加。椎间盘轴向应变随纤维体积比的增加而降低，随髓核横截面积的增加而略有降低。部分模型的轴向应变（椎间盘高度变化量除以初始高度）在文献的结果范围内（0.08 ± 0.01）[11, 28]，说明部分参数的组合模型可以模拟正常椎间盘的力学行为。
2.3 椎间盘转动角度
本部分结果（见图4）为所有组合模型在载荷2（720 N压力与5.5 Nm扭矩）下的转动角度。椎间盘轴向转动受到材料和纤维体积比例的影响很大。然而，髓核的横截面积对椎间盘轴向转动的影响相对较小。部分参数组合模型的椎间盘旋转在体内测量结果范围内[27, 29]，M2S2P4、M2S3P4和M2S4P4模型的椎间盘旋转在体外测量结果范围内[27]，这进一步说明部分参数的组合模型可以模拟正常椎间盘的力学行为。
通过拟合结果可以看出，当纤维体积比例较小时，其变化对椎间盘的旋转角影响很大；纤维材料为M1时，当体积比例超过60%后，其变化对椎间盘的旋转角影响较小；纤维材料为M2时，当体积比例超过20%后，其变化对椎间盘的旋转角影响很小。此外，可以观察到髓核面积的比例对椎间盘旋转角度的影响相对较小。
2.4 椎间盘应力、应变与突出变形
本部分结果（见图5）为模型M2S2P4在载荷1（520 N压力与±5°转动角）和载荷2（720 N压力与5.5 Nm扭矩）下的最大主应变和应力。在520 N压力下纤维最大主应变发生在纤维环内侧靠近终板处，数值为0.07，在520 N压力与5°转动角下，纤维最大主应变在纤维环后外侧与端面连接区域，数值为0.16，结果与文献[11]中的数据一致。纤维环基质在单独720 N压力以及720 N压力与5.5 Nm扭矩的组合下应力分布是相同的，最大应力都在纤维环内部和靠近端板处，与文献研究结果也相同[30]，进一步说明模型可以模拟正常椎间盘的力学行为。
图5
模型M2S2P4中纤维最大主应变、纤维环基质最大应力分布与椎间盘中间横截面最外侧变形
Figure5.
Maximum principal strain of fibers, maximum stress distribution of fiber ring matrix and the bulge of the outermost region of the IVD in the mid-transverse plane of model M2S2P4
在本文建立的椎间盘有限元模型中，纤维环后部较薄，而且材料的弹性模量较小，因此，在轴向压力作用下，纤维环后部的凸起最大。所以，建立患者特定结构特性的有限元模型，对于准确研究椎间盘的生物力学特性具有重要意义[17]。
3 讨论与结论
通过上述模拟结果分析与文献对比，可以认为本研究创建的腰椎间盘有限元模型有效性得到了验证，可以用于进一步研究。在本研究中，我们发现纤维的材料性质和体积比例对椎间盘的生物力学特性有很大的影响，特别是对椎间盘的旋转角度（见图4）。例如，当纤维采用材料M1、纤维体积比为S1时，计算的椎间盘旋转角是实验结果[27, 29]的10倍。但是，它们对椎间盘的压力和高度的影响相对较小（见图2～3）。因此，对椎间盘有限元模型进行验证时，应观察椎间盘的旋转角度是否在实验范围内。此外，我们发现，髓核的面积比例对椎间盘的生物力学影响不大（见图2～4），因此，在建立椎间盘的有限元模型时，将髓核的面积比例控制在25%～50%[3, 31]是合理的。
虽然文献中使用的纤维材料性能有很大的不同[1, 9, 11, 32]，但通过调整纤维的体积比例（或横截面积），模拟结果也可以在实验范围内（见图4）。一方面这说明不同文献中纤维的体积可能不同，因为只有很少的文献[3, 7]说明有限元模型中纤维的体积比例或横截面积，如果模型中的纤维体积较小，则使用的弹性模量较大。这可以解释不同文献中纤维弹性模量的差异，特别是实验研究与有限元研究中的差异。因为在纤维力学实验中，样品是单一的片层，是纤维环基体和纤维束的组合，而不是单一的纤维[1-2, 4]。因此，有限元模型中纤维的体积可能小于实验体积，则模型中纤维的弹性模量将大于实验弹性模量。另一方面实验中纤维弹性模量的差异可能来源于椎间盘不同程度的退化，椎间盘的退化可能导致其弹性模量增大，但这需要进一步实验验证。
在轴向压缩和轴向扭矩（载荷2）共同作用下，所有模型计算的平均髓核压力均在体外测量结果[27]范围内（见图2），但是只有部分模型计算的椎间盘高度变化和旋转角度在体外模拟[11, 27]和实验测量结果范围内[28]（见图3～4）。这说明两点：第一，验证了本文模型的有效性；第二，可以通过调整纤维的弹性模量或体积比，使模型更加真实合理。为了给模型标定提供指导或参考，本文总结了两种材料特性下纤维体积比与椎间盘旋转角度之间的关系，如图4所示。结果表明，当纤维体积比例较小时，其变化对椎间盘的旋转角影响很大；而当纤维体积比例增大时（M1纤维材料体积比例超过60%，M2纤维材料体积比例超过20%），其变化则对椎间盘的旋转角影响很小。此外，可以观察到髓核面积的比例对椎间盘旋转角度的影响相对较小（见图4）。
不同方向载荷对椎间盘应力或应变分布的影响不同。例如在轴向压力或在轴向压力与较小扭矩的组合作用下，纤维环后部鼓起最大，说明纤维环后部在轴向压缩载荷下易受损伤，在轴向转动载荷下纤维环后外侧的位移或变形最大，说明后外侧纤维环在旋转载荷下易受损伤（见图5）。因此，可以通过患者椎间盘损伤位置判断损伤原因。需要注意的是在相同载荷下纤维和纤维环基质的应力和应变分布是相同的，所以在研究椎间盘运动范围和应力分布时可以用基质代替纤维（适当增大基质弹性模量）简化模型。但是由于纤维和纤维环基质材料特性不同，在研究椎间盘断裂时需要创建详细的纤维模型。另外，加载速率可影响椎间盘的生物力学[3, 33-34]，但因为加载速率通常不被视为模型参数[35]，因此本研究未考虑加载速率。如果进行更广泛的详细研究，则应考虑加载速率[35]。
重要声明
利益冲突声明：本文全体作者均声明不存在利益冲突。
作者贡献声明：董瑞春负责方案设计、模拟计算、数据整理分析和正文撰写；刘忠和时明负责方案设计、结果分析与讨论；郭云强负责模拟计算、数据整理分析和撰写草稿；安钰坤和石舟负责模拟计算和结果分析。
