-
摘要: 为了模拟和预测肌体组织复杂的再生修复过程，提出基于力学环境和血液供给条件建立骨折愈合仿真模型.针对骨折固定的力学条件和生物学因素，用一种时间动态模型模拟二期骨折愈合阶段的机械稳定性、血管再生和组织分化之间复杂的作用关系.与以往模型不同，本研究建立骨折的三维几何模型，通过有限元法计算骨痂局部力刺激，并与模糊逻辑相结合，将血液浓度作为时-空状态变量引入到模型中，描述骨痂力学及组织分化过程.通过前进欧拉法进行组织浓度等时间步长的迭代更新，在Visual Studio 2012环境下实现愈合进程模拟.最后，利用仿真模型预测稳定与不稳定环境下骨间动度随骨折愈合时间的变化情况，并将仿真结果数据与实验数据进行对比，结果表明，仿真结果与实验数据在趋势和数值上都有较好的吻合，仿真结果数据全部分布在实验数据平均偏差范围内.该结果验证了骨折愈合模型的精确性以及在模拟骨折愈合过程方面的优势.Abstract: A fracture healing model considering mechanical environment and blood supply conditions is proposed to simulate and predict the complex regenerative repair process for tissue. For mechanics condition of fracture fixation and biological factors, a dynamic spatio-temporal model is developed to simulate the complex interactions of mechanical stability, revascularization and tissue differentiation in secondary fracture healing. Unlike previous study, a three-dimensional finite element model is established. The blood perfusion regarded as a spatio-temporal state variable is included into the model to simulate the revascularization process. With finite element method and fuzzy logic, the dynamic model can describe the callus mechanics and biological processes of tissue differentiation. The callus healing process is simulated in Visual Studio 2012 by forward Euler integration method over equidistant time steps iterative loop. Finally, the model predicates the course of interfragmentary movement of two groups with a different axial stability. Through the comparison with the experiment data, it turned out that the simulation result distributed within the scope of the average deviation of the experimental data, which corresponds well to the experiment curve trend and value. The agreement of simulation results with experiment results verifies the accuracy of the fracture healing model and the advantage of the simulating healing process.
-
Key words:
- Fracture healing /
- three-dimensional model /
- fuzzy rules /
- tissue differentiation /
- finite element
-
骨折是一种常见的创伤, 骨折的高发性使得对骨折机理及促进愈合的研究尤为迫切.与其他组织损伤修复不同的是, 骨折不是靠纤维结缔组织连接, 而是骨组织的完全再生[1].骨组织再生是一个极其复杂的生物学修复过程, 受到多种因素的影响, 这些影响因素大致可分为生物学因素和生物力学因素.事实上, 不是所有的骨折都可以被修复, 有时会有不修复或拖延修复, 由此引发患肢疼痛和功能障碍.关于骨折愈合过程及影响骨折愈合速度和质量的研究已经取得了一些进展, 但是受到研究手段和骨折愈合过程复杂又不能直接观察的限制, 仍有约5$%$\$\sim$ $10\ \%$的骨折因各种原因发生延迟愈合或不愈合.如何利用有效方法搭建一个平台使其能够预测骨折愈合这个复杂的过程、可以探寻骨折愈合的机理、可以无限次地通过试验去找到提高手术成功率、降低不愈合风险、提高愈合质量的影响因素和有效办法.建立骨愈合过程的数字化仿真系统为以上问题的解决提供了最佳途径.
Claes和Heigele[2]分别建立了骨折后1周、4周和8周三个愈合阶段的二维轴对称绵羊截骨有限元模型.该模型依靠单元特性的改变来模拟组织分化和骨痂硬化的过程建立了局部力刺激与新组织形成的连接关系提出了定量的组织分化理论.但是该模型只能描述三个特定时期的骨折愈合过程. Ament和Hofer[3]于2000年提出了一种骨折愈合过程中组织分化的模型该模型基于一系列从实验中推导出的模糊逻辑规则以应变能密度作为控制组织分化的力学刺激给出了膜内骨化、软骨骨化和萎缩的组织转化过程.然而该模型没有考虑血液因素对愈合的影响在愈合过程的模糊规则描述方面不够全面.爱尔兰学者Lacroix和Prendergast[4]基于量化组织分化理论建立了两相介质有限元模型在骨折愈合中组织分化的时间进程.将应变和流速作为组织分化的两个力学调控因素成功预测了不同骨折间隙下从炎症后的肉芽组织到骨塑型阶段的愈合过程[5]. Bailón-Plaza和Meulen[6]利用有限差分方法研究了生长因子作用下的骨折愈合进程模拟了连续的组织分化和骨痂中不同细胞演变的过程但没有考虑力学因素对愈合的影响. 2008年Geris等[7]提出了一种基于Bailón-Plaza模型的生物进程调节模型该模型在考虑力学因素的同时将血液浓度等生物学变量引入模型中模拟骨折愈合生物进程. García等[8]基于Prendergast理论用数学模型阐明了力学环境和组织形成的复杂关系预测了愈合初期骨痂生长的几何形态虽然该模型考虑了细胞转化、迁移、增殖以及细胞间质的连续的变化过程描述了组织形成和力学环境的复杂关系但是在愈合初期忽略了断端周围未分化的骨痂组织. Kuiper和Doblaré等[9-12]在2000~2008年相继提出用热弹性模量分析来改变有限元网格几何图形从而实现预测骨组织生长的动态愈合过程.以上研究都没有考虑动态血管重建在愈合中的作用.直到2010年德国乌尔姆大学Simon等[13]提出了一种模拟二期骨折愈合过程中的机械稳定性、血管重建和组织分化的动态模型模型将血供情况和局部力刺激作为变量通过模糊规则模拟组织分化的动态过程预测了骨折断端不同的力学条件下的组织分化和血管再生. 2014年Steiner等[14-15]成功预测了多种形式的作用力对骨折愈合进程的影响但是该研究仅对二维轴对称形式的骨折模型进行了仿真.
本文以三维理想横向骨折几何模型为基础用一种动态迭代模型来仿真二次骨折愈合中随时间变化的血管再生和组织分化过程.主要以骨重建理论和模糊理论为基础基于影响骨折愈合的力学环境和组织分化血液供给条件进行骨折愈合的力学调控模型以及包括21条模糊规则的组织分化的模糊模型的建立在Visual Studio 2012环境下实现组织成分更新过程.
通过构建骨折愈合仿真模型能够针对不同的患者进行个体化骨折愈合过程的仿真与预测.通过复杂骨折愈合过程的预测对医生制订正确的手术方案提供指导进而提高手术成功率提高骨折愈合质量可以有效减少骨折不愈合和延迟愈合的情况帮助骨折患者通过有效的治疗手段恢复健康生活减轻由此带来的社会经济负担.
1. 方法
1.1 骨折三维几何模型
将二维断层DICOM图像导入到MIMICS软件10.01中经过图像分割和三维重建等图像处理过程建立真实的长骨断骨三维几何模型并生成STL文件.由于在MIMICS软件中生成的三维几何模型是壳体不是实体模型因此需要将MIMICS生成的STL文件导入GEOMAGIC进行实体化并进行模型的平滑处理与优化输出IGES文件.将生成的IGES文件导入PROE软件中去除多余的部分建立预定义骨痂区域的骨折愈合几何模型该模型就是有限元网格划分的输入模型.
1.2 骨折力学模型及有限元计算
长骨骨折后骨痂的力学模型被看成是各种纯组织混合的各向同性的弹性体而皮质骨本身应该是非均匀各向异性的.然而在骨折愈合早期由于力是传递给骨痂的因此皮质骨力学性能的相关性低于肉芽组织故假设皮质骨的属性为各向同性.基于弹性力学的基本方程建立骨折愈合单元的应力应变关系: $\pmb\sigma={D}\pmb\varepsilon$其中${D}$为弹性矩阵.弹性矩阵的参数由两个力学参数弹性模量${E}$和泊松比$\nu$确定.
Carter和Hayes[16]通过大量实验建立了弹性模量与表观密度的关系式$E=3\ 790\ \rho^3$. Shefelbine[17]等针对骨痂组织分化基本成分和含量进行研究得出骨痂单元的弹性模量和各组织的弹性模量以及各组织浓度三次方有关单元泊松比与各组织浓度和各组织泊松比有关.
$ E=\sum\limits_{\tau \in T}E_\tau c_\tau^3 $
(1) $ \nu=\sum\limits_{\tau \in T}c_\tau \nu_\tau $
(2) 其中${T}$是骨、软骨和软组织集合$c_\tau$是集合中局部各组织浓度$E_\tau$是混合物对应的弹性模量$\nu_\tau$是混合物对应的泊松比.混合物各组织成分的力学特性见表 1.
表 1 组织成分的力学特性Table 1 Material properties of tissue types组织成分 弹性模量(MPa) 泊松比 皮质骨 10 000 0.36 编织骨 4 000 0.36 纤维软骨 200 0.45 结缔组织 3 0.30 根据虚功原理得
$ \begin{align} \delta \pmb d^{e{\rm T}}{F}^e=&\ \delta \pmb d^{{e}{\rm T}}\int_{V^{e}}{B}^{\rm T}\pmb\sigma \rm dV=\nonumber\\ &\ \delta \pmb d^{{e}{\rm T}}\int_{V^{e}}{B}^{\rm T}{D}{B}\pmb d^{e}\rm dV \end{align} $
(3) 其中${F}^e$是作用在单元上的等效节点力.由于等式两边的虚位移$\delta \pmb d^e$是任意的因此得到平衡方程
$ \begin{align} {F}^e={K}^e\pmb d^e \end{align} $
(4) 其中${K}^e$称为单元刚度矩阵.
在线性四面体有限元方法求解中求解平衡方程的技术问题常常最终归结为线性代数方程组的求解问题.对于$n$阶线性方程组${A}\pmb x=\pmb b$的数值解法分为直接法和迭代法.在骨痂愈合模型中有限元方程的解算在经过共轭梯度法和LU分解法等多种方法的尝试后确定采用LU分解法进行求解.
在有限元求解程序中创建两个子函数进行有限元方程的求解在void Deal::GetShift($\cdot$)中进行调用获取位移两个子函数分别为:
1) bool CMathEasyFun::lu(double *$a$ int *$pivot$ int $n$) //LU分解.
2) bool CMathEasyFun::guass(double const *$lu$ int const *$p$ double *$b$ int $n$) //由LU分解计算线性方程的解.
1.3 调控组织分化的力学激励
$ \begin{align} \pmb S:=\left[ \begin{array}{c} \varepsilon_0 \\ \gamma_0 \end{array} \right] \end{align} $
(5) 其中，
$ \varepsilon_0:=\frac{1}{3}(\varepsilon_1+\varepsilon_2+\varepsilon_3) $
(6) $ \gamma_0:=\frac{1}{\sqrt{2}}\sqrt{(\varepsilon_1-\varepsilon_2)^2+(\varepsilon_1-\varepsilon_3)^2+(\varepsilon_2-\varepsilon_3)^2} $
(7) 式(6)用应变第一不变量表示膨胀应变又称为静力学应变代表体积的改变. $\varepsilon_1$ $\varepsilon_2$ $\varepsilon_3$是主应变.式(7)是应变偏张量表示的畸变应变代表形状的改变.
1.4 外固定器力学模型
骨折愈合模拟仿真中采用的外部固定器是1997年Claes等进行的骨折间隙和骨端稳定性对骨折愈合影响实验中设计的如图 1所示.
1.5 骨间动度的计算
外固定下的骨折愈合系统在外力的作用下骨折断端会产生微小移动即骨间动度.骨折断端的微动不断刺激着骨折的愈合同时骨痂组织不断分化骨痂的硬度随之增大会使骨间动度变小.当骨痂逐渐分化成骨时骨间动度减小为零骨折达到愈合状态.因此骨间动度既是判断骨折是否愈合的指标也是影响骨折愈合进程的必要力学因素.
假设绵羊在正常行走过程中承受的最大轴向作用力为500 N即$F=500~{\rm N}$在皮质骨末端2 mm处施加轴向固定约束将力全部作用于骨折愈合系统中则整体作用力可分为两部分:固定力和愈合力[18].
$ \begin{align} F(\mu)=F_{\rm fix}(\mu)+F_{\rm cal}(\mu)=F_{\rm fix}(\mu)+k_{\rm cal}\mu \end{align} $
(8) 其中$F_{\rm fix}(\mu)$是固定力$F_{\rm cal}(\mu)$是愈合力$k_{\rm cal}$是骨痂刚度.
为了更好地确定位移与固定力和愈合力的关系将骨痂与固定架弹簧系统分离开在每次时间步长迭代都以骨痂的线性部分计算为程序的初始阶段在顶部节点处设定1 mm的单位位移$\mu_{\rm unit}$那么可以知道轴向反作用力$F_N$根据胡克定律在迭代步i时骨痂的硬度$k_{\rm cal}$可表示为
$ \begin{align} k_{\rm cal}=\frac{F_N}{\mu_{\rm unit}} \end{align} $
(9) 如图 2所示整个愈合系统的非线性力-位移曲线可以看成是骨痂的线性曲线与固定架弹簧系统的非线性曲线的叠加所得的非线性分段函数.通过反函数即可求得全部外力作用下实际骨间动度$\mu$实现骨痂组织材料的迭代更新.每一次的组织分化的材料属性的更新都会伴随着一次骨间动度的计算.
1.6 生物组织浓度的状态变量
在整个骨折模拟仿真中愈合区域被离散成有限的组织单元.每个单元都是由不同的组织形式混合而成具有各自的生物组织特性.本研究假设各时段的骨痂是由骨组织、软骨组织和纤维结缔组织以一定的含量比混合而成.用空间位置和时间确定了描述生物组织浓度状态的数字模型.
$ \begin{align} {C}(\pmb x,t):=\left[ \begin{array}{c} C_{\rm perf}(\pmb x,t) \\ C_{\rm bone}(\pmb x,t) \\ C_{\rm cart}(\pmb x,t) \end{array} \right] \end{align} $
(10) 其中$x\in \Omega\subset {\bf R}^3$代表恒定的预定义愈合范围. $t$ $\in$ $[0+\infty)$. $C_{\rm perf}(\pmb x,t)$是血流浓度是仿真血管再生和组织分化之间复杂相互作用过程的状态变量. $C_{\rm bone}(\pmb x,t)$是编织骨的浓度$C_{\rm cart}(\pmb x,t)$是纤维软骨的浓度.
1.7 模糊控制组织分化模型
用充血(血供)、软骨密度和骨密度三个组织浓度变量相邻单元的充血和骨密度两个相邻组织状态变量膨胀应变和畸变应变两个力学刺激作为模糊输入经过21条模糊规则描述组织分化过程.充血、骨密度和软骨密度的改变量作为输出.首先建立7个输入3个输出语言变量的语言值然后以所计算的应变的组织学结果及细胞培养实验[19-20]为基础建立如图 3所示的隶属度函数.
1.8 组织分化模糊规则
模糊规则主要描述愈合进程中组织转变过程规则的推导源于大量的医学知识这些规则是模糊模型的基础本研究采用"if $A$ and $B$ then $C$"的语句描述其主要形式如表 2所示.
表 2 二期骨折愈合的模糊规则Table 2 The fuzzy rules of the fracture healing分化过程 血供 相邻血供 骨密度 邻骨密度 软骨密度 膨胀应变 畸变应变 $\Delta$血供 $\Delta$骨密度 $\Delta$软骨密度 1 低 中或高 $- $ $- $ $- $ 非负中 非中 增高 $- $ $- $ 2 中 高 $- $ $- $ $- $ 非负中 非中 增高 $- $ $- $ 3 高 $- $ $- $ $- $ $- $ 非负中 非中 增高 $- $ $- $ 4 $- $ $- $ $- $ $- $ $- $ 负中 $- $ 降低 $- $ $- $ 5 高 $- $ $- $ 中或高 低 负低 低 $- $ 增高 $- $ 6 高 $- $ $- $ 中或高 低 正低 低 $- $ 增高 $- $ 7 $- $ $- $ 低或中 $- $ 低 负中 非过载 $- $ $- $ 增高 8 $- $ $- $ $- $ $- $ 中或高 负中 非过载 $- $ $- $ 增高 9 $- $ $- $ $- $ $- $ 高 负低 非过载 $- $ $- $ 增高 10 不低 $- $ $- $ 中或高 中或高 正低 零 $- $ 增高 降低 11 不低 $- $ $- $ 中或高 中或高 负中 零 $- $ 增高 降低 12 不低 $- $ $- $ 中或高 中或高 负中 低 $- $ 增高 降低 13 不低 $- $ $- $ 中或高 中或高 负低 低 $- $ 增高 降低 14 不低 高 高 $- $ 低 负低 低 $- $ 增高 降低 15 不低 高 高 $- $ 低 负低 零 $- $ 增高 降低 16 高 $- $ $- $ $- $ $- $ 零 零 $- $ 降低 $- $ 17 $- $ $- $ $- $ $- $ $- $ 负过载 $- $ 降低 降低 降低 18 $- $ $- $ $- $ $- $ $- $ 正过载 $- $ 降低 降低 降低 19 $- $ $- $ $- $ $- $ $- $ $- $ 过载 降低 降低 降低 20 $- $ $- $ $- $ $- $ $- $ 非负过载 $- $ 不变 不变 不变 21 $- $ $- $ $- $ $- $ $- $ 非正过载 $- $ 不变 不变 不变 规则1~4描述了血管生成血肿的刺激使得骨折处的毛细血管再生血管的再生是从已存在的血管中产生新血管的过程并在力学的调控下产生生物学效应而诱导血管重建.以血管的力学生物学为基础研究者开展了应力诱导血管重建的差异蛋白质组学的研究得出低应力可能通过抑制蛋白的表达从而诱导血管重建[21]的结论.因此血管生成与相邻单元的血液浓度及局部力刺激大小有关.故在适当的应变条件下至少一个单元有相同或更高充血量血液浓度才会上升血管生成.过高的膨胀应变则会阻碍血管再生.
规则5和规则6描述了膜内骨化.膜内骨化是一种骨生长的成骨方式尤其是在成骨过程中没有软骨生成时膜内骨化是结缔组织转变成骨组织的直接转化形式.初期骨髓间充质干细胞开始增殖聚集血管生成活性增加毛细血管开始增长骨髓间充质干细胞进而发展成为成骨细胞开始产生类骨质.与此同时类骨质开始以相同的速率矿化.在骨表面由成骨细胞促使生成新的类骨质而成骨细胞逐渐变成骨细胞[22-23].这种成骨方式仅仅在稳定的力学环境和充足的血液供应条件下并且相邻单元具有较高骨密度时才会发生.
规则7~9描述了软骨的形成.在固定不良骨折断端活动较大同时伴随血液供给不足的情况下断端周围组织内的间充质细胞才会分化成软骨细胞继而产生软骨[24].
规则10~13描述了软骨钙化过程软骨骨痂形成后成骨细胞在高应力刺激下合成并分泌出溶胶原、蛋白多糖和糖蛋白构成骨的有机质.溶胶原逐渐聚合成胶原纤维并和成骨细胞一起被基质包围逐渐发生钙化.这个过程的发生与血液供给无关.实验证明在骨折端持续的压缩即低力刺激下可抑制软骨内钙化而在间歇性高力刺激下会促进软骨钙化的发生[25].因此软骨钙化而引起的硬度的增加可以看成是在较高的力刺激下骨密度的增加和软骨密度的减少.
规则14和规则15描述软骨骨化过程在软骨钙化区内营养物质逐渐弥散发生障碍软骨细胞逐渐凋亡钙化基质逐渐降解释放血管生长因子此时大量的毛细血管和成骨细胞侵入矿化的软骨基质逐渐产生类骨质发生骨化现象[26].因此在血液供给良好的条件下完全钙化的软骨才能发生骨化.在本模型中完全钙化的软骨是以高的骨密度和低的软骨密度表达的.规则不仅预测了骨密度的增加而且也预测了相同数量软骨密度的减少.
除了正常状态下的骨痂分化过程规则16~21描述在过度力加载条件下引起的组织萎缩或者没产生诱发组织分化条件时保持的原有参数不变.模糊规则主要描述的组织转变过程如图 4所示.
在模糊输入输出隶属度函数和模糊规则确立以后建立FIS文件.为了减少计算时间方便编程本文将模糊控制运算所表达的输入输出对应关系建立成模糊查询表的形式用Matlab下的系统测试工具收集构造查询表所需的数据资料.查询表与模糊逻辑推理系统的效果几乎一样它们都是输入变量空间到输出变量空间的一种映射查询表法节省了大量重复计算的时间.基于建立的FIS数据结构通过系统测试生成mat文件.用C++编程语言对此生成的数据进行匹配查询.将应变和各组织浓度匹配最符合的输出值提取出来即为模糊逻辑计算值.
2. 结果和讨论
建立一个如图 5(a)所示的横向截骨三维骨折愈合模型.中间部分为预定义的骨痂范围尺寸的设定与对比试验的尺寸一致.上下两端为截取的原始骨内外直径分别为12 mm和16 mm且横截面积为28$\pi\ \rm{mm^2}$.
将建立的三维几何模型导入到网格划分软件中(自主研发)进行体网格的划分网格数为14\448如图 5(b)所示.由于网格划分软件会生成许多数据而本研究只需要节点坐标和单元编号将文本文件导入到Matlab中进行预处理只提取出目标数据生成后续有限元计算所需要的单元编号和节点坐标两个txt文本格式的文件.节点坐标文件中三列数据分别代表每个节点的空间坐标值.而单元编号文件中的四列分别为每个单元的四个节点的节点序号.
随着模型的复杂程度的增高和网格单元个数的增加计算也会相应复杂.特别是对模型的整体刚度矩阵等大型稀疏矩阵以及各单元组织浓度的迭代更新的计算和存储不仅要对其算法进行优化而且需要有较高配置的计算机硬件的支持.本研究采用的运行系统类型为64位的Windows7操作系统处理器参数为: Inter(R) Xeon(R) E5606@2.13 GHz四核安装内存为8 GB RAM显卡1 GB.在Visual Studio 2012环境下首先进行模块和类成员及成员函数的基本设定[27-29]. Deal类中主要实现文件的信息加载、数据的预处理、有限元的解算、模糊运算查询以及数据的更新和输出.信息加载及模型的预处理部分通过Deal类中成员函数Start($\cdot$)和Init($\cdot$)进行初始化读取网格划分后的单元编号和节点坐标扫描所有单元获取单元的数目.然后通过成员函数PreCompute($\cdot$)进行预计算并基于有限元理论运用GetB(double* fB int iElementIndex); GetD($\cdot$); GetAll_K($\cdot$)等功能函数获取几何矩阵、弹性矩阵、单元刚度矩阵进行总体刚度矩阵的封装继而解算力学有限元模型得到调控组织分化的力学刺激.然后通过功能函数Search($\cdot$)进行模糊运算的查询匹配获得组织浓度的改变量进行材料属性的更新计算同时更新数据类.将计算出的各组织单元的材料属性赋予有限单元中.当迭代达到平衡且收敛时输出平均骨密度、血液浓度等结果.
图 6为三维外固定横向截骨模型的外加载荷和固定约束的示意图.在模型顶端施加载荷大小为$F$ $=$ $500~{\rm N}$假设这是在绵羊正常行走中足够满足跖骨可能承受的最大的力.底部为固定约束部分其自由度均为0包括三个方向的位移和三个方向的旋转.外部固定架的安放使得模型只能发生轴向位移其他方向自由度均为0.
为了验证仿真模型的精确性方便仿真结果与实验结果进行对比初始条件和边界条件都是根据Claes的实验条件来设定如表 3所示.
表 3 初始条件和边界条件(%)Table 3 Initial condition and boundary condition (%)血液浓度 骨密度 软骨密度 结缔组织 骨皮质 100 100 0 0 断端间骨痂 0 0 0 100 骨痂外边缘 30 0 0 100 2.1 骨折愈合仿真结果的精确性验证
应用骨折愈合模拟程序本研究在2 mm和3 mm两种不同骨折间隙环境下针对骨间轴向稳定性对骨折愈合的影响进行分析.两种环境下外部固定器在相同力加载条件下允许的最大折断间位移即骨间动度分别为0.25 mm和1.25 mm即2 mm为相对稳定情况3 mm为不稳定情况.根据固定弹簧系统的非线性曲线形式确定两种不同间隙条件下外部固定器的力学非线性行为如图 7所示.将模拟程序计算所得的骨间动度、血液浓度以及弹性模量等变量在等时间步长的迭代更新数据以.txt文档的形式输出.程序运行结束后共约150次迭代步.
本研究预测了稳定与不稳定两种环境下骨间动度$\mu$随骨折愈合时间的变化.将该结果与Claes等于1997年做的经典绵羊跖骨横向截骨骨折愈合过程中各组羊的骨间动度随时间变化的实验数据进行对比分析如图 8所示. 图 8(a)和图 8(b)分别为2 mm间隙稳定情况和3 mm间隙不稳定情况的模拟结果.
图 8 骨间动度仿真与实验数据的对比Fig. 8 Comparison between interfragmentary movement over time and test results从图 8可以看出随着愈合时间(横坐标)的不断增多仿真曲线和实验曲线的骨间动度(纵坐标)都呈现逐渐减小的趋势这符合骨折愈合的特征.由于力刺激和血液供应的共同作用骨痂开始组织分化骨折愈合初始时刻骨痂几乎都被看成由结缔组织组成的刚度较小骨间动度相对较大随着分化的组织转变开始出现软骨和骨骨痂的刚度逐渐提高承受外力的量级也逐渐增大骨间动度随之变小直至为零.
德国乌尔姆大学的Simon等于2010年建立的骨折愈合仿真模型是基于二维骨折模型设计的依据实验数据建立了8条模糊规则.本研究建立的仿真模型是基于三维骨折模型设计的建立7个输入3个输出语言变量的语言值然后以所计算的应变的组织学结果及细胞培养实验建立隶属度函数依据生物学和医学基础理论建立所有模糊规则共21条.将本研究的愈合模拟数据与Simon等2010年的模拟结果以及各组绵羊骨折愈合实验数据平均值进行对比如图 9所示.其中图 9(a)和图 9(b)分别为2 mm间隙稳定组和3 mm间隙不稳定组的对比数据.
图 9 骨间动度仿真与Simon仿真结果的对比Fig. 9 Comparison between simulated interfragmentary movement and Simon results如图 9所示本研究建立的骨折愈合模型的仿真结果与Simon模拟获得的数据在数值上非常接近在大部分区域本研究建立的骨折愈合模型的仿真结果比Simon等得到的仿真结果更贴近实验数据且分布在实验数据平均偏差范围内.这与本研究利用三维骨折模型进行有限元仿真和设定了21条模糊规则有关.在相同的力学条件和血液供给条件下本研究建立的骨折愈合仿真模型得到的骨间动度数值与文献中实验测得的骨间动度数值之间的平均误差分别为0.029 (稳定组)和0.075 (不稳定组).产生误差的原因可能来源于本文对骨的生物力学特性做了一定的简化并且骨折的愈合是由多种因素共同调控本研究仅考虑了影响骨折愈合的力学刺激和血液供给条件而忽略了生长因子等因素的作用.
2.2 骨间稳定性对愈合的影响分析
从图 9两种不同的骨折稳定性情况模拟的位移曲线可以看出对稳定的2 mm间隙情况在46天左右完成桥接愈合而相对不稳定的3 mm间隙情况在65天左右完成愈合桥接即稳定固定能加快骨桥接过程这个预测结果与实验观测情况一致.
随着骨折的愈合骨痂组织分化过程中毛细血管逐渐生成骨断端间血液浓度不断提高直至达到理想的100%血液浓度而收敛如图 10所示. 图 11为标准血供下稳定组和不稳定组愈合阶段骨痂随着组织分化过程中组织成分的不断改变刚度随时间变化的关系图.
从模拟程序仿真的血液浓度随时间的变化图 10可以看出随着愈合时间的发展局部血液浓度值呈上升趋势两组曲线的总体趋势相同.相对稳定的固定组血液浓度增长的速率即图中曲线的斜率大于相对不稳定组.血液浓度达到100%理想充血的时间相对较早且在同一愈合时刻稳定组的血液浓度高于不稳定组血管再生的速度快.产生这一现象的原因在于骨折后血管的再生是在已存在的血管上生成新血管过程断端的运动会破坏新生的血管使血管形成受阻.模拟程序预测的结果与实验观测一致且验证了骨折断端的稳定性与血管生成的关系.随着组织的不断分化结缔组织从软骨转变成编织骨从而使骨痂的刚度不断增高逐渐趋于原始骨.对比图 11的两组曲线稳定组刚度的变化率要快于相对不稳定组愈合速度快.
3. 结论
本研究首次在骨折愈合过程仿真中引入三维骨折模型和21条模糊规则并在Visual Studio 2012环境下进行C++编程结合Matlab软件对组织分化模糊模型的预计算实现了基于力学环境和血液供给条件的骨折愈合仿真对于同尺寸的羊拓骨利用本研究建立的骨折愈合模型对骨间动度随愈合时间的变化情况进行仿真通过仿真结果与实验结果的对比验证了骨折模型的精确性.由于本研究建立的骨折模型是在Visual Studio 2012平台上实现的因此可以很容易地将未来的研究成果以模块替换的方式加进仿真平台中包括复杂几何模型、更精确的网格划分算法、非线性力学模型、多尺度的血液供给模型等.
骨折愈合仿真的成功实施将产生以下几点影响和作用:
1) 该系统能作为一个仿真试验平台用于各种适应人体生物力学要求的内固定器材的评价和优化设计; 用于非生物学因素对骨组织愈合速度和质量影响的测试分析; 由于该仿真模型的可扩展性该平台亦可实现对部分影响骨愈合进程的生物学因素进行分析和研究.
2) 该系统能够针对不同的患者进行个体化骨折愈合过程的仿真与预测.
3) 该系统能够实现复杂骨折愈合过程的预测对医生制订正确的手术方案提供指导.
4) 该系统可以利用计算机中建立的仿真模型进行多次重复试验研究不需要真实的生物学试验节省时间提高效率节省费用避免人道主义的争议.
目前的研究在力学模型、生物学参数和实验研究方面还有一定的局限性在未来的研究中我们将做多组不同力学环境下的骨折愈合实验通过对比来进一步讨论所建立的骨折愈合仿真模型的精确性.我们也将在力学模型和生物学模型两个方面开展更加深入的研究工作在保证模型精确性的前提下进一步实现通过模型的仿真能够提取更多中间参数为探寻骨折愈合机理研究提供更多有价值的数据信息.
-
图 8 骨间动度仿真与实验数据的对比
Fig. 8 Comparison between interfragmentary movement over time and test results
图 9 骨间动度仿真与Simon仿真结果的对比
Fig. 9 Comparison between simulated interfragmentary movement and Simon results
表 1 组织成分的力学特性
Table 1 Material properties of tissue types
组织成分 弹性模量(MPa) 泊松比 皮质骨 10 000 0.36 编织骨 4 000 0.36 纤维软骨 200 0.45 结缔组织 3 0.30 
下载: 导出CSV
表 2 二期骨折愈合的模糊规则
Table 2 The fuzzy rules of the fracture healing
分化过程 血供 相邻血供 骨密度 邻骨密度 软骨密度 膨胀应变 畸变应变 $\Delta$血供 $\Delta$骨密度 $\Delta$软骨密度 1 低 中或高 $- $ $- $ $- $ 非负中 非中 增高 $- $ $- $ 2 中 高 $- $ $- $ $- $ 非负中 非中 增高 $- $ $- $ 3 高 $- $ $- $ $- $ $- $ 非负中 非中 增高 $- $ $- $ 4 $- $ $- $ $- $ $- $ $- $ 负中 $- $ 降低 $- $ $- $ 5 高 $- $ $- $ 中或高 低 负低 低 $- $ 增高 $- $ 6 高 $- $ $- $ 中或高 低 正低 低 $- $ 增高 $- $ 7 $- $ $- $ 低或中 $- $ 低 负中 非过载 $- $ $- $ 增高 8 $- $ $- $ $- $ $- $ 中或高 负中 非过载 $- $ $- $ 增高 9 $- $ $- $ $- $ $- $ 高 负低 非过载 $- $ $- $ 增高 10 不低 $- $ $- $ 中或高 中或高 正低 零 $- $ 增高 降低 11 不低 $- $ $- $ 中或高 中或高 负中 零 $- $ 增高 降低 12 不低 $- $ $- $ 中或高 中或高 负中 低 $- $ 增高 降低 13 不低 $- $ $- $ 中或高 中或高 负低 低 $- $ 增高 降低 14 不低 高 高 $- $ 低 负低 低 $- $ 增高 降低 15 不低 高 高 $- $ 低 负低 零 $- $ 增高 降低 16 高 $- $ $- $ $- $ $- $ 零 零 $- $ 降低 $- $ 17 $- $ $- $ $- $ $- $ $- $ 负过载 $- $ 降低 降低 降低 18 $- $ $- $ $- $ $- $ $- $ 正过载 $- $ 降低 降低 降低 19 $- $ $- $ $- $ $- $ $- $ $- $ 过载 降低 降低 降低 20 $- $ $- $ $- $ $- $ $- $ 非负过载 $- $ 不变 不变 不变 21 $- $ $- $ $- $ $- $ $- $ 非正过载 $- $ 不变 不变 不变
下载: 导出CSV
表 3 初始条件和边界条件(%)
Table 3 Initial condition and boundary condition (%)
血液浓度 骨密度 软骨密度 结缔组织 骨皮质 100 100 0 0 断端间骨痂 0 0 0 100 骨痂外边缘 30 0 0 100
下载: 导出CSV
-
[1] Kayabasi O, Erzineanli F. Finite dement modeling and analysis of a new cemented hip prosthesis. Advances in Engineering Software, 2006, 37 (7):477-483 doi: 10.1016/j.advengsoft.2005.09.003 [2] Claes L E, Heigele C A. Magnitudes of local stress and strain along bony surfaces predict the course and type of fracture healing. Journal of Biomechanics, 1999, 32 (3):255-266 doi: 10.1016/S0021-9290(98)00153-5 [3] Ament C, Hofer E P. A fuzzy logic model of fracture healing. Journal of Biomechanics, 2000, 33 (8):961-968 doi: 10.1016/S0021-9290(00)00049-X [4] Lacroix D, Prendergast P J. A mechano-regulation model for tissue differentiation during fracture healing:analysis of gap size and loading. Journal of Biomechanics, 2002, 35 (9):1163-1171 doi: 10.1016/S0021-9290(02)00086-6 [5] Prendergast P J, Huiskes R, Soballe K. Biophysical stimuli on cells during tissue differentiation at implant interfaces. Journal of Biomechanics, 1997, 30 (6):539-548 doi: 10.1016/S0021-9290(96)00140-6 [6] Bailón-Plaza A, van Der Meulen M C H. A mathematical framework to study the effects of growth factor influences on fracture healing. Journal of Theoretical Biology, 2001, 212 (2):191-209 doi: 10.1006/jtbi.2001.2372 [7] Geris L, Gerisch A, Vander Sloten J, Weiner R, Van Oosterwyck H. Angiogenesis in bone fracture healing:a bioregulatory model. Journal of Theoretical Biology, 2008, 251 (1):137-158 doi: 10.1016/j.jtbi.2007.11.008 [8] García-Aznar J M, Kuiper J H, Gómez-Benito M J, Doblaré M, Richardson J B. Computational simulation of fracture healing:influence of interfragmentary movement on the callus growth. Journal of Biomechanics, 2007, 40 (7):1467-1476 doi: 10.1016/j.jbiomech.2006.06.013 [9] Kuiper J H, Ashton B A, Richardson J B. Computer simulation of fracture callus formation and stiffness restauration. In: Proceedings of the 12th Conference of the European Society of Biomechanics. Dublin, Ireland: Royal Academy of Medicine, 2000. 61-61 [10] Doblaré M, García J M, Gómez M J. Modelling bone tissue fracture and healing:a review. Engineering Fracture Mechanics, 2004, 71 (13-14):1809-1840 doi: 10.1016/j.engfracmech.2003.08.003 [11] Gómez-Benito M J, García-Aznar J M, Kuiper J H, Doblaré M. A 3D computational simulation of fracture callus formation:influence of the stiffness of the external fixator. Journal of Biomechanical Engineering, 2006, 128 (3):290-299 http://www.ncbi.nlm.nih.gov/pubmed/16706578 [12] Isaksson H, van Donkelaar C C, Huiskes R, Ito K. A mechano-regulatory bone-healing model incorporating cell-phenotype specific activity. Journal of Theoretical Biology, 2008, 252 (2):230-246 doi: 10.1016/j.jtbi.2008.01.030 [13] Simon U, Augat P, Utz M, Claes L. A numerical model of the fracture healing process that describes tissue development and revascularisation. Computer Methods in Biomechanics and Biomedical Engineering, 2011, 14 (1):79-93 doi: 10.1080/10255842.2010.499865 [14] Steiner M, Claes L, Ignatius A, Simon U, Wehner T. Disadvantages of interfragmentary shear on fracture healing-mechanical insights through numerical simulation. Journal of Orthopaedic Research, 2014, 32 (7):865-872 doi: 10.1002/jor.22617 [15] Steiner M, Claes L, Ignatius A, Niemeyer F, Simon U, Wehner T. Prediction of fracture healing under axial loading, shear loading and bending is possible using distortional and dilatational strains as determining mechanical stimuli. Journal of the Royal Society Interface, 2013, 10 (86):1123-1126 http://europepmc.org/abstract/med/23825112 [16] Carter D R, Hayes W C. The compressive behavior of bone as a two-phase porous structure. The Journal of Bone & Joint Surgery, 1977, 59 (7):954-962 doi: 10.1007/978-1-4471-5451-8_116 [17] Shefelbine S J, Augat P, Claes L, Simon U. Trabecular bone fracture healing simulation with finite element analysis and fuzzy logic. Journal of Biomechanics, 2005, 38 (12):2440-2450 doi: 10.1016/j.jbiomech.2004.10.019 [18] Claes L, Augat P, Suger G, Wilke H J. Influence of size and stability of the osteotomy gap on the success of fracture healing. Journal of Orthopaedic Research, 1997, 15 (4):577-584 doi: 10.1002/(ISSN)1554-527X [19] Duda G N, Eckert-Hübner K, Sokiranski R, Kreutner A, Miller R, Claes L. Analysis of interfragmentary movement as a function of musculoskeletal loading conditions in sheep. Journal of Biomechanics, 1998, 31 (3):201-210 https://www.sciencedirect.com/science/article/pii/S0021929097001279 [20] Kaspar D, Seidl W, Neidlinger-Wilke C, Ignatius A, Claes L. Dynamic cell stretching increases human osteoblast proliferation and cicp synthesis but decreases osteocalcin synthesis and alkaline phosphatase activity. Journal of Biomechanics, 2000, 33 (1):45-51 doi: 10.1016/S0021-9290(99)00171-2 [21] Guo Z Y, Yan Z Q, Bai L, Zhang M L, Jiang Z L. Flow shear stress affects macromolecular accumulation through modulation of internal elastic lamina fenestrae. Clinical Biomechanics, 2008, 23 (1):S104-S111 https://www.sciencedirect.com/journal/clinical-biomechanics/vol/23/suppl/S1 [22] 金大地, 史占军.骨与骨折愈合.成都:电子科技大学出版社, 1994. 15-16Jin Da-Di, Shi Zhan-Jun. Bone and Fracture Healing. Chengdu:University of Electronic Science and Technology Press, 1994. 15-16 [23] 李瑞宗.骨肿瘤诊断图谱.天津:天津科技翻译出版公司, 1993. 3-5Li Rui-Zong. Atlas of Bone Tumor Diagnosis. Tianjin:Tianjin Science and Technology Translation and Publishing Company, 1993. 3-5 [24] 金岩.组织工程学原理与技术.西安:第四军医大学出版社, 2004. 8-10Jin Yan. Principle and Technology of Tissue Engineering. Xi'an:The Fourth Military Medical University Press, 2004. 8-10 [25] 周树权.应力与微动对骨折愈合的影响.医药前沿, 2012, 2 (8):328-329 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=hainanyx201001051Zhou Shu-Quan. The effect of stress and micro motion on fracture healing. Frontiers of Medicine, 2012, 2 (8):328-329 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=hainanyx201001051 [26] 王志鑫, 陈浩宏.骨折愈合学.武汉:湖北科学技术出版社, 1995. 16-19Wang Zhi-Xin, Chen Hao-Hong. Fracture Healing. Wuhan:Hubei Science and Technology Press, 1995. 16-19 [27] 郑阿奇, 丁有和. Visual C++教程.北京:机械工业出版社, 2015. 38-54Zheng A-Qi, Ding You-He. Visual C++ Course. Beijing:Machinery Industry Press, 2015. 38-54 [28] 侯增广, 赵新刚, 程龙, 王启宁, 王卫群.康复机器人与智能辅助系统的研究进展.自动化学报, 2016, 42 (12):1765-1779 http://www.aas.net.cn/CN/abstract/abstract18966.shtmlHou Zeng-Guang, Zhao Xin-Gang, Cheng Long, Wang Qi-Ning, Wang Wei-Qun. Recent advances in rehabilitation robots and intelligent assistance systems. Acta Automatica Sinica, 2016, 42 (12):1765-1779 http://www.aas.net.cn/CN/abstract/abstract18966.shtml [29] 宋亚飞, 王晓丹, 雷蕾.基于直觉模糊集的时域证据组合方法研究.自动化学报, 2016, 42 (9):1322-1338 http://www.aas.net.cn/CN/abstract/abstract18921.shtmlSong Ya-Fei, Wang Xiao-Dan, Lei Lei. Combination of temporal evidence sources based on intuitionistic fuzzy sets. Acta Automatica Sinica, 2016, 42 (9):1322-1338 http://www.aas.net.cn/CN/abstract/abstract18921.shtml 期刊类型引用(7)
1. 李志云，董佩龙. 严重创伤后骨折卧床制动时间对骨折愈合时间的影响. 海南医学. 2023(06): 795-800 . 
百度学术2. 王沫楠. 力和组织内氧气调控的骨折愈合过程仿真. 自动化学报. 2021(01): 79-91 . 本站查看3. 崔勇，王琛，宋晓. 基于驻极体材料的机械天线式低频通信系统仿真研究. 自动化学报. 2021(06): 1335-1342 . 本站查看4. 吴佳，魏成建. 有限元分析在骨折愈合模拟中的研究进展. 医学食疗与健康. 2020(12): 214-215 . 百度学术5. 付瑞森，杨海胜. 骨折愈合的数值模拟方法. 生物医学工程学杂志. 2020(05): 930-935 . 百度学术6. 林国钦，肖丹，牟朋林. 自制骨折闭合复位装置在山羊下肢骨折中的效果研究. 中国实用医药. 2020(35): 206-208 . 百度学术7. 张艳，李莹，张雪松. 血管内皮生长因子和骨形态发生蛋白-2在糖尿病病人骨组织内表达水平的研究. 骨科. 2019(06): 529-533 . 百度学术其他类型引用(1)
-
下载:
计量
- 文章访问数: 1849
- HTML全文浏览量: 279
- PDF下载量: 697
- 被引次数: 8
