一种图像相似块评估方法与流程

本发明设计计算机领域，特别涉及一种冷冻电镜图像相似块评估的方法。
背景技术：
：随着科学技术的迅速发展，人们对生物的认识，不再仅仅局限于宏观阶段，人们想要更加深刻了解认识生物的微观世界，根据结构决定功能，要了解生物大分子的功能和内在机制，如复制、繁殖，我们必须首先知道其结构。因为生物样品内在结构的复杂性和分子量巨大，解析其三维结构一直是影响是后续研究其功能的瓶颈，直到上世纪中后期，伴随着材料科学的飞速发展，越来越多的精密仪器设备逐渐出现，冷冻电子显微镜等精密仪器设备渐渐用于生物样品研究。经过几十年的发展，测定生物大分子结构的技术以及方法，都获得了很大的发展，主要有以下三种：X射线晶体学，NMR(核磁共振波谱学)，冷冻电子显微镜技术(简称冷冻电镜技术，Cryo-EM技术)，X射线晶体学要求得到待测生物样品的结晶，NMR需要对生物样品进行提纯，这两种方法都有局限性，因为生物样品内在结构的复杂性和分子量巨大，解析其三维结构一直是影响着其后续功能研究的瓶颈，直到上个世纪中后期，伴随着材料科学的飞速发展，越来越多的精密仪器设备逐渐出现，电子显微镜等渐渐用于生物样品的研究。冷冻电镜成像的一个缺点是，许多的生物样品都对辐射敏感，电子束可能会损伤样品，为了防止这种情况发生，通常将样品放在很低的电子束下成像，然而，这容易导致图像的信噪比非常低，图像的噪声增强，得到的图像常被泊松噪声和高斯噪声污染，除此之外，环境因素也会造成图像噪声，如磁场变化，机械振动，声振动，热不稳定性，电磁透镜等。另外，噪声还和检流器有关，近似表现为泊松分布，在将模拟投影图像转换成数字图像的时候也会带入高斯噪声，因此，需要针对噪声的不同表现形式去消除冷冻电镜图像中的噪声。另外，图像中的噪声也会影响冷冻电镜图像三维重构过程中的单颗粒挑选和校准步骤，因此，冷冻电镜图像去噪的主要目的是在减少噪声的同时，尽可能多保留图像的细节信息，提高图像质量。常用的电镜图像去噪方法，例如高斯滤波器技术，能够消除噪声，但是该方法也使得图像的边缘变得模糊，从整体上降低了图像质量，而且，图像边缘不够清晰，严重影响后续的颗粒挑选。要获得生物大分子的结构，我们首先需要得到其投影图像，但是，低温电镜技术拍摄的冷冻电镜图像有着弱相位、低信噪比、低对比度、背景强度不均衡、颗粒内部纹理不规则等特点，这给准确的测定生物大分子的立体结构带来很大的困难，为了解决这个问题，需要对投影图像进行去噪，提升图像的视觉效果，最大限度的恢复图像质量。在分析原始电子冷冻电镜图像的基础上，我们发现基于图像块的去噪方法能够有效消除原始图像中的噪声成分，但是，同时也发现，该方法首先需要找到与参考块相似的图像块，而如何评价两个图像块是否相似，成为我们解决这个问题的关键所在。在许多经典的基于图像块的图像处理算法中，经常使用欧拉距离来评价两个图像块是否相似，但是这种方式存在一定的局限性，这是因为，欧拉距离只考虑了图像的灰度值，同时，图像块子空间并不完全是欧拉空间，使用欧拉距离并不是一个很好的相似性衡量标准，因此，在选择相似块时，需要考虑图像子空间的结构和一些结构特征。因此，找到一种能够准确评估相似块的简单且高效的方法，成为当务之急。技术实现要素：图像块匹配(Block-Matching)算法是建立在图像信息的冗余性和相关性基础上，通过计算候选块集合X与参考块集合R之间的距离，找到候选块所属的参考块类，常用的块匹配算法有K近邻搜索法，该方法属于局部搜索，搜索速度相对较快，但是只能得到局部最优解，另一种匹配方法是全局搜索方法，虽然可以得到全局最优解，但是该方法耗时久，匹配算法的优劣已经严重影响着其他的后续处理。图像块子空间并不完全是欧拉空间，在判断相似块时，需要考虑图像子空间的结构和图像块的结构特征。本发明提出一种图像相似块评估方法，为一种基于测地距离的相似块匹配算法，使用测地距离代替欧拉距离来判断两个图像块是否相似。本发明提出的基于测地距离的相似块衡量方法，用于基于图像块的去噪算法中，能有效提高相似块的准确率，能提高去噪效果。一种图像相似块评估方法，包括以下步骤：(1)取冷冻电子显微镜图像中两个大小为K*K的图像块SA，SB；(2)计算图像块SA中每个像素点对应的梯度值tAi，，并保存在矩阵DSA中；计算图像块SB中每个像素点对应的梯度值tBi，并保存在矩阵DSB中；(3)计算第i个点的权重weight1，weight2：weight1＝0.5*(value[Ai]-value[Bi])2其中：value[Ai]、value[Bi]分别代表图像块SA和SB中对应相同位置的Ai、Bi两点的灰度值；weight2＝0.5*(tAi+tBi+tan|α-β|)tAi，tBi分别代表图像块SA和SB中第i个像素点的梯度值；α，β分别代表第i个像素点像素值变化最大方向和最小变化方向的夹角；(4)计算图像块中相同位置的Ai，Bi两点的测地距离：d(SAi,SBi)＝weight1+weight2；(5)计算SA，SB之间的测地距离：d(SA,SB)=1K*KΣi=1k*kd(SAi,SBi)]]>其中：i为图像块中的第i个像素点，i＝1,2……K*K；(6)比较d(SA,SB)与预定义阈值T的大小，判断图像块SA，SB是否相似。在本发明中，步骤(1)所述的两个图像块是通过在整个冷冻电镜图像中随机抽取的。在本发明中，所述K表示图像块的大小。大小为K的图像块中共有K*K个像素点。K为2-50，优选为5-20，更优选为6-10。在本发明中，步骤(2)中梯度值通过matlab中的gradient函数求得。在本发明中，步骤(3)中所述的α为[0,π]。在本发明中，β为[0,π]。在本发明中，步骤(4)中所述灰度值通过以下方法获得：冷冻电子显微镜图像本身就是灰度图像，通过分析头文件格式，用matlab读取后，存在矩阵中的数值就是图像的灰度值。在本发明中，步骤(5)中K的值与步骤(1)中K的值相同。在本发明中，步骤(6)中所述T的值与图像块大小K有关。在本发明中，T小于2K*K，T为2-70，优选为5-60，更优选为10-50。在本发明中，步骤(6)中：1)如果d(SA,SB)≤T，两个图像块相似；2)如果d(SA,SB)＞T，两个图像块不相似。在本发明中，在大小为N*N的冷冻电镜图像中随机抽取两个大小为K*K的图像块SA和SB，使用以上方法判断两个图像块是否相似。在本发明中，所述N为2m，m为4-50，优选为5-20，更优选为6-15。在本发明中，N*N表示冷冻电子显微镜图像大小。在本发明中，通过计算图像中随机抽取的两个图像块之间的测地距离来衡量他们之间的相似性，处理对象是图像块。在本发明中，测地距离计算方法除了考虑图像的灰度值外，还考虑了图像的梯度值对测地距离的影响，具体的计算公式是含有图像灰度和梯度的联合表达式。在本发明中，测地距离计算公式同时考虑像素值和灰度值对测地距离的影响，各占1/2的比例，同时，考虑灰度值最大变化跟最小变化方向，采用在[0,π]区间递增的正弦函数来描述两点梯度变化。在本发明中，计算的是图像块中的每个像素点到另一图像块相同位置像素点之间的测地距离，迭代多次，直到所有的像素点都参与计算，进行迭代但是不更新图像块。在本发明中，先计算图像块中K*K个点之间的测地距离，然后通过加权平均得到两个图像块之间的测地距离，通过测地距离来判断图像块是否相似，主要用于基于图像块的去噪算法设计中。在本发明中，块匹配是最简单和有效的寻找相似块的方法，但是该方法效率较低，在实际的图像处理应用中，通常使用一个比参考块稍微大一点的局部窗口来代替全局以寻找相似块，使用局部滑动窗口的方法的基本思想是，假设能在一个比较小的区域能找到参考块的许多相似块，但是，在实际中，图像的某些显著特征，如角，圆边，并不会在某一邻域出现，在非重复模式中使用局部相似性可能会出现较大的误差，在这种情况下，在整个图像域全局搜索相似块会更加合适。为了解决这两个问题，本发明提出了一种新的解决方法，该方法中使用测地距离代替欧拉距离来判断两个图像块是否相似，然后将所有的相似块聚集成相似块组，求相似块组的均值，然后把每一个图像块减去相似块组的均值，通过减法操作能够移除相似块的直流成分，但是并不会改变图像块的重要的结构特征，然后对所有减去直流成分的图像块进行先验学习。在本发明中，在现有的相似块匹配算法中，通常使用候选块Sx和参考块之间的欧拉距离来衡量其相似性。然而，当两点之间存在大量不在计算考虑区域的点时，用欧拉距离来计算两点之间的距离通常是无效的，因为它没有考虑到局部连通性，所以欧拉距离有着空间上的局限性，为了克服这种局限，在本发明中，我们将采用测地距离代替欧式距离，通过计算两个图像块之间的测地距离来衡量相似性。测地距离用于数据处理，通常表现在分类和相似性比较，本发明通过计算两个图像块之间的测地距离来衡量相似性。在图像域中能够用二维离散函数来表示图像，图像的梯度由二维离散函数得到，梯度方向就是灰度值的最大变化方向，因此，可以使用两点的梯度值来描述测地距离。α，β是两个梯度方向的夹角，取值范围均为[0,π]，灰度变化最大方向与灰度变化最小方向的夹角。在本发明中，小波变换虽然有效的利用了自然图像的稀疏属性，但是在小波分解时，依赖于所采用的小波基函数，缺乏平移不变性，在去除噪声的时候会丢失大量的细节信息，而现有的基于块的BM3D算法，使用硬阈值和维纳滤波，该方法依赖于阈值的选择，而且该方法中，使用欧式距离来衡量相似块，有一定的局限性，基于全局字典来对图像进行去噪，没有充分考虑到自然图像的图像块之间的非局部自相似先验知识。当先验学习的思想与图像的稀疏和冗余表示相结合，基于图像相似块先验学习的去噪方法，使用字典对图像进行去噪。字典去噪的理论依据是，理想图像在适当的过完备字典下存在稀疏表示，噪声破坏了这种稀疏表示，通过选择或设计适当的字典,求出自然图像在该字典下的稀疏表示就可以达到减弱或消除噪声的目的。使用字典去噪的一个先验知识就是该信号的稀疏性，基于图像块先验知识的去噪方法能够有效的保留图像的局部信息以达到较好的实验效果，通过学习得到的字典比使用固定稀疏基具有更好的去噪性能。MRC图像中存在大量的全同颗粒，对图像进行分块可以有效的利用这些全同颗粒的属性，能取得更好的实验效果。因此，本发明结合基于测地距离的相似块匹配方法和图像块的非局部自相似先验知识(NSS)，对NSS进行先验学习，构造相似块组的字典，通过求解加权稀疏编码模型的最优解，得到相似块组的稀疏编码，利用图像的稀疏表示对图像块进行去噪，最后重构所有的图像块得到去噪后的图像。该方法充分考虑到流形空间中的距离计算方法，使用测地距离准确的选择相似块，同时，兼顾相似块之间的内部先验知识(NSS)和外部先验知识(稀疏性)。MRC图像去噪的目的是消除图像中的噪声，提高图像的衬度和信噪比，为后续的单颗粒挑选和二维投影图像分类提供足够的有效信息。假设观察图像y，无噪图像x，噪声v，则有y＝x+v，v～N(0,σ2)，因此，图像去噪问题转化为通过观察图像y求得图像x的估计值使得最小，即均方误差MSE最小，从而能得到最大PSNR值，达到最优的去噪效果。非局部均值能够在消除噪声的同时有效的保留图像的细节信息，能用于多种不同的图像。为了有效利用图像的结构信息和非局部均值的能有效消除噪声的性质，本发明结合冷冻电镜图像的特点，对图像进行分组分块处理，从一幅MRC图像中抽取若干图像块，选取N个参考块，将所有的相似块聚集成相似块组，这里使用上文提到的测地距离来衡量相似块，然后将所有的相似块聚集成一个相似块组，共N个相似块组，每个相似块组包含M个相似块，用ym表示图像y中的图像块，用xm表示图像x中的图像块，因此，图像去噪问题可以转换成求最小MSE问题，即其中因为非局部均值能够抑制噪声，字典表示能够有效表示图像中的非噪声信号，因此，结合非局部均值和字典表示，可知因此通过求解字典D和稀疏编码系数达到对图像块进行去噪的效果，最后通过聚合所有图像块得到去噪之后的MRC图像。对图像进行分块分组处理，主要有以下优点：能够充分利用图像块之间的先验信息，对其进行建模，更适合并行运算，提高效率。图4是使用本发明提出的方法对MRC图像进行去噪的流程图。用表示观察图像y中的M个大小为p*p的相似块，其中用表示这M个图像块的均值，下一步的工作就是对进行先验学习得到K个高斯分布。在基于图像块的处理方法中，通常认为所有的图像块是独立采样。每个GMM模型由K个高斯分布组成，通过对所有的无噪图像块进行学习得到K个高斯成分。在本发明中，考虑到对图像块的均值操作能够抑制噪声，最大程度上保留图像块的结构特征，本发明采用对进行先验学习的方法得到K个高斯成分。因此，在本发明中的似然函数表示为这里的πk是权值因子，表示第k个高斯成分被选中的概率，在本发明中，N个相似块组相互独立，因此全局目标似然函数可以表示为为了下文计算方便，对目标似然函数取对数通常的方法是对目标似然函数求导，通过令其导数为0，求解方程。但是，lnL表达式中含有K个高斯成分的累加操作，对数函数不能进行化简，因此，不能通过求导直接得到最大值。本发明采用从GMM模型中随机选点，通过EM算法求解。通过GMM学习，可以得到能够描述图像块结构特征的K个高斯分布，每一个高斯分布称为一个高斯成分，下文通过贝叶斯方法，为每一个相似块组选择合适的高斯成分，并对相似块组去噪。对每个相似块组计算它由第k个高斯成分生成的概率。因为相似块组有着相同的高斯分布。假设第k个高斯成分能有效的描述根据贝叶斯概率公式有取对数得这里对不同的k取值，C都是相同的。对每一个k值，k＝1,2,......K，分别计算最大后验概率的对数值比较k取不同值时，的大小，选择能使取得最大值的k值，对应的高斯成分被选择用来对进行后续处理，此时Σk表示第k个高斯成分的协方差矩阵，对∑k进行奇异值分解(SVD)，即Σk＝DΛDT，其中D是正交特征向量矩阵，Λ是特征值对角矩阵，特征向量D表征非局部自相似统计结构，因此，可以用D作为稀疏编码的字典，用α表示稀疏编码系数，因此本发明中稀疏编码模型约束条件可以表示为用α表示稀疏编码的系数，v是噪声，w是α的权值向量，表示2范数，||wTα||1表示1范数，D是正交矩阵，DDT＝I，|D|＝±1，I是单位矩阵。根据最大后验概率MAP，可知根据贝叶斯公式，可得噪声v满足v～N(0,σ2)分布，因此，根据高斯分布概率密度函数有疏编码系数α满足拉普拉斯分布，因此，其中c是常数，所以可得其中ε是接近0的正数。因为D是正交矩阵，DDT＝I，|D|＝±1，因此取可得又因为所以令因此稀疏编码约束条件也可以表示为对αi进行求导，可得即通常写成下列形式其中(a)+＝max(a,0)，sgn(zi)是符号函数。定义函数SoftMAP(gi,|τi|＝sgn(gi)(|gi|-|τi|)+，因此所以图像块中的噪声和不准确的相似块集合会影响GMM学习，反过来影响字典和加权稀疏编码矩阵，本发明提出的测地距离能提高相似块组的准确性，从而提高字典的的正确率，最后使用字典D和加权稀疏编码的组合得到去噪后的图像块在本发明中，先使用GMM模型对进行先验学习，得到K个高斯成分，当第k个高斯成分被选择作为的最合适的高斯成分时，我们对其协方差Σk进行SVD分解，得到字典D，根据稀疏编码模型求解对应的加权稀疏编码，然后对图像块组分别进行去噪。使用本发明方法求得MRC图像中各图像块的估计值，最后重构所有的图像块得到去噪后的图像当图像的某一位置出现多个估计值时，通过加权平均计算最终估计。通过更新噪声方差进行多次迭代，η是常量，从而提高去噪效果。与现有技术相比较，本发明的技术方案具有以下有一技术效果：1、计算简单，快速。由于图像块是由若干像素点构成，通过计算点与点之间的测地距离进行累加评价得到图像块间的测地距离，在运算过程中，直接调用matlab中编译的梯度函数，节约了整体计算时间。2、准确性更高。传统的评估方法，如欧式距离，但是由于图像块子空间并不完全是欧拉空间，因此，使用欧拉距离并不是一个很好的相似性衡量标准，而基于流形空间的测地距离则能更好的描述子空间，使用测地距离作为相似性衡量标准，能够更准确的搜索到所有相似块。3、欧拉空间是流形空间的一个特例，如果是原始图像中的相邻点，就直接用欧式距离代替测地距离，否则用上文提到的距离公式计算测地距离。4、本发明的技术方案的去噪算法在本质上表现为能在最大程度上去除噪声，同时又能保持原始图像有效信息的完整性，并且具有相对较低的计算时间复杂度和空间复杂度。5、本发明的技术方案基于图像块的思想有效利用了冷冻电子显微镜图像中全同颗粒的属性，测地距离同时兼顾图像的灰度值和梯度值，表述更准确。附图说明图1为本发明的两个图像块相似性评估流程图。图2为本发明一个具体实施例中的相似块示意图。图3为本发明操作示意图。图4为使用本发明提出的方法对MRC图像进行去噪的流程图。具体实施方式为了使本领域的技术人员更好的理解本申请的技术方案，下面将结合本申请实施例中的附图，对本申请实施例中的技术方案进行清楚、完整的描述。应当理解，此处所描述的具体实施用例仅仅用以解释本发明，并不用于限定本发明。一种图像块相似性评估方法，包括以下步骤：步骤A，在一幅冷冻电镜图像中随机抽取两个大小相同的图像块。A1、分别用矩阵存储两个图像块的像素值。A2、分别用矩阵存储两个图像块的梯度值步骤B，计算测地距离：B1、计算权值weight1，weight2B2、计算图像块中对应点的测地距离；B3、计算图像块之间的测地距离步骤C，比较测地距离和固定的阈值T的大小：步骤D，判断两个图像块是否相似：更具体地：一种图像相似块评估方法，包括以下步骤：(1)取冷冻电子显微镜图像中两个大小为K*K的图像块SA，SB；(2)计算图像块SA中每个像素点对应的梯度值tAi，并保存在矩阵DSA中；计算图像块SB中每个像素点对应的梯度值tBi，并保存在矩阵DSB中；((3)计算第i个点的权重weight1，weight2：weight1＝0.5*(value[Ai]-value[Bi])2其中：value[Ai]、value[Bi]分别代表图像块SA和SB中对应相同位置的Ai、Bi两点的灰度值；weight2＝0.5*(tAi+tBi+tan|α-β|)，tAi，tBi分别代表图像块SA和SB中第i个像素点的梯度值；α，β分别代表第i个像素点像素值变化最大方向和最小变化方向的夹角；(4)计算图像块中相同位置的Ai，Bi两点的测地距离：d(SAi,SBi)＝weight1+weight2；(5)计算SA，SB之间的测地距离：d(SA,SB)=1K*KΣi=1k*kd(SAi,SBi)]]>其中：i为图像块中的第i个像素点，i＝1,2……K*K；(6)比较d(SA,SB)与预定义阈值T的大小，判断图像块SA，SB是否相似。在本发明中，步骤(1)所述的两个图像块是通过在整个冷冻电镜图像中随机抽取的。在本发明中，所述K表示图像块的大小。大小为K的图像块中共有K*K个像素点。K为2-50，优选为5-20，更优选为6-10。在本发明中，步骤(2)中梯度值通过matlab中的gradient函数求得。在本发明中，步骤(3)中所述的α为[0,π]。在本发明中，β为[0,π]。在本发明中，步骤(4)中所述灰度值通过以下方法获得：冷冻电子显微镜图像本身就是灰度图像，通过分析头文件格式，用matlab读取后，存在矩阵中的数值就是图像的灰度值。在本发明中，步骤(5)中K的值与步骤(1)中K的值相同。在本发明中，步骤(6)中所述T的值与图像块大小K有关。在本发明中，T小于2K*K，T为2-70，优选为5-60，更优选为10-50。在本发明中，步骤(6)中：1)如果d(SA,SB)≤T，两个图像块相似；2)如果d(SA,SB)＞T，两个图像块不相似。在本发明中，在大小为N*N的冷冻电镜图像中随机抽取两个大小为K*K的图像块SA和SB，使用以上方法判断两个图像块是否相似。在本发明中，所述N为2m，m为4-50，优选为5-20，更有选为6-15。在本发明中，N*N表示冷冻电子显微镜图像大小。实施例1从图3可知，假设两个夹角分别是α和β，α表示点Ai灰度值变化最大方向与灰度值变化最小方向的夹角，β表示点Bi灰度值变化最大方向与灰度值变化最小方向的夹角,定义两点之间的测地距离d(SAi,SBi)。定义权重weight1，weight2：weight1＝0.5*(value[Ai]-value[Bi])2其中：value[Ai]、value[Bi]分别代表图像块中相同位置的Ai、Bi两点的灰度值；weight2＝0.5*(tAi+tBi+tan|α-β|)tAi，tBi分别代表图像块SA和SB中第i个像素点的梯度值；α，β分别代表第i个像素点像素值变化最大方向和最小变化方向的夹角；计算图像块中相同位置的Ai，Bi两点的测地距离：d(SAi,SBi)＝weight1+weight2图像块是由若干像素点构成，给定两个图像块SA，SB，大小为p*p，图像块SA，SB之间的测地距离d(SA,SB)=1p*pΣi=1p*pd(SA1,SB1)---(3.5)]]>其中，i为图像块中的第i个像素点，i＝1,2.....p*p，比较d(SA,SB)与固定的阈值T的大小，如果d(SA,SB)＜T，认为图像块SA，SB是相似块，图1是使用本发明提出的测地距离来衡量相似块的流程图。以冷冻电镜图像N＝1280，K＝6，T＝50为例，举例说明本发明工作流程。如图1所示，计算图像块之间的测地距离的流程图：包括如下步骤：(1)在大小为N*N的冷冻电镜图像中随机抽取两个大小为K*K的图像块SA和SB。(2)计算SA和SB之间的测地距离：SA=1041731451601101412031251457615258235255178161811171251692201167390156152244233129144171162178215157151]]>SB=120169152154112139200128140851486524025017016087601201682211107084150158229234134140178170175200150157]]>调用matlab中的gradient函数分别求SA，SB的梯度值，分别存储在矩阵DSA，DSB中。DSA=69.020.5-6.5-17.5-9.531-78-29-24.53.5-9.0-9420.0-28.5-47-48.5-2236.044.047.5-26.5-73.5-1317-4.044.040.5-57.5-44.515.0-9.03.526.5-10.5-32-6.0]]>DSB=49.016.0-7.5-20.0-7.527.0-72.0-30.0-21.54.0-10.0-83.010.0-35.0-45.0-41.5-50.0-27.048.050.5-29.0-75.5-13.014.08.039.538.0-47.5-476.0-8.0-1.515.0-12.5-21.57.0]]>计算位于图像块中相同位置的两点的距离di，存储在矩阵R中R=163.50-52.2523.5-34.2565-219.5-303524.5-56.75-71.5152.50-182.5-203.75101396.522157-35-25.510.7572.5022.50-19.50-13259.25-5.7577-36.75434-115.5018.5058.2514-39.25177]]>计算得到SA，SB间测地距离d(SA,SB)＝32.07(3)比较d(SA,SB)和T的大小，其中T为50：d(SA,SB)＜T(4)判断是否相似：图像块SA和SB是相似块。当前第1页1 2 3