CN102024076B - 仿生骨支架内部微观孔结构的建模方法 - Google Patents

Abstract

本发明涉及一种仿生骨支架内部微观孔结构的建模方法。本方法基于多约束背包问题模型，采用椭球体作为基本单元体，骨支架模型的最小包容盒作为约束空间，运用混合遗传算法来求解支架内部微观孔结构的负模型；然后，通过不含微观孔的支架模型与微观孔结构的负模型之间的布尔减运算，来构造包含微观孔结构的仿生骨支架模型。本发明在保证仿生骨支架内部微观孔结构具有较高的孔隙率和表面积的情况下，改善其内部微观孔形状单一、分布不均匀的情况，增大其内部连通孔道的扭曲度，从而提高仿生骨支架模型的性能，使其具有较强的生物活性。

Description

仿生骨支架内部微观孔结构的建模方法

技术领域

本发明涉及到骨组织工程领域的一种仿生骨支架的建模方法，特别是仿生骨支架内部微观孔结构的建模方法。

背景技术

在骨组织工程研究领域中，成骨细胞与仿生骨支架的复合体在体外的共同培养，是构建组织工程化骨的起始和关键，也是目前国内外研究的热点。从骨组织工程角度来说，仿生骨支架起到了骨细胞外基质的作用，是对骨细胞外基质结构和功能的仿生。仿生骨支架是一种由非均质多微孔连通的三维支架结构，它创造一种有利于骨细胞的粘附、生长、增殖和功能发挥的微环境。医学研究表明，要使仿生骨支架具有生物活性，它必须是具备合适的孔隙率和多孔贯通的空间网状结构，以利于组织液的渗透、骨细胞和血管的生长。因此，仿生骨支架模型应该满足以下的要求：

(1)具有很大的表面积，以利于细胞的黏附、增殖和分化，以及生长因子等生物信号分子的负载；

(2)具有合适的孔尺寸、高的孔隙率，既要保持一定的力学性能，又要为细胞、血管等组织的生长提供合理的空间；

(3)微观孔呈多样性，以利于骨单位的形成，且满足不同组织生长的需要；

(4)微观孔分布相对均匀，且相互导通，以利于细胞生长、营养物质和代谢产物的流动。

关于仿生骨支架内部微观孔结构的建模，国内外已经报道了一些相关的研究成果。S.J.Hollister等提出了基于自组织和拓扑优化的均匀化理论来构建骨支架的微观孔结构，以解决支架的力学性能和孔隙率之间的冲突；杨楠和饶高等提出了采用蒙特卡洛法来构建骨支架的微观孔结构，该方法先采用蒙特卡洛法生成骨支架微观孔结构的二维模型，然后通过拉伸和叠加构建其的三维模型；吴懋亮和李莉敏等提出了采用CAD方法来构建骨支架的微观孔结构，该方法通过单元模型的规则空间组合构建微观孔结构的负型，采用布尔运算获得骨支架微观孔结构的三维模型；刘丰等提出了材料点单元堆积造型方法来构建骨支架的微观孔结构，以建立具有材料梯度和结构梯度的仿生骨支架模型。以上这些建模方法均可以生成仿生骨支架的微观孔结构，但是仍然存在以下的问题：

(1)内部微观孔形状比较单一，微观孔分布比较规则，不利于骨细胞的黏附和增殖；

(2)有些建模算法比较复杂，人工干预较多，建模效率相对较低；

(3)特别是采用蒙特卡洛法构建的支架模型，层与层之间容易产生盲孔，微观孔的连通性比较差。

发明内容

本发明的目的在于针对已有技术存在的缺陷，提供一种仿生骨支架内部微观孔结构的建模方法。在保证仿生骨支架具有较高的孔隙率和表面积的情况下，改善其内部微观孔形状单一、分布不均匀的情况，增大其内部连通孔道的扭曲度，从而提高仿生骨支架模型的性能，使其具有较强的生物活性。

为了实现上述的目的，本发明采用下述的技术方案。一种仿生骨支架内部微观孔结构的建模方法，其特征在于建模步骤如下：

第一步：计算仿生骨支架模型的最小包容盒

输入仿生骨支架的实体模型，基于OBB原理计算支架的最小包容盒，作为微观孔结构负模型的空间约束边界。

第二步：构建支架微观孔结构的单元体模型

如图1所示，通过统计分析发现，人体骨骼内部的微观孔是不规则的分形体，这些微观孔的表面积和体积之比更接近于椭球体。因此，选择椭球体作为负模型的单元体，其模型定义为一个八元组：

ELL＝(a，β，x_c，y_c，z_c，θ_x，θ_y，θ_z)    (1)

式中，ELL表示椭球单元体模型，α是椭球体的长轴参数，β是长短轴比例因子，(x_c，y_c，z_c)是椭球体中心坐标，(θ_x，θ_y，θ_z)是长轴与三个坐标平面的夹角。

第三步：构建支架微观孔结构的负模型

根据多约束背包问题模型来构建支架微观孔结构的负模型，以支架的最小包容盒约束空间作为“背包”，椭球单元体作为“物品”，不断地“填充”支架的约束空间，直到所有单元体构成的实体模型满足孔隙率、连通性等约束条件为止。单元体在约束空间中无序均匀地排列，就构成了内部微观孔结构的负模型，该模型描述了仿生骨支架内部微观孔的个体特征和空间排列特征。支架内部微观孔结构负模型的数学描述如下：

$\left\{\begin{array}{c}\mathrm{Maxf}(v_1,v_2,\mathrm{\&Lambda;},v_n)\&GreaterEqual;\mathrm{\&eta;}{V}_z\\ \mathrm{st}.|b_i-b_j|<2(c_i-c_j)<a_i+a_j,L_i>0.5b_i\\ i,j=\mathrm{1,2},\mathrm{\&Lambda;},n,i\&NotEqual;j\end{array}\right.---\left(2\right)$

式中，Maxf(*)表示求(v₁，v₂，...，v_n)体积和的最大值，v_i表示第i个椭球体在模型中的实际体积，V_z表示不含微观孔的支架体积，η表示孔隙率，a_i、b_i表示第i个椭球体的长、短轴参数，c_i-c_j表示相交的椭球体i、j之间的中心距，L_i表示第i个椭球体的中心到支架边界的距离。

第四步：求解支架微观孔结构的负模型

负模型的求解过程就是计算出负模型中所有椭球体的形状参数(α，β)、位置参数(x_c，y_c，z_c)和姿态参数(θ_x，θ_y，θ_z)，这属于典型的多约束背包问题的求解。采用混合遗传算法来解决该类问题是非常有效的，该算法可定义为一个九元组：

HGA＝(C，E，P₀，M，Φ，Γ，Ψ，D，T)      (3)

式中，HGA表示混合遗传算法模型，C表示个体的编码形式，E表示个体的适应度函数，P₀表示初始群体，M表示群体大小，Φ是选择算子，Γ是交叉算子，Ψ是变异算子，D是扰动算子，T是算法的终止条件。

1、浮点数编码结构

采用二维浮点数组来组织个体基因，每个单元体占据一个基因位。数组第一维维数等于单元体的总数n，第二维维数为8，分别存放单元体的8个参数。编码方案如下表1所示：

表1浮点数编码方案

2、适应度函数

在遗传算法中，适应度较高的个体遗传到下一代的概率较大，适应度较低的个体遗传到下一代的概率较小。个体适应度函数定义为：

$\mathrm{Fit}\left(v\right)=\frac{\mathrm{\&Sigma;}{v}_i}{\mathrm{\&eta;}{V}_z}---\left(4\right)$

式中，Fit(*)表示求*的个体适应度值，v_i表示第i个椭球体在模型中的实际体积，η表示孔隙率，V_z表示不含微孔的支架体积。

3、初始群体

初始群体的采样过程就是随机产生n个椭球体模型，模型参数包括形状参数(α，β)、位置参数(x_c，y_c，z_c)和姿态参数(θ_x，θ_y，θ_z)。根据骨组织工程的要求，支架内部微观孔的尺寸通常取100～400μm。a在[100，400]范围内均匀产生随机数进行赋值，β在[0.5，1]范围内均匀产生随机数进行赋值，(x_c，y_c，z_c)在支架的约束空间内均匀产生随机数进行赋值，(θ_x，θ_y，θ_z)在[0⁰，90⁰]范围内均匀产生随机数进行赋值。编码串长度8n决定了构成微观孔结构负模型的单元体数量，n值由下列公式计算：

式中，n表示负模型中单元体的数量，V_c是长方体包围盒的体积，

是所有椭球体长轴参数的平均值，β是长短轴比例因子。

群体大小M表示群体所含个体的数量，对算法的效率和群体的多样性产生重要的影响。一般建议M的取值范围为20～100，本算法中取M＝80。

4、遗传算子

(1)选择算子Φ：采用比例选择与最优保存策略相结合的混合选择算子，进行优胜劣汰操作，即当前群体中适应度最高的个体不参与交叉运算和变异运算，用它来替换本代群体中经过交叉、变异操作后所产生的适应度最低的个体。

(2)交叉算子Γ：为了避免近亲交叉操作而降低种群的多样性，采用均匀交叉算子进行个体的交叉操作，交叉概率P_c＝0.6。

(3)变异算子Ψ：以概率P_m＝0.09产生一个均匀随机数{1，2，...，n}，来选取变异基因位。根据该基因位的基因值范围均匀产生一个随机数，来替换该基因位的值，实现均匀变异操作。

(4)扰动算子D：以扰动概率P_d＝0.05产生一个满足约束的随机个体，替代群体中适应度最高的个体。对解空间的局部进行扰动，提高全局的搜索能力。

第五步：支架内部微观孔结构的生成

将不含微观孔的支架模型与微观孔结构的负模型导入相同的坐标空间，进行两个模型的布尔减运算，得到包含微观孔结构的仿生骨支架模型。

本发明与现有技术相比较，具有如下显而易见的突出的实质性特点和显著优点：本发明基于多约束背包问题模型，采用椭球体作为基本单元体，骨支架模型的最小包容盒作为约束空间，运用混合遗传算法来求解支架内部微观孔结构的负模型；然后，通过不含微观孔的支架模型与微观孔结构的负模型之间的布尔减运算，来构造包含微观孔结构的仿生骨支架模型。本方法在保证仿生骨支架内部微观孔结构具有较高的孔隙率和表面积的情况下，改善其内部微观孔形状单一、分布不均匀的情况，增大其内部连通孔道的扭曲度，从而提高仿生骨支架模型的性能，使其具有较强的生物活性。

附图说明

图1是人体骨骼的内部微观孔结构；

图2是支架微观孔结构建模的具体实施流程；

图3是某缺损骨的仿生骨支架模型；

图4是微观孔结构的负模型；

图5是含微观孔结构的仿生骨支架模型。

具体实施方式

本发明的一个优选实施例结合附图详述如下：

仿生骨支架内部微观孔结构建模方法的具体实施流程如图2所示，主要包括以下五个步骤：

步骤1计算骨支架模型的最小包容盒

输入仿生骨支架外观模型M_sld，基于OBB原理计算M_sld的最小矩形包容盒(X_b，Y_b，Z_b)，作为微观孔结构负模型M_sub的空间约束边界；

步骤2构建支架微观孔结构的单元体模型

根据公式(1)所描述的模型参数组设计椭球单元体的数据结构，定义单元体模型类：

ELL＝(a，β，x_c，y_c，z_c，θ_x，θ_y，θ_z)

式中，α为椭球单元体的长轴参数，β为椭球单元体的长短轴比例因子，x_c为椭球体单元体几何中心点的X坐标值，y_c为椭球体单元体几何中心点的Y坐标值，z_c为椭球体单元体几何中心点的Z坐标值，θ_x为椭球单元体的长轴与YOZ坐标平面的夹角，θ_y为椭球单元体的长轴与XOZ坐标平面的夹角，θ_z为椭球单元体的长轴与XOY坐标平面的夹角。

步骤3构建支架微观孔结构的负模型

根据微观孔结构负模型的数学描述公式(2)，设计构成负模型的单元体参数及数量的约束条件：

|α_iβ_i-α_jβ_j|＜2(c_i-c_j)＜α_i+α_j    (6)

L_i＞0.5α_iβ_i                         (7)

$n\&GreaterEqual;\frac{\mathrm{\&Sigma;}{v}_i}{\mathrm{\&eta;}{V}_z}---\left(8\right)$

根据公式(5)和公式(8)计算负模型中椭球体的数量n，确定个体编码串的长度8n；根据公式(6)和公式(7)的约束条件，采用随机函数均匀随机产生M*n个椭球单元体模型，从而构成M组个体；每个椭球单元体模型的参数取值范围如下：

α∈[100，400]，β∈[0.5，1]，xc∈[0，X_b]，yc∈[0，Y_b]，zc∈[0，Z_b]，(θx，θy，θz)∈[0，π/2]

步骤4求解支架微观孔结构的负模型

(a)根据表1的编码方案，对M组个体分别进行浮点数编码，从而创建初始群体P(t)；

(b)根据选择算子Φ、交叉算子Γ、变异算子Ψ、扰动算子D对群体P(t)进行选择、交叉、变异和扰动运算，得到新的群体P(t+1)；根据公式(4)计算每个个体的适应度Fit(v)_j，j＝1，2，...，M；

(c)搜索群体P(t+1)中适应度最大的个体，如果第j个个体适应度最大，且Fit(v)_j≥η，则转入(d)；否则，P(t+1)→P(t)，转入(b)；

(d)对第j个个体进行解码，得到n个椭球体的参数组，输出微观孔的负模型M_sub；

步骤5支架内部微观孔结构的生成

将模型M_sld与M_sub进行布尔减操作，得到含微观孔结构的仿生骨支架M_pore，算法结束。

下面以实施例来说明本发明所提出的算法：如图3所示，利用医学图像处理系统构建某缺损骨的仿生骨支架模型；并基于OBB原理计算其最小包容盒(见图3中的虚线部分)，作为微观孔结构负模型的约束空间；然后，基于背包问题模型的构造思想，以椭球作为单元体不断地“填充”约束空间，并运用混合遗传算法进行求解，得到微观孔结构的负模型，如图4所示；最后，通过不含微观孔的仿生骨支架模型和负模型之间的布尔减运算，构造出含微观孔结构的仿生骨支架模型，如图5所示。

以上所述仅为本发明的较佳实施而已，并不用以限制本发明。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

本说明书中未作详细描述的内容属于本领域专业技术人员公知的现有技术。

Claims (3)

1.一种仿生骨支架内部微观孔结构的建模方法，其特征在于建模步骤为：

（1）基于不含微观孔的仿生骨支架模型，运用OBB原理计算仿生骨支架模型的最小包容盒；

（2）在最小包容盒的约束空间内，以椭球体作为负模型的单元体，构建支架微观孔结构的单元体模型；具体步骤为：以椭球作为基本单元体，来构建支架微观孔结构的单元体模型，该模型定义为一个八元组：

式中，ELL表示椭球单元体模型，α是椭球的长轴参数，β是长短轴比例因子，(x _c , y _c , z _c )是椭球中心坐标，(θ _x , θ _y , θ _z )是长轴与三个坐标平面的夹角

（3）以椭球单元体模型作为负模型的构建单元，基于多约束背包问题模型，构建支架微观孔结构的负模型；具体步骤为：基于多约束背包问题模型，来构建支架内部微观孔结构的负模型，该模型的数学表达如下：

式中，Maxf(*)表示求(v ₁ ,v ₂ ,…，v _n )体积和的最大值，v _i 表示第i个椭球体在模型中的实际体积，η表示孔隙率，V _z 表示不含微观孔的支架体积，a _i 、b _i 表示第i个椭球体的长、短轴参数，c _i -c _j 表示相交的椭球体i、j之间的中心距，L _i 表示第i个椭球体的几何中心到支架边界的距离；

（4）针对所构建的支架微观孔结构的负模型，利用混合遗传算法求解支架微观孔结构负模型的最优解；具体步骤为：利用混合遗传算法求解微观孔结构负模型的最优解，该算法定义为一个九元组：

式中，HGA表示混合遗传算法模型，C表示个体的编码形式，E表示个体的适应度函数，P ₀ 表示初始群体，M表示群体大小，Φ是选择算子，Г是交叉算子，Ψ是变异算子，D是扰动算子，T是算法的终止条件

（5）将不含微观孔的仿生骨支架模型与最优化的负模型进行布尔减运算，以生成支架的内部微观孔结构。

2.根据权利要求1所述的仿生骨支架内部微观孔结构的建模方法，其特征在于：所述混合遗传算法中，以孔隙率作为个体适应度的评价指标，个体适应度函数定义为：

式中，Fit(*)表示求*的个体适应度值，v _i 表示第i个椭球体在模型中的实际体积，η表示孔隙率，V _z 表示不含微观孔的支架体积。

3.根据权利要求1所述的仿生骨支架内部微观孔结构的建模方法，其特征在于：所述混合遗传算法中，所设计的遗传算子包括选择算子Φ、交叉算子Г、变异算子Ψ、扰动算子D；Φ采用比例选择与最优保存策略相结合的混合选择算子；Г采用均匀交叉算子；Ψ采用均匀变异算子；D设计为以扰动概率P _d =0.05产生一个满足约束的随机个体，替代群体中适应度最高的个体。

Images