2. 常州市中医医院骨伤科, 常州, 213003
2. Department of Orthopedics, Changzhou TCM Hospital, Changzhou, 213003, China
了解航天飞行中航天员的体能变化并设法维持其正常体能, 是航天医学研究的重点课题。以往研究表明, 航天员在长期失重状态下, 骨质不断流失,导致椎体失重性骨质疏松,严重危害宇航员的身体健康。因此, 探索骨质流失对腰椎承载能力的影响已经成为亟待解决的关键课题。
早期的医学研究者主要通过体外实验对腰椎承载能力进行研究,但这种方法无法得到腰椎内在应力、应变的结果。有限元方法(Finite element algorithm, FEA)由于具有成本低、计算周期短,能得到传感器难于测量的内部应力应变等优点,被广泛应用到生物力学的研究中,通过对骨骼几何模型相似性[1-2]、材料性质相似性[3-4]以及边界条件[5]和载荷相似性[6]的深入研究,得到的结果也越来越准确。Rohlmann等[7]利用有限元分析法建立了人体脊椎L1~L3段有限元模型,模拟计算了人体在前屈和后伸两种运动状态下的应力应变分布;Tsouknidas等[8]基于CT数据建立了腰椎有限元模型,通过该模型确定了腰椎节段在复杂载荷下的机械响应;Zhu等[9]整合了不同研究中建立的脊椎局部模型,组成一个完整的脊柱有限元模型。以前的大部分研究都是关于正常脊椎在不同载荷条件下的内部应力应变情况, 而尚未有关于脊椎失重性骨质疏松情况下生物力学性质变化的探究工作。
本文主要根据脊椎的解剖结构, 建立了包括L3~L5运动节段的正常腰椎三维模型, 并通过模型处理构建了失重性骨质疏松的腰椎模型, 然后利用有限元方法, 对正常腰椎以及失重性骨质疏松腰椎进行生物力学分析, 对比了正常腰椎以及失重性骨质疏松腰椎在同等条件下的内部应变以及应力。
1 材料和方法 1.1 材料获取采用由常州市中医医院提供的正常男性脊椎CT数据,脊椎范围的CT数据包括379张CT序列,将获得的CT断层图像以DICOM格式储存。
1.2 L3~L5腰椎模型的建立以及分析 1.2.1 CT图像的处理将DICOM格式的CT断层图像直接导入mimics19.0,然后经过图像处理生成三维腰椎骨骼模型(下面以L5腰椎建立过程为例),具体步骤如下:(1)利用阈值分割(Thresholding)工具将骨骼与周边的组织肌肉区分出来,设定的灰度值为Bone阈值选项,范围为226~2 944;(2)使用编辑蒙板(Editing masks)工具中的Draw和Erase命令逐层对图像进行处理,为后续有限元操作对图像的边缘进行修改以及对图像中的空洞进行填充; (3)通过区域增长(Region growing)工具将L5腰椎单独提取出来,和其他椎体分离。依次通过上述步骤对L3和L4以及两块椎间盘进行图像处理,从而获得L3~L5运动节段的三维模型(图 1)。
![]() |
图 1 L3~L5腰椎三维模型 Figure 1 L3~L5 lumbar three-dimensional model |
1.2.2 体网格的建立
将初步重建的3D模型在mimics自带的网格划分程序3-matic Research中划分面网格,具体流程如下:(1)执行Quality Preserving Reduce Triangles命令缩减三角形网格数量;(2)在原网格的基础上优化网格,对三角网格进行自动优化(Adaptive remesh);(3)对非均匀或出现局部过密或过疏的网格进行手动修改,最终得到优化后的全腰椎面网格模型。在Research中将面网格转化为体网格,建立全腰椎体网格单元模型。最终的模型共有143 298个单元。
1.2.3 材料的赋值本文将骨骼当作各向同性材料来处理。骨骼的表观密度以及材料属性可以通过CT图像反映出来。研究成果表明,骨骼的表观密度与CT灰度值存在近似的线性关系[10-11],表达式为
$ \mathit{\rho }{\rm{ = }}\frac{{{\mathit{\rho }_{{\rm{max}}}}\mathit{HU}}}{{\mathit{H}{\mathit{U}_{{\rm{max}}}}}} $ | (1) |
式中:ρmax为骨骼最大密度值;ρ为表观密度;HUmax为骨骼皮质骨的最大灰度值。正常腰椎最大骨骼密度为2 g/cm3,失重性骨质疏松腰椎最大骨骼密度为1.6 g/cm3。
在各项同性的情况下,骨骼材料属性只须用杨氏模量E和泊松比λ两个弹性常数来表示。而CT扫描是轴向进行的,故可采用在轴向上的材料属性与表观密度之间的经验公式。表达式为
$ \mathit{E}{\rm{ = 1}}\;{\rm{904}}{\mathit{\rho }^{{\rm{1}}{\rm{.64}}}} $ | (2) |
$ \mathit{\lambda }={\rm{0}}{\rm{.3}} $ | (3) |
采用均匀赋值法给L3, L4和L5腰椎进行赋值,将CT灰度值等间隔取样,以间隔中点的CT灰度值代表间隔中所有体单元的CT值,最后用经验公式计算出腰椎的表观密度值,并通过表观密度和杨氏模量的经验公式为腰椎赋值。张国栋等[12]通过对比分析不同材料数目的股骨有限元模型,得到赋予骨骼10种材料属性即能满足有限元需求。
根据正常腰椎和骨质疏松腰椎的经验公式,通过mimics为所有体单元添加材质。表 1为正常腰椎和失重性骨质疏松腰椎的表观密度、杨氏模量以及泊松比的赋值情况。
![]() |
表 1 正常以及骨质疏松腰椎材料特性 Table 1 Normal and osteoporotic lumbar material properties |
椎间盘则主要分为纤维环以及髓核两部分,本文不考虑椎间盘的退变,故设定正常腰椎以及失重性骨质疏松腰椎的椎间盘材料特性相同,如表 2所示。
![]() |
表 2 腰椎间盘材料特性 Table 2 Properties of lumbar intervertebral disc |
材料赋值结束后,将已经赋值材料属性的三维模型以inp格式输出至abaqus中进行有限元分析。
1.2.4 定义接触以及相互作用人体腰椎在直立负重条件下,腰椎终板与椎间盘始终贴合在一起,腰椎与椎间盘间不存在相对位移,而各个腰椎间的相对位移很小,因此设定上下椎体终板和椎间盘的接触为无相对位移的绑定来模拟在直立负重条件下脊椎的形变情况。
1.2.5 边界条件设定L5腰椎下表面所有节点以及下关节突在各个方向上自由度为0;人体直立负重时,通过手臂将重量传递给肩部,然后通过脊柱承重,向下传递,故采取在L3上表面所有节点分别施加竖直向下0.4 MPa的载荷和10 N·m的弯矩来模拟椎体正常站立以及身体后屈的受力情况。
2 数值分析(1) 腰椎和椎间盘在受轴向载荷条件下,椎体和椎间盘会产生明显的应力应变:
① 正常腰椎和椎间盘在受轴向载荷条件下,力主要通过椎体外层皮质骨以及椎弓传递,由于L3上表面直接受载荷冲击,L3椎体应力大于其他椎体,椎体应力随着椎体序号增加而减小(图 2(a))。而椎弓处的应力随着椎体序号增加而增大(图 2(b)),正常腰椎最大应力为4.537 MPa。
![]() |
图 2 正常腰椎有限元分析结果 Figure 2 Results of finite element analysis of normal lumbar |
② 失重性骨质疏松腰椎和椎间盘在受轴向载荷条件下,力依旧主要通过椎体外层皮质骨以及椎弓传递,椎体应力随着椎体序号增加而减小(图 3(a)),椎弓处的应力随着椎体序号增加而增大(图 3(b)),其最大应力为3.529 MPa,较正常腰椎最大应力减小22%。
![]() |
图 3 骨质疏松腰椎有限元分析结果 Figure 3 Results of finite element analysis of lumbar spine |
③ 正常腰椎以及失重性骨质疏松腰椎在受相同载荷条件下,正常腰椎应力范围为4.183×10-5~4.537 MPa,失重性骨质疏松腰椎应力范围为2.797×10-8~3.529 MPa。
④ 正常腰椎和失重性骨质疏松腰椎都会发生微小的形变,L3腰椎的形变最为明显,随着腰椎序号的增加其形变量减小,而在相同载荷情况下,正常腰椎最大形变量为0.759 8 mm(图 4(a)),失重性骨质疏松腰椎最大形变量为0.823 8 mm(图 4(b)), 失重性骨质疏松腰椎最大形变量较正常腰椎最大形变量增大了8.4%。
![]() |
图 4 椎体形变 Figure 4 Vertebral deformation |
(2) 腰椎和椎间盘在受弯矩条件下:
① 当腰椎受弯矩向后弯曲时,椎弓部分整体的应力远大于椎体部分的应力,第4和第5腰椎的椎体部位所受应力很小。
② 在受相同弯矩条件下,正常腰椎最大应力为22.62 MPa(图 5(a)), 相对应的失重性骨质疏松腰椎最大应力为22.46 MPa(图 5(b)),其应力大小及分布差别不大。
![]() |
图 5 腰椎有限元分析结果 Figure 5 Results of lumbar finite element analysis |
③ 相同载荷条件下,正常腰椎最大形变量为1.365 mm(图 6(a)),而失重性骨质疏松腰椎最大形变量为1.609 mm(图 6(b)),失重性骨质疏松腰椎最大形变量较正常腰椎最大形变量增大了15.2%。
![]() |
图 6 椎体形变 Figure 6 Vertebral deformation |
3 结束语
腰椎受轴向载荷条件下应力传导的主要承载为外侧的皮质骨以及椎弓部分,腰椎结构承载能力主要由皮质骨以及椎弓承担,在受弯矩条件下,应力传导的主要承载为椎弓部分。
有限元结果显示:当正常腰椎发生骨质流失并成长为骨质疏松腰椎时, 对比正常腰椎在相同轴向载荷情况下最大应力减小22%,最大形变量增大了8.4%。当腰椎在受弯矩情况下,骨质疏松腰椎较正常腰椎应力分布及大小相近,而最大形变量增大了15.2%。这一结果可以看出,骨质结构承载能力不仅依赖于应力的最大值, 而且也依赖于材料强度, 当骨质流失发生时, 腰椎骨骼材料性能会发生变化, 这包括了强度极限和杨氏模量。当骨质材料的杨氏模量减小时,腰椎的结构变得更加柔软,从而导致最大应力随着骨质流失而减小。另一方面,虽然最大应力减小了, 但随着骨质流失, 材料强度会相应减小, 从而导致结构的承载能力也在下降。
[1] |
陈中, 邢跃刚, 罗绍华.
利用数字几何技术重建个性化骨骼模型[J]. 中国组织工程研究, 2016, 20(39): 5846–5851.
DOI:10.3969/j.issn.2095-4344.2016.39.011 CHEN Zhong, XING Yuegang, LUO Shaohua. Reconstruction for patient-specific bone model based on digital geometry processing technique[J]. Chinese Journal of Tissue Engineering Research, 2016, 20(39): 5846–5851. DOI:10.3969/j.issn.2095-4344.2016.39.011 |
[2] |
MEIJER G J, HOMMINGA J, VELDHUIZEN A G, et al.
Influence of interpersonal geometrical variation on spinal motion segment stiffness:Implications for patient-specific modeling[J]. Spine, 2011, 36(14): 929–935.
DOI:10.1097/BRS.0b013e3181fd7f7f
|
[3] |
TRAVERT C, JOLIVET E, SAPIN-DE B E, et al.
Sensitivity of patient-specific vertebral finite element model from low dose imaging to material properties and loading conditions[J]. Medical & Biological Engineering & Computing, 2011, 49(12): 1355–1361.
|
[4] |
蔡康健, 王丽珍, 姚杰, 等.
腰椎椎体有限元建模的最优单元尺寸和材料属性分布及建模方法[J]. 医用生物力学, 2016, 31(2): 135–141.
DOI:10.3871/j.1004-7220.2016.02.135 CAI Kangjian, WANG Lizhen, YAO Jie, et al. The optimal element size, material property distributions and modeling methods for finite element modeling of lumbar vertebra[J]. Journal of Medical Biomechanics, 2016, 31(2): 135–141. DOI:10.3871/j.1004-7220.2016.02.135 |
[5] |
白石柱, 李涤尘, 赵铱民, 等.
上颌骨有限元分析中边界约束条件的研究[J]. 临床口腔医学杂志, 2006, 22(12): 720–723.
DOI:10.3969/j.issn.1003-1634.2006.12.006 BAI Shizhu, LI Dichen, ZHAO Yimin, et al. The boundary constraint study in finite element analysis of maxillary[J]. Journal of Clinical Stomatology, 2006, 22(12): 720–723. DOI:10.3969/j.issn.1003-1634.2006.12.006 |
[6] |
赵峰, 韩彦峰, 胡江峰.
弹性模量和初始应力对种植体骨界面应力分布影响的三维有限元分析[J]. 中国口腔种植学杂志, 2006, 11(2): 55–58.
ZHAO Feng, HAN Yanfeng, HU Jiangfeng. Three-dimensional finite element method analysis of relation of implant elastic modulus and initial stress and bone-implant surface stress distribution[J]. Chinese Journal of Oral Implantology, 2006, 11(2): 55–58. |
[7] |
ROHLMANN A, BAUER L T, BERGMANN G, et al.
Determination of trunk muscle forces for flexion and extension by using a validated finite element model of the lumbar spine and measured in vivo data[J]. Journal of Biomechanics, 2006, 39(6): 981–989.
DOI:10.1016/j.jbiomech.2005.02.019
|
[8] |
TSOUKNIDAS A, MICHAILIDIS N, SAVVAKIS S, et al.
A finite element model technique to determine the mechanical response of a lumbar spine segment under complex loads[J]. Journal of Applied Biomechanics, 2012, 28(4): 448–456.
DOI:10.1123/jab.28.4.448
|
[9] |
ZHU R, ZANDER T, DREISCHARF M, et al.
Considerations when loading spinal finite element models with predicted muscle forces from inverse static analyses[J]. Journal of Biomechanics, 2013, 46(7): 1376–1378.
DOI:10.1016/j.jbiomech.2013.03.003
|
[10] |
RHO J Y, HOBATHO M C, ASHMAN R B.
Relations of mechanical-properties to density and CT numbers in human bone[J]. Med Eng Phys, 1995, 17(5): 347–355.
DOI:10.1016/1350-4533(95)97314-F
|
[11] |
WIRTZ D C, SCHIFFERS N, PANDORF T, et al.
Critical evaluation of known bone material properties to realize anisotropic FE-simulation of the proximal femur[J]. Journal of Biomechanics, 2000, 33(10): 1325–1330.
DOI:10.1016/S0021-9290(00)00069-5
|
[12] |
张国栋, 廖维靖, 陶圣祥, 等.
股骨有限元分析赋材料属性的方法[J]. 中国组织工程研究, 2009, 13(43): 8436–8441.
DOI:10.3969/j.issn.1673-8225.2009.43.006 ZHANG Guodong, LIAO Weijing, TAO Shengxiang, et al. Methods for material assignment of finite element analysis with femurs[J]. Chinese Journal of Tissue Engineering Research, 2009, 13(43): 8436–8441. DOI:10.3969/j.issn.1673-8225.2009.43.006 |