CN105167849B - 一种基于蚁群算法的血管三维路径规划方法 - Google Patents

Abstract

本发明公开了一种基于蚁群算法的血管三维路径规划方法，本方法的操作步骤如下：1数据导入、2血管建模、3血管中心线提取、4建立血管中心线网络拓扑结构、5蚂蚁侦察算法、6参数初始化、7启发式信息计算、8概率选择、9信息素动态挥发、10信息素增量计算、11信息素更新、12规划结束判断、13结果输出。该方法在血管中心线提取的基础上采用蚂蚁侦察算法，通过改进蚁群算法，并综合考虑导管直径、血管长度、最小直径、最大曲率和最大挠率辅助外科医生规划手术的最优路径。该方法提血管介入手术术前路径规划的可靠性，保证了导管的通过性，能够为外科医生提供一种新的手术路径参考标准。

Description

一种基于蚁群算法的血管三维路径规划方法

技术领域

本发明涉及的是一种基于蚁群算法的血管三维路径规划方法，属于血管三维路径规划技术领域。

背景技术

血管介入手术具有出血少，创伤小，恢复快的优点，因此在血管疾病的治疗上得到广泛应用。目前全国大部分医院的血管介入手术还是依靠外科医生徒手操作导管导丝并借助X射线等医学成像技术和虚拟现实等计算机辅助技术完成手术。在有限的技术条件下，从医生的操作水平、手术的危险性、术后并发症等因素考虑，术中通常选用较粗的主动脉作为手术路径，不考虑其他可能的最优路径。

血管三维路径规划是在三维几何空间内，寻找一条从手术切口到病灶部位的一条最优路径，是一种全局规划问题，是术前规划的重要手段。血管内三维路径规划的关键技术是血管中心线的提取和导航路径的规划。提取中心线的方法有拓扑细化法、距离变换法、势能场法、基于水平集法和基于分割的方法等。导航路径规划方面，随着智能算法的兴起，有学者采用Dijkstra算法和A星算法在提取的血管中心线上搜索到两点之间的最短路径。但这些方法在规划过程中并没有考虑到血管的其他特性。

蚁群算法是仿生学中群体智能算法之一，它是1991年由意大利学者M.Dorigo等人受到真实世界中蚂蚁觅食的行为的启发后而提出。算法核心即蚂蚁能以更高的概率选中残留信息素浓度最高的路径，越来越多的蚂蚁会被吸引到这条路径上，形成一种正反馈原理来找到一条离食物源最短的路径。目前国内外基于蚁群算法的三维路径规划主要集中在无人机、潜水器、三维管道、机器人等方面，暂没有文献将蚁群算法应用在血管三维路径规划方面。

发明内容

本发明的目的在于针对现有血管路径规划技术存在的问题和不足，提供一种基于蚁群算法的三维路径规划方法，该方法在血管中心线提取的基础上采用蚂蚁侦察算法，通过改进蚁群算法，并综合考虑导管直径、血管长度、直径、曲率和挠率辅助外科医生规划手术的最优路径。

为达到上述目的，本发明的构思是：对CTA医学影像数据进行血管建模并提取中心线，在给定起始点的血管中心线网络中在此基础上采用蚂蚁侦察算法，通过改进蚁群算法，综合导管直径和血管其他特性，如血管长度、最小直径、最大曲率和最大挠率因素进行全局路径规划。

根据本发明的构思，本发明采用以下技术方案实现：

一种基于蚁群算法的血管三维路径规划方法，其特征在于操作步骤如下：1数据导入、2血管建模、3血管中心线提取、4建立血管中心线网络拓扑结构、5蚂蚁侦察算法、6参数初始化、7启发式信息计算、8概率选择、9信息素动态挥发、10信息素增量计算、11信息素更新、12规划结束判断、13结果输出。

所述步骤1数据导入：导入一组完整的医学CTA图像数据(DICOM)；

所述步骤2血管建模：将导入的CTA图像数据采用基于3D水平集的方法进行血管建模；

所述步骤3血管中心线提取：将建立的血管模型采用基于Voronoi图和Eikonal方程的方法提取血管中心线，并求得最大内切球半径R；

所述步骤4建立血管中心线网络拓扑结构：根据中心线提取结果建立网络拓扑结构，确定结点和边的数量，并计算各段路径的距离、最大曲率、最大挠率和最小直径；

4-1距离计算，采用第一类曲线积分来计算空间曲线的长度。假设参数化形式的空间曲线方程为

因此可得空间曲线的距离公式为：

4-2最大曲率计算，曲率(curvature)就是针对曲线上某个点的切线方向角对弧长的转动率算出中心线上每个点曲线的曲率，中心线的曲率反映了血管弯曲的程度，而曲率越大则血管弯曲程度越大，必将增加导管通过的难度，降低血管的可通过性。对于参数化形式空间曲线方程的曲率计算为：

根据此公式计算该路径上每个点的曲率并求得该路径的最大曲率C_max；

4-3最大挠率计算，挠率(torsion)的绝对值度量了曲线上邻近两点的次法向量之间的夹角对弧长的变化率,中心线上的挠率可以反映出血管扭曲的程度，对于盆腔内血管的扭曲程度是比较大的，而挠率越大血管扭曲程度越大，也将增加导管通过的难度，降低血管可通过性，降低安全系数。对于参数化形式空间曲线方程的挠率计算为：

根据此公式计算该路径上每个点的挠率并求得该路径的最大挠率T_max；

4-4最小直径计算，根据最大内切球半径R可计算中心线上每个点的直径D，由此计算各路径上的最小直径D_min，假设手术导管的直径为D_c，要求这段血管路径的最小直径D_min＞D_c以保证导管能够穿过血管；

所述步骤5蚂蚁侦察算法：确定起始点，执行具有自我复制功能的蚂蚁侦察算法。

由于血管网络拓扑结构不是一个完全图，有的结点之间不存在回路，必然存在一些端结点。若蚂蚁在搜寻过程中，走到了一个端结点而并非终点时，要么就直接退出寻路，要么就退回到端结点的上一结点继续寻路。若选择前者，虽然会加快收敛，但是欠缺全局考虑，过早退出；若选择后者，虽然从全局考虑，但会降低算法效率。因此，本发明在确定起始点后，让蚂蚁先对整个网络拓扑结构进行侦察，区分并去除这类非终点端结点后再执行蚁群算法。这类侦察蚂蚁具备自我复制功能，在侦察过程中可以根据顶点V的度数N，即边数，相应自我复制(除了已访问的边之外)，再分别往其他顶点继续侦察。

5-1根据起点v和终点u的度数分别放置对应数量的蚂蚁Mv＝Nv，Mu＝Nu；

5-2每只蚂蚁侦察顶点V_i时先将其保存至各自的路由表Ta_i中，做侦察标记；

5-3每只蚂蚁判断侦察的顶点V_i的度数N_i，

5-3-1若n≠1，则复制n-2只蚂蚁副本，继续侦察；

5-3-2若n＝1，则报告该顶点V_i为异常点，剔除该点，修改其上一顶点V_i-1的度数N_i-1，更新整体网络拓扑结构，进入步骤5-4；

5-4根据路由表Ta_i沿原路退回到上一顶点V_i-1，若n＝1，则执行步骤5-3-2；若n≠1，则进入步骤5-5；

5-5若所有顶点都已侦察完，则终止侦察，进入参数初始化6；否则进入步骤5-3；

所述步骤6参数初始化：初始化各参数，如设置蚂蚁数量M，信息素强度常量Q，迭代总数，导管直径D_c，将蚂蚁置于初始位置；

所述步骤7启发式信息计算：启发式函数对算法的收敛性和稳定性有着重要的影响，在蚁群算法中，启发式信息η一般是两点之间距离的反比，即1/L_I，J，本发明规划路径的标准是手术安全和路径最优，因此路径规划中不仅要考虑到血管的长度L_I，J、最大曲率C_max和最大挠率T_max，还要考虑到血管的最小直径D_min，启发式信息定义为

η_I，J(t)表示t时刻结点I(x，y，z)到结点J(x′，y′，z′)的启发信息；x，y，z和x′，y′，z′为点的三维坐标；D_min表示这条路径的最小直径，D_c表示导管直径，要求D_min＞D_c；L_I，J表示这条的路径长度；C_max表示这段路径的最大曲率，T_max，表示这段路径的最大挠率。

所述步骤8概率选择：设有M只蚂蚁，n个结点，a_i(t)表示t时刻位于某结点I的蚂蚁数量，则在蚁群算法的每一步路径选择中，蚂蚁m按照概率决定下一步往哪条路上移动。

表示第t时刻蚂蚁m由点I移动到点J的概率，allowed_m＝{1，2，...，n}-t_m表示蚂蚁下一步可以选择的节点，t_m表示蚂蚁曾经访问过的节点；α表示信息素启发因子，表示残留信息素量在蚂蚁运动时所起的作用，值越大说明蚂蚁更倾向于选择其它蚂蚁经过的路径，蚂蚁之间的协作能力越强；β表示期望值启发因子，表明启发信息在蚂蚁选择路径时的受重视程度，值越大，状态转移规律越接近贪心规则。τ_I，J(t)表示t时刻点I到点J路径上的信息素。

所述步骤9信息素动态挥发：信息素挥发速度会随着时间的推移受到温度、湿度等因素而变化，是一个动态的变化的过程，越复杂的路径，启发信息越小，挥发速度越快，信息素残留越少，因此信息素的挥发系数ρ随着启发式信息η_I，J(t)的变化而动态挥发，所有蚂蚁到达终点后，按如下关系计算，

ρ_I，J(t)表示t时刻在I到J路径上的挥发系数，η_I，J(t)表示t时刻I到J路径上的信息素，表示t时刻所有路径上的启发信息，k表示某只蚂蚁在t时刻所经过的路径数。

所述步骤10信息素增量计算：信息素更新模型是基本蚁群算法的随机搜索与快速收敛重要环节。本发明针对问题本身的全局性最优要求，采用蚁周模型；考虑到血管最小直径和导管直径在手术中的影响，修改蚁周模型如下：

Q为信息素强度常量，是蚂蚁在一个循环中在所经过的路径上释放的信息素的总量，在一定程度上影响着算法的收敛速度；要求血管路径的最小直径D_min＞D_c，L_m表示第m只蚂蚁在此次循环中所经过的路径长度，根据每只蚂蚁在其路由表Ta_m中的结果计算。

所述步骤11信息素更新：初始时刻各路径上的信息素量是相等的，当蚂蚁完成一次循环后，信息素会随着时间的推移逐渐挥发，因此要对信息素浓度要进行更新，在蚂蚁进入下一个循环之前对相应路径上的信息素做如下更新：

τ_I，J(t+1)＝(1-ρ)τ_I，J(t)+Δτ_I，J(t，t+1)

ρ(O＜ρ＜1)是信息素动态挥发系数，1-ρ是信息素动态残留因子，表示第m只蚂蚁在本次循环中在路径(I,J)上留下的信息素量，Δτ_I，J(t，t+1)表示本次循环当中所有经过路径(I,J)的蚂蚁留下的信息素增量。

所述步骤12规划结束判断：当达到最大迭代次数则退出循环，否则进入启发式信息计算7，继续循环；

所述步骤13结果输出：整理路由表输出结果。

本发明与现有技术相比较，具有如下显而易见的突出实质性特点和显著优点：

(1)本发明将蚁群算法应用到血管三维路径规划技术领域。

(2)本发明改进了蚁群算法中的启发式信息模型，综合考虑了导管直径、血管长度、最小直径、最大曲率和最大挠率，提高了术前规划路径的可靠性，保证了导管的通过性。

(3)本发明改进了蚁群算法中的信息素增量模型，结合了血管最小直径和导管直径在规划中的影响，加快了算法的收敛。

附图说明

图1为本发明一种基于蚁群算法的血管三维路径规划方法的流程图。

具体实施方式

下面结合附图对本发明的实施方式作进一步详细的说明。

实施例一：

参见图1，一种基于蚁群算法的血管三维路径规划方法，其特征在于操作步骤如下,1数据导入、2血管建模、3血管中心线提取、4建立血管中心线网络拓扑结构、5蚂蚁侦察算法、6参数初始化、7启发式信息计算、8概率选择、9信息素动态挥发、10信息素增量计算、11信息素更新、12规划结束判断、13结果输出。

所述步骤1数据导入：导入一组完整的医学CTA图像数据(DICOM)；

所述步骤2血管建模：将导入的CTA图像数据采用基于3D水平集的方法进行血管建模；

所述步骤3血管中心线提取：将建立的血管模型采用基于Voronoi图和Eikonal方程的方法提取血管中心线，并求得最大内切球半径R；

所述步骤4建立血管中心线网络拓扑结构：根据中心线提取结果建立网络拓扑结构，确定结点和边的数量，并计算各段路径的距离、最大曲率、最大挠率和最小直径；

4-1距离计算，采用第一类曲线积分来计算空间曲线的长度；

4-2最大曲率计算，根据参数化形式空间曲线方程的曲率计算公式计算该路径上每个点的曲率并求得该路径的最大曲率C_max；

4-3最大挠率计算，根据参数化形式空间曲线方程的挠率计算公式计算该路径上每个点的挠率并求得该路径的最大挠率T_max；

4-4最小直径计算，根据最大内切球半径R计算中心线上每个点的直径D，由此计算各路径上的最小直径D_min，假设手术导管的直径为D_c，则要求这段血管路径的最小直径D_min＞D_c以保证导管能够穿过血管。

所述步骤5蚂蚁侦察算法：在确定起始点后，让蚂蚁先对整个网络拓扑结构进行侦察，区分并去除这类非终点端结点后再执行蚁群算法。这类侦察蚂蚁具备自我复制功能，在侦察过程中可以根据顶点V的度数N，即边数，相应自我复制(除了已访问的边之外)，再分别往其他顶点继续侦察；

5-1根据起点v和终点u的度数分别放置对应数量的蚂蚁Mv＝Nv，Mu＝Nu；

5-2从顶点集V中选择一个顶点V_i，判断顶点V_i是否已侦察；

5-2-1若V_i没有侦察，则每只蚂蚁侦察顶点V_i时先将其保存至各自的路由表Ta_j中，做侦察标记；

5-2-2若V_i已经侦察，则返回到5-2；

5-3每只蚂蚁判断侦察的顶点V_i的度数N_i，

5-3-1若n≠1，则复制n-2只蚂蚁副本，继续侦察；

5-3-2若n＝1，则报告该顶点V_i为异常点，剔除该点，修改其上一顶点V_i-1的度数N_i-1，更新整体网络拓扑结构，进入步骤5-4；

5-4根据路由表Ta_i沿原路退回到上一顶点V_i-1，若n＝1，则执行步骤5-3-2；若n≠1，则进入步骤5-5；

5-5若所有顶点都已侦察完，则终止侦察，进入参数初始化6；否则进入步骤5-3；

所述步骤6参数初始化：初始化各参数，如设置蚂蚁数量M，信息素强度常量Q，迭代总数，导管直径D_c，将蚂蚁置于初始位置。

实施例二：本实施例与实施例一基本相同，特别之处如下：

所述步骤7启发式信息计算：本发明规划路径的标准是手术安全和路径最优，因此路径规划中不仅要考虑到血管的长度L_I，J、最大曲率C_max和最大挠率T_max，还要考虑到血管的最小直径D_min；

所述步骤8概率选择：在蚁群算法的每一步路径选择中，蚂蚁m按照概率公式的比率决定下一步往哪条路上移动；

所述步骤9信息素动态挥发：信息素挥发速度会随着时间的推移受到温度、湿度等因素而变化，是一个动态的变化的过程，越复杂的路径，启发信息越小，挥发速度越快，信息素残留越少，因此信息素的挥发系数ρ随着启发式信息η_I，J(t)的变化而动态挥发。

实施例三：本实施例与实施例二基本相同，特别之处如下：

所述步骤10信息素增量计算：信息素更新模型是基本蚁群算法的随机搜索与快速收敛重要环节。本发明针对问题本身的全局性最优要求，采用蚁周模型；考虑到血管最小直径和导管直径在手术中的影响，按照改进后信息素增量模型进行计算；

所述步骤11信息素更新：初始时刻各路径上的信息素量是相等的，当蚂蚁完成一次循环后，信息素会随着时间的推移逐渐挥发，因此要对信息素浓度要进行更新，在蚂蚁进入下一个循环之前对相应路径上的信息素做相应更新；

所述步骤12规划结束判断：当达到最大迭代次数则退出循环，否则进入启发式信息计算7，继续循环；

所述步骤13结果输出：整理路由表输出结果。

Claims (4)

1.一种基于蚁群算法的血管三维路径规划方法，其特征在于操作步骤如下：1数据导入、2血管建模、3血管中心线提取、4建立血管中心线网络拓扑结构、5蚂蚁侦察算法、6参数初始化、7启发式信息计算、8概率选择、9信息素动态挥发、10信息素增量计算、11信息素更新、12规划结束判断、13结果输出；

所述步骤1数据导入：导入一组完整的医学CTA图像数据；

所述步骤2血管建模：将导入的CTA图像数据采用基于3D水平集的方法进行血管建模；

所述步骤3血管中心线提取：将建立的血管模型采用基于Voronoi图和Eikonal方程的方法提取血管中心线，并求得最大内切球半径R；

所述步骤4建立血管中心线网络拓扑结构：根据中心线提取结果建立网络拓扑结构，确定结点和边的数量，并计算各段路径的距离、最大曲率、最大挠率和最小直径；

4-1距离计算，采用第一类曲线积分来计算空间曲线的长度；

4-2最大曲率计算，根据参数化形式空间曲线方程的曲率计算公式计算该路径上每个点的曲率并求得该路径的最大曲率C_max；

4-3最大挠率计算，根据参数化形式空间曲线方程的挠率计算公式计算该路径上每个点的挠率并求得该路径的最大挠率T_max；

4-4最小直径计算，根据最大内切球半径R计算中心线上每个点的直径D，由此计算各路径上的最小直径D_min，假设手术导管的直径为D_c，则要求这段血管路径的最小直径D_min＞D_c以保证导管能够穿过血管；

所述步骤5蚂蚁侦察算法：在确定起始点后，让蚂蚁先对整个网络拓扑结构进行侦察，区分并去除这类非终点端结点后再执行蚁群算法；这类侦察蚂蚁具备自我复制功能，在侦察过程中能够根据顶点V的度数N，即边数，相应自我复制，再分别往其他顶点继续侦察；

5-1根据起点v和终点u的度数分别放置对应数量的蚂蚁Mv＝Nv，Mu＝Nu；

5-2从顶点集V中选择一个顶点V_i，判断顶点V_i是否已侦察；

5-2-1若V_i没有侦察，则每只蚂蚁侦察顶点V_i时先将其保存至各自的路由表Ta_i中，做侦察标记；

5-2-2若V_i已经侦察，则返回到5-2；

5-3每只蚂蚁判断侦察的顶点V_i的度数N_i，

5-3-1若n≠1，则复制n-2只蚂蚁副本，继续侦察；

5-3-2若n＝1，则报告该顶点V_i为异常点，剔除该点，修改其上一顶点V_i-1的度数N_i-1，更新整体网络拓扑结构，进入步骤5-4；

5-4根据路由表Ta_i沿原路退回到上一顶点V_i-1，若n＝1，则执行步骤5-3-2；若n≠1，则进入步骤5-5；

5-5若所有顶点都已侦察完，则终止侦察，进入参数初始化6；否则进入步骤5-3；

所述步骤6参数初始化：初始化各参数，如设置蚂蚁数量M，信息素强度常量Q，迭代总数，导管直径D_c，将蚂蚁置于初始位置；

所述步骤7启发式信息计算：本规划路径的标准是手术安全和路径最优，因此路径规划中不仅要考虑到血管的长度L_I，J、最大曲率C_max和最大挠率T_max，还要考虑到血管的最小直径D_min；

所述步骤10信息素增量计算：信息素更新模型是基本蚁群算法的随机搜索与快速收敛重要环节；针对问题本身的全局性最优要求，采用蚁周模型；考虑到血管最小直径和导管直径在手术中的影响，按照改进后信息素增量模型进行计算，修改蚁周模型如下：

Q为信息素强度常量，是蚂蚁在一个循环中在所经过的路径上释放的信息素的总量，在一定程度上影响着算法的收敛速度；要求血管路径的最小直径D_min＞D_c，L_m表示第m只蚂蚁在此次循环中所经过的路径长度，根据每只蚂蚁在其路由表Ta_m中的结果计算。

2.根据权利要求1所述的基于蚁群算法的血管三维路径规划方法，其特征在于所述步骤8概率选择：在蚁群算法的每一步路径选择中，蚂蚁m按照概率公式的比率决定下一步往哪条路上移动。

3.根据权利要求1所述的基于蚁群算法的血管三维路径规划方法，其特征在于所述步骤9信息素动态挥发：信息素挥发速度会随着时间的推移受到温度、湿度因素而变化，是一个动态的变化的过程，越复杂的路径，启发信息越小，挥发速度越快，信息素残留越少，因此信息素的挥发系数ρ随着启发式信息η_I，J(t)的变化而动态挥发。

4.根据权利要求1所述的基于蚁群算法的血管三维路径规划方法，其特征在于所述步骤11信息素更新：初始时刻各路径上的信息素量是相等的，当蚂蚁完成一次循环后，信息素会随着时间的推移逐渐挥发，因此要对信息素浓度要进行更新，在蚂蚁进入下一个循环之前对相应路径上的信息素做相应更新；

所述步骤12规划结束判断：当达到最大迭代次数则退出循环，否则进入启发式信息计算7，继续循环；

所述步骤13结果输出：整理路由表输出结果。