基于双边随机投影的光学相干层析成像去噪方法

摘要

基于双边随机投影的光学相干层析成像去噪方法，属于图像去噪处理技术领域。采用该发明的方法不需要大量训练数据，利用少量相邻帧即可在去除散斑的同时保留原有生物结构信息，运算效率较高。基于双边随机投影的去噪算法能有效的抑制散斑噪声，获得了较好的视觉效果，且计算复杂度较低。该方法的提出为医生进行视网膜疾病诊断以及改善计算机自动辅助视网膜疾病分析的性能提供了有效帮助，且计算时间可以满足临床的实时性要求。

说明书

技术领域

本发明涉及一种基于双边随机投影(BRP)的光学相干层析成像去噪算法，属于图像去噪处理技术领域。

背景技术

随着商业化频域OCT技术的问世，实现了三维高分辨率的成像，为青光眼、年龄相关性黄斑变性和糖尿病黄斑水肿等视网膜疾病的定量分析和研究提供了重要帮助。然而，由于光学相干层析成像(OCT)成像系统不可避免地引入散斑噪声，导致了图像的质量退化，严重影响了人们对图像信息的准确获取，给医生进行疾病诊断以及计算机辅助疾病诊断分析带来了极大挑战。如何抑制散斑噪声一直是国内外学者致力研究的重要课题。

目前针对散斑去噪主要有硬件和软件两种处理方法。其中，硬件方法主要依靠改进光学系统架构，该方法能够有效抑制散斑噪声，但也会使得OCT成像系统结构更加复杂，造价更高。与硬件方法相比，软件方法无需改变系统结构，更加灵活、易实现。常用的软件方法能在一定程度上有效去除散斑噪声，但减小噪声的同时也会破坏图像的细节和纹理信息，且部分方法需要设定不同的阈值收缩以达到最佳去噪效果。也有用深度学习和基于稳健性主成分分析算法(RPCA)等方法来解决散斑问题的研究，然而也存在需要大量图像样本进行训练、训练时间较长、噪声的适应性较差及分解效率低等问题，导致得到的图像去噪模型还有待于提升。

发明内容

本发明需解决的技术问题是，提供一种能有效地抑制散斑噪声且计算复杂度较低的基于双边随机投影的光学相干层析成像去噪算法。

为了解决上述问题，本发明提出的基于双边随机投影的光学相干层析成像去噪方法，包括以下步骤：

S01：图像预处理，对视网膜图像进行对齐预处理，增强相邻帧之间的相关性；

其中，对于每个OCT图像，在慢扫描方向，即Y方向(y-z视图)，相邻的两个B扫描图像之间存在很大失真，B扫描图像对齐步骤如下：首先分割出OCT图像中视网膜的上边界，采用多分辨率表面检测方法进行分割。分割后计算所有B扫描图像中视网膜的上边界的平均位置，并以此为基准，估计各B扫描图像的位移值。在计算上边界的平均位置时只使用除中心凹外左侧部分和右侧部分的像素点。各B扫描图像根据视网膜的上边界的平均位置向上或向下进行平移，从而对齐各B扫描图像。

S02:将一帧包含n个像素的B扫描图像看作是矩阵中的列向量，相邻的组织结构相似的m帧B扫描图像则可用矩阵表示为

X＝L+S+N (1)

其中，L表示无噪生物结构图像组成的低秩矩阵，S表示由运动伪差引起的稀疏矩阵,N表示噪声矩阵。要从扫描的OCT图像中提取无噪图像即等同于从观测信号X中恢复信号L。S03:设低秩矩阵L有rank(L)≤r，且S矩阵中的非零元素个数满足card(S)≤k，则去噪问题可表示为最小化代价函数求解问题：

其中，r表示低秩矩阵的秩，k表示矩阵稀疏度参数，k越小，矩阵越稀疏；S04:由于公式(2)中矩阵L和S均为未知矩阵，给定r和k，公式(2)的优化问题可转化为交替求解公式(3)及公式(4)直至收敛。

其中，L

S05:根据双边随机投影算法，对于矩阵

Y

其中，Y

S06:X的r秩快速低秩逼近可由公式(6)计算。分别对Y

Y

其中，Q

对于公式(4)中S

其中，

S07:经过多次迭代，直到

本发明与现有技术相比，具有如下优点和有益效果：采用该发明的方法不需要大量训练数据，利用少量相邻帧即可在去除散斑的同时保留原有生物结构信息，运算效率较高。基于双边随机投影的去噪算法能有效的抑制散斑噪声，获得了较好的视觉效果，且计算复杂度较低。该方法的提出为医生进行视网膜疾病诊断以及改善计算机自动辅助视网膜疾病分析的性能提供了有效帮助，且计算时间可以满足临床的实时性要求。

附图说明

图1是OCT图像预处理前后的三维渲染图；其中，(a)为预处理前上表面三维渲染图，(b)为预处理后上表面三维渲染图。

图2是原始图像与本发明提供的去噪方法对图像处理后的对比图；其中，(a)为原始图像，(b)为去噪图像。

图3是一维原始图像深度信息信号图与去噪图像深度信息信号图的对比图；其中，(a)为一维原始图像深度信息信号图，(b)为去噪图像深度信息信号图。

图4是去噪前的原始图像。

图5是RPCA算法去噪结果图。

图6是本专利算法去噪结果图。

具体实施方式

以下将结合附图对本发明作进一步详细说明。

S01：图像预处理，对视网膜图像进行对齐预处理如图1所示，对于每个OCT图像，在慢扫描方向，即Y方向(y-z视图)，相邻的两个B扫描图像之间存在很大失真，B扫描图像对齐步骤如下：首先分割出OCT图像中视网膜的上边界，采用多分辨率表面检测方法进行分割。分割后计算所有B扫描图像中视网膜的上边界的平均位置，并以此为基准，估计各B扫描图像的位移值。在计算上边界的平均位置时只使用除中心凹外左侧部分和右侧部分的像素点。各B扫描图像根据视网膜的上边界的平均位置向上或向下进行平移，从而对齐各B扫描图像。图1显示了视网膜上表面在校正前后的效果。如结果图所示，对齐之前相邻两个B扫描图之间有很多跳变，视网膜上表面高低起伏。对齐之后相邻两个B扫描图之间平滑连续，视网膜上表面也趋于平滑。对齐之后大大减少了运动伪差的影响，视网膜各层趋于平滑，视网膜成像的形状更接近真实的视网膜形状。

S02:将一帧包含n个像素的B扫描图像看作是矩阵中的列向量，相邻的组织结构相似的m帧B扫描图像则可用矩阵表示为

X＝L+S+N (1)

其中，L表示无噪生物结构图像组成的低秩矩阵，S表示由运动伪差引起的稀疏矩阵,N表示噪声矩阵。要从扫描的OCT图像中提取无噪图像即等同于从观测信号X中恢复信号L。

S03:设低秩矩阵L有rank(L)≤r，且S矩阵中的非零元素个数满足card(S)≤k，则去噪问题可表示为最小化代价函数求解问题：

S04:由于公式(2)中矩阵L和S均为未知矩阵，给定r和k，公式(2)的优化问题可转化为交替求解公式(3)及公式(4)直至收敛。

S05:根据双边随机投影算法，对于矩阵

Y

其中，

S06:X的r秩快速低秩逼近可由公式(6)计算。分别对Y

Y

对于公式(4)中S

其中，

S07:经过多次迭代，直到

将双边随机投影散斑去噪算法的求解步骤总结如下：

步骤1，初始化高斯随机矩阵A

步骤2，利用公式(5)计算右随机投影矩阵Y

步骤3，更新A

步骤4，更新A

步骤5，根据公式(7)对Y

步骤6，根据公式(9)求解S

经过多次迭代，直到

在原始图像和去噪图像中选取对应位置的A扫描图像(图2中竖线位置)分别绘制一维深度信息信号图(图3)。随着深度逐渐增大，视网膜从上而下对应的结构组织分别为神经纤维层、神经节细胞层、内从状层、内核层、外从状层、外核层(outer nuclear layer，ONL)、内节/外节、视网膜色素上皮层以及脉络膜等。对比图3中去噪处理前后的深度信息信号图可见，去噪处理后信号曲线的光滑度以及视网膜生物结构信息的对比度都有明显的提升。

为了更好地评价本发明方法的客观性能，将原始图像与本文提出的双边随机投影去噪方法以及RPCA算法进行对比分析。图5和图6为不同去噪方法处理所得的去噪效果对比图。原始图像以及各算法处理得到的去噪图像的量化分析结果如表1所示。

由表1可知，采用双边随机投影去噪方法处理的图像整体更加光滑，对比度更高，图像质量更好。与RPCA算法相比，改进后的BRP方法在信噪比、对比度信噪比以及等效视数指标上分别提高了1.22dB、0.84dB和59.5，去噪效果更好。除此之外，对各方法的计算复杂度也进行了评估。在Intel(R)Core(TM)2Duo CPU处理器和8GB RAM的PC机，MATLAB 2017a软件上进行了测试，用RPCA算法和BRP算法对实验数据集中单幅B扫描图像的处理时间分别为4s和0.3s。BRP算法大大降低了计算复杂度。

以上，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何不经过创造性劳动想到的变化或替换，都应该属于本发明的保护范围之内。

表1