本发明涉及医学图像处理的技术领域，尤其涉及一种基于蒙特卡洛树搜索的血管3d/2d配准方法，以及基于蒙特卡洛树搜索的血管3d/2d配准装置。
背景技术：
目前，微创介入治疗是血管疾病的主要治疗手段，在介入治疗过程中，手术器械的操作以x射线血管造影(xra)图像为指导。造影剂通过导管注入感兴趣的动脉并进行成像，这种成像模式在显示血管管腔方面有令人满意的表现。在xra中，血管内导航的手术器械也可以清楚地显示出来。然而，由于xra缺乏空间信息，在单视图二维投影的指导下进行准确的介入操作对介入医师来说是困难的。因此，在介入手术中，医生经常使用旋转c型臂获得的多视图血管造影图像，然而这样会加大造影剂的注入且会给患者带来负担。为了解决这个问题，可以将术前计算机断层扫描血管造影(cta)图像与术中xra成像相结合使用。通过将三维血管投影覆盖在二维实时图像上来增强介入图像，医生和病人均可以从多模态数据融合和可视化中受益。为了实现这一目标，3d/2d配准技术是其中获得良好的对齐和对应关系的关键。
3d/2d配准方法通常使用一张术前的3d图像和作为配准源的多张术中2dx射线造影图来实现配准。3d和2d图像数据以及c型臂的成像几何参数需要作为配准模型的输入。对于x射线图像序列，通常认为3d/2d配准是3d/2d t的。旋转成像是临床常用的一种多平面x射线成像技术，在这种情况下，多平面3d/2d配准可视为3d/2.5d。由于3d/2d配准是3d/2d t和3d/2.5d的基础，因此本文只讨论3d图像到单帧单平面2d图像的配准。根据配准技术的性质，将3d/2d配准方法分为基于灰度的和基于特征的配准方法。
基于灰度的方法通常是通过优化3d术前图像和2d术中图像投影的相似性测度来实现的。数字重建放射成像(drr)和最大强度投影(mip)是两种常见的通过计算机断层成像(ct)图像生成模拟x射线投影的方法。hipwell分析了drr/mip与x射线图像的六种相似度指标，其中模式强度和梯度差异表现最好。考虑到基于灰度的方法使用整个图像灰度信息进行配准，因此它们对背景异常值非常敏感。此外，当配准具有大尺度变换的数据时，基于优化的方法可能产生较小的捕获范围。
基于特征的3d/2d配准依赖于从两种模态的图像中提取的具有一致性特征。中心线是血管配准中最常用的特征表示。迭代最近点(icp)方法将点云配准分解为包含匹配和配准阶段的交替连续的过程。对于匹配阶段，可以通过查找点的最小欧氏距离来分配点对应。baka介绍了一种将二维点反投影到三维空间，然后执行经典icp过程的3d/2d配准方法。rivest-henault利用血管中心预先计算距离变换，并构造目标函数，可以加速配准过程。benseghir提出了迭代最近邻曲线(icc)方法，该方法利用血管分支作为基于最近邻关系的配对元素，然后估计出最小化配对分支之间距离和的变换。这些类icp方法对噪声和异常值很敏感，因为它们将对应关系限制为一对一的匹配。基于最近邻关系的匹配也会导致这些方法严重依赖于初始姿态。
由于硬分配策略对噪声和异常值的敏感性，在概率分配框架中，软分配策略将一对一对应关系松弛为一对多。在高斯混合模型(gmm)和期望最大化(em)算法的基础上，myronenko提出了相干点漂移(cpd)方法，该方法强制gmm中心体进行相干移动，以保持点集的全局拓扑结构。kang使用与cpd相同的框架进行3d/2d点云配准。在考虑透视投影的非线性特性时，采用粒子群优化算法求解配准参数的最优估计。baka提出了一种从jianandvemuri的方法扩展而来的方向约束ogmm方法，利用方向和位置估计两个点集的l2距离，并对其进行优化，实现3d/2d冠状动脉配准。由于ogmm利用了中心线的定位，因此对噪声数据具有较高的精度和鲁棒性。
由于血管拓扑是跨形态和维度的不变属性，因此图匹配成为一种有效的血管配准方法。血管图的匹配可以描述为先估计血管分叉点的对应关系，然后将分叉点视为顶点来匹配连接它们之间的曲线。serradell将血管配准描述为寻找最大可能对应关系的搜索过程，并使用先验搜索来加速该过程。
pinheiro基于待匹配图的拓扑一致性，将血管匹配描述为一种树搜索方法，并采用蒙特卡罗树搜索来解决该问题。moriconi通过设置节点和边属性来定义图匹配的亲密度函数，然后将二次分配问题的函数最大化，得到节点匹配。
上述几个方法均与3d/3d或2d/2d血管系统匹配/配准相关。对于3d/2d血管匹配，投影中的重叠问题决定了该方法需要强调异常值和噪声。icc方法通过两个节点之间的配对曲线来使用血管拓扑。对于有噪声的二维图，利用邻域关系可以约束候选对象。通过保证分叉点的连通性，benseghir在icc方法的基础上提出了一种分而治之的树形中心线匹配方法。liu将3d和2d血管视为树形拓扑，将树形拓扑表示为序列。然后，通过序列对序列的匹配，实现节点的血管匹配。通过保持前后关系，利用拓扑排序算法提取序列，然后按顺序遍历。
由于透视投影的非线性性质，基于数值优化的3d/2d配准容易陷入局部极值，使得这些方法对初始配准位姿敏感。初始化操作是3d/2d配准的关键。一般来说，术前和术中图像来自不同的设备。大多数配准方法的捕获范围不足以覆盖两种采集设备的坐标系统之间的转换。markelj介绍了一些初始化方法，如患者位置和方向的对齐、相应标记对的注册以及手动初始化。然而，一种利用内在特征的自动初始化方法可能更适用于术中3d/2d的配准。varnavas预先计算了大范围3d位姿的二维投影模板，并使用术中图像评估模板与二维透视的广义霍夫变换之间的相似性，以获得初始对准。miao建立了从金属植入体轮廓中提取的二维轮廓的形状上下文编码库，并采用jensen-shannondivergence算法作为快速库匹配的匹配度量。gouveia提出了一种基于回归的cta和xra的初始配准方法，该方法将二维投影图像特征与刚性变换参数联系起来。
为了获得较大的捕获范围或对初始位姿不敏感的配准方法，基于“先匹配再求变换”的框架比基于优化的方法更合适。假定已知3d和2d点的对应关系，3d/2d配准的目标类似于计算机视觉领域的pnp(perspective-n-point)问题。它的目的是通过将对点变换的估计问题简化为对四个控制点坐标的估计问题来确定摄像机的姿态。
技术实现要素：
为克服现有技术的缺陷，本发明要解决的技术问题是提供了一种基于蒙特卡洛树搜索的血管3d/2d配准方法，其能够不依赖配准的初始位姿，3d图像和2d图像的血管图匹配精确度高。
本发明的技术方案是：这种基于蒙特卡洛树搜索的血管3d/2d配准方法，其包括以下步骤：
(1)利用3d和2d血管拓扑一致性来实现血管图匹配，将3d血管图和2d血管图的匹配表示为顶点对的集合和边对的集合；
(2)将血管匹配过程看作求解3d和2d血管点的稠密匹配，先估计顶点之间的稀疏匹配，再匹配连接两个顶点之间的边，从而得到两个血管图的稠密匹配；
(3)根据血管边匹配的可以分解为连续状态的特性，构建一个搜索树；
(4)在搜索树的每个节点上，使用一个封闭解来计算配准结果，设计一个评价分数；
(5)使用一个蒙特卡洛树搜索mcts的变种用于遍历该搜索树，评价分数作为树中每个节点的奖励值，mcts的目标是在搜索树空间中找到最具有最高奖励的节点。
本发明利用3d和2d血管拓扑一致性来实现血管图匹配，将3d血管图和2d血管图的匹配表示为顶点对的集合和边对的集合；将血管匹配过程看作求解3d和2d血管点的稠密匹配，先估计顶点之间的稀疏匹配，再匹配连接两个顶点之间的边，从而得到两个血管图的稠密匹配；根据血管边匹配的可以分解为连续状态的特性，构建一个搜索树；在搜索树的每个节点上，使用一个封闭解来计算配准结果，设计一个评价分数；使用一个蒙特卡洛树搜索mcts的变种用于遍历该搜索树，评价分数作为树中每个节点的奖励值，mcts的目标是在搜索树空间中找到最具有最高奖励的节点，从而能够不依赖配准的初始位姿，3d图像和2d图像的血管图匹配精确度高。
还提供了基于蒙特卡洛树搜索的血管3d/2d配准装置，其包括：
匹配结果表示模块，其利用3d和2d血管拓扑一致性来实现血管图匹配，将3d血管图和2d血管图的匹配表示为顶点对的集合和边对的集合；
稠密匹配模块，其将血管匹配过程看作求解3d和2d血管点的稠密匹配，先估计顶点之间的稀疏匹配，再匹配连接两个顶点之间的边，从而得到两个血管图的稠密匹配；
搜索树的构造模块，其根据血管边匹配的可以分解为连续状态的特性，构建一个搜索树；
计算结果及节点评分模块，其在搜索树的每个节点上，使用一个封闭解来计算配准结果，设计一个评价分数；
最优结果寻找模块，其使用一个蒙特卡洛树搜索mcts的变种用于遍历该搜索树，评价分数作为树中每个节点的奖励值，mcts的目标是在搜索树空间中找到最具有最高奖励的节点。
附图说明
图1a是待配准的3d血管模型；图1b是3d血管中心线；图1c是2d造影图像；图1d是2d血管中心线。
图2是基于mcts的血管匹配的一次迭代。在根节点上选择一对血管顶点，然后通过选择、扩展、仿真和反向传播四个步骤对树进行扩展。
图3示出了根据本发明的基于蒙特卡洛树搜索的血管3d/2d配准方法的流程图。
具体实施方式
如图3所示，这种基于蒙特卡洛树搜索的血管3d/2d配准方法，其包括以下步骤：
(1)利用3d和2d血管拓扑一致性来实现血管图匹配，将3d血管图和2d血管图的匹配表示为顶点对的集合和边对的集合；
(2)将血管匹配过程看作求解3d和2d血管点的稠密匹配，先估计顶点之间的稀疏匹配，再匹配连接两个顶点之间的边，从而得到两个血管图的稠密匹配；
(3)根据血管边匹配的可以分解为连续状态的特性，构建一个搜索树；
(4)在搜索树的每个节点上，使用一个封闭解来计算配准结果，设计一个评价分数；
(5)使用一个蒙特卡洛树搜索mcts的变种用于遍历该搜索树，评价分数作为树中每个节点的奖励值，mcts的目标是在搜索树空间中找到最具有最高奖励的节点。
本发明利用3d和2d血管拓扑一致性来实现血管图匹配，将3d血管图和2d血管图的匹配表示为顶点对的集合和边对的集合；将血管匹配过程看作求解3d和2d血管点的稠密匹配，先估计顶点之间的稀疏匹配，再匹配连接两个顶点之间的边，从而得到两个血管图的稠密匹配；根据血管边匹配的可以分解为连续状态的特性，构建一个搜索树；在搜索树的每个节点上，使用一个封闭解来计算配准结果，设计一个评价分数；使用一个蒙特卡洛树搜索mcts的变种用于遍历该搜索树，评价分数作为树中每个节点的奖励值，mcts的目标是在搜索树空间中找到最具有最高奖励的节点，从而能够不依赖配准的初始位姿，3d图像和2d图像的血管图匹配精确度高。
优选地，所述步骤(1)中，血管图表示为
优选地，所述步骤(1)中，根据公式(1),3d和2d血管结构的配准表示为找到3d血管的最优变换t，在该变换t作用后的3d血管的投影与2d血管达到最佳对齐，
其中
3d点集
所述步骤(2)中，3d血管通常被认为是一个无环图，2d血管则用有环图来表示，因为在2d平面上有多个投影的三维分支重叠，从而形成伪分叉和连接边。因此，3d和2d的血管拓扑通常是不同的
优选地，所述步骤(2)中，两个边上点的稠密匹配通过均匀插值实现，将该过程描述为从
优选地，所述步骤(3)中，根据
优选地，所述步骤(4)中，评价分数q为公式(2)
对于给定的
本方法使用一个蒙特卡洛树搜索(mcts)的变种用于遍历该搜索树，公式(2)中定义的分数作为树中每个节点的奖励值。mcts的目标是在搜索树空间
优选地，所述步骤(5)包括以下分步骤：
(5.1)选择：在当前搜索树中选择最紧急的可拓展节点，从根节点开始自顶向下，在每一层通过贪婪策略实现；可拓展节点意味着根据当前节点的匹配状态找到至少一对新的匹配边；节点的紧急度通过公式(3)计算
其中qsim表示以当前结点作为根节点的子树的可能最大的奖励值，n表示一共迭代的次数，n表示当前结点的访问次数；因此，第二项用于惩罚当前结点过高访问次数，以鼓励对未访问或者较少访问分支的探索；γ是用来平衡对已知的节点的开发和对未访问结点的探索；
(5.2)拓展：对于选中的节点，对其添加至多nexp个子节点；当存在大于nexp的子节点，计算他们的奖励，并选择较高的nexp个子节点；
(5.3)模拟：将每个新扩展的子节点作为根节点，以深度优先的方式向下访问节点；在每一层中，随机选择子节点，直到到达最深处的叶子节点；这个过程重复nsim次，然后计算这些访问叶节点的奖励，将最高奖励作为扩展节点的qsim；
(5.4)反向传播：沿着被选中节点到根节点的路径，更新qsim值和节点访问次数n；
这四个步骤迭代地进行直至达到最大的计算量nmax或者理论上最高的奖励qmax，迭代结束后，具有最高的奖励q的节点从搜索中获得，且节点中对应的匹配
本领域普通技术人员可以理解，实现上述实施例方法中的全部或部分步骤是可以通过程序来指令相关的硬件来完成，所述的程序可以存储于一计算机可读取存储介质中，该程序在执行时，包括上述实施例方法的各步骤，而所述的存储介质可以是：rom/ram、磁碟、光盘、存储卡等。因此，与本发明的方法相对应的，本发明还同时包括一种基于蒙特卡洛树搜索的血管3d/2d配准装置，该装置通常以与方法各步骤相对应的功能模块的形式表示。该装置包括：
匹配结果表示模块，其利用3d和2d血管拓扑一致性来实现血管图匹配，将3d血管图和2d血管图的匹配表示为顶点对的集合和边对的集合；
稠密匹配模块，其将血管匹配过程看作求解3d和2d血管点的稠密匹配，先估计顶点之间的稀疏匹配，再匹配连接两个顶点之间的边，从而得到两个血管图的稠密匹配；
搜索树的构造模块，其根据血管边匹配的可以分解为连续状态的特性，构建一个搜索树；
计算结果及节点评分模块，其在搜索树的每个节点上，使用一个封闭解来计算配准结果，设计一个评价分数；
最优结果寻找模块，其使用一个蒙特卡洛树搜索mcts的变种用于遍历该搜索树，评价分数作为树中每个节点的奖励值，mcts的目标是在搜索树空间中找到最具有最高奖励的节点。
以上所述，仅是本发明的较佳实施例，并非对本发明作任何形式上的限制，凡是依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与修饰，均仍属本发明技术方案的保护范围。
1.基于蒙特卡洛树搜索的血管3d/2d配准方法，其特征在于：其包括以下步骤：
(1)利用3d和2d血管拓扑一致性来实现血管图匹配，将3d血管图和2d血管图的匹配表示为顶点对的集合和边对的集合；
(2)将血管匹配过程看作求解3d和2d血管点的稠密匹配，先估计顶点之间的稀疏匹配，再匹配连接两个顶点之间的边，从而得到两个血管图的稠密匹配；
(3)根据血管边匹配的可以分解为连续状态的特性，构建一个搜索树；
(4)在搜索树的每个节点上，使用一个封闭解来计算配准结果，设计一个评价分数；
(5)使用一个蒙特卡洛树搜索mcts的变种用于遍历该搜索树，评价分数作为树中每个节点的奖励值，mcts的目标是在搜索树空间中找到最具有最高奖励的节点。
2.根据权利要求1所述的基于蒙特卡洛树搜索的血管3d/2d配准方法，其特征在于：所述步骤(1)中，血管图表示为
3.根据权利要求2所述的基于蒙特卡洛树搜索的血管3d/2d配准方法，其特征在于：所述步骤(1)中，根据公式(1),3d和2d血管结构的配准表示为找到3d血管的最优变换t，在该变换t作用后的3d血管的投影与2d血管达到最佳对齐，
其中
4.根据权利要求3所述的基于蒙特卡洛树搜索的血管3d/2d配准方法，其特征在于：所述步骤(2)中，两个边上点的稠密匹配通过均匀插值实现，将该过程描述为从πε到π的满射
5.根据权利要求4所述的基于蒙特卡洛树搜索的血管3d/2d配准方法，其特征在于：所述步骤(3)中，根据πε的定义以及血管的拓扑连续性，当给定初始的顶点匹配对，沿着血管图逐步地增加新的匹配边
6.根据权利要求5所述的基于蒙特卡洛树搜索的血管3d/2d配准方法，其特征在于：所述步骤(4)中，评价分数q为公式(2)
对于给定的πε，刚性变换经过两个满射得到t＝f1(f2(πε))，q的第一项与3d和2d点集的重叠度相关，3d点的投影
7.根据权利要求6所述的基于蒙特卡洛树搜索的血管3d/2d配准方法，其特征在于：所述步骤(5)包括以下分步骤：
(5.1)选择：在当前搜索树中选择最紧急的可拓展节点，从根节点开始自顶向下，在每一层通过贪婪策略实现；可拓展节点意味着根据当前节点的匹配状态找到至少一对新的匹配边；节点的紧急度通过公式(3)计算
其中qsim表示以当前结点作为根节点的子树的可能最大的奖励值，n表示一共迭代的次数，n表示当前结点的访问次数；因此，第二项用于惩罚当前结点过高访问次数，以鼓励对未访问或者较少访问分支的探索；γ是用来平衡对已知的节点的开发和对未访问结点的探索；
(5.2)拓展：对于选中的节点，对其添加至多nexp个子节点；当存在大于nexp的子节点，计算他们的奖励，并选择较高的nexp个子节点；
(5.3)模拟：将每个新扩展的子节点作为根节点，以深度优先的方式向下访问节点；在每一层中，随机选择子节点，直到到达最深处的叶子节点；这个过程重复nsim次，然后计算这些访问叶节点的奖励，将最高奖励作为扩展节点的qsim；
(5.4)反向传播：沿着被选中节点到根节点的路径，更新qsim值和节点访问次数n；
这四个步骤迭代地进行直至达到最大的计算量nmax或者理论上最高的奖励qmax，迭代结束后，具有最高的奖励q的节点从搜索中获得，且节点中对应的匹配πε和变换t是本方法的最终结果。
8.基于蒙特卡洛树搜索的血管3d/2d配准装置，其特征在于：其包括：
匹配结果表示模块，其利用3d和2d血管拓扑一致性来实现血管图匹配，将3d血管图和2d血管图的匹配表示为顶点对的集合和边对的集合；
稠密匹配模块，其将血管匹配过程看作求解3d和2d血管点的稠密匹配，先估计顶点之间的稀疏匹配，再匹配连接两个顶点之间的边，从而得到两个血管图的稠密匹配；
搜索树的构造模块，其根据血管边匹配的可以分解为连续状态的特性，构建一个搜索树；
计算结果及节点评分模块，其在搜索树的每个节点上，使用一个封闭解来计算配准结果，设计一个评价分数；
最优结果寻找模块，其使用一个蒙特卡洛树搜索mcts的变种用于遍历该搜索树，评价分数作为树中每个节点的奖励值，mcts的目标是在搜索树空间中找到最具有最高奖励的节点。
技术总结