輕度青少年特發性脊柱側凸有限元模型構建及分析

吳曉薇，李鑒軼,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胸椎的應力變化趨勢與靜態時的應力呈現相反的趨勢，提示該患者在進行后伸及旋轉活動時，椎體的應力變化可能有利于改變畸形。