一种基于先验标记特征的精准图像配准算法

刘天弼, 冯 瑞

（复旦大学 计算机科学技术学院, 上海 201203）

0 引 言

基因芯片与第二代DNA测序是两种重要的高通量基因组学研究技术[1]. 要想准确快速地读取高通量基因芯片的数据, 需要准确地拍摄到芯片图像. 第二代DNA基因测序是通过对基因芯片拍摄图片, 再对图像进行扫描读取碱基序列. 因此, 测序仪对基因芯片拍摄出高质量的图片是测序过程中碱基识别的重要前提. 在拍摄图像的分辨率已经确定的情况下, 单个芯片单元所占的像素数越少, 则整张图像包含的单元越多, 测序的通量就越高. 然而, 要准确地读取基因芯片上的信息，需要拍摄的芯片单元处于非常准确的位置, 对于仅占少量像素的芯片单元而言, 少量的位置偏差就会对信息的读取造成很大的影响. 因此, 将拍摄的视场与芯片位置精准对齐, 对提高第二代基因测序技术的准确度和通量都具有十分重要的意义.

要对基因芯片拍摄出理想的图片, 需要每个芯片单元都有固定的位置坐标. 摄像仪拍摄的图片,反映的是当前视场的位置信息；而计算当前视场与芯片位置之间的偏差, 可视为当前拍摄图片与理想情况下拍摄的图片之间的图像配准问题.

1 图像配准技术研究

图像配准技术种类繁多, 通常每一种配准技术都是针对某一具体场景而设计的, 不同的应用环境要综合各方面的因素来选取相应的图像配准技术.

传统的图像配准算法主要分为基于灰度的图像配准和基于特征的图像配准两大类[2−3]. 基于灰度的图像配准需要遍历计算图像中的灰度、梯度[4]、互信息[5−8]等, 从而确定两幅图像间的空间变换关系,因此运行速度较慢, 计算量较大. 自21世纪以来, 图像配准主要采用基于特征的方法, 利用图像的线特征、面特征、块特征、边缘特征等确定空间变换关系. 从图像中提取关键点的技术历史悠久且种类繁多, 如SIFT (Scale Invariant Feature Transform)[9]、SURF (Speeded Up Robust Features)[10]、SUSAN(Smallest Univalue Segment Assimilating Nucleus)角点[11]、MIC (Minimum Intensity Change)角点[12]、BRISK (Binary Robust Invariant Scalable Keypoints)检测[13]等, 近些年来广泛使用的特征点均有良好的尺度不变特征和抗噪声性能. 此外, 可使用图像中的直线段、边缘曲线等线特征[14−16], 或闭合区域[17−18]等面特征进行配准运算. 基于特征的配准算法自由度较高, 但是快速的特征提取一般仅利用图像中的小部分信息, 会导致精度降低; 具有良好区分度的特征描述往往需要复杂的计算, 反而导致效率降低.

目前，大多数关于图像配准的研究都涉及深度学习. 在过去的几年中, 深度学习使计算机视觉任务具有先进的性能, 如图像分类、目标检测和语义分割等. 使用深度学习做图像配准需要准备合适的数据集进行训练, 以期得到优秀的配准算法模型. 从功能上看, 基于深度学习的图像配准可以分为两大类[19]:①利用深度卷积神经网络估计两幅图像的相似性度量; ②直接利用深度回归网络预测转换参数. 前者进行相似性度量, 仍然需要传统数字图像处理进行迭代做配准优化, 运算耗时较久; 后者是通过大量学习获得配准映射能力, 但只能用于非刚性配准, 无法达到高精度要求.

根据深度学习的种类划分, 此类配准方法也可以分为基于监督学习的配准与基于非监督学习的配准两大类. 基于监督学习的配准, 是通过两张作为配准对的图片输入得到映射向量, 即真实变形场.其关键在于构建数据集时需提供与配准对相对应的真实变形场作为Ground Truth, 作为数据集标签[20−22]. 基于非监督学习的配准方法在训练时只需要提供配准对, 不需要标签[23−25], 完成训练后可通过一张图像计算获得配准图像.

然而, 基于深度学习本身的特点, 其用于高精度配准的场合尚存在诸多问题: ①深度学习依赖于庞大的训练数据集; ②算法训练成果是基于机器学习经验, 无法满足高精准的配准; ③单一场景下的图像配准对泛化能力要求不高, 而深度神经网络计算量庞大, 面对高分辨率图像的复杂计算反而成为负担.

本文通过在基因芯片上的特定位置设置标记, 并确定在理想状态下拍摄的图片中所有标记的位置坐标, 从而达到图像配准的目标；在得到实时拍摄的图像后, 通过对标记的位置进行捕捉可初步将视场中心与标准图像的中心对齐; 然后对标记上的关键点进行捕捉, 得到一系列特征点的坐标; 再对关键点的映射关系进行拟合, 计算出高精度的坐标和角度偏差结果.

本文设计的用于准确捕捉图像特征的新颖卷积方法, 具有良好的抗干扰能力, 其算法简捷, 且易于移植使用硬件加速; 同时使用图像全范围多处采样进行拟合的计算方法, 保证了图像配准达到较高的精确度. 经实践验证和仿真实验分析, 本文算法能取得精确的偏差计算结果, 并且在应用场景下具有很好的鲁棒性.

2 基于先验标记特征的精准图像配准算法

2.1 先验标记设计

基因芯片上布满了携带DNA信息的探针对, 在布满探针对的区域内, 开辟一块不放置探针对的独立区域, 即可看作芯片内的标记. 标记的作用在拍摄的芯片图像中体现出来, 每个探针均携带荧光信号, 在图像上可看作是1个带有灰度信息的单元（cell）, 而标记处不放置探针, 所以始终呈现黑色状态.

在芯片的水平或竖直方向上, 设置1条直线型的标记, 称为track标记; 1条水平track标记与1条竖直track标记相交于一点, 称为cross标记.

先验标记的作用就是人为地在图像中设置关键特征, 作为图像配准过程中实现映射的特征采样点. 几种先验track标记设计示例如图1所示.

图 1 3种track标记设计示意图Fig. 1 Schematic diagram of three track mark designs

图1中黑色线段即track标记, 而cross标记是依赖于track标记的相交呈现的. 设计track标记遵循以下准则.

(1) 所有track标记宽度统一.

(2) 单条track标记必须为水平方向或者竖直方向.

(3) track标记在视场范围内的分布应呈左右、上下轴对称.

(4) cross标记的分布能够体现出对图像的采样比较均衡.

在芯片上设定好先验目标之后, 根据芯片成像的分辨率确定全部标记的位置、宽度及长度的标准数据, 作为标准位置图像的特征. 这样就将视场对芯片位置做校准的问题, 转化为当前视场拍摄到的图片与标准图像的配准问题.

2.2 先验标记的捕捉

2.2.1 track标记的捕捉

视场中, 带有亮度信息的cell随机分布, 而track标记则全部呈现为黑色, 可通过提取其灰度特征进行捕捉.

以竖直的track标记为例, 设其成像宽度为w像素, 长度为L像素. 若拍摄图像的高为H像素, 宽为W像素, 则对于任意一种track标记设计方式, 都有L≤H. 如图2所示, 在竖直方向上截取一段含有track标记的图像, 其中track标记的数量为N个.

图 2 于竖直方向截取图像做水平卷积示意图Fig. 2 Schematic diagram of the image in the vertical direction for horizontal convolution

取卷积核大小为l×w, 其中,l为竖直方向长度,w为水平方向宽度, 且l＜L. 在水平方向上对当前截取的l×W范围上做卷积操作. 卷积核有以下两种设计方式.

(1)均型: 卷积核所有元素值等于1, 相当于卷积范围内覆盖的像素响应均等.

(2)谷型: 卷积核所有元素值为正, 均值等于1. 卷积核行向量的左右边缘元素值最大, 中央最小,列方向上是行向量的复制，相当于卷积覆盖的像素以中间列的响应最小.

均型卷积适用于对抗拍摄场景中噪声的情况; 谷型卷积更适合对抗cell间光线存在串扰的情况.关于卷积核设计对算法鲁棒性的影响, 本文将在实验部分做进一步分析.

设卷积核为k, 截取l×W范围的图像为I′, 执行卷积运算

即

执行卷积运算时, 为捕捉track标记的中心位置, 图像需要在左右方向做像素填充(padding), 其上下方向不做. padding像素的灰度为纯白色. padding数量为

卷积结果为W个元素的一维张量, 即

卷积结果r的索引j与截取图像的x坐标相对应. 将r的元素做升序(Ascend, Asc)排列, 前N个元素对应的索引即为搜索到的track标记的x坐标

再将此N个坐标重新排序, 则截取图像范围内从左至右的track标记的x坐标就全部搜索完成, 即

竖直方向上捕捉track标记的位置, 使用上述卷积核的转置kT, 在竖直方向上执行卷积运算, 即可通过类似的算法得到Ytrack.

2.2.2 cross标记的捕捉

如图3所示, 在cross标记范围的横纵方向上执行track标记捕捉运算, 即可得到cross标记上的交点坐标

其中,X是在cross标记位置上捕捉的x坐标集合,Y是在cross标记位置上捕捉的y坐标集合.

图 3 卷积运算捕捉cross标记示意图Fig. 3 Schematic diagram of the convolution operation to capture cross marks

2.3 基于先验标记特征的图像配准算法

根据当前视场拍摄到的图片与标准图像的配准, 即可计算出拍摄的图片与标准位置之间的坐标偏差和角度偏差. 提高配准精度的基本依据是准确的关键点位置坐标, 提高配准精度的基本方法是通过多采样做最优拟合.

2.3.1 坐标原点对齐

在做位置偏差的精确计算之前, 需要图像的中心位置与标准位置对齐, 将中心位置视作坐标原点.其原因有如下两点.

(1)便于分析特征点配准的映射关系.

(2)寻找cross标记既定范围, 便于捕捉cross标记.

基于2.1节先验track标记的设计准则(3), 视场范围内的track标记呈轴对称分布, 则分别计算相互对称的2个track标记的坐标偏差, 其均值能够消除旋转对坐标偏差的影响, 获得中心点的坐标偏差.

竖直track标记有u行, 每行v个; 水平track标记有k列, 每列l个. 特别是, 对于贯穿型的track标记, 可将竖直track线合理地划分为u段, 使用相应的卷积核进行捕捉, 仍看作u×v的track分布;同理, 将水平track线合理地划分为k段, 仍看作k×l的track分布.

通过坐标变换将当前拍摄图像的中心与标准位置对齐, 即旋转中心与标准中心位置对齐. 因位置坐标是以整数表示, 因此中心位置对齐的误差在1像素左右.

2.3.2 旋转角度与错位拟合

根据先验设计, 视场内存在M×N个cross标记, 通过对cross标记的捕捉, 可得到M×N个交点的位置坐标.

假设当前拍摄图像与标准位置的角度偏差为φ, 因图像具有旋转不变特征, 所以图像坐标系原点位置可任意指定且不影响旋转角度的计算. 取图像中任意一点研究其旋转的位置特性, 如图4所示,像素点A经旋转一定角度之后, 处于A′的位置.

图 4 单个像素错位映射Fig. 4 Mapping of single pixel misalignment

前文根据cross标记已捕捉到M×N个交点坐标, 是符合式(15)的解. 对点坐标进行处理, 得

根据式(15)可得

3 实验分析

两幅图像的配准结果是否优秀, 目前并没有一个很明确的标准. 但是本文具体讨论的配准是一幅图像向标准位置的配准, 针对位置校准问题的评价标准较为明确, 即坐标与角度偏差的准确度.

3.1 位置与角度配准精度分析

实验中, 使用不同型号的基因测序仪进行拍摄: 一种拍摄的图像为矩形, 分辨率为(2 560 × 2 160)像素; 另一种拍摄的图像为正方形, 成像分辨率为(5 012 × 5 012)像素. 使用不同型号的基因芯片, 芯片上的探针数量不同, track标记均使用非均匀贯穿型设计.

基因芯片上的探针经拍摄后形成单元cell, 其对应的分辨率为(a×a)像素. 可认为位置校准的结果误差符合正态分布, 通过多次校准, 取cell位置处灰度最大值的校准结果为标准位, 则可估算出坐标校准的误差期望, 并以图像边缘位置的多个cell位置误差推测出角度配准的准确度. 不同芯片在不同仪器上进行拍摄和配准, 位置与角度配准误差如表1所示, 其中, 仪器型号用分辨率表示. 基于实际拍摄情况, cell周围的像素会受到光线串扰的影响, 因此在估算时对cell周围的像素保留了灰度串扰值. 实际使用中, 位置误差在0.5像素以下即可视为镜头与基因芯片精准对齐.

表 1 位置与角度配准误差Tab. 1 Position and angle deviation in registration

从实验结果可以看到, 在位置校准能够保障精度的情况下, 可以提高研究对象的拍摄分辨率, 还可以减少单个信息单元占用的像素数, 能够大大提升计算通量.

3.2 卷积核设计分析

一般来说, 当出现track标记不在视场范围内、track标记错位对齐的现象, 已属于设备故障问题,本算法不展开故障检测问题的讨论. 对于角度偏差的适用范围, 与标记捕捉时的卷积核设计有关.

3.2.1 卷积核长度与容错角度

无论是对track标记还是cross标记的捕捉, 在搜索x坐标或y坐标时, 卷积运算的作用是完全一致的, 即通过灰度特征提取得到track标记位置. 因此, 卷积核应设计成宽度与track线的宽度一致.

视场中cell的边长等于track标记的宽度, 以均值为最高亮度的25%的灰度随机填充cell. 以track标记中心位置为圆心逐步增加图像的旋转角度, 当通过卷积捕捉到track标记位置坐标偏差时,记当前角度为卷积核的容错角度.

越大的卷积核其抗干扰能力越强, 所以设计卷积核时应当使用较大的长度. 然而图5的结果表明,越长的卷积核的容错角度越小. 在实际应用中, 拍摄系统配准前的角度偏差很容易控制在1°以下, 因此本文讨论的卷积运算能够很好地应用在拍摄系统中.

图 5 卷积核长度与容错角度的关系Fig. 5 Relationship between the length of the convolution kernel and the angle of fault tolerance

3.2.2 卷积切面与干扰

图像配准的干扰主要有两个: 噪声和光线串扰. 噪声在拍摄环境中不可避免, 而拍摄对象中的带有亮度的信息单元在成像时会对周围的像素产生光线串扰.

对先验标记进行捕捉的卷积核, 其宽度等于track标记的宽度. 为方便起见, 取竖直方向track标记对应的卷积核为研究对象, 卷积核内元素是相同的行向量的堆叠. 因此, 取行向量作为卷积核的切面, 研究其形状与鲁棒性的关系. 同样地, 对应于水平方向的track标记, 将列视为切面.

选取如图6所示的几种切面形状的卷积核, 分别是均型、峰型、谷型. 图6 中, 纵坐标表示卷积元素的值(切面不同位置的数值), 横坐标表示卷积切面元素的位置.

在没有噪声和串扰的理想情况下, track线上的卷积结果应当为0, 其他位置的卷积期望为一个正值. 这种卷积特性可以很好地区别track标记与业务区域, 并准确定位track的中心位置.

为了研究卷积核设计对track的捕捉能力, 需要研究track附近的卷积结果之间的区分能力. 在卷积运算接近track标记并经过整个标记直到离开, 这段范围的卷积值应当是先递减再递增, 出现一个极小值, 极小值出现的位置就是捕捉到的track标记位置.

图7为通过卷积捕捉track标注的示意图. 从图7中可以看到, 普通场景下的业务区域卷积结果是明显高于track标记上的, 在track标记位置附近呈现凸函数形状, 因此track捕捉能力仅需要考察track线及其左右相邻范围的卷积值的梯度绝对值, 其值越大表明对位置的区分能力越强.下面通过仿真, 分析理想状态下的卷积运算的区分能力，并进一步分析不同卷积核面对噪声和光线串扰的抗干扰能力.

图 6 不同的卷积核切面Fig. 6 Cross sections of various convolution kernels

图 7 通过卷积捕捉track标记Fig. 7 Capturing a track mark by convolution

(1) 理想状态下卷积核的区分能力

取卷积核宽度为9像素, 在没有任何干扰的理想状态下, 业务区域以25%的亮度单元随机分布, 各卷积核在此状态下的区分度如图8所示. 关注区域在[–5, 5]上, 此区域之外的梯度会逐步趋向0附近.

理想状态下的卷积以均值型的最为理想, 在任何位置处都有稳定的区分能力. 峰型卷积核在track标记左右边缘处区分能力最强. 谷型卷积核在track中心位置处区分能力最强. track捕捉能力只体现在关注区域中, 因此红色标识的峰型卷积核的鲁棒性最差.

(2) 噪声环境下卷积核的区分能力

对理想环境增加不同程度的白噪声, 以噪声均值N0区别噪声的严重程度, 同样对上文列举的5种卷积核进行对比. 对比后发现5种卷积核对track的区分能力有所变化, 结果如图9所示.

图 8 不同卷积核在理想状态下的区分度Fig. 8 Differentiation of various convolution kernels under ideal conditions

图 9 不同卷积核在噪声环境下的区分度Fig. 9 Differentiation of various convolution kernels in a noisy environment

从图9中可以看出, 随着噪声能量的增加, 所有的卷积核对track标记的区分能力都在下降. 在[–5, 5]之间, 梯度最低值代表该卷积核区分能力最差值. 当最差区分能力低到一定程度时, 应对光线的极端分布情况的能力就会变差. 因此，从图9中还可以看出, 峰型卷积核并不适用于track标记的捕捉, 而均型卷积核保持着稳定的区分能力.

(3)光线串扰环境下卷积核的区分能力

增加了40%的光线串扰之后, 5种卷积核对track的区分能力有所变化, 结果如图10所示. 图10表明, 红色曲线代表的峰型卷积核区分能力变得更差, 因此峰型卷积核完全不适合在捕捉track标记时使用； 在有光线串扰的情况下, 均值型卷积核的区分能力也有所下降, 而谷型卷积核则表现出很好的区分能力, 并且比均值型卷积核具有更好的稳定性.

图 10 不同卷积核在串扰环境下的区分度Fig. 10 Differentiation of various convolution kernels in a light crosstalk environment

鉴于谷型卷积核对光线串扰具有优秀的抗干扰能力, 因此本文对谷型形状与串扰之间的相应关系做了进一步观察. 如图11所示, 设计5种谷型深度不同的卷积切面, 其卷积核各元素之和与均值型卷积核保持一致. 图11中, 纵坐标表示卷积元素的值(切面不同位置的数值), 横坐标表示卷积切面元素的位置. 选取10%、20%、30%、40% 这4个等级的串扰, 5种卷积核对track标记的区分能力如图12所示.

图 11 5种抗串扰卷积核切面Fig. 11 Cross sections of five anti−crosstalk convolution kernels

图 12 5种抗串扰效果Fig. 12 Anti−crosstalk effects of five convolution kernels

从图12可看到, 10%串扰时卷积核(2)(图12( a))的抗串扰性能最好, 20%串扰时卷积核(3)(图12 (b))的抗串扰性能最好, 30%串扰时卷积核(5)(图12( c))的抗串扰性能最好, 串扰超过40%之后, 上述5种谷型卷积核性能均变得很差. 因此可以得出结论： 谷型卷积核的“谷深”越深, 应对光线串扰的能力越强； 谷型卷积核的设计需对应具体的串扰强度, 从频域角度计算其滤波性能.

综上所述, 从实验结果可知，卷积操作对捕捉track标记具有很高的准确性和很强的鲁棒性, 在白噪声的环境下, 均型卷积核始终保持了较好的track标记区分能力, 而在光线串扰的环境中, 需要使用谷型卷积核以保证其区分能力.

4 结束语

本文针对基因测序仪需要将镜头与基因芯片精准对齐的问题, 提出了一种基于先验标记的图像配准方法: 计算当前视场与理想位置的偏差. 该方法通过在基因芯片上人工设计标记, 再分析所拍摄图像的灰度特征, 利用卷积运算对图像特征的良好提取能力捕捉标记, 使用拟合的方式计算偏差. 从实际应用和仿真分析可知, 本文算法对位置偏差的计算结果具有良好的精确度, 误差限制在0.5像素以下,并且有良好的抗干扰能力, 在第二代基因测序技术提高测序准确度和通量的需求下具有很好的应用前景.