[发明专利]基于XCT图像非刚度配准的目标鼠解剖结构图谱获取方法

摘要

本发明公开了一种基于XCT图像非刚度配准的目标鼠解剖结构图谱获取方法，其基本步骤是首先，将Digimouse模型设定为参考鼠解剖结构图谱，同时将Digimouse对应的XCT图像设定为参考图像；其次，对目标鼠进行XCT成像获取目标图像并进行预处理；然后，利用非刚度图像配准技术构建参考图像至目标图像的配准映射矩阵；最终，将配准映射矩阵作用于参考鼠解剖结构图谱，构建出目标鼠解剖结构图谱，从而实现对目标鼠各个组织器官的标定。本发明所使用配准方法精度高，能够通过图像配准的方法简单易行的实现小鼠组织器官在小鼠XCT图像上的标识；本发明提出方法同样适用于其他医学应用领域如人脑结构学研究，即可以利用标准人脑解剖学图像实现目标人脑的解剖结构的获取。

主权项

一种基于XCT图像非刚度配准的目标鼠解剖结构图谱获取方法，其特征在于，其基本步骤是：首先，将Digimouse模型设定为参考鼠解剖结构图谱，同时将Digimouse对应的XCT图像设定为参考图像；其次，对目标鼠进行XCT成像获取目标图像并进行预处理；然后，利用非刚度图像配准技术构建参考图像至目标图像的配准映射矩阵；最终，将配准映射矩阵作用于参考鼠解剖结构图谱，构建出目标鼠解剖结构图谱；具体步骤如下：步骤一、设定参考鼠解剖结构图谱和参考图像：将Digimouse结构设定为参考鼠解剖结构图谱Ar，将所述Digimouse的XCT图像设定为参考图像Ir；步骤二、目标鼠XCT图像的获取与预处理：利用XCT设备对目标鼠进行全身成像，获取目标鼠的XCT图像；对目标鼠XCT图像进行仿射变换，将目标鼠XCT图像中小鼠的头部与背部的朝向变换至与参考图像Ir相同；将变换后的XCT图像作为目标图像It；步骤三、构建初步配准映射矩阵Mc和初步配准图像Ic：利用图像分割算法分别分割出参考图像Ir中的小鼠骨骼区域与表皮区域，及目标图像It中的小鼠骨骼区域与表皮区域；利用边缘检测算法分别提取上述参考图像Ir中小鼠骨骼区域和表皮区域及目标图像It中的小鼠骨骼区域和表皮区域共计四个区域的边界轮廓，对所述四个区域的边界轮廓进行等概率采样，从而计算出参考鼠骨骼特征点集Lrb、参考鼠表皮特征点集Lrs、目标鼠骨骼特征点集Ltb以及目标鼠表皮特征点集Lts；利用TPS‑RPM算法的开源程序，计算参考鼠骨骼特征点集Lrb与目标鼠骨骼特征点集Ltb之间的骨骼特征点对应矩阵Cb；在计算时，将目标鼠骨骼特征点集Ltb设定为目标点集，将参考鼠骨骼特征点集Lrb设定为待配准点集，并设定TPS‑RPM所用模拟退火算法中的初始温度系数及迭代终止条件，通过迭代计算出骨骼特征点对应矩阵Cb；利用TPS‑RPM算法的开源程序，计算参考鼠表皮特征点集Lrs与目标鼠表皮特征点集Lts之间的表皮特征点对应矩阵Cs；在计算时，将目标鼠表皮特征点集Lts设定为目标点集，将参考鼠表皮特征点集Lrs设定为待配准点集，并设定TPS‑RPM所用模拟退火算法中的初始温度系数及迭代终止条件，通过迭代计算出表皮特征点对应矩阵Cs；利用TPS‑RPM算法的开源程序求取骨骼特征点对应矩阵Cb与表皮特征点对应矩阵Cs的具体过程如下：所述参考鼠骨骼特征点集Lrb与目标鼠骨骼特征点集Ltb之间的骨骼特征点对应矩阵Cb及参考鼠表皮特征点集Lrs与目标鼠表皮特征点集Lts之间的表皮特征点对应矩阵Cs均为模糊对应矩阵，在模糊对应矩阵中各个元素均为[0,1]之间的浮点数，用以描述两点之间对应程度的强弱；所述TPS‑RPM算法在计算参考鼠骨骼特征点集Lrb与目标鼠骨骼特征点集Ltb之间的骨骼特征点对应矩阵Cb时，利用模拟退火算法计算如下能量方程的极小值点：E(Cb,f)=Σi=1NrbΣj=1Ntbcb,ij||lrb,i-f(ltb,j)||2+λ||Lf||2+TΣi=1NrbΣj=1Ntbcb,ijlogcb,ij-γΣi=1NrbΣj=1Ntbcb,ij---(1)]]>式(1)中，cb,ij为对应矩阵Cb中的元素；Nrb为参考鼠骨骼特征点集Lrb所包含特征点的个数，Ntb为目标鼠骨骼特征点集Ltb所包含点特征的个数；lrb,i为参考鼠骨骼特征点集Lrb中的第i个特征点的坐标，ltb,j为目标鼠骨骼特征点集Ltb中的第j个特征点的坐标；f为薄板样条函数；||Lf||2为对f的平滑性约束，其中L为微分算子；T为模拟退火过程中用于控制模糊对应程度的温度系数；λ为先验平滑度权值；γ为鲁棒性控制权值；随着温度系数T逐渐减小，对应矩阵Cb与薄板样条函数f交替迭代；对应矩阵Cb的最终迭代结果即为所求骨骼特征点对应矩阵Cb；所述TPS‑RPM算法在计算参考鼠表皮特征点集Lrs与目标鼠表皮特征点集Lts之间的表皮特征点对应矩阵Cs时，利用模拟退火算法计算如下能量方程的极小值点：E(Cs,f)=Σi=1NrsΣj=1Ntscs,ij||lrs,i-f(lts,j)||2+λ||Lf||2+TΣi=1NrsΣj=1Ntscs,ijlogcs,ij-γΣi=1NrsΣj=1Ntscs,ij---(2)]]>式(2)中，cs,ij为对应矩阵Cs中的元素；Nrs为参考鼠表皮特征点集Lrs所包含特征点的个数，Nts为目标鼠表皮特征点集Lts所包含点特征的个数；lrs,i为参考鼠表皮特征点集Lrs中的第i个点的坐标，lts,j为目标鼠表皮特征点集Lts中的第j个点的坐标；其余变量含义与式(1)中相同；随着温度系数T逐渐减小，对应矩阵Cs与薄板样条函数f交替迭代；对应矩阵Cs的最终迭代结果即为所求表皮特征点对应矩阵Cs；分别将求得的骨骼特征点对应矩阵Cb及表皮特征点对应矩阵Cs作用于参考鼠骨骼特征点集Lrb及参考鼠表皮特征点集Lrs，得到初步配准后的骨骼特征点集Lcb及初步配准后的表皮特征点集Lcs；利用上述求得的初步特征点配准的结果，构建初步配准局部位移矩阵P：P=LrbLrb-LcbLrsLrs-Lcs=[x,y,z,Δx,Δy,Δz]---(3)]]>将上述矩阵P转化为三组四维数据点集Px＝{(x,y,z,Δx)}，Py＝{(x,y,z,Δy)}，Pz＝{(x,y,z,Δz)}，将Δx、Δy及Δz分别视为点(x,y,z)的函数值，即Δx＝G1(x,y,z)，Δy＝G2(x,y,z)，Δz＝G3(x,y,z)；采用多层次B样条拟合算法分别对所述三组四维数据点集Px＝{(x,y,z,Δx)}，Py＝{(x,y,z,Δy)}，Pz＝{(x,y,z,Δz)}进行拟合，用于计算所述参考图像Ir中各个像素沿着x、y、z三个方向上的位移量；利用所述多层次B样条拟合方法拟合数据点集Px＝{(x,y,z,Δx)}的过程详细描述如下：所述多层次是指，利用覆于参考图像Ir上的一组逐渐加密的立方控制网格Φ0,Φ1,...,Φk,...,Φh依次对迭代更新的四维数据点集进行B样条拟合，并将所求h层拟合函数之和作为最终多层次B样条拟合函数；其中，Δ0ξ＝Δx，Δk+1ξ＝Δkξ‑gk(x,y,z)，gk(x,y,z)为第k层B样条拟合结果；所述第k层B样条拟合过程描述如下：假设第k层控制网格Φk尺寸为Kx×Ky×Kz，则第k层B样条拟合函数如下式所示：gk(x,y,z)=Σl=03Σm=03Σn=03Bl(dx)Bm(dy)Bn(dz)φk,(l+i,m+j,n+k)---(4)]]>式(4)中，φk,(l+i,m+j,n+k)为位于控制网格Φk中坐标为(l+i,m+j,n+k)的控制节点值；l,m,n∈{0,1,2,3}；Bl、Bm及Bn分别为第l,m,n阶B样条基函数，其中0至3阶B样条函数的表达式描述如下:B0(δ)=(1-δ)3/6B1(δ)=(3δ3-6δ2+4)/6B2(δ)=(-3δ3+3δ2+3δ+1)/6B3(δ)=δ3/6---(5)]]>所述式(4)中，控制网格Φk中各个控制节点的取值通过以下两步计算：(a)计算中每个数据点对控制网格Φk中各个控制节点取值的影响量：以中的一个数据点p＝(xp,yp,zp,Δkξp)为例进行说明如下：数据点p＝(xp,yp,zp,Δkξp)对控制网格Φk中各个控制节点的影响矩阵表现为一个尺寸为Kx×Ky×Kz的矩阵Ψp；为计算简便，定义与矩阵Ψp尺寸相同的两个矩阵Γp与Ωp；所述矩阵Ψp、Γp与Ωp中坐标为(l+i,m+j,n+k)的元素分别通过式(6)计算：ψp,(l+i,m+j,n+k)=γp,(1+i,m+j,n+k)Δkξpωp,(l+i,m+j,n+k)γp,(l+i,m+j,n+k)=Bl(dxp)Bm(dyp)Bn(dzp)ωp,(l+i,m+j,n+k)=Σl=03Σm=03Σn=03[Bl(dxp)Bm(dyp)Bn(dzp)]2---(6)]]>式(6)中，l,m,n∈{0,1,2,3}；在矩阵Ψp、Γp与Ωp中，除所述坐标为(l+i,m+j,n+k)共计64个元素以外其余位置，ψp、γp与ωp均为0；(b)求取控制网格Φk中每个控制节点的值综合中每个数据点对控制网格Φk中各个控制节点取值的影响，求取制网格Φk中每个控制节点的值；控制网格Φk中坐标为(a,b,c)的控制节点φk,(a,b,c)取值为：式(7)中，γp,(a,b,c)、ωp,(a,b,c)、ψp,(a,b,c)分别为Ψp、Γp与Ωp中坐标为(a,b,c)的元素值；至此，第k层B样条拟合函数gk(x,y,z)确立；综合各层次拟合函数，计算出多层次B样条拟合函数利用所述求取的多层次B样条拟合函数g(x,y,z)，计算出参考图像Ir中各个像素沿x轴的位移；同理，利用多层次B样条拟合四维数据点集Py＝{(x,y,z,Δy)}和Pz＝{(x,y,z,Δz)}，从而计算出参考图像Ir中各个像素沿y、z轴向中的位移量，由此构建参考图像Ir的初步配准映射矩阵Mc；利用所述初步配准映射矩阵Mc，反向求解出初步配准图像Ic；在构建初步配准图像时，灰度的赋值采用三线性插值方法；步骤四、构建精细配准映射矩阵Mf：利用图像分割算法，提取出初步配准图像Ic中的小鼠表皮区域及骨骼区域，并利用边缘检测算法分别提取出所述小鼠表皮区域及骨骼区域的轮廓；将小鼠表皮及骨骼的轮廓相互叠加，并使用矩形网格进行采样，由此得到一组初步配准图像特征点集L'c；利用块匹配方法求取初步配准图像特征点集L'c中每个特征点在目标图像上的对应位置，以L'c中的任意一点pl＝(xp,yp,zp)为例，将所述过程描述如下：在初步配准图像Ic中以坐标(xp,yp,zp)为中心选取尺寸为N1×N1×N1的立方邻域T，在目标图像It中以坐标(xp,yp,zp)为中心选取N2×N2×N2的立方邻域S，其中N2>N1；将T作为模板，S作为搜索域，在S区域中搜索与T具有最大相似度的子区域s1，并将s1的中心点pl′作为点pl的对应位置；以此类推，依次寻找出初步配准图像特征点集L'c中各点在目标图像It上的对应位置，由此构建出精细配准特征点集Lf；利用初步配准图像特征点集L'c及精细配准特征点集Lf构建精细配准局部位移矩阵Q：Q＝[L'c,L'c‑Lf]＝[x,y,z,Δx,Δy,Δz]  (8)将上述矩阵Q转化为三组四维数据点集Qx＝{(x,y,z,Δx)}，Qy＝{(x,y,z,Δy)}，Qz＝{(x,y,z,Δz)}；与步骤三中多层次B样条拟合过程一致，分别对所述三组四维数据点集Qx＝{(x,y,z,Δx)}，Qy＝{(x,y,z,Δy)}进行多层次B样条拟合，从而分别计算出初步配准图像Ic中各个像素沿x、y、z三个轴向的位移量，由此构建初步配准图像Ic的精细配准映射矩阵Mf；步骤五、构建目标鼠解剖结构图谱：将所述参考鼠解剖结构图谱Ar投影至像素坐标系，得到像素坐标系下的参考鼠解剖结构图谱按顺序依次将初步配准映射矩阵Mc及精细配准映射矩阵Mf作用于所述像素坐标系下的参考鼠解剖结构图谱使所述像素坐标系下的参考鼠解剖结构图谱产生与步骤三初步配准及步骤四精细配准过程相同的变形，获得像素坐标系下的配准鼠解剖结构图谱将所述像素坐标系下的配准鼠解剖结构图谱投影至物理坐标系下，得到物理坐标系下的配准鼠解剖结构图谱Af，所述物理坐标系下的配准鼠解剖结构图谱Af即为目标鼠解剖结构图谱。