锥束计算机断层扫描图像与X光图像配准方法

摘要

本发明公布了一种锥束计算机断层扫描图像与X光图像配准方法，利用基于混合残差卷积神经网络的回归模型建立二维X光图像与三维锥束CT图像非刚性变形参数之间的关联，基于混合残差卷积神经网络与变形参数迭代优化方法，实现可靠的在线二维三维图像配准；包括：提取得到图像通道；训练基于混合残差卷积神经网络的回归模型；进行基于回归的二维三维非刚性配准；迭代优化变形参数；由此得到由变形参数确定的最终的体图像，实现基于迭代回归的二维三维图像配准。本发明能够实现可靠的在线二维三维图像配准，可用于口腔临床基于二维与三维图像评价治疗及分析颅面生长。

说明书

技术领域

本发明涉及口腔临床医学和计算机视觉技术领域，具体涉及锥束计算机断层扫描图像与X光图像配准方法。

背景技术

口腔正畸治疗通常会持续数年，在治疗过程中采集的多张图像将用于评价颌面结构的刚性与非刚性形态差异，特别是对于未成年人患者采集的图像之间的差异还涉及到结构生长的非刚性形态变化。在临床口腔正畸中广泛使用锥束计算机断层扫描(锥束CT)图像之前，二维X光图像是唯一记录颌面形态变化的介质。非刚性的二维三维图像配准是从治疗中的X光图像中获取三维图像的关键步骤。传统的二维三维图像配准方法依赖迭代优化以最小化数字重建X光片(digitally-reconstructed-radiographs，DRR)与目标X光片之间的差异，并求解三维图像的空间刚性与非刚性变换。基于迭代优化的配准过程由于需要大量在线DRR计算而耗时。相对于刚性的图像配准，非刚性的二维三维图像配准由于涉及到大量参数的求解通常需要更大的时间代价。传统三维图像配准的测度被用于进行二维三维图像配准，例如离散直方图测度与变分互信息测度等。为了缓解测度计算中的时间代价，基于特征的二维三维图像配准使用基于诸如角点，线段等特征的测度。但是医学图像中准确检测结构的几何特征本身也是一个具有挑战性的任务，特别是考虑到颅面X光图像中由于结构重合造成图像模糊，不易于自动获取精确的结构几何特征。此外，基于特征的二维三维图像配准方法仍然使用迭代优化，在迭代优化的过程中会扩散特征检测误差。基于统计的表面网格与灰度模型可以降低二维三维图像配准中参数空间，但是对于复杂结构，仅仅调节统计模型中的有限参数难于获取统计模型与目标结构的紧密拟合。基于统计模型的二维三维图像配准仍然依赖于大量的在线DRR计算，也会进一步增加配准的时间代价。利用回归模型可以在二维图像与三维图像变形参数之间建立联系,不需要进行在线的DRR。基于部分最小二乘回归已被用于从二维X光图像重建三维图像，但是其中仍需要额外的三维代理模型用于模型学习。回归森林也被用于构造二维X光图像特征空间与三维体图像特征空间之间的映射，但是该映射依赖于二维与三维图像的特征提取。卷积神经网络在大量的计算机视觉与医学图像处理任务中具有很好的性能。基于卷积神经网络的方法用于进行二维三维图像配准，但是现有基于卷积神经网络回归模型的配准工作仅仅处理刚性二维三维图像配准中少量参数，并且其中的单次回归也难于处理颅面三维图像中的所涉及的非刚性复杂结构变化。

发明内容

为了克服上述现有技术的不足，本发明提供一种锥束计算机断层扫描图像与X光图像配准方法，基于卷积神经网络回归方法，进行二维X光图像与三维锥束CT图像之间的非刚性配准，可有效地改进二维三维图像配准，并实现可靠的在线二维三维图像配准，可用于口腔临床基于二维与三维图像评价治疗及分析颅面生长。

本发明的原理是：利用基于卷积神经网络的回归模型建立二维图像与三维图像非刚性变形参数之间的关联。在基于卷积神经网络的回归模型中，引入残差连接构造混合残差网络用于增强图像中的信息的传播以及回归模型的训练。该方法从一组锥束CT图像中提取锥束CT图像相对于参照锥束CT图像的非刚性变换参数，并对锥束CT图像的非刚性变形参数构造子空间，子空间中每一个点都对应三维锥束CT图像的非刚性变形参数，可用于对参照锥束CT图像进行非刚性的变形。在子空间中随机采样的变形点对，并构造该点对对应的三维锥束CT图像，通过DRR技术获取三维体图像的二维投影。以子空间中的变形点对的差异与对应的DRR图像对训练回归模型。不同与以往基于单步回归的配准，本发明引入回馈机制，通过迭代不断优化变形参数。在迭代过程中，基于卷积神经网络的回归模型不断更新三维图像变形参数，对参照图像施加获取的变形参数得到新的三维锥束CT图像，该图像的DRR投影如果与目标X光图像一致则迭代结束。否则，以当前新的三维锥束CT图像的DRR投影更新回归模型的输入图像点对，利用基于卷积神经网络的回归模型更新三维图像变形参数。该方法结合基于混合残差卷积神经网络的回归模型与变形参数迭代优化机制，利用混合残差卷积神经网络中的多尺度特征融合与迭代优化可以有效地改进二维三维图像配准，实现可靠的在线二维三维图像配准。

本发明提供的技术方案是：

一种锥束计算机断层扫描图像与X光图像配准方法，利用基于混合残差卷积神经网络的回归模型建立二维X光图像与三维锥束CT图像非刚性变形参数之间的关联，基于混合残差卷积神经网络与变形参数迭代优化方法，实现可靠的在线二维三维图像配准；包括如下步骤(流程参见图1)：

1)提取得到图像通道，作为基于混合残差卷积神经网络的回归模型的输入图像：

本发明具体采用各项异性扩散滤波通道及相关的梯度图像来增强X光图像，得到图像通道；

各项异性扩散滤波过程由如下热方程(式1)确定：

式1中，通量函数c控制扩散的速率，为图像I关于时刻t的偏导数，Δ为图像拉普拉斯算子，▽为梯度算子。热方程的离散数值解为式2：

其中，常数α与数值解的稳定有关，I^(t)与I^(t-1)分别对应t与t-1时刻的各项异性的扩散滤波图像，▽为梯度算子。基于梯度值定义系数β_r为式3：

β＝exp(-||▽I||²/κ)(式3)

参数r取值为1～4，从1变化到4对应梯度估计过程中的上下左右四个方向。各项异性的扩散滤波图像与梯度图像被用于回归模型的输入，其中t∈{0,…,8}。因而回归模型的输入图像(包括原始的X光图像、各项异性的扩散滤波图像及梯度图像；原始的X光图像指对X光图像不需要进行诸如角点、分割线段等图像特征提取，而将整张X光图像作为卷积神经网络回归模型的输入)有18个通道其中I^(t)表示t时刻各项异性扩散滤波通道，表示t时刻各项异性扩散滤波通道对应的梯度图像。

2)训练基于混合残差卷积神经网络的回归模型；

混合残差卷积神经网络的回归模型包含：卷积网络部分与全连接网络部分。其中，卷积网络部分是对图像进行自动的特征提取与编码。卷积网络部分的基本单位是卷积块。卷积块由常见的残差卷积网络构成，每个卷积块中包含两个卷积层，每个卷积层后都有一个修正线性单元(RELU)。感知域的大小是3，步长为1。卷积块中还包含一个2×2×2的池化层。在卷积块中的第一个卷积层的输入与第二个卷积层的输出之间添加残差连接以增强信息扩散。本发明中，卷积网络部分一共包含四个卷积块。全连接网络部分包含三个全连接层用于从图像特征中估计三维图像的变形参数。在本发明中引入长的残差连接，该连接被加在多尺度的卷积块网络与第一个全连接层之间，利用长残差连接可以有效辅助网络的训练，并融合多尺度的特征以估计精确的三维图像的变形参数。回归模型中对于来自卷积块1、2、3的特征图进行降采样，其降采样率分别为1/8、1/4、1/2，降采样后与卷积块4的特征图融合,作为全连接网络的输入。在图1中，长的残差连接用侧面箭头表示，其对应的特征图的降采样利用步长为8、4与2的卷积计算完成。与四个卷积块特征融合对应的全连接层后，网络中还有两个全连接层。

21)将基于混合残差卷积神经网络的回归模型损失函数L定义为预测的三维图像变形参数与真实的三维图像变形参数之间的欧式距离，如式4：

L＝||δT^p-δT^g||²(式4)

其中，δT^p为预测的三维图像变形参数，δT^g为真实的三维图像变形参数；

22)从一组锥束CT图像构造混合残差卷积神经网络的训练数据，训练基于混合残差卷积神经网络的回归模型；

221)从一组锥束CT图像中提取锥束CT图像相对于参照锥束CT图像的非刚性变换参数，并对锥束CT图像的非刚性变形参数构造子空间，子空间中每一个点都对应三维锥束CT图像的非刚性变形参数，可用于对参照锥束CT图像进行非刚性的变形；

222)在子空间中随机采样变形点对，并构造该点对对应的三维锥束CT图像，计算三维体图像的二维投影，得到对应的DRR图像对；

223)以子空间中的变形点对的差异与对应的DRR图像对训练基于混合残差卷积神经网络的回归模型；

3)基于回归的二维三维非刚性配准，利用基于混合残差卷积神经网络的回归模型建立二维X光图像与三维锥束CT图像非刚性变形参数之间的关联，获得参照图像V_ref对应的三维非刚性变形参数，使得变形之后得到的锥束CT图像的二维DRR投影与目标X光图像I_tar一致；

31)将目标X光图像记作I_tar，参照锥束CT图像记作V_ref；

32)将图像对(I_tar,I_ref)作为卷积神经网络的变形参数回归模型的输入，其中，I_ref表示参照锥束CT图像V_ref的DRR投影；

33)卷积神经网络的变形参数回归模型的输出即为图像对所对应的三维图像之间的非刚性变形参数；

具体地，本发明利用基于B样条非刚性变形模型对三维锥束CT图像进行非刚性的图像变形，其中在B样条控制网格节点上的偏移参数决定三维图像的非刚性变形；三维B样条的张量积被用于定义非刚性变形，其中在参照图像V_ref上定义控制网格G；二维三维图像配准的目标函数g定义为式5：

其中，P表示DRR投影函数，用于从变形后的三维锥束CT图像中生成二维X光图像；d表示基于互信息的二维图像测度，用于度量变形后的三维锥束CT图像的DRR投影与目标X光图像之间的差异；三维空间变换δT＝f(I_tar，I_ref)，其中f表示基于卷积神经网络的回归函数；M表示基于B-样条张量积的三维非刚性变形，用式6表示：

其中，B表示三维B-样条变形函数的基函数；

4)迭代优化变形参数；

41)首先，将目标X光图像与参照三维锥束CT图像的DRR投影图像组成的图像对作为回归模型的输入；

42)利用回归模型输出的变形参数，对参照三维锥束CT图像进行变形；

43)将变形后的三维锥束CT图像的DRR投影与目标X光图像进行比较，如果变形后的三维锥束CT图像的DRR投影与目标X光图像的图像均方差距离小于预先定义的阈值，则二者相似，即当前变形后的参照体图像与目标X光图像一致，返回当前变形后的三维锥束CT图像；否则，将变形后的三维锥束CT图像的DRR投影图像更新输入图像对中的参考图像，再次将更新图像对作为回归模型的输入，进一步估计该图像对所对应的三维图像变形参数；

44)建立基于迭代回归的二维三维图像配准的目标函数为式7：

其中，M表示基于B-样条张量积的三维非刚性变形，P表示DRR投影函数，I_tar为目标X光图像，为(t-1)步迭代中得到锥束CT图像。在第i次迭代中，控制网格G^(t)由之前迭代中得到的变形参数确定，如式8：

其中，G⁽⁰⁾为初始时刻的控制网格；t时刻的变形参数

45)经过t次迭代后得到的最终的体图像V^(t)由之前所有迭代得到的变形参数确定，表示为式9：

其中，M表示基于B-样条张量积的三维非刚性变形，为(t-1)步迭代中得到锥束CT图像，G^(t)为t时刻的控制网格，δT^(t)为t时刻的变形参数；

由此得到由变形参数确定的最终的体图像V^(t)，实现基于迭代回归的二维三维图像配准。

由于包含迭代优化，仍然有在线的DRR计算与二维的图像测度计算，但是，本发明中，为了实现二维三维图像配准仅仅需要有限次数的迭代，其效率远优于传统的利用迭代优化的方法。

与现有技术相比，本发明的有益效果是：

利用本发明提供的方法，可以有效地进行在线的二维三维图像配准，进行二维X光图像与三维锥束CT图像之间配准，从X光图像估计对应的三维锥束CT图像的非刚性变形参数，其中基于卷积神经网络的回归模型可以描述二维X光图像与对应的三维图像非刚性变形之间的非线性映射。迭代优化机制可有效获取与目标X光图像一致的三维锥束图像。本发明可用于口腔临床基于二维与三维图像评价治疗以及分析颅面生长。

在本发明具体实施中，计算了配准后的三维锥束CT图像的DRR投影图像与目标X光图像之间的差异，度量了前颅底、下颌外缘曲线距离，其误差都小于0.5毫米。因此，采用本发明基于迭代卷积神经网络回归模型的二维三维图像配准的精度，可满足临床口腔对二维三维图像配准的精度要求。

附图说明

图1为本发明方法的流程框图。

具体实施方式

下面结合附图，通过实施例进一步描述本发明，但不以任何方式限制本发明的范围。

本发明提供的基于卷积神经网络回归的方法进行二维X光图像与三维锥束CT图像之间的非刚性配准。利用基于卷积神经网络的回归模型建立二维X光图像与三维锥束CT图像非刚性变形参数之间的关联。结合基于混合残差卷积神经网络的回归模型与变形参数迭代优化机制，实现可靠的在线二维三维图像配准。

图1为本发明方法的流程框图。下面参照附图对本发明做进一步说明。

步骤一：提取图像通道；

不同于现有的传统基于特征的二维三维图像配准方法，本发明方法接收整张X光图像作为卷积神经网络回归模型的输入，不需要进行诸如角点、分割线段等特征提取。由于X光图像是单通道的灰度图像，本发明使用一系列各项异性扩散滤波通道以及相关的梯度图像增强输入X光图像。各项异性的扩散滤波过程由如下热方程(式1)确定：

式1中，通量函数c控制扩散的速率，为图像I关于时刻t的偏导数，Δ为图像拉普拉斯算子，▽为梯度算子。热方程的离散数值解为式2：

其中，常数α与数值解的稳定有关，I^(t)与I^(t-1)分别对应t与t-1时刻的各项异性的扩散滤波图像，▽为梯度算子。基于梯度值定义系数β_r为式3：

β＝exp(-||▽I||²/κ)(式3)

参数r取值为1～4，从1变化到4对应梯度估计过程中的上下左右四个方向，κ为常数。各项异性的扩散滤波图像与梯度图像被用于回归模型的输入，其中t∈{0,…,8}。因而回归模型的输入图像有18个通道其中I^(t)表示t时刻各项异性扩散滤波通道，表示t时刻各项异性扩散滤波通道对应的梯度图像。

步骤二：训练混合残差卷积神经网络；

本发明利用基于混合残差卷积神经网络的回归模型建立二维图像与三维图像非刚性变形参数之间的关联。基于卷积神经网络的回归模型包含两个部分：卷积网络部分与全连接网络部分。其中卷积网络部分是对图像进行自动的特征提取与编码。卷积网络部分的基本单位是卷积块。卷积块网络由常见的残差卷积网络构成，每个卷积块中包含两个卷积层，每个卷积层后都有一个修正线性单元(RELU)。所有卷积滤波器的感知域的大小都是3，步长为1。卷积块网络中还包含一个2×2×2的池化层。在卷积块中的第一个卷积层的输入与第二个卷积层的输出之间添加残差连接以增强信息扩散。本发明中，卷积网络部分一共包含四个卷积块网络。

全连接网络部分包含三个全连接层用于从图像特征中估计三维图像的变形参数。在本发明中引入长的残差连接，该连接被加在多尺度的卷积块网络与第一个全连接层之间，利用长残差连接可以有效辅助网络的训练，并融合多尺度的特征以估计精确的三维图像的变形参数。回归模型中对于来自卷积块1、2、3的特征图进行降采样，其降采样率分别为1/8、1/4、1/2，降采样后与卷积块4的特征图融合,作为全连接网络的输入。在图1中长的残差连接用侧面箭头表示，其对应的特征图的降采样利用步长为8、4与2的卷积计算完成。与四个卷积块特征融合对应的全连接层后，网络中还有两个全连接层。基于卷积神经网络的回归模型损失函数L定义为预测的三维图像变形参数与真实的三维图像变形参数之间的欧式距离：

L＝||δT^p-δT^g||²(式4)

其中，δT^p为预测的三维图像变形参数，δT^g为真实的三维图像变形参数。

本发明从一组锥束CT图像构造混合残差卷积神经网络的训练数据。从一组锥束CT图像中提取锥束CT图像相对于参照锥束CT图像的非刚性变换参数，并对锥束CT图像的非刚性变形参数构造子空间，子空间中每一个点都对应三维锥束CT图像的非刚性变形参数，可用于对参照锥束CT图像进行非刚性的变形。在子空间中随机采样的变形点对，并构造该点对对应的三维锥束CT图像，通过DRR技术获取三维体图像的二维投影。以子空间中的变形点对的差异与对应的DRR图像对训练基于混合残差卷积神经网络的回归模型。

步骤三：基于回归的二维三维配准；

利用基于混合残差卷积神经网络的回归模型建立二维X光图像与三维锥束CT图像非刚性变形参数之间的关联。将目标X光图像记作I_tar，参照锥束CT图像记作V_ref。二维三维图像非刚性配准的目的是推导参照图像V_ref对应的三维非刚性变形参数，使得变形之后得到的锥束CT图像的二维DRR投影与目标X光图像I_tar一致。基于卷积神经网络的变形参数回归模型的输入为图像对(I_tar,I_ref),其中I_ref表示参照锥束CT图像V_ref的DRR投影，输出为图像对所对应的三维图像之间的非刚性变形参数。本发明利用基于B样条非刚性变形模型，其中在B样条控制网格节点上的偏移参数决定三维图像的非刚性变形。三维B样条的张量积被用于定义非刚性变形，其中在参照图像V_ref上定义控制网格G。二维三维图像配准的目标函数g定义为式5：

其中，P表示DRR投影函数，用于从变形后的三维锥束CT图像中生成二维X光图像；d表示基于互信息的二维图像测度，用于度量变形后的三维锥束CT图像的DRR投影与目标X光图像之间的差异；三维空间变换δT＝f(I_tar，I_ref)，其中f表示基于卷积神经网络的回归函数；M表示基于B-样条张量积的三维非刚性变形，用式6表示：

其中，B表示三维B-样条变形函数的基函数。

步骤四：迭代优化变形参数。

不同与基于单步回归的二维三维图像配准，本发明引入迭代优化机制处理颌面部结构的精细形状。首先，目标X光图像与参照三维锥束CT图像的DRR投影图像组成的图像对作为回归模型的输入。利用回归模型输出的变形参数对参照三维锥束图像进行变形，将变形后的三维锥束图像的DRR投影与目标X光图像进行比较，如果二者相似，即当前变形后的参照体图像与目标X光图像一致，则返回当前变形后的三维锥束CT图像。否则，利用将变形后的三维锥束图像的DRR投影图像更新输入图像对中的参考图像。再次将更新图像对作为回归模型的输入，进一步估计该图像对所对应的三维图像变形参数。

基于迭代回归的二维三维图像配准的目标函数为式7：

其中，M表示基于B-样条张量积的三维非刚性变形，P表示DRR投影函数，I_tar为目标X光图像，为(t-1)步迭代中得到锥束CT图像。在第i次迭代中，控制网格G^(t)由之前迭代中得到的变形参数确定，如式8：

其中，G⁽⁰⁾为初始时刻的控制网格。t时刻的变形参数

经过t次迭代后得到的最终的体图像V^(t)由之前所有迭代得到的变形参数确定，表示为式9：

其中，M表示基于B-样条张量积的三维非刚性变形，为(t-1)步迭代中得到锥束CT图像，G^(t)为t时刻的控制网格，δT^(t)为t时刻的变形参数。由于包含迭代优化，仍然有在线的DRR计算与二维的图像测度计算，但是，本发明中，为了实现二维三维图像配准仅仅需要有限次数的迭代，其效率远优于传统的利用迭代优化的方法。

为了验证基于迭代卷积神经网络回归模型的二维三维图像配准的精度，本发明具体实施中，计算了配准后的三维锥束CT图像的DRR投影图像与目标X光图像之间的差异，度量了前颅底，下颌外缘曲线距离，其误差都小于0.5毫米，可满足临床口腔对二维三维图像配准的精度要求。

利用本发明的方法，可以有效地进行在线的二维三维图像配准，从X光图像估计对应的三维锥束CT图像的非刚性变形参数，其中基于卷积神经网络的回归模型可以描述二维X光图像与对应的三维图像非刚性变形之间的非线性映射。迭代优化机制可有效获取与目标X光图像一致的三维锥束图像。

需要注意的是，公布实施例的目的在于帮助进一步理解本发明，但是本领域的技术人员可以理解：在不脱离本发明及所附权利要求的精神和范围内，各种替换和修改都是可能的。因此，本发明不应局限于实施例所公开的内容，本发明要求保护的范围以权利要求书界定的范围为准。