基于层次式优化的多幅眼底图像融合方法

魏丽芳 余 轮

1(福建农林大学计算机与信息学院，福州 350002)2(福州大学物理与信息工程学院,福州 350002)

基于层次式优化的多幅眼底图像融合方法

魏丽芳1余 轮2*

1(福建农林大学计算机与信息学院，福州 350002)2(福州大学物理与信息工程学院,福州 350002)

针对眼底图像配准后直接叠加产生的接缝及如何保证多幅图像融合后细节信息不丢失的问题，提出层次式优化的多频带眼底图像融合方法。该方法通过多阈值分割获取掩模图像并计算其欧式距离得到各层次图像蒙版；根据欧氏距离值及拉普拉斯能量和设计每层蒙版图像的改进加权系数。构建基于信息熵、空间频率和清晰度的图像融合联合客观评价方法进行分析。最后，利用配准误差及重叠率对图像序列进行分组融合后进入下一层，对于奇数幅图像序列中没有组合的图像直接放入下一层。再根据分组规则重新分组进行优化融合，实现层次式的多幅眼底图像优化融合。通过对75 组取自福建省附属第一医院眼科及眼底相机实验系统的图像序列，涉及4 898 组图像对(正常眼底图像2 952 对，病变眼底图像1 946 对)进行测试与验证，结果表明，提出方法在有效除去拼接缝的同时，融合图像在重叠区域的RMSE值约为(0.1±0.05)像素。提出的融合方法在客观评价和主观视觉效果之间取得了较好的平衡。

眼底图像； 图像融合； 蒙版图像； 加权系数； 层次式融合

引言

图像融合是将两幅图像中的信息综合到一幅图像中并以可视化方法显示的技术[1-3]。图像融合的前提是两幅图像已经实现较好的配准。图像融合效果的好坏将直接影响拼接图像的视觉效果。对于待融合的两幅图像，若在重叠部分，简单的取第一幅图像和第二幅图像进行叠加会产生图像模糊和拼接痕迹，无法满足人的视觉要求。图像融合的主要思想是实现图像在衔接处平滑过渡。

在糖尿病早期诊断和治疗过程中，利用眼底图像配准、融合技术等可从多个方面辅助眼科医生进行病变检测和治疗跟踪[4-8]。以消除拼接缝为目的的眼底图像融合，主要考虑在实现图像间平滑过渡的同时保证图像信息不损失。而多幅眼底图像融合，需要考虑目标图像生成过程中产生的误差累积及目标图像越来越大而导致计算量越来越大等问题。文献[2]以充分保留原图像信息考虑，提出基于金字塔分解的多层次局部极值的医学图像融合方法，该方法在保持原图细节信息同时有效实现了多模态医学图像融合，但该方法对噪声影响严重。Cattin等采用多频带方法对眼底图像融合，通过距离值蒙版设计权值函数，实现任意重叠图像之间的多频带融合，完成眼底图像的无缝拼接，证明了多频带融合方法在眼底图像融合中的适用性[4]；Poletti 等在重叠区域设计基于距离幂次的加权系数，实现多幅眼底图像的融合，不必进行图像分解而有效去除灰度变化不大图像之间形成的拼接缝[6]；但该方法中幂次值需依据经验值获取，且对于灰度变换较大的眼底图像边缘拼接处仍可见接缝。

针对以上分析和存在的问题，提出基于加权蒙版的层次式优化多频带眼底图像融合方法。该方法在所有图像完成良好全局配准基础上[7]，获取融合蒙版及蒙版加权系数，实现基于配准误差和重叠率的层次式多幅眼底图像的多频带融合，在保证图像质量的同时降低运算量。

1 材料与方法

所提出的方法在图像两两配准的基础上，通过距离图像获取蒙版图像，并通过距离图像在重叠区域内的距离值及拉普拉斯能量和来确定蒙版图像加权系数，根据最小误差原则设计层次式多频带融合的多图像融合方法。实现流程如图1所示，首先要实现两两图像的良好融合，进而构建层次式的多幅图像融合方法。

图1 层次式优化融合流程图Fig.1 The flow chart of hierarchical optimization fusion

1.1掩模图像与距离图像

提取掩模图像实际上是将图像分割为背景和前景，为准确提取掩模图像，对眼底图像的B通道图像进行多阈值分割。设眼底图像直方图的NT个波峰为Peak(1)，Peak(2)，…,Peak(NT)，其中，

Peak(1)

(1)

由于眼底图像背景为较暗的灰度区域，像素值近似为0，因此选择第一个波谷作为区分背景和前景的阈值。设T_MIN为Peak(1)和Peak(2)峰值间的最小波谷值，Mf为提取的背景掩模图像，即目标区域，R为原图像

(2)

对于掩模图像Mf目标区域内的点(即值为1的点)计算其到目标区域边界的最短距离(采用欧氏距离)的正值d(i,j)，对目标区域以外的点及边界轮廓上的点到边界轮廓的距离值设为0[8-10]，表示为

(3)

(4)

式中,D1和D2分别为两幅原图像的距离变换矩阵，若(i,j)属于目标区域，而在重叠区域内

(5)

蒙版图像获取过程见图2。通过精确配准获得两幅图像之间的坐标变换关系，即空间位置关系；利用图像分割技术获得掩模图像；计算掩模的图像 距离矩阵；由距离矩阵计算得到蒙版图像。

图2 蒙版图像获取过程Fig.2 The process of mask image acquisition

1.2蒙版加权系数的设计

在各金字塔层的高斯金字塔图像是原图分辨率变化，故仍保留着原图像的结构信息，存在幅值的变化，因此利用EOG(the energy of image gradient)等清晰度指标来选择图像的变换域系数是合理的。而清晰度指标中拉普拉斯能量和具有最好的区分能力[12-15]。而基于距离值和拉普拉斯能量和的加权系数可以使图像之间的灰度平稳过渡。

因此，通过距离图像在重叠区域内的距离值和拉普拉斯能量和相结合来确定加权系数。

在重叠区域内，某个点(i,j)的基于距离值的权重w1(i,j)是通过原图像的距离图像中的值来确定的

(6)

式中，F代表重叠区域。

高斯金字塔图像在第l层的拉普拉斯能量ML及拉普拉斯能量和为SML[12-13,15]：

(7)

(8)

式中，取step=1，表示系数间的可变间距，Mnum和Nnum表示窗口的大小，Mnum=64/2i-1，GI为金字塔图像,Nnum=64/2i-1。

关于拉普拉斯能量和的权重w2(i,j)表示为

(9)

将表示距离值的权重w1(i,j)和拉普拉斯能量和的权重w2(i,j)相结合形成重叠区域内新的权重系数，该权重系数的设计主要是以提高具有高清晰度图像的权重值为目的。

(10)

1.3多频带眼底图像融合

以两幅待融合图像实现Ii和Ij融合过程为例，Mmask为图像Ii和Ij通过距离变换和加权得到的融合蒙版，两幅图像多频带融合过程的具体实现描述如下：

(1)将待拼接图像Ii和Ij映射到基准图像、所在的坐标系下。

(2)获取融合蒙版图像Mmask，计算蒙版图像的加权系数。

(4)获得合成金字塔C，在C的金字塔第l层

LCl(i,j)=GRl(i,j)*LIi,l(i,j)+

GR′*LIj,l(i,j)

(11)

(5)将C从金字塔高层到底层进行扩张和叠加得到最终的融合图像。

1.4层次式多幅图像融合

层次式融合一方面可以减少误差累积，另一方面层次式方法与依次融合方法相比，两种方法执行相同的融合操作次数，但每次融合得到融合图像的尺寸增加速度较慢，从而降低运算量[7,13]。层次式融合的基本思想描述如下：

(1)根据全局配准变换关系矩阵，将N幅图像映射到统一坐标系下，同时记录两两图像之间的配准误差和所有图像掩模图像提取的阈值T_MIN(i)，其中，i=1,2,…,N；第一层融合则根据配准误差最小原则对所有图像进行分组，即具有最小误差的两个图像归为一组进行融合，而对于奇数幅图像序列，将没有组合的一幅图像放入下一层。

(2)在第二层融合时，首先将Tf=min(T_MIN)作为阈值对经过第一层融合获得图像进行阈值分割获取任两幅图像之间的重叠区域，对具有重叠且重叠最小的图像进行重新分组融合，同样，奇数幅图像的处理同第一层融合处理方法相同；然后在下一层重复这样的组合进行融合，直到生成目标图像。图3给出了多幅图像层次方式融合示意图，假设图像序列数N为偶数，根据最小配准误差原则组合，如图像i,j,t∈{1,2,…,N-1}、Ii与Ij、In与Im、It与IN组合，根据层次融合思想逐层实现多图融合。

图3 层次式融合示意图 Fig.3 The schemes of hierarchical fusion

1.5融合效果评价

源图像数据一部分由福建省附属第一医院眼科提供，另一部分源图像来自眼底相机实验系统，数据分辨率主要有3种，分别为2 048像素×1 360像素、1 920像素×1 296像素和13 072像素×2048 像素，图像视角包括60°、45°和20°，图像获得方式有两种：散瞳方式和免散瞳方式。图像数据分别来自400 多位眼底受检者(包括左、右眼)，其中有358 位受检者的眼底发生病变。测试中，通过对其中75 组图像序列(每组2 到18 幅图像不等)，涉及4 898 组图像对(正常眼底图像2 952对，病变眼底图像1 946对)进行实验。

图像融合效果的评价大体上可分为主观评价方法和客观评价方法两种。主观评价方法主要是通过评价者的目视效果进行判定，即凭借主观感知进行评价、判定融合结果的优劣。客观评价方法通常是利用图像的统计参数进行判定，即给出量化指标，模拟人类视觉系统感知机制进行衡量图像融合效果。往往对融合效果的评价需要主观评价与客观评价相结合进行综合评价，即对融合图像质量在主观评价基础上，进行客观定量评价[14,16-17]。

实验时，通过各种眼底图像融合结果示例体现提出方法在主观评价方面的性能，同时考虑单模态眼底图像融合的目的，通过融合图像自身性能的评价指标：熵、空间频率和清晰度来进行讨论和分析，并根据RMSE值和SNR值的对比，进行重叠区域图像融合后信息保持性能的分析。

所进行的眼底图像融合，是以图像拼合消除拼接接缝为目的的，没有涉及参考图像，故选择不需要考虑参考图像的图像熵(H)、空间频率(SF)和清晰度(G)进行融合性能的客观评价。

图像的熵是包含平均信息量多少的度量

(12)

式中，p(i)为灰度i的分布概率，L的范围是0～255。融合后图像的熵值大小反映了图像包含信息量的多少，熵值越大说明融合的效果相对越好。

已知融合图像I大小为M×N，M为图像行数，N为图像列数。

行频率为

(13)

列频率为

(14)

则空间频率为

(15)

空间频率反映图像的全面活跃水平，值越大表示图像越清晰，清晰度(G)采用梯度法衡量

(16)

为了强调融合后眼底图像的清晰度指标，设计基于熵、空间频率和清晰度的、用不同权值的联合客观评价指标(J)，根据各指标意义可知，该值越大表示性能越好：

J=(2G+2SF+H)/5

(17)

在验证方法保持原图像信息的性能测试中，采用均方根误差 (root mean squared error，RMSE)和信噪比 (signal-to-noise ratio，SNR)进行验证[17]

(18)

(19)

式中，In(n=1,2)为建立拼接图像的原始图像，Iover为拼接图像中I1和I2重叠的部分，P为重叠区域的面积。

2 结果

为从主观评价角度进行分析，首先给出两组多幅眼底图像融合结果示例，如图4所示。其中，(b)为(a)原图序列1散瞳正常眼底图像的融合结果，(d)为(c)原图序列2免散瞳正常眼底图像融合结果。从融合结果图可以看出，融合图像无明显的接缝，保留了图像细节结构信息，实现了图像间的平滑过渡，即在主观评价方面获得较好结果。

图4 眼底图像融合结果示例。(a) 原图序列1；(b) 原图序列1融合结果；(c) 原图序列2；(d) 原图序列2融合结果Fig.4 The examples of fundus image fusion.(a) Original image array1；(b) The fusion results original image array1；(c) Original image array2；(d) The fusion results original image array2

为了分析所提出融合方法的性能，对实验图像序列分别采用加权平均法、文献[6]距离幂次(幂次经验值n=2)变换获得蒙版、文献[8]根据距离值确定蒙版加权系数和提出方法进行实验。图5～图7为3组示例图像数据的结果，其中图5为正常散瞳眼底图像，图6为正常免散瞳眼底图像，图7为发生增殖性糖尿病视网膜五期玻璃体病变图像,图5～图7中的(a)和(b)为原图像对，(c)为简单加权平均法融合结果，(d)为文献[6]方法融合结果，(e)为文献[8]方法融合结果，(f)为本方法的融合结果。

由融合结果可看出，在图5(c)、 图6(c)和图7(c)中的加权平均法的拼接痕迹有所改善，但是拼接痕迹仍可见。在图5(d)、 图6(d)和图7(d)中采用文献[6]方法的接缝几乎不可见，但在重叠边缘与图像背景交界处仍存在细微接缝(接缝见箭头指示位置)。而在图5(e) (f)、 图6(e) (f)和图7(e) (f)中采用文献[8]方法和本方法的融合结果中均无明显的拼接痕迹，得到了图像间的平滑过渡，并保留了图像的有效结构信息，主观上观察差别不大。

表1 融合结果客观评价性能对比Tab.1 The performance comparison of objective evaluation

图5 正常散瞳图像融合结果对比。(a) 原图1；(b) 原图2；(c) 加权平均法；(d)文献[6]方法；(e)文献[8]方法；(f)本方法Fig.5 The contrast of normal mydriatic fundus image fusion result.(a) Original image 1； (b) Original image 2；(c) Weighted average method；(d) Literature [6] method；(e) Literature [8] method；(f) The proposed method

图6 正常免散瞳眼底图像融合结果对比。(a)原图1； (b)原图2；(c)加权平均法；(d)文献[6]方法；(e)文献[8]方法；(f)本方法Fig.6 The contrast of normal nonmydriatic fundus image fusion result.(a) Original image 1； (b) Original image 2；(c) Weighted average method；(d) Literature [6] method；(e) Literature [8] method；(f) The proposed method

表1给出了图5～图7中三组图像数据应用不同方法计算融合图像的熵、空间频率和清晰度及联合评价指标的值。由表1可知，提出方法具有较高的熵值，而对发生病变的眼底图像的空间频率略小于未进行滤波的加权平均法，这是由于发生病变的眼底图像有时使得图像混沌模糊，经过多频带方法使得图像经过多次高斯滤波处理，使已经混沌的图像信息减弱，但这对于高分辨率的眼底图像来说是可以接受的。通过突出强调空间频率和清晰度的联合评价指标可以看出，提出方法的优越性明显。

表2为不同方法在重叠区域融合图像的RMSE和SNR值的对比。根据表2各组数据的RMSE和SNR值可以看出，提出方法在重叠区域内，对原图像信息的保持比文献[6]和文献[8]方法具有一定的优越性，且RMSE保持在(0.10±0.05)像素，虽与加权平均方法相比对原信息的保持上略差一些，但加权平均法在融合结果中却存在明显接缝。因此，综合评价主观和客观评价可知，提出方法可以在保证图像细节结构信息良好保持的基础上，有效消除拼接接缝，符合人的视觉要求，满足医生临床诊断要求，在客观指标和主观视觉效果之间取得一个很好的平衡。

图7 病变散瞳眼底图像融合结果对比1。(a)原图1； (b)原图2；(c)加权平均法；(d)文献[6]方法；(e)文献[8]方法；(f)本文方法Fig.7 The contrast of disease mydriatic fundus image fusion result1.(a) Original image 1；(b) Original image 2；(c) Weighted average method；(d) Literature [6] method； (e) Literature [8] method；(f) The proposed method

表2 原图像信息保持性能测试Tab.2 The performance of image information maintained

3 讨论

针对亮度不均引起的接缝及如何保证多幅图像融合后细节信息不丢失的问题，同时需要考虑由于眼底图像的曲面特性，使得配准变换后的图像融合要考虑这种空间差异而引起的接缝问题，本文提出了加权蒙版与层次式融合相结合的多频带眼底图像融合方法。从主观评价实验结果的图4(b)和(d)可知，重叠区域内图像细节信息完整，拼接处过度自然，多幅图像变换后无明显变形，有效解决了拼接接缝及多图拼接误差累积的问题，说明本文图像融合框架合理，通过层次式融合思想减少误差累积确实可行。在图5、图6和图7的对比结果中可以看出，所提出的方法对几种不同类型眼底图像的融合均是有效的，与文献[6]距离幂次(幂次经验值n=2)变换获得蒙版相比，可获得较好的无缝拼接结果，和文献[8]根据距离值确定蒙版加权系数获得的融合结果相比，在主观评价中差异不大。

而从表1的客观评价分析可知，采用基于距离图像在重叠区域内的距离值与拉普拉斯能量和相结合来确定加权系数的方法有一定优越性。在表2原图像信息保持性能测试对比中可以看出，所提出的方法在主观评价和客观评价获得较好平衡的同时，可较好地保持原图的信息，其性能与已有方法相比，有相近或更好的性能。同时，从理论分析可知，采用层次式多幅图像融合方法减少运算量的方法是有效的。

另外，在考虑如何增强多幅图像融合体系的性能方面，将进一步研究多图融合时大范围重复区域多次出现时，融合图像的计算方法包括多次出现区域的处理和中心定位图像的选取。

4 结论

针对眼底图像的特点及对眼底图像融合的要求，提出基于多频带的层次式优化多幅眼底图像融合方法，在保证图像质量的同时，实现图像之间的平滑过渡。提出方法考虑到眼底图像之间重叠区域形状的多变性，获取眼底掩模图像及其距离图像，通过比较距离图像重叠区域内任一点在各自距离图像中的值而产生蒙版图像，利用距离值与各层的高斯金字塔图像的拉普拉斯能量和来确定蒙版加权系数，根据配准误差和重叠率实现层次式融合，降低多幅图像融合的计算量。实验分析表明，提出方法在实现任意重叠形状眼底图像之间平滑过渡的同时，保证了图像的质量。提出方法不仅适用于眼底图像，对于其他图像数据的融合也具有一定参考价值。

[1] James AP,Dasarathy BV.Medical image fusion:a survey of the state of the art [J].Information Fusion,2014,19:4-19.

[2] Xu Zhiping.Medical image fusion using multi-level local extrema [J].Information Fusion,2013,9(2014):38-48.

[3] 张泾周,李婷,吴疆.医学图像的小波变换融合算法研究 [J].中国生物医学工程学报,2008,27(4):521-525.

[4] Cattin P,Bay H,Gool VL,etal．Retina mosaicing using local features [C] // Larsen R,Nielsen M,Sporring J,eds.Medical Image Computing and Computer-Assisted Intervention,Copenhagen:Springer-Verlag,2006,4191:185-192．

[5] Daftari IK,Mishra KK,O’Brien JM,etal．Fundus image fusion in eyeplay software:an evaluation of a novel technique for ocular melanoma radiation treatment planning [J].Medical Physics,2010,37(10):5199-5207

[6] Poletti E,Benedetti G,Ruggeri A.Super-image mosaic of infant retinal fundus:selection and registration of the best-quality frames from videos [C]// Sergio Cerutti,Riccardo Barbieri,eds.Conference Proceedings IEEE EMBS.Osaka:IEEE,2013:5883-5886.

[7] 魏丽芳,潘林,郑炳锟,等.优化的多幅眼底图像拼接方法 [J].中国图象图形学报,2011,16(12):2175-2182.

[8] Choe TE,Cohen I,Lee M,etal.Optimal global mosaic generation from retinal images [C] //Rajapakse,eds.The 18th International Conference on Pattern Recognition,HongKong:IEEE Computer Society,2006：6812-6841．

[9] Szeliski R.Computer Vision:Algorithms and Applications [M].New York:Springer,2010:27-59.

[10] Bailey DG．An efficient euclidean distance transform [C] // Reinhard K,Jovisa Z,eds.International Workshop on Combinatorial Image Analysis.Auckland:Springer-Verlag,2004:94-408．

[11] Wang Jun,Tan Ying.Efficient euclidean distance transform algorithm of binary images in arbitrary dimensions [J].Pattern Recognition,2013,46(1):230-242.

[12] 屈小波,闫敬文,杨贵德．改进拉普拉斯能量和的尖锐频率局部化Contourlet域多聚焦图像融合方法 [J].光学精密工程,2009，(5)：1023-1033．

[13] Ishita D,Bhabatosh C.Multi-focus image fusion using a morphology-based focus measure in a quad-tree structure [J].Information Fusion,2013,14(2):136-146.

[14] Padma G,Vinod K.CT and MR image fusion scheme in nonsubsampled contourlet transform domain [J].Journal of Digital Imaging,2014,27(3):407-418.

[15] Wang Zhaobin,Ma Yide,Gu J.Multi-focus image fusion using PCNN [J].Pattern Recognition,2010,43(6):2003-2016．

[16] 阿都建华,王邦平,王珂,等.基于剪切波变换的医学图像融合算法[J].中国生物医学工程学报,2013,32(3):284-291.

[17] Piccinini F,Bevilacqua A,Lucarelli E.Automated image mosaics by non-automated light microscopes:the MicroMos software tool [J].J Microsc,2013,252(3):226-250.

HierarchicalOptimizationforMulti-FundusImageFusion

Wei Lifang1Yu Lun2*

1(CollegeofComputerandInformationFujianAgricultureandForestryUniversity，Fuzhou350002,China)2(CollegeofPhysicsandInformationEngineeringFuzhouUniversity，Fuzhou350002,China)

In order to remove the seams and ensure information details are not missed in the image fusion,the multi-band blending image fusion method for multi-fundus image was proposed in this paper based on hierarchical optimization fusion.The mask image is achieved by multi-threshold segmentation and Euclidean distance transform.And the improved weighted coefficients are designed based on the distance value in the overlap region and the layers of Gaussian image sum of Laplacian energy for the mask image.The combined objective evaluation was proposed that it is based on information entropy,spatial frequency and definition.The registration error and overlap rate are used to realize hierarchical optimization multi-band image fusion method for grouping fusion.And the image alone is put into next level directly.75 group image sequences are applied involving 4 898 set of images (normal fundus image 2 952 pairs and disease fundus image 1 946 pairs) for testing and validation.The results show that the proposed method is effective that removes the seam and obtains the RMSE value of about (0.1±0.05) pixels in overlap region.The proposed method can make equilibrium between objective evaluation and subjective visual effect.

fundus image; image fusion; mask image; weighted coefficient; hierarchical fusion

10.3969/j.issn.0258-8021.2015.02.006

2014-05-28，录用日期:2014-11-11

国家自然科学基金(60827002)

TP391.09

A

0258-8021(2015) 02-0168-08

*通信作者(Corresponding author)，E-mail:yulun@fzu.edu.cn