多谱段短波红外小鼠静脉荧光观测成像研究

张 瑞，汤心溢，朱雯青

1. 中国科学院大学，北京 100049 2. 中国科学院上海技术物理研究所，上海 200083 3. 中国科学院红外探测与成像技术重点实验室，北京 100049

引 言

生物光学窗口是我们根据生物组织对不同波段的光穿透能力的不同而定义的有利于生物组织进行光学成像的光学窗口。 根据光在生物组织中的衰减长度与波长之间的关系，一般地，我们定义第一生物光学窗口为700～900 nm的光波段，定义第二生物光学窗口为1 000～1 400 nm的光波段，定义第三生物光学窗口为1 500～1 700 nm的光波段。 在这三个光学窗口中，光在生物组织中的衰减长度都比较大，穿透生物组织的能力也比较强。

InGaAs短波红外探测器响应范围很好的覆盖了从900～1 700 nm的波段，精准的包含了第二和第三生物光学窗口。 在生物组织成像应用中短波红外具有以下的成像优势： 生物组织光学损伤小，成像深度大，成像信噪比高，空间和时间成像分辨率高。 根据这些特点配合适当的光学系统，可以实现针对小鼠静脉的显微短波红外光谱观测。

由于显微成像的对焦区域集中的特点，使得小鼠静脉荧光光谱反馈在单张图像信息较少，因此我们设计了一种新颖的融合算法。 融合算法在图像领域有了很好的发展。 主要分为两大类，一类是基于深度学习的方法，另一类也是本工作引入的传统算法。 基于深度学习的方法，需要大量的数据集，这一点在文献11中有了很好的表述，而基于传统算法的融合则可以盲设计。 在传统算法中GF为基于梯度滤波的融合算法，其鲁棒性差，且对细节保持较弱。 NSCT[1]算法则会在平滑区域引入噪声，MFMGRW[2]、PCA[3]以及基于小波变换的融合算法fusedDctVar[4]、fusedDctVarCv[5]除了会在平滑区域引入噪声之外，其在对焦区域检测也缺失鲁棒性。

针对短波红外生物观测的特点，设计了一种基于InGaAs探测器的红外小鼠静脉成像实验，设计了一种高时空分辨率的短波红外光谱生物探测系统，同时针对显微成像的特点，提出了一种多尺度梯度域引导滤波的多焦距融合算法，很好的实现了对小鼠静脉的观测。

1 短波红外光谱成像

设计的InGaAs短波红外探测系统样机如图1(a)所示，其结构原理如图1(b)所示。 我们设计的短波相机最大分辨率达到640×512，其积分时间最大可以达到5 000 ms，可以最大程度接收弱光成像，以减少观测中不必要的强激光对生物组织造成损伤。 配合不同的光学系统，可以实现对生物组织的宏观和显微成像。

图1 InGaAs短波红外小鼠静脉生物成像(a)： 相机样机; (b): 相机原理; (c): 现场荧光实验图; (d): 观察所用的显微镜头Fig.1 Vein imaging of mice with InGaAs camera(a)： The demo of camera; (b): Design framework of camera; (c): Fluorescence experiment; (d): Lens

实验采用了显微光学成像系统如图1(d)所示，可以对小鼠的静脉精准成像。 我们实现了1 100，1 250和1 350 nm等多个波段，不同积分时间下的小鼠静脉数据采集，实验现场如图1(c)所示。 注射荧光显影液后，激光照射在小鼠表面。 在显微光学镜头下，由于对小鼠静脉对焦区域较小，我们设计了针对静脉图像特点的多焦距融合算法。

2 单光谱多焦距融合算法

小鼠静脉交织重叠，构成了一张静脉网络。 单张图像获得的静脉成像信息少，多张图像又缺乏信息的连贯性。 有鉴于此，我们设计了一种基于梯度域引导滤波的对焦像素检测的多焦距融合算法。 通过梯度域引导滤波最大程度的获取原始图像中的高频静脉成分，针对交织网络场景的特点设计融合规则，来获取信息量更加丰富的静脉图像。

2.1 梯度域引导滤波金字塔对焦像素检测

我们设计的融合算法框架如图2所示，通过多尺度梯度域引导滤波提取图像对焦像素区域，然后根据获得的对焦区域求解融合决策函数，最后进一步改进了融合函数，可以获得更好的融合鲁棒性。 下面具体介绍算法的内容。

图2 融合算法总体框架流程Fig.2 Framework of algorithm

2.1.1 梯度域引导滤波

梯度域引导滤波GDGF(gradient domain guided image filter)由引导滤波[6]改进而来。 通过对边缘感知加权函数的改进，可以更大程度的提取图像的边缘细节。 原理如下：

GDGF的损失函数如式(1)所示

(1)

式(1)中，γk为边缘感知滤波因子，ΓI, k为边缘感知加权滤波因子，根据引导图像I计算像素k的重要性。γk与ΓI, k分别定义如式(2)和式(3)

(2)

(3)

其中，N是引导图像的总像素数，ε=(0.001×L)2，L是输入图像的动态范围。Sk为引导图像I像素k周围两个窗口像素的方差的积，即窗口因子Ωr(k)的方差，Ωr(k)表示以像素k为中心半径为r的正方形窗口，S表示所有的Si(i=1，…，N)，μS表示Si的均值。

Sk=σI, 1(k)σI, r(k)

(4)

Sk同时计算了两个窗口函数的方差，这两个窗口函数分别是3×3和(2r+1)×(2r+1)。

最优化值ak和bk分别由式(5)和式(6)得到

(5)

bk=μp, r(k)-akμI, r(k)

(6)

最后输出值为

qi=μa(i)Ii+μb(i)

(7)

其中μa(i)和μb(i)分别表示ak与bk在窗口处的均值，按如式(8)和式(9)计算

(8)

(9)

梯度域引导滤波相对引导滤波可以更好的提取图像的细节，针对窗口大小的控制，可以获取不同尺度的图像频率信息，基于此特点，我们设计了多尺度梯度域金字塔对焦像素提取算法，同时可以利用梯度域引导滤波对边缘的感知特点，针对小鼠静脉的边缘信息可以更好的保存。

2.1.2 单光谱对焦像素检测

在一般的光电成像系统中，系统函数是一个低通滤波器。 模糊区域的系统函数比清晰区域的高频分量更低[7-8]，因此在对焦区域有更多的高频信息。 在多焦距融合算法中首要的就是对焦像素区域检测，我们可以通过对高频分量的检测实现对焦区像素区域的检测[9]。

针对静脉图像条状曲线纹理的特点，设计的基于GDGF金字塔对焦像素检测[10-12]，主要分为以下两个步骤。 首先是多尺度高频分量提取，基于GDGF保持高频细节的特点，通过多尺度的形式，增加算法对图像细节提取的鲁棒性，其原理如图2黄色与蓝色框所示。

其中I是输入图像，初始的滤波半径是rinitial=1，eps为滤波器的参数，设定为0.05，于是得到图像细节的分解函数如式(10)所示

(10)

式(10)中，ri作为迭代因子，其迭代方式为ri=2×r(i-1)，于是根据多尺度分解函数，得到图像细节金字塔Di，其中(i=0, 1, 2, …，n)。 最后通过式(11)获取n的数值

(11)

于是得到IA,IB两张待融合图像的对焦像素细节分量(focus pixel decision map，FPDM)为

(12)

2.2 融合规则

在多焦距融合中，清晰的区域一定含有最多的对焦像素，而模糊的区域则含有最少的对焦像素。 FPDMA与FPDMB分别表示待融合图像IA,IB的对焦像素区域，初步的粗融合规则CF(coarse-fusion)定义如式(13)

(13)

由于生物组织中静脉的交织重叠，图像中所反映的单一焦距下静脉成像会交错出现在不同的焦距的图像中。 这些图像中难免会出现小区域孔洞，而围成这些孔洞的静脉图像是不可以被抛弃的，如图3所示。 再者，如果直接使用这个CF粗融合规则作为最后的融合规则的话，最后得到的融合图像又会不可避免的引入融合噪声，所以对这个CF融合规则做一个精细化处理，如式(14)—式(16)所示

图3 静脉纹理细节的网孔说明Fig.3 Illustration of texture detail

ICF(x,y)=CF(x,y)IA(x,y)+

(1-CF(x,y))IB(x,y)

(14)

FFM(x,y)=GDGF(ICF(x,y), CF(x,y))

(15)

Ioutput=FFM(x,y)IA(x,y)+

(1-FFM(x,y))IB(x,y)

(16)

其中FFM(x,y)表示最终的决策融合映射函数(final fusion map，FFM)。

3 实验与结果讨论

3.1 成像

通过不同波段的短波红外激光照射，实现1 100，1 250和1 350 nm等不同波段的小鼠静脉观测，InGaAs短波探测器成像结果如图4所示。

图4 采集到的不同波段的静脉成像图Fig.4 Vein images of different wavelengths

从图4中我们可以清楚的看到小鼠的静脉短波红外光谱纹理(高亮细节)，四周深色的须状为小鼠的毛发，但是由于光学对焦在不同的区域，单张图像中无法全面的显示小鼠静脉信息，并且不同的波段的红外光在小鼠机体组织中穿透力也有不同。

3.2 单光谱多焦距图像融合

将获得的不同谱段InGaAs短波红外小鼠静脉图像部署到多个算法平台上包括GF，NSCT，MFMGRW，PCA，fusedDctVar以及fusedDCTVarCv等6个算法，并与我们的算法结果进行比较。 同时采用一些评价指标包括Qabf，SCD，MS_SSIM以及SSIM等4个指标[13]，来综合评价我们的算法。

我们设计的基于梯度域引导滤波的金字塔图像融合算法结果如图5—图7所示。 在视觉效果上本算法两焦距融合的边界平滑度更好，且针对原始图起到了一定的降噪作用。 此外，从1 100，1 250和1 350 nm三个波段的图像融合结果中可以看出，其他的几种算法包括GF，NSCT，MFMGRW，PCA，fusedDctVar以及FusedDctVarCv等在不同程度上引入额外的噪声，同时针对对焦区域检测上也缺乏鲁棒性。 在GF算法中，1 250 nm的的融合图像中直接丢失了A图像的对焦区域，其他几种算法也有不同程度的细节丢失。 这些问题在我们设计的算法上得到了很好的解决。 采用了Qabf，SCD，MS_SSIM, SSIM等四个融合图像评价指标，其结果如表1—表3所示，我们的算法综合评价最高，其中在MS_SSIM这个指标上我们得到了全优。 我们设计的算法成功的将两个光谱图像的对焦区域融合到了一张图中，使得我们在单张光谱图像中获取更多的小鼠静脉信息。

表1 1 100 nm 5 000 ms各个融合算法评价指标得分Table 1 Evaluation index scores of each fusion algorithm at 1 100 nm 5 000 ms

表2 1 250 nm 5 000 ms各个融合算法评价指标得分Table 2 Evaluation index scores of each fusion algorithm at 1 250 nm 5 000 ms

表3 1 350 nm 5 000 ms各个融合算法评价指标得分Table 3 Evaluation index scores of each fusion algorithm at 1 350 nm 5 000 ms

4 结 论

设计了InGaAs短波红外小鼠静脉荧光光谱观测系统，通过短波红外光谱照射，实现了针对小鼠静脉短波红外光谱荧光观测成像。 小鼠静脉交织重叠，为了在单张图像中呈现更多小鼠静脉光谱的细节，提出的一种基于梯度域引导滤波金字塔的图像融合算法，先通过梯度域引导滤波金字塔算法检测出图像的对焦像素区域，在基于梯度域引导滤波设计新的融合规则。 最终在短波红外光谱的照射下，实现了针对小鼠静脉的不同波段红外光谱下的广域观测。 由于多焦距成像是因为景深的因素，即短波红外光谱对生物组织的穿透特性，未来可以开展针对生物组织静脉的深度短波红外光谱成像研究工作，进一步探索短波红外光谱对生物组织观测的作用。

图6 1 250 nm波长5 000 ms积分时间融合结果Fig.6 Fusion result of 5 000 ms integration time at 1 250 nm wavelength

图7 1 350 nm 5 000 ms各个融合算法结果Fig.7 Fusion result of 5 000 ms integration time at 1 350 nm wavelength