本研究基于计算机断层扫描(CT)图像数据建立并验证正常人颈椎C4~7三维有限元模型, 为研究中医手法治疗颈椎慢性疾病的生物力学机制提供模型平台。基于受试者颈部CT图像, 依次运用Mimics 17.0、Geomagic 12.0及Abaqus 6.13等软件创建正常人C4~7节段有限元模型。在模型上分别模拟前屈、后伸、左右侧弯和左右旋转工况, 计算椎体间相对动度(ROM), 将计算的结果与文献结果进行对比分析, 并观察模型在1 Nm载荷下6种工况下模型的主要应力分布情况。本研究成功建立了正常人颈椎C4~7三维有限元模型, 共包含591 459单元、121 446节点, 模拟了椎体、椎间盘、韧带、关节等几何结构与材料特性。模型在前后屈伸、左右侧弯和左右旋转6种工况下的ROM与实验研究数据基本一致, 在1 Nm扭矩或弯矩载荷下, 模型主要应力分布基本反映了正常人颈椎生理活动时的主要应力分布情况。本研究建立的正常人颈椎C4~7三维有限元模型精确逼真, 符合颈椎的生物力学特性, 可用于研究中医手法治疗颈椎慢性疾病的生物力学分析。
引用本文： 邓真, 王辉昊, 牛文鑫, 兰天鹰, 王宽, 詹红生. 正常人下颈椎C4~7节段三维有限元模型的建立与验证. 生物医学工程学杂志, 2016, 33(4): 652-658. doi: 10.7507/1001-5515.20160108 复制
引言
颈椎病是目前临床最常见的疾病之一。近年来，随着移动便携式设备的快速发展和工作生活方式的改变，颈椎病的发病越来越年轻化，并且没有减缓的趋势[1]。下颈椎(C4~7节段)上下连接着颈椎活动幅度和角度最大的上颈椎以及作为颈椎活动固定基础的胸椎，其稳定性相对较差，生理功能比较特殊，同时，下颈椎发生退行性病变及创伤的概率较高。为了了解颈椎的生物力学特性和损伤机制，指导临床对颈椎病的预防、诊断和治疗，早期研究者常用离体实验和在体实验的方式进行，但这些研究方法实验步骤复杂、周期长、效率低并且不能获得脊柱内部结构的应力特点，因此，多种生物力学模型应运而生。随着数字化技术的不断提高，越来越多的学者采用数值模型计算的方法进行应力和应变分析[2]。其中，三维有限元模型可以克服离体、在体实验或模型的缺点，近年来被广泛应用于脊柱生物力学方向的研究[3]。
早在1994年，Bozic等[4]就依据计算机断层扫描(computed tomography, CT)数据建立了单个C4椎体的三维有限元模型，并对其进行了简单的应力分析；随后Yoganandan等[5]建立了相对完整的C4~6颈椎节段有限元模型，并与离体实验的结果进行了比较。此后，国外涌现大量不同节段、不同材料属性和不同有限元单元类型的颈椎三维有限元模型[6]，并且将模型结果与离体实验结果对比验证。国内在1999年，李雪迎等[7]建立了颈椎C1~T1节段的三维有限元模型，通过模拟不同力学条件下颈椎牵引的生物力学过程，得出最有效的颈椎牵引条件。随着计算机技术的不断提高，越来越多的研究人员加入和参与到该领域，构建了不同特点、应用于不同研究领域的各式颈椎三维有限元模型。2015年，王辉昊等[8]基于颈椎CT血管造影图像建立了人体全颈椎及椎动脉流固耦合模型，从解剖形态和活动功能对模型进行了验证，该模型为深入了解颈椎及其附属结构的生物力学机制提供了新思路。
中医骨伤科手法是中医治疗颈椎病的一大特色疗法，其有效性得到了全球广泛的认同[9-10]。然而，针对手法的生物力学机制研究相对有限，为了科学、合理地解释和阐述中医手法治疗颈椎病的机制，结合目前脊柱有限元领域已有的研究基础，选择运用有限元模型来模拟和观察中医手法作用于颈椎的结构应力特点等生物力学机制成为了必然的趋势。
基于以上背景和研究基础，本研究建立了正常人颈椎C4~7节段有限元模型，通过静态加载生理载荷，计算相邻椎体间的相对活动度(range of motion，ROM)，并将结果与离体实验及有限元模型结果进行比较验证，为后续研究手法的生物力学作用奠定实验基础并提供模型平台。
1 材料和方法
1.1 数据获取
招募1名健康女性志愿者(30岁，身高163 cm，体重62 kg)，既往无颈椎病病史，无颈椎外伤史及手术史，没有颈痛及上肢症状。采用64层螺旋CT机进行颈椎扫描；志愿者取仰卧位，颈肩背放松，保持颈椎中立位，调整扫描断面与身体垂直。扫描条件：140 kV，200 mA，层厚0.624 mm。分辨率512×512像素，共获得162幅二维断层图像，以DICOM格式保存。
1.2 有限元建模
将志愿者CT数据导入比利时Materialise公司的交互式医学影像控制系统(Mimics17.0软件)，通过界定阈值、形成蒙板、区域增长、编辑等过程，提取颈椎C4~7各个节段，并生成三维模型，使用Remesh功能对三维模型进行初步光顺处理。将结果以STL格式导入逆向工程软件Geomagc12.0。在Geomagic12.0中通过点云、多边形、精确曲面等模块，对模型进行降噪、填充、去除特征及NURBS曲面片拟合处理，将结果以IGES格式保存。
将IGES格式椎体导入有限元分析软件Abaqus6.13中，在其前处理Part和Assembly模块将相邻两个椎体装配并进行布尔运算，得到两个椎体之间与终板紧密贴合的椎间盘模型，并对其进行相应区域划分，划分为髓核区域和纤维环区域，结果以IGES格式保存。
采用只受拉不受压的truss单元对韧带进行模拟，共建立5种韧带，包括前纵韧带(anterior longitudinal ligament, ALL)、后纵韧带(posterior longitudinal ligament, PLL)、黄韧带(ligamentum flavum, LF)、关节囊韧带(capsular ligament, CL)、棘间韧带(interspinous ligament, ISL)。韧带的连接位置以及相应韧带truss单元的数量均参考文献[11-12]。本模型采用非线性面-面无摩擦接触关系模拟关节间的相互作用。椎体使用三维10节点实体单元(C3D10)进行划分，椎间盘(包括髓核和纤维环)均采用增强沙漏的C3D8R单元进行划分。椎间盘与椎体之间通过共用终板面网格实现共节点，从而完成两者之间的关系。韧带使用truss单元进行定义和网格划分。单元类型及材料属性参数参考文献[13-14]，详细数据见表 1。
表1
颈椎有限元模型各部件单元类型及材料属性
Table1.
Mechanical properties of the material and element types of different parts of the cervical spine model
1.3 有限元分析
约束C7下终板所有节点的全部6个自由度作为边界条件。在C4上终板建立一个参考点，该点位于颈椎旋转中心上，将其与C4上终板所有单元节点耦合。根据右手定律，对参考点分别施加±0.5、±1、±1.5、±2 Nm纯扭矩，模拟颈椎前屈、后伸、左轴向旋转、右轴向旋转、左侧弯、右侧弯6种工况下的活动。前屈后伸的方向参考X、Y、Z全局坐标(X-Y平面为水平面、X-Z平面为冠状面、Y-Z平面为矢状面)，扭转方向参考颈曲切线方向，侧弯时垂直于颈曲方向并与矢状面平行。通过边界条件和加载条件的设置，对模型进行6种工况加载，分析不同条件下模型各部分的应力应变和位移。
2 结果
2.1 有限元模型
建立的正常人下颈椎C4~7节段三维有限元模型模拟了4个椎体、3个椎间盘、5种韧带，共计591 459单元和121 446节点，清晰完整地模拟了颈椎的几何特性和内在材料属性，见图 1。
图1
有限元模型
Figure1.
Three-dimensional finite element model
2.2 有限元模型的验证
本研究共计算了模型在前屈、后伸、左轴向旋转、右轴向旋转、左侧弯、右侧弯6种工况下的关节ROM。为了保证本研究验证结果的可靠性，课题组通过大量文献回顾，纳入了Moroney等[15]、Wheeldon等[16]、Nightingale等[17-18]、Panjabi等[19]和Lysell[20]等共6篇经典且公认较好的离体实验文献，将本研究结果与它们进行比较。颈椎相应节段ROM如图 2~4所示。
图2
不同颈椎节段在不同前屈(正)和后伸(负)弯矩下的关节活动度比较
Figure2.
Comparison of ROM of the present model response against in vitro studies in different segments under flexion (positive) and extension (negative)
图3
不同颈椎节段在1 Nm扭矩下轴向旋转的关节活动度比较
Figure3.
Comparison of ROM of the present model response against in vitro studies in different segments under axial rotation with the load of 1 Nm
图4
不同颈椎节段在1 Nm扭矩下侧向弯曲的关节活动度比较
Figure4.
Comparison of ROM of the present model response against in vitro studies in different segments under lateral bending with the load of 1 Nm
在不同大小和方向的弯矩作用下，本模型表现出非对称性，即前屈时ROM略大于后伸时ROM，模型与Wheeldon等[16]及Nightingale等[17-18]离体实验结果一致性较好，前屈与后伸工况下，最大的ROM均发生在C4、5节段，最小的ROM为C6、7节段。在前屈时(弯矩为正)，本模型C4~6节段ROM与Wheeldon等[16]离体实验数据更接近，后伸时(弯矩为负)本模型C4~6节段则与Nightingale等[17-18]离体实验数据更接近。C6、7节段则与两个离体实验均有较好的一致性，在数据上与Nightingale等[17-18]离体实验数据更接近。
在轴向旋转工况时，本模型最大的ROM同样发生在C4、5节段，然后依次递减。C4~7的ROM均比Moroney等[15]离体测量结果要小，C4、5的ROM比Panjabi等[19]离体实验数据要小，而C5~7的ROM则与Panjabi等[19]离体实验的数据一致。
在左右侧弯工况时，最大ROM依然发生在C4、5节段，最小为C6、7节段。模型结果整体上比Lysell[20]的测量结果要大，比Moroney等[15]实验结果要小；而与Panjabi等[19]的离体实验相比，C4、5节段结果偏小，C5~7节段结果一致。
2.3 颈椎C4~7模型的应力分布
对模型6种工况下施加的扭矩或弯矩均为1 Nm，以von Mises等效应力作为主要指标，对模型主要部件(椎体和椎间盘)进行应力分布情况观察。结果如图 5所示。
图5
不同工况1 Nm扭矩或弯矩作用下，颈椎C4~7模型椎体和椎间盘应力分布情况
Figure5.
Stress nephogram of C4-7 cervical spinal segments and intervertebral disc under different status
前屈时，椎体最大应力集中在C7皮质骨上侧、椎弓根、椎体后侧以及后纵韧带和棘间韧带连接等区域，主要是由于颈椎前屈活动时，C7作为基准面承受由上而下传递的载荷，同时后纵韧带和棘间韧带相对于前纵韧带而言，能更有效地制约颈椎的前屈；椎间盘主要应力集中在终板前侧和后侧，前侧由于受到压力，应力由上到下逐次增加，而后侧由于受到拉力，应力由上到下逐渐减少。后伸时，椎体前侧应力相对于后侧更为集中和明显，这是由于前纵韧带对抗椎体后伸运动产生的效果；同时，椎间盘主要应力集中在终板后侧，这是由于后伸时，前纵韧带分担了大多数拉力，椎间盘前侧受到的拉力较小的缘故。左右侧弯时，最大应力分布在两侧椎体、椎弓根和椎间盘连接处以及小关节处；椎间盘受力也集中在两侧，同侧更明显。左右旋转时，椎体应力主要集中在旋转同侧椎体、小关节及椎弓根区域；椎间盘应力呈近似“U”形分布，髓核区域受力最小，前侧纤维环受力明显小于其他三侧纤维环区域，但最大应力集中在椎间盘的侧后方。
3 讨论
颈椎是连接大脑与躯干的唯一桥梁，颈部创伤和损伤会给人体带来巨大的损害和生活质量的下降。下颈椎(C4~7节段)作为整个颈椎序列的根基，上与颈椎活动幅度和角度最大的上颈椎相连，下与相对基准的胸椎相连，其稳定性相对较差，生理功能比较特殊，临床上发生退行性病变、劳损性损伤及创伤的概率较高。目前，颈部劳损性损伤多采用中医手法治疗[21]，临床疗效明确，但其所利用的生物力学机制和作用机制研究却仍显匮乏。
本研究的目的是通过构建正常人C4~7节段三维有限元模型，为研究中医手法治疗颈椎慢性疾病提供模型平台。数字虚拟模型的构建不能完全模拟和替代离体实验的结果，但如若构建的有限元模型能够反映大部分离体实验测量得到的相关信息，那么运用该模型进行相关生物力学分析就能在不使用尸体标本进行实验的前提下，反映和预测相关生物力学问题。本研究详细地构建了正常人C4~7节段三维有限元模型，充分模拟了椎体、椎间盘、韧带、钩椎关节的生理功能和材料属性。目前，多数研究[22-24]建立的颈椎有限元模型忽略了钩椎关节的影响，对其不进行任何设定，事实上，中医骨伤科所强调的颈椎“骨错缝筋出槽”现象常发生在颈椎钩椎关节部位，该部位的结构和功能是不能被忽略的，很有可能是引起颈椎症状的关键所在，中医正骨手法也擅长对该部位进行整复治疗[25]。本研究认为不能完全忽视钩椎关节对有限元模型造成的影响，因此，我们将其设定为无摩擦面-面接触来模拟颈椎钩椎关节的状态。
由于本研究采取的验证方法为文献结果对比，课题组更注重本模型是否能反映离体实验测量的结果，因此我们纳入的对比文献均为经典且公认较好的离体实验。从实验总体结果来看，本研究所建立的颈椎模型能够有效地反映人体颈椎的各项活动，在不同载荷作用下，本模型计算的关节ROM的结果与文献报道一致，在不同大小弯矩下，本模型结果与Wheeldon等[16]和Nightingale等[17-18]的离体实验均有较好的一致性。
本模型在部分情况下计算结果与离体实验并不一致。我们分析导致该结果发生的原因可能是：离体实验数据的测量采用的是不同年龄、不同椎间盘退行性变的尸体标本，计算的是多个尸体标本不同节段的平均值，而本研究所选取的是一名30岁健康女性志愿者的CT扫描图像数据，椎间盘、椎体无任何退变。同时由于本模型在构建过程中存在一定的模型简化步骤，所选取的材料属性来源于有限的文献报道，加之不同实验所选取的标本不同，几何特性和材料属性在本质上是有差异的。就本文所选取的三篇离体实验文献([16-18])而言，其结果也并非完全一致。Zhang等[14]对三项离体实验[19, 26-27]结果进行比较分析发现，三项研究在所得数值上存在较大的差异，但是这三项研究数据所反映的颈椎椎体之间ROM的总体趋势是一致的。因此，总地来说本模型所反映的结果是真实、可以接受的。
在模型的应力方面，本模型尽量接近真实地体现了正常人颈椎在进行6个方向活动时主要的应力区域，不同的活动姿势时，颈椎应力分布有其各自的特点，也体现了不同解剖结构在人体运动时所承担的主要功能。
总而言之，通过验证，本模型所构建的C4~7节段三维有限元模型是可靠、有效的，能够反映正常人颈椎的功能活动，其所展现的6种工况下的ROM与离体实验数据结果相比，总体趋势是一致的，在C6、7节段则完全反映了离体实验的结果。该模型可以用于分析椎体、椎间盘、小关节和颈部主要韧带的生物力学特性。而中医手法治疗颈椎病正是将治疗手法的力作用在椎体上，通过正骨手法改善椎间盘、韧带以及小关节与周围肌肉组织的动力平衡与静力平衡，从而达到治疗目的[28]。以往的研究多数是从临床治疗效果反推其机制，尤其是作为核心特色治疗技术的中医骨伤手法，其作用的部位清晰，疗效明确，但是其具体作用过程和作用机制却鲜有研究。三维有限元模型正好补充了这块短板，通过有限元模型可以完整地模拟手法治疗的过程，观察内部结构应力、应变和解剖结构位置改变，从而为完整地阐述中医手法疗效提供了必不可少的实验基础。
人体颈椎结构很复杂，本研究尽可能真实地完成了颈椎几何模型和材料属性的模拟，但在模型构建过程中不可避免地使用了一些简化步骤，如部分几何特征去除、几何三角面片数量缩减等。目前，对于钩椎关节的材料属性的定义以及关节面的设定尚存在争议，也有研究[29]使用特殊类型的设定来模拟钩椎关节，而本研究以无摩擦面-面接触关系来定义钩椎关节之间的相互作用也只能部分模拟钩椎关节之间的生理功能，这样设定是否完全正确尚需进行进一步研究和讨论。此外，由于肌肉组织参数的设定非常复杂，其作用也不容易控制[30-31]，因此本研究忽略了肌肉对颈椎的影响，并未对肌肉组织进行有限元模型重建。颈椎主要相关的血管(如椎动脉)的模拟从几何上没有困难，但在进行分析时涉及流固耦合分析领域，目前该方向的研究较少[8]，其模型的构建以及与椎体之间的关系设定等均是难点，本研究同样未对血管进行模拟和分析。因此，为了更好、更完整地模拟颈椎真实的生理和功能状况，肌肉和血管的纳入以及颈椎解剖上更细微结构的模拟将成为日后研究的重心。
引言
颈椎病是目前临床最常见的疾病之一。近年来，随着移动便携式设备的快速发展和工作生活方式的改变，颈椎病的发病越来越年轻化，并且没有减缓的趋势[1]。下颈椎(C4~7节段)上下连接着颈椎活动幅度和角度最大的上颈椎以及作为颈椎活动固定基础的胸椎，其稳定性相对较差，生理功能比较特殊，同时，下颈椎发生退行性病变及创伤的概率较高。为了了解颈椎的生物力学特性和损伤机制，指导临床对颈椎病的预防、诊断和治疗，早期研究者常用离体实验和在体实验的方式进行，但这些研究方法实验步骤复杂、周期长、效率低并且不能获得脊柱内部结构的应力特点，因此，多种生物力学模型应运而生。随着数字化技术的不断提高，越来越多的学者采用数值模型计算的方法进行应力和应变分析[2]。其中，三维有限元模型可以克服离体、在体实验或模型的缺点，近年来被广泛应用于脊柱生物力学方向的研究[3]。
早在1994年，Bozic等[4]就依据计算机断层扫描(computed tomography, CT)数据建立了单个C4椎体的三维有限元模型，并对其进行了简单的应力分析；随后Yoganandan等[5]建立了相对完整的C4~6颈椎节段有限元模型，并与离体实验的结果进行了比较。此后，国外涌现大量不同节段、不同材料属性和不同有限元单元类型的颈椎三维有限元模型[6]，并且将模型结果与离体实验结果对比验证。国内在1999年，李雪迎等[7]建立了颈椎C1~T1节段的三维有限元模型，通过模拟不同力学条件下颈椎牵引的生物力学过程，得出最有效的颈椎牵引条件。随着计算机技术的不断提高，越来越多的研究人员加入和参与到该领域，构建了不同特点、应用于不同研究领域的各式颈椎三维有限元模型。2015年，王辉昊等[8]基于颈椎CT血管造影图像建立了人体全颈椎及椎动脉流固耦合模型，从解剖形态和活动功能对模型进行了验证，该模型为深入了解颈椎及其附属结构的生物力学机制提供了新思路。
中医骨伤科手法是中医治疗颈椎病的一大特色疗法，其有效性得到了全球广泛的认同[9-10]。然而，针对手法的生物力学机制研究相对有限，为了科学、合理地解释和阐述中医手法治疗颈椎病的机制，结合目前脊柱有限元领域已有的研究基础，选择运用有限元模型来模拟和观察中医手法作用于颈椎的结构应力特点等生物力学机制成为了必然的趋势。
基于以上背景和研究基础，本研究建立了正常人颈椎C4~7节段有限元模型，通过静态加载生理载荷，计算相邻椎体间的相对活动度(range of motion，ROM)，并将结果与离体实验及有限元模型结果进行比较验证，为后续研究手法的生物力学作用奠定实验基础并提供模型平台。
1 材料和方法
1.1 数据获取
招募1名健康女性志愿者(30岁，身高163 cm，体重62 kg)，既往无颈椎病病史，无颈椎外伤史及手术史，没有颈痛及上肢症状。采用64层螺旋CT机进行颈椎扫描；志愿者取仰卧位，颈肩背放松，保持颈椎中立位，调整扫描断面与身体垂直。扫描条件：140 kV，200 mA，层厚0.624 mm。分辨率512×512像素，共获得162幅二维断层图像，以DICOM格式保存。
1.2 有限元建模
将志愿者CT数据导入比利时Materialise公司的交互式医学影像控制系统(Mimics17.0软件)，通过界定阈值、形成蒙板、区域增长、编辑等过程，提取颈椎C4~7各个节段，并生成三维模型，使用Remesh功能对三维模型进行初步光顺处理。将结果以STL格式导入逆向工程软件Geomagc12.0。在Geomagic12.0中通过点云、多边形、精确曲面等模块，对模型进行降噪、填充、去除特征及NURBS曲面片拟合处理，将结果以IGES格式保存。
将IGES格式椎体导入有限元分析软件Abaqus6.13中，在其前处理Part和Assembly模块将相邻两个椎体装配并进行布尔运算，得到两个椎体之间与终板紧密贴合的椎间盘模型，并对其进行相应区域划分，划分为髓核区域和纤维环区域，结果以IGES格式保存。
采用只受拉不受压的truss单元对韧带进行模拟，共建立5种韧带，包括前纵韧带(anterior longitudinal ligament, ALL)、后纵韧带(posterior longitudinal ligament, PLL)、黄韧带(ligamentum flavum, LF)、关节囊韧带(capsular ligament, CL)、棘间韧带(interspinous ligament, ISL)。韧带的连接位置以及相应韧带truss单元的数量均参考文献[11-12]。本模型采用非线性面-面无摩擦接触关系模拟关节间的相互作用。椎体使用三维10节点实体单元(C3D10)进行划分，椎间盘(包括髓核和纤维环)均采用增强沙漏的C3D8R单元进行划分。椎间盘与椎体之间通过共用终板面网格实现共节点，从而完成两者之间的关系。韧带使用truss单元进行定义和网格划分。单元类型及材料属性参数参考文献[13-14]，详细数据见表 1。
表1
颈椎有限元模型各部件单元类型及材料属性
Table1.
Mechanical properties of the material and element types of different parts of the cervical spine model
1.3 有限元分析
约束C7下终板所有节点的全部6个自由度作为边界条件。在C4上终板建立一个参考点，该点位于颈椎旋转中心上，将其与C4上终板所有单元节点耦合。根据右手定律，对参考点分别施加±0.5、±1、±1.5、±2 Nm纯扭矩，模拟颈椎前屈、后伸、左轴向旋转、右轴向旋转、左侧弯、右侧弯6种工况下的活动。前屈后伸的方向参考X、Y、Z全局坐标(X-Y平面为水平面、X-Z平面为冠状面、Y-Z平面为矢状面)，扭转方向参考颈曲切线方向，侧弯时垂直于颈曲方向并与矢状面平行。通过边界条件和加载条件的设置，对模型进行6种工况加载，分析不同条件下模型各部分的应力应变和位移。
2 结果
2.1 有限元模型
建立的正常人下颈椎C4~7节段三维有限元模型模拟了4个椎体、3个椎间盘、5种韧带，共计591 459单元和121 446节点，清晰完整地模拟了颈椎的几何特性和内在材料属性，见图 1。
图1
有限元模型
Figure1.
Three-dimensional finite element model
2.2 有限元模型的验证
本研究共计算了模型在前屈、后伸、左轴向旋转、右轴向旋转、左侧弯、右侧弯6种工况下的关节ROM。为了保证本研究验证结果的可靠性，课题组通过大量文献回顾，纳入了Moroney等[15]、Wheeldon等[16]、Nightingale等[17-18]、Panjabi等[19]和Lysell[20]等共6篇经典且公认较好的离体实验文献，将本研究结果与它们进行比较。颈椎相应节段ROM如图 2~4所示。
图2
不同颈椎节段在不同前屈(正)和后伸(负)弯矩下的关节活动度比较
Figure2.
Comparison of ROM of the present model response against in vitro studies in different segments under flexion (positive) and extension (negative)
图3
不同颈椎节段在1 Nm扭矩下轴向旋转的关节活动度比较
Figure3.
Comparison of ROM of the present model response against in vitro studies in different segments under axial rotation with the load of 1 Nm
图4
不同颈椎节段在1 Nm扭矩下侧向弯曲的关节活动度比较
Figure4.
Comparison of ROM of the present model response against in vitro studies in different segments under lateral bending with the load of 1 Nm
在不同大小和方向的弯矩作用下，本模型表现出非对称性，即前屈时ROM略大于后伸时ROM，模型与Wheeldon等[16]及Nightingale等[17-18]离体实验结果一致性较好，前屈与后伸工况下，最大的ROM均发生在C4、5节段，最小的ROM为C6、7节段。在前屈时(弯矩为正)，本模型C4~6节段ROM与Wheeldon等[16]离体实验数据更接近，后伸时(弯矩为负)本模型C4~6节段则与Nightingale等[17-18]离体实验数据更接近。C6、7节段则与两个离体实验均有较好的一致性，在数据上与Nightingale等[17-18]离体实验数据更接近。
在轴向旋转工况时，本模型最大的ROM同样发生在C4、5节段，然后依次递减。C4~7的ROM均比Moroney等[15]离体测量结果要小，C4、5的ROM比Panjabi等[19]离体实验数据要小，而C5~7的ROM则与Panjabi等[19]离体实验的数据一致。
在左右侧弯工况时，最大ROM依然发生在C4、5节段，最小为C6、7节段。模型结果整体上比Lysell[20]的测量结果要大，比Moroney等[15]实验结果要小；而与Panjabi等[19]的离体实验相比，C4、5节段结果偏小，C5~7节段结果一致。
2.3 颈椎C4~7模型的应力分布
对模型6种工况下施加的扭矩或弯矩均为1 Nm，以von Mises等效应力作为主要指标，对模型主要部件(椎体和椎间盘)进行应力分布情况观察。结果如图 5所示。
图5
不同工况1 Nm扭矩或弯矩作用下，颈椎C4~7模型椎体和椎间盘应力分布情况
Figure5.
Stress nephogram of C4-7 cervical spinal segments and intervertebral disc under different status
前屈时，椎体最大应力集中在C7皮质骨上侧、椎弓根、椎体后侧以及后纵韧带和棘间韧带连接等区域，主要是由于颈椎前屈活动时，C7作为基准面承受由上而下传递的载荷，同时后纵韧带和棘间韧带相对于前纵韧带而言，能更有效地制约颈椎的前屈；椎间盘主要应力集中在终板前侧和后侧，前侧由于受到压力，应力由上到下逐次增加，而后侧由于受到拉力，应力由上到下逐渐减少。后伸时，椎体前侧应力相对于后侧更为集中和明显，这是由于前纵韧带对抗椎体后伸运动产生的效果；同时，椎间盘主要应力集中在终板后侧，这是由于后伸时，前纵韧带分担了大多数拉力，椎间盘前侧受到的拉力较小的缘故。左右侧弯时，最大应力分布在两侧椎体、椎弓根和椎间盘连接处以及小关节处；椎间盘受力也集中在两侧，同侧更明显。左右旋转时，椎体应力主要集中在旋转同侧椎体、小关节及椎弓根区域；椎间盘应力呈近似“U”形分布，髓核区域受力最小，前侧纤维环受力明显小于其他三侧纤维环区域，但最大应力集中在椎间盘的侧后方。
3 讨论
颈椎是连接大脑与躯干的唯一桥梁，颈部创伤和损伤会给人体带来巨大的损害和生活质量的下降。下颈椎(C4~7节段)作为整个颈椎序列的根基，上与颈椎活动幅度和角度最大的上颈椎相连，下与相对基准的胸椎相连，其稳定性相对较差，生理功能比较特殊，临床上发生退行性病变、劳损性损伤及创伤的概率较高。目前，颈部劳损性损伤多采用中医手法治疗[21]，临床疗效明确，但其所利用的生物力学机制和作用机制研究却仍显匮乏。
本研究的目的是通过构建正常人C4~7节段三维有限元模型，为研究中医手法治疗颈椎慢性疾病提供模型平台。数字虚拟模型的构建不能完全模拟和替代离体实验的结果，但如若构建的有限元模型能够反映大部分离体实验测量得到的相关信息，那么运用该模型进行相关生物力学分析就能在不使用尸体标本进行实验的前提下，反映和预测相关生物力学问题。本研究详细地构建了正常人C4~7节段三维有限元模型，充分模拟了椎体、椎间盘、韧带、钩椎关节的生理功能和材料属性。目前，多数研究[22-24]建立的颈椎有限元模型忽略了钩椎关节的影响，对其不进行任何设定，事实上，中医骨伤科所强调的颈椎“骨错缝筋出槽”现象常发生在颈椎钩椎关节部位，该部位的结构和功能是不能被忽略的，很有可能是引起颈椎症状的关键所在，中医正骨手法也擅长对该部位进行整复治疗[25]。本研究认为不能完全忽视钩椎关节对有限元模型造成的影响，因此，我们将其设定为无摩擦面-面接触来模拟颈椎钩椎关节的状态。
由于本研究采取的验证方法为文献结果对比，课题组更注重本模型是否能反映离体实验测量的结果，因此我们纳入的对比文献均为经典且公认较好的离体实验。从实验总体结果来看，本研究所建立的颈椎模型能够有效地反映人体颈椎的各项活动，在不同载荷作用下，本模型计算的关节ROM的结果与文献报道一致，在不同大小弯矩下，本模型结果与Wheeldon等[16]和Nightingale等[17-18]的离体实验均有较好的一致性。
本模型在部分情况下计算结果与离体实验并不一致。我们分析导致该结果发生的原因可能是：离体实验数据的测量采用的是不同年龄、不同椎间盘退行性变的尸体标本，计算的是多个尸体标本不同节段的平均值，而本研究所选取的是一名30岁健康女性志愿者的CT扫描图像数据，椎间盘、椎体无任何退变。同时由于本模型在构建过程中存在一定的模型简化步骤，所选取的材料属性来源于有限的文献报道，加之不同实验所选取的标本不同，几何特性和材料属性在本质上是有差异的。就本文所选取的三篇离体实验文献([16-18])而言，其结果也并非完全一致。Zhang等[14]对三项离体实验[19, 26-27]结果进行比较分析发现，三项研究在所得数值上存在较大的差异，但是这三项研究数据所反映的颈椎椎体之间ROM的总体趋势是一致的。因此，总地来说本模型所反映的结果是真实、可以接受的。
在模型的应力方面，本模型尽量接近真实地体现了正常人颈椎在进行6个方向活动时主要的应力区域，不同的活动姿势时，颈椎应力分布有其各自的特点，也体现了不同解剖结构在人体运动时所承担的主要功能。
总而言之，通过验证，本模型所构建的C4~7节段三维有限元模型是可靠、有效的，能够反映正常人颈椎的功能活动，其所展现的6种工况下的ROM与离体实验数据结果相比，总体趋势是一致的，在C6、7节段则完全反映了离体实验的结果。该模型可以用于分析椎体、椎间盘、小关节和颈部主要韧带的生物力学特性。而中医手法治疗颈椎病正是将治疗手法的力作用在椎体上，通过正骨手法改善椎间盘、韧带以及小关节与周围肌肉组织的动力平衡与静力平衡，从而达到治疗目的[28]。以往的研究多数是从临床治疗效果反推其机制，尤其是作为核心特色治疗技术的中医骨伤手法，其作用的部位清晰，疗效明确，但是其具体作用过程和作用机制却鲜有研究。三维有限元模型正好补充了这块短板，通过有限元模型可以完整地模拟手法治疗的过程，观察内部结构应力、应变和解剖结构位置改变，从而为完整地阐述中医手法疗效提供了必不可少的实验基础。
人体颈椎结构很复杂，本研究尽可能真实地完成了颈椎几何模型和材料属性的模拟，但在模型构建过程中不可避免地使用了一些简化步骤，如部分几何特征去除、几何三角面片数量缩减等。目前，对于钩椎关节的材料属性的定义以及关节面的设定尚存在争议，也有研究[29]使用特殊类型的设定来模拟钩椎关节，而本研究以无摩擦面-面接触关系来定义钩椎关节之间的相互作用也只能部分模拟钩椎关节之间的生理功能，这样设定是否完全正确尚需进行进一步研究和讨论。此外，由于肌肉组织参数的设定非常复杂，其作用也不容易控制[30-31]，因此本研究忽略了肌肉对颈椎的影响，并未对肌肉组织进行有限元模型重建。颈椎主要相关的血管(如椎动脉)的模拟从几何上没有困难，但在进行分析时涉及流固耦合分析领域，目前该方向的研究较少[8]，其模型的构建以及与椎体之间的关系设定等均是难点，本研究同样未对血管进行模拟和分析。因此，为了更好、更完整地模拟颈椎真实的生理和功能状况，肌肉和血管的纳入以及颈椎解剖上更细微结构的模拟将成为日后研究的重心。
