人骨三維有限元分析的研究進(jìn)展

鄭 鋒（綜述），林海濱（審校）

（莆田學(xué)院附屬醫(yī)院 1急診科，2骨科，福建 莆田 351100）

有限元分析的基本原理是根據(jù)幾何外形、材料特性以及受力條件等因素將彈性物體離散為有限的體單元，這些體單元只在有限個(gè)節(jié)點(diǎn)上相交接，力通過結(jié)點(diǎn)傳遞，導(dǎo)致每個(gè)體單元的變形、任意體單元或節(jié)點(diǎn)的應(yīng)力分布可通過多種簡單的方程式來求解。常用Ansys及Abaqus等商業(yè)軟件對骨骼進(jìn)行有限元分析，具有高精度、低成本、零風(fēng)險(xiǎn)、可重復(fù)的優(yōu)點(diǎn)，迄今國內(nèi)外諸多學(xué)者從不同應(yīng)用目的對此進(jìn)行了大量研究，其有效性、優(yōu)越性已在基礎(chǔ)試驗(yàn)及臨床應(yīng)用中得到充分證明。

去軟組織的人骨有限元分析是其他復(fù)雜類型有限元分析的基礎(chǔ)，然而，即使是試驗(yàn)條件頗為單一的去軟組織人骨有限元分析，其實(shí)踐過程、條件設(shè)置、結(jié)果評判及臨床應(yīng)用等方面仍未能形成統(tǒng)一標(biāo)準(zhǔn)。該文針對去軟組織人骨的有限元分析進(jìn)行綜述，旨在對三維重建技術(shù)及有限元分析在骨骼力學(xué)研究的一般過程及實(shí)際應(yīng)用進(jìn)行歸納總結(jié)。

1 去軟組織人骨三維有限元分析的一般過程

1．1 三維重建 CT/磁共振成像（magnatic resonace imaging，MRI）以無創(chuàng)手段對人體進(jìn)行三維重建，針對人骨的有限元分析大體上包括CT/MRI掃描-三維重建-體網(wǎng)格劃分-賦予材料屬性-有限元分析等基本步驟。醫(yī)學(xué)三維重建軟件的基本原理，是以疊加的方式對掃描之Dicom格式圖像進(jìn)行三維重建，三維模型由多層面曲線構(gòu)成，每一層面曲線都是由大量規(guī)則、不規(guī)則的自由曲線構(gòu)成。重建質(zhì)量差的三維模型在有限元分析中缺乏相應(yīng)的方程式描述，任何細(xì)微的重建錯(cuò)誤都會(huì)使運(yùn)算量呈幾何級(jí)增加，甚至導(dǎo)致運(yùn)算失敗。三維重建的每一步驟，圖像編輯-重建三維模型-表面光滑處理等，都會(huì)影響三維模型質(zhì)量［1］，在很大程度上決定有限元分析與真實(shí)的接近程度。

受限于CT/MRI圖像的信噪比和空間分辨率，軟組織三維重建效果仍難以令人滿意［2］，仍然無法對肌腱、肌肉等軟組織進(jìn)行良好的三維重建，此外，軟組織的力學(xué)性能存在人種差異、個(gè)體差異、生理狀態(tài)及運(yùn)動(dòng)狀態(tài)的差異，其材料屬性往往難以通過試驗(yàn)手段進(jìn)行界定，涵括復(fù)合組織在內(nèi)的人骨有限元分析迄今仍難以實(shí)現(xiàn)。

1．2 體單元?jiǎng)澐?任何復(fù)雜的彈性物體均由無限個(gè)質(zhì)點(diǎn)組成，具有無限個(gè)自由度，從而不能求解，而有限元分析法可根據(jù)試驗(yàn)要求不同，將彈性物體分解成不同類型、一定數(shù)量的體單元，實(shí)現(xiàn)從無限到有限的轉(zhuǎn)變，從而達(dá)到求解的目的。基于三維重建后模型的形狀，Ansys可自動(dòng)化對骨骼三維模型進(jìn)行體網(wǎng)格劃分，一般所得體單元數(shù)量多在幾萬、十幾萬個(gè)甚至更多，一般可以達(dá)到有限元分析的需要。特定目的之研究如骨小梁有限元分析的高分辨率有限元，1個(gè)5 mm×5 mM×5 mm小梁骨微有限元模型可有數(shù)十萬個(gè)單元，其求解需要專用高效的方法及多臺(tái)計(jì)算機(jī)并行處理來分析［3］。有限元分析可根據(jù)彈性物體的幾何材料特性以及受力條件采用塊單元、殼單元、面單元、纜式單元等不同種類，亦有文獻(xiàn)對此進(jìn)行了描述［4］，不同單元類型并無本質(zhì)區(qū)別。

1．3 賦予材料屬性 從本質(zhì)來說，骨骼是由多種（趨于無限種）材料構(gòu)成的復(fù)合體，這些材料各具不同力學(xué)特性，并具有明顯個(gè)體差異，整骨有限元模型的材料屬性呈異質(zhì)性分布［5］。研究表明，從CT數(shù)據(jù)中提取邊緣光滑的骨表面三維幾何形態(tài)，經(jīng)軟件自動(dòng)網(wǎng)格化，由CT值得到的材料屬性被賦予到這些有限元網(wǎng)格上，這種方法比基于體素的方法能更加準(zhǔn)確地評估骨表面應(yīng)變［6］。

進(jìn)行骨骼有限元分析的重要條件之一是確定骨骼內(nèi)不同部位的材料屬性，其基本原理可簡單概括如下。①將骨小梁視為連續(xù)體。利用三維定量計(jì)算斷層掃描將每個(gè)體素直接轉(zhuǎn)換為立方體有限單元［7］，相應(yīng)的，骨小梁被視為一個(gè)連續(xù)體，通過給定的力學(xué)性能與密度之間的回歸關(guān)系，賦予代表小梁骨的單元的材料屬性［8］。②簡化骨骼材料特性。骨骼是一種各向異性的生物材料，但其各向異性彈性常數(shù)太多，故一般簡化為正交各向異性，同時(shí)由于骨骼各向異性較弱，故尚可進(jìn)一步簡化為各向同性［9，10］。③骨骼表觀密度與 CT值的關(guān)系。研究表明，骨骼表觀密度ρ與CT值（Hu）具有近似的線性關(guān)系［11］（灰度值 Gv=CT 值 Hu+1024）;骨骼材料特性與骨骼表觀密度存在冪指數(shù)關(guān)系的經(jīng)驗(yàn)公式［12］。綜合上述，可以用CT值-材料屬性經(jīng)驗(yàn)公式來進(jìn)行賦值后進(jìn)行有限元分析，目前尚無基于國人測量數(shù)據(jù)的此類經(jīng)驗(yàn)公式報(bào)道。

有學(xué)者認(rèn)為構(gòu)成骨骼的材料呈均勻分布，將骨骼簡單區(qū)分成密質(zhì)骨和松質(zhì)骨，分別賦予彈性模量和泊松比等材料屬性后進(jìn)行有限元分析［13］，其結(jié)果與真實(shí)生理狀態(tài)的骨骼相異。亦有學(xué)者將股骨材料屬性分為10種，據(jù)此進(jìn)行有限元分析［14］，但未進(jìn)行其有效性驗(yàn)證。張國棟等［1］通過理論推算結(jié)合體外力學(xué)驗(yàn)證的方式，論證了將骨骼三維模型材料屬性等分成10種可以達(dá)到有限元分析的目的。

1．4 三維有限元分析 根據(jù)試驗(yàn)?zāi)康?，在Ansys中可以對骨骼實(shí)施幾乎任意形式的力學(xué)分析，如加壓、扭曲、拉伸等靜態(tài)力學(xué)分析，評估在日常活動(dòng)的動(dòng)態(tài)力學(xué)分析，骨科手術(shù)模擬如內(nèi)固定植入物及假體力學(xué)分析，包括韌帶等軟組織在內(nèi)的關(guān)節(jié)力學(xué)分析等等。

2 骨強(qiáng)度有限元分析及骨折預(yù)測

骨骼材料在外力作用下抵抗永久變形和斷裂的能力稱為骨強(qiáng)度，骨強(qiáng)度概念與骨折密不可分，按外力作用的性質(zhì)不同，骨強(qiáng)度主要包括屈服強(qiáng)度、抗拉強(qiáng)度、抗壓強(qiáng)度、抗彎強(qiáng)度等。研究表明，當(dāng)外力超過骨強(qiáng)度時(shí)，就會(huì)發(fā)生骨折［15-17］。骨折的發(fā)生除了與骨骼材料分布、材料力學(xué)性能有關(guān)之外，與骨骼形狀密不可分［18］。對于去軟組織人骨的有限元分析，最值得關(guān)注的是骨強(qiáng)度及骨折預(yù)測。Keyak等［19］評估了基于定量CT的股骨近端有限元模型預(yù)測站立及側(cè)身跌倒時(shí)骨折部位及類型的能力，認(rèn)為有限元分析預(yù)測骨折部位的準(zhǔn)確率約為60%～70%。國外學(xué)者就分析和優(yōu)化椎體成形術(shù)進(jìn)行了大量研究，在提高脊柱骨折的治療效果方面，有限元分析技術(shù)具有一定優(yōu)勢［20］。

雖然有限元分析可以獲得應(yīng)力及應(yīng)變、節(jié)點(diǎn)位移、物體的整體剛度及應(yīng)變能密度等指標(biāo)，但骨折預(yù)測需要指定某些材料特性，如破壞特性，才有可能計(jì)算出骨折在什么部位、什么時(shí)候以及如何發(fā)生。骨折部位與破壞負(fù)載的確定取決于對骨組織及整骨破壞標(biāo)準(zhǔn)的選擇［19-22］。目前對這些情況的估計(jì)沒有標(biāo)準(zhǔn)化的算法。

3 骨密度-骨強(qiáng)度-骨折風(fēng)險(xiǎn)預(yù)測

有限元分析的類型因不同研究目的而異，其本質(zhì)原因仍然是骨強(qiáng)度，針對人骨來說，骨密度是反映骨強(qiáng)度較為理想的指標(biāo)，骨強(qiáng)度及骨折預(yù)測涉及非常龐大的計(jì)算量，其實(shí)踐過程需要進(jìn)行規(guī)范的培訓(xùn)，在當(dāng)前尚無標(biāo)準(zhǔn)化的應(yīng)用軟件的情況下，以骨密度指標(biāo)代替骨強(qiáng)度及骨折預(yù)測，似乎更具有臨床應(yīng)用的可行性。

一般認(rèn)為，骨密度是反映骨強(qiáng)度的有效指標(biāo)，甚至Johnell等［23］對12個(gè)人群，39 000個(gè)人進(jìn)行研究后，認(rèn)為骨骼強(qiáng)度80%取決于骨密度。有限元分析已表明在關(guān)鍵區(qū)域骨密度的略微增加就會(huì)提高側(cè)身跌倒姿勢下預(yù)期的股骨破壞負(fù)載，與更大范圍非特異的骨密度增加的效應(yīng)類似［22］。Van Rietbergen等［24］以微有限元模型來評估正常和疏松股骨的應(yīng)力及應(yīng)變分布。在行走時(shí)，與正常骨相比，疏松骨上的應(yīng)變量更大，其分布也更不均勻。因此，在疏松股骨中有較大比例的骨組織存在被破壞的危險(xiǎn)。

國外學(xué)者就骨質(zhì)疏松癥的骨密度-骨強(qiáng)度進(jìn)行了有益的探討［25］，致力于創(chuàng)建新的診斷標(biāo)準(zhǔn)來避免單純依靠骨密度測量所帶來的偏差。有限元分析已用于分析應(yīng)用甲狀旁腺激素和阿侖膦酸鈉治療絕經(jīng)后骨質(zhì)疏松女性患者骨強(qiáng)度的變化。研究表明，治療1年后、2年后，股骨強(qiáng)度有了明顯的提高，骨強(qiáng)度的提高主要源于骨小梁密度增加［26］。張國棟等［27，28］以股骨、脊柱進(jìn)行了基于三維重建技術(shù)及有限元分析的骨密度測量的探討，從方法學(xué)上實(shí)現(xiàn)了對骨骼表觀密度的測量，并為實(shí)現(xiàn)骨密度-骨強(qiáng)度-骨折風(fēng)險(xiǎn)預(yù)測的統(tǒng)一分析作初步準(zhǔn)備。骨密度-骨強(qiáng)度-骨折風(fēng)險(xiǎn)預(yù)測的統(tǒng)一分析可在數(shù)字化植入物等骨骼相關(guān)領(lǐng)域研究中應(yīng)用，在某些病理狀態(tài)下，如骨質(zhì)疏松癥、股骨頭壞死骨密度診斷，可作為輔助的量化指標(biāo)。

4 不足及展望

去軟組織人骨有限元分析是大多數(shù)人骨有限元分析的基礎(chǔ)以及最終體現(xiàn)，實(shí)現(xiàn)骨密度-骨強(qiáng)度-骨折預(yù)測的統(tǒng)一分析是數(shù)字骨科學(xué)領(lǐng)域值得關(guān)注的熱點(diǎn)，然而迄今仍無相應(yīng)的指標(biāo)及可靠結(jié)論。在一定荷載下，整骨骨折是否發(fā)生、最大荷載、骨折部位、應(yīng)變變化及骨折發(fā)生的過程，尚未實(shí)現(xiàn)統(tǒng)一分析，這與如下因素有關(guān)。①整骨的“破壞負(fù)載”及“骨強(qiáng)度”的計(jì)算界定仍然存在巨大爭議:骨折是否首先發(fā)生在應(yīng)力最大的位置，對于骨骼中具有不同力學(xué)性能的構(gòu)成材料的破壞載荷如何界定，最大應(yīng)力變化與骨折的關(guān)系，導(dǎo)致骨折的靜態(tài)載荷與動(dòng)態(tài)載荷之間的關(guān)系等等，其中整骨的“破壞負(fù)載”及“骨強(qiáng)度”的計(jì)算界定是關(guān)鍵因素。②骨折的位置:發(fā)生骨折的部位與承受不同形式、大小的荷載之間的關(guān)系，仍無明確結(jié)論。③應(yīng)變與骨折關(guān)系:骨骼在承受荷載同時(shí)發(fā)生應(yīng)變（壓縮、扭曲、延伸等），一般認(rèn)為，不同年齡段人群、不同個(gè)體及不同生理狀態(tài)的骨骼的載荷-應(yīng)變規(guī)律不同，以什么指標(biāo)判斷應(yīng)變與骨折的關(guān)系，目前尚無定論。④骨密度與骨折的關(guān)系。大量研究表明，骨密度與骨折有著密切的聯(lián)系，由于個(gè)體差異的原因，以及目前骨密度測量方式存在一定不足，導(dǎo)致骨密度-骨強(qiáng)度的關(guān)系未能確定，骨密度有可能是判斷骨折風(fēng)險(xiǎn)較好的指標(biāo)之一，但需要大量的基礎(chǔ)試驗(yàn)及臨床調(diào)查研究。⑤活體骨折的判斷:對于活體而言，骨折的判斷尤為困難，由于骨折的發(fā)生除了與骨強(qiáng)度下降有關(guān)，同時(shí)與肌肉、韌帶等軟組織密切相關(guān)，由于這些軟組織的力學(xué)性能較骨骼更難以明確，尤其目前軟組織的三維重建仍存在技術(shù)上的困難，故軟組織在防止骨折發(fā)生中起到什么樣的作用仍不能明確。

從骨強(qiáng)度這一骨折發(fā)生最本質(zhì)的原因著手，探討有限元分析預(yù)測骨折及骨折風(fēng)險(xiǎn)的有效指標(biāo)，進(jìn)行大樣本的有限元分析與體外力學(xué)試驗(yàn)的比較，分析其回歸規(guī)律，加上更精確地輸入及更合理的建模技術(shù)，同時(shí)要整合各個(gè)水平（如骨組織及整骨）上的研究成果并以之分析現(xiàn)有的大型臨床數(shù)據(jù)等等，有望實(shí)現(xiàn)人骨有限元分析的標(biāo)準(zhǔn)化。

［1］張國棟，廖維靖，陶圣祥，等．股骨有限元分析賦材料屬性的方法［J］．中國組織工程研究與臨床康復(fù)，2009，13（43）:8436-8441．

［2］Chappell KE，Stonebridge-Foster A，Glover A，etal．Magic angle effects in MR Neurography［J］．AJNR AMJNeuroradiol，2004，25（3）:431-440．

［3］Morgan EF，Yeh OC，Keaveny TM．Damage in trabecular bone at small strains［J］．Eur JMorphol，2005，42（1/2）:13-21．

［4］蘇佳燦，張春才，陳維華，等．骨盆及髖臼三維有限元模型單元選擇及構(gòu)建生物力學(xué)意義［J］．中國臨床康復(fù)，2004，8（26）:5550-5552．

［5］Lengsfeld M，Schmitt J，Alter P，etal．Comparison of geometrybased and CT voxel-based finite elementmodelling and experimental validation［J］．Med Eng Phys，1998，20（7）:515-522．

［6］Taddei F，Cristofolini L，MartelliS，etal．Subject-specific finite elementmodels of long bones:An in vitro evaluation of the overall accuracy［J］．JBiomech，2006，39（13）:2457-2467．

［7］Keyak JH，F(xiàn)ourkas MG，Meagher JM，etal．Validation of an automatedmethod of three-dimensional finite elementmodelling ofbone［J］．JBiomed Eng，1993，15（6）:505-509．

［8］Crawford RP，Cann CE，Keaveny TM．Finite element models predict in vitro vertebral body compressive strength better than quantitative computed tomography［J］．Bone，2003，33（4）:744-750．

［9］Peng L，Bai J，Zeng XL，etal．Comparison of isotropic and orthotropicmaterial property assignments on femoral finite elementmodels under two loading conditions［J］．Med Eng Phys，2006，28（3）:227-233．

［10］彭亮，曾小麗，白凈．基于CT數(shù)據(jù)的股骨三維有限元建模方法［J］．清華大學(xué)學(xué)報(bào)（自然科學(xué)版），2007，47（3）:416-419．

［11］Rho JY，Hobatho MC，Ashman RB．Relationsofmechanical properties to density and CT numbers in human bone［J］．Med Eng Phys，1995，17（5）:347-355．

［12］Wirtz DC，Schiffers N，Pandorf T，etal．Criticalevaluation of known bonematerial properties to realize anisotropic FE-simulation of the proximal femur［J］．JBiomech，2000，33（10）:1325-1330．

［13］姜海波，葛世榮．基于CT掃描人體股骨的有限元分析［J］．工程力學(xué)，2007，24（10）:156-159．

［14］馬劍雄，馬信龍，張清功，等．有限元模型評價(jià)股骨正常站立位的生物力學(xué)特性［J］．中國組織工程研究與臨床康復(fù)，2007，12（35）:6823-6826．

［15］Bessho M，Ohnishi I，Matsumoto T，etal．Prediction of proximal femur strength using a CT-based nonlinear finite element method:differences in predicted fracture load and site with changing load and boundary conditions［J］．Bone，2009，45（2）:226-231．

［16］Bessho M，Ohnishi I，Matsuyama J，etal．Prediction of strength and strain of the proximal femur by a CT-based finite elementmethod［J］．JBiomech，2007，40（8）:1745-1753．

［17］Yu XZ，Guo YM，Li J，etal．Finite elementanalysisof impact loads on the femur［J］．Chin JTraumatol，2007，10（1）:44-48．

［18］Wu SY，Yang L，Qi J，etal．Quantitative study on proximal femur of osteoporotic elderly women with 16-slice CT［J］．China JRadiol，2005，39（11）:1165-1170．

［19］Keyak JH，Rossi SA，Jones KA，etal．Prediction of fracture location in the proximal femur using finite element models［J］．Med Eng Phys，2001，23（9）:657-664．

［20］Keller TS，Kosmopoulos V，Lieberman IH．Vertebroplasty and kyphoplasty affect vertebralmotion segment stiffness and stress distributions:amicrostructural finite-element study［J］．Spine，2005，30（11）:1258-1265．

［21］Keyak JH，Rossi SA．Prediction of femoral fracture load using finite elementmodels:an examination of stress-and strain-based failure theories［J］．JBiomech，2000，33（2）:209-214．

［22］Oden ZM，SelvitelliDM，Bouxsein ML．Effectof local density changes on the failure load of the proximal femur［J］．JOrthop Res，1999，17（5）:661-667．

［23］Johnell O，Kanis JA，Oden A，etal．Predictive value of BMD for hip and other fractures［J］．J Bone Miner Res，2005，20（7）:1185-1194．

［24］Van Rietbergen B，Huiskes R，Eckstein F，etal．Trabecular bone tissue strains in the healthy and osteoporotic human femur［J］．J Bone Miner Res，2003，18（10）:1781-1788．

［25］Kanis JA，Black D，Cooper C，etal．International Osteoporosis Foundation;National Osteoporosis Foundation．A new approach to the development of assessment guidelines for osteoporosis［J］．Osteoporos Int，2002，13（7）:527-536．

［26］Keaveny TM，Hoffmann PF，Singh M，etal．Femoral bone strength and its relation to cortical and trabecular changes after treatment with PTH，alendronate，and their combination as assessed by finite element analysis of quantitative CT scans［J］．JBone Miner Res，2008，23（12）:1974-1982．

［27］張國棟，陶圣祥，毛文玉，等．基于三維重建技術(shù)及有限元分析的骨密度測量［J］．中國組織工程研究與臨床康復(fù)，2010，14（9）:1539-1544．

［28］張國棟，毛文玉，廖維靖，等．基于三維重建技術(shù)及有限元分析的脊柱骨密度測量及其意義［J］．中國臨床解剖學(xué)雜志，2010，28（1）:78-84．