轻度青少年特发性脊柱侧凸有限元模型构建及分析
2021-07-30吴晓薇李鉴轶杜冰冉杨巧华范天成李路韬
吴晓薇,李鉴轶,2,杜冰冉,杨巧华,范天成,李路韬
1.南方医科大学基础医学院人体解剖学教研室,广东省医学生物力学重点实验室,广州 510515;2.南方医科大学第七附属医院,广东 佛山 528244
青少年特发性脊柱侧凸(adolescent idiopathic scoliosis,AIS)是发生于青春期最常见的脊柱畸形,2013年9月至2014年7月期间,中山大学附属第一医院脊柱侧弯中心筛查了99695名广东省初、高中学生,发现AIS的发病率为5.14%[1]。AIS会导致患者的脊柱扭曲、畸形,胸廓变形,不仅影响患者的体型外貌,严重时还会限制心肺功能,给患者带来身心双重痛苦。关于青少年特发性脊柱侧凸发生的原因有许多学说,包括遗传、激素水平、生物力学因素、肌骨系统及神经系统等因素[2~10],但确切的发病机制尚不清楚。随着计算机技术的突飞猛进,有限元分析技术作为生物医学领域仿真人体结构力学研究的一个重要实验手段,被广泛应用于AIS的病因研究中。但目前AIS的三维有限元分析主要针对重度患者,尚未有针对轻度AIS患者的三维有限元分析。本研究利用有限元分析技术构建轻度AIS患者的有限元模型并验证模型的有效性,以期将轻度AIS患者的有限元模型用于该疾病的相关研究中。
1 材料和方法
1.1 临床材料
1.1.1 研究对象 选取中山大学附属第一医院脊柱侧弯中心女性志愿者1名,12岁,诊断为AIS,三弯畸形,顶椎分别为第3胸椎、第11胸椎、第4腰椎,Cobb角分别为20°、15°、18°。该志愿者既往无脊柱外伤史及手术史,未合并其他脊柱损伤。收集该志愿者螺旋CT断层图片489张,扫描层厚1.3 mm,扫描范围从寰椎至骶髂关节平面,包含完整的椎骨、横突、棘突等结构,获得的原始CT影像数据以DICOM格式保存输出。
1.1.2 实验软件 医学影像三维重建软件Mimics 21.0(Materialise公司,比利时),逆向工程软件Geomagic Wrap 2017(Raindrop Geomagic公司,美国),计算机辅助设计(computer-aided design,CAD)软件Solidworks 2018(达索系统,法国),有限元力学分析计算机软件Ansys 17.0(ANSYS公司,美国)。
1.2 研究方法
1.2.1 构建有限元模型(1)首先建立三维模型。将收集到的DICOM格式原始CT影像数据导入Mimics 21.0软件,新建蒙版(Mask)调节灰度值(单位:HU),以获得最接近实际需求的骨骼状态,通过蒙版编辑(Edit Mask)对所需提取的骨与其他骨进行分割,利用三维计算功能(Calculate)对需要提取模型的蒙版进行计算导出脊柱的三维模型(图1)。利用三维工具(3D tools)的光滑功能(smooth)对三维模型进行处理,选择对应的三维模型后再以STL格式将其导出。(2)网格划分、优化平滑。将所得的STL文件导入Geomagic软件中,默认单位为mm,网格大小确定为0.4。对骨骼进行整体的“删除钉状物”、“去除特征”、“砂纸”、“减少噪音”、“快速光滑”、“松弛”等操作,使网格更加平滑(图2)。之后通过“网格医生”对骨骼整体进行检查,当网格医生“确诊”无误后,便可以使用“精确曲面”功能进行曲面的操作,主要包括轮廓线的编辑和曲面片的基本编辑。生成曲面片之后,进行“构造格栅”操作,当格栅的构造不存在问题以后,就可以直接进行拟合曲面了。将上述完成的模型分别以WRAP格式及STP格式文件导出进行保存。(3)装配、组合、绘制、调整。在Solidworks软件中创建装配体,把所有的骨骼零件导入装配体中,使用“配合”功能来固定零件。当装配件不存在干涉后,到“零件”格式下对皮质骨和松质骨进行组合,绘制软骨、椎间盘(含纤维环和髓核)、终板,其中髓核体积约占椎间盘的50%~60%。(4)材料赋值、添加韧带。根据不同材料的弹性模量与泊松比,对皮质骨、松质骨、椎间盘、关节软骨等不同部位的实体进行材料赋值,并使用Spring弹簧单元替代韧带。该模型设置了前纵韧带、后纵韧带、黄韧带、横突间韧带、棘上韧带和棘间韧带,通过定义弹簧的刚度值来代表各韧带的约束力大小[11],见表1~2。
图2 椎体三维有限元模型a:未经平滑处理的椎体 b:经Geomagic Wrap平滑处理后的椎体Fig.2 3D finite element model of vertebral bodya:the vertebra before smoothing;b:the vertebra after smoothing by Geomagic Wrap
表1 AIS患者三维有限元模型不同部位的材料属性Tab.1 Material properties of different parts of the AIS3D finite
表2 不同韧带的弹簧刚度(N/mm)Tab.2 Spring stiffness of different ligaments(N/mm)
1.2.2 有效性验证(1)几何形态验证。将构建好的AIS患者脊柱有限元模型与该患者原始的CT数据图像进行重叠比对,从几何形态上比较模型与实体CT图像的符合情况。(2)力学数据验证。结合该患者畸形顶椎所处的节段,从有限元模型中将T1~4、T9~12、L1~4分段提取出来,包括其间的椎间盘及韧带,将下端椎体的下表面各向活动进行约束,然后在上端椎体的上表面施加不同方向的1.5 N·m纯力矩,模拟屈伸、侧屈、旋转运动,分别测量上述3个节段在施加载荷后各方向上的活动度(range of motion,ROM),将活动度情况与Busscher等[12]的研究结果进行比较。同理进行L1~5节段验证,从有限元模型中将L1~5节段提取出来,在L1椎体的上表面向下施加均匀分布的400 N载荷,约束L5椎体下表面各向活动,施加不同方向的1.5 N·m纯力矩,模拟侧屈、旋转运动,测量该节段各方向的ROM,将活动度情况与Zander等[13]的实验值进行比较。同理,从有限元模型中将C7~L1节段提取出来,在C7椎体的上表面向下施加均匀分布的40 N载荷,约束L1椎体下表面各向活动,施加2 N·m的纯力矩模拟旋转运动,测量该节段的ROM,将活动度情况与Watkins等[14]的研究结果进行比较。
1.2.3 应力分布 约束有限元模型的骶髂关节下表面,在C6的上表面向下施加均匀分布的40 N载荷,并在该表面施加不同方向的1.5 N·m纯力矩,模拟前屈、后伸、左侧屈、右侧屈、左旋转、右旋转6种运动,提取该模型各椎体的应力分布结果。
2 结果
2.1 有限元模型
本实验成功建立了轻度AIS第6颈椎至第5腰椎的三维有限元模型,模型的总节点数为2561811个、总单元数为1547806个(图3~4)。
图3 AIS患者的三维有限元模型(正面观)图4 AIS患者的三维有限元模型(网格划分)Fig.3 3D finite element model of the AIS(Anterior view)Fig.4 3D finite element model of the AIS(Grid division)
2.2 模型有效性
2.2.1 几何形态的验证 将AIS患者脊柱有限元模型与其CT数据图像进行重叠比对,可见二者具有良好的一致性(图5)。
图5 AIS患者脊柱三维有限元模型与CT图像的比对情况Fig.5 Comparison between CTimages and 3D finite element model of the AIS
2.2.2 力学数据验证(1)T1~4、T9~12、L1~4分节段验证。该有限元模型各节段在加载1.5 N·m纯力矩后,模拟屈伸、侧屈、旋转运动,活动度与Busscher等[12]的研究结果对比,胸椎上段的旋转活动和腰椎的侧屈、旋转活动稍有差异(图6)。(2)提取该有限元模型L1~5腰椎节段,加载1.5 N·m纯力矩后,模拟侧屈、旋转运动,活动度与Zander等[13]的实验基本一致(图7)。(3)提取该有限元模型C7~L1节段,加载2 N·m纯力矩后,模拟旋转运动,活动度与Watkins等[14]的研究结果基本一致(图8)。
图6 不同节段活动度 本研究与Busscher等[12]对比Fig.6 Comparison of the range of motion in different segment(thisstudy versus Busscher)
图7 L1~~L5节段活动度 本研究与Zander等[13]对比Fig.7 Comparison of the range of motion in L1~L5 segment(this study versus Rohlmann)
图8 C7~L1节段旋转活动度 本研究与Watkins等[14]对比Fig.8 Comparison of the range of motion of axial rotation in C7~L1 segment(this study versus Watkins)
2.3 应力分布
2.3.1 全节段椎体应力变化 模拟该有限元模型进行屈伸、侧屈、旋转运动后,各椎体应力与静态时应力进行对比的结果见图9。
图9 应力分布对比Fig.9 Comparison of stress distribution
2.3.2 特征椎体应力变化 提取畸形最显著的第3胸椎处的应力分布云图(图10、11),可见各向活动时的应力变化均较静态时有较大差异,其中后伸、旋转活动时,第2~4胸椎的应力变化趋势与静态时相反。
图10 T3节段的应力分布云图(a:静态b:后伸c:左旋d:右旋)Fig.10 Stress distribution in T3 segment(a:static state;b:backward extension;c:left rotation;d:right rotation)
3 讨论
利用有限元分析技术进行AIS相关研究的基础是构建一个有效的有限元模型。对脊柱有限元模型的有效性进行验证的方法主要有4种[15],包括将有限元模型进行相同条件下的尸体标本力学实验、对比力学实验数据,与有限元模型相关文献所报道的实验结果进行比较,模拟手术操作并与术后相关数据进行比较,以不同体位下的脊柱正侧位片相关数据做评估。由于AIS患者的畸形形态各不相同,无论是以往的尸体标本或者目前的真实案例,都难以找到身体条件、疾病特点几近类似的样本;加之轻度患者无手术指征,因此难以通过常规的方法去严格验证轻度AIS患者有限元模型的有效性。结合病例特点和研究目的,本研究通过将有限元模型与该患者原始CT图像作比较及观察模型在不同工况下不同节段的脊柱活动度是否与既往文献的报道相仿两方面来进行模型的有效性验证。
从几何形态的对比图可看出,本研究构建的有限元模型与原始CT影像数据所得的骨骼形态具有良好的一致性,体现出了患者的三弯畸形且侧凸角度均较小,直观表明该模型有效。
本研究的对象是一名三弯畸形的AIS患者,因此采取分节段的方式进行力学数据的验证。结合研究对象的顶椎位置和既往研究报道,首先提取T1~4、T9~12、L1~43个节段进行验证,结果发现胸椎上段的旋转活动和腰椎的侧屈、旋转活动情况与既往Busscher实验结果存在1°~2°的差异。由于脊柱的活动度主要存在于颈椎和腰椎,为探究上述差异的原因,本实验提取该模型的L1~5、C7~L1节段,以上述节段的活动度代表腰段和颈腰段的活动度。
图11 T3节段的应力分布对比(后伸、旋转与静态)Fig.11 Comparison of stress distribution in T3 segment(backward extension and axial rotation versus static state)
既往对脊柱的研究多聚焦于特定节段,如颈椎、腰椎,而全脊柱的力学实验较少,因此本研究首先选取2009年Busscher等[12]利用尸体T1~L4节段进行的生物力学测试结果进行比较。但是,Busscher等[12]采用的标本为老年男性尸体(死亡时的平均年龄为72岁),而本研究的志愿者是一名12岁的少女,为消除年龄导致的脊柱活动度差异,本研究将实验力矩进行减半处理,其余实验条件与Busscher等[12]保持一致。从结果部分可以看出,本研究构建模型的活动度情况与Busscher等[12]实验结果差异不大,较明显的差异存在于旋转运动上。因此,本实验又提取模型的L1~5、C7~L1节段进行模拟运动,并将活动度情况与既往类似研究继续进行对比。
通过与Busscher等[12]全脊柱分节段的实验结果进行对比,可以发现,本研究构建的有限元模型中T1~4节段的屈伸和侧屈活动、T9~12节段的各向活动及L1~4节段的屈伸活动情况与Busscher的结果一致;进一步将该模型的L1~5节段的侧屈和旋转活动、C7~L1节段的旋转活动情况与Zander等[13]、Watkins等[14]的研究结果进行对比,发现上述活动情况高度一致。
综合分析以上3种运动的模拟结果,认为本研究构建的轻度AIS患者的有限元模型是有效的,可以将该模型用于轻度AIS的进一步研究中。
从该有限元模型全节段椎体应力及特征椎体应力与静态时应力进行对比的结果,发现畸形最显著的第3胸椎(Cobb角20°),在模拟后伸、旋转活动时,第2~4胸椎的应力变化趋势与静态时的应力呈现相反的趋势,提示该患者在进行后伸及旋转活动时,椎体的应力变化可能有利于改变畸形。