基于模板匹配和不变矩的浮游动物快速测量方法

徐 帅，杨俊毅，郑旻辉，谢尚微

(1.国家海洋局 第二海洋研究所，浙江 杭州 310012；2.国家海洋局 海洋生态系统与生物地球化学重点实验室，浙江 杭州 310012)

0 引言

浮游生物在海洋生态圈中扮演着重要角色，是海洋生态学中重要的研究对象。目前海洋生物学家对浮游生物的研究手段主要以现场采样、实验室人工鉴定与计数的传统观测手段为主，人工识别费时费力，远不能适应和满足现代海洋科学发展的需要。

随着图像处理技术的不断发展，使得基于图像处理的浮游生物自动识别成为可能，目前已有研究者取得瞩目的成果。例如，CULVERHOUSE et al[1]基于细胞形状和表面纹理特征开发了甲藻分类软件DiCANN；LOKE et al[2]基于浮游动物轮廓特征利用最近邻分类规则进行甲藻识别；HU et al[3]以形态特征为基础建立神经网络分类器，然后以表面结构特征为基础进行浮游生物图像识别；TANG et al[4]结合经典的几何外型特征以及不变矩等特征对海洋浮游生物图像实行自动识别分类；YANG et al[5]提出基于内容特征的浮游生物图像识别算法。

随着近年来水下观测技术快速发展，浮游生物观测方法逐渐从传统方法向原位观测分析方向发展，受原位观测平台计算资源的限制，上述根据生物特征的智能识别方法存在一些不足，如神经网络方法样本训练量大且效率低下，支持向量机方法虽然识别率较高但需提取大量高维特征，前期数据准备计算量大等，这些都制约着其在原位分析中的应用。针对这些不足，本文提出将计算机视觉中的模板匹配方法与不变中心矩描述方法相结合，将模板匹配作为浮游动物识别的初筛环节，在模板匹配限定的子空间内采用比较不变矩描述方法判定目标动物，然后测量目标动物的尺寸信息，统计动物个数。这种方法与支持向量机等识别方法相比速度快，计算量小，不仅能够验证显微图像中是否存在模板表示的目标动物，还能快速找到模板对应的动物目标的质心位置，即可方便地得出图像中动物的长度、宽度等尺寸信息，适合应用在浮游动物的原位观测分析中。

1 图像分析方法

1.1 方法步骤

本文所提出的方法步骤如图1所示，先对图像进行预处理，然后确定感兴趣区域，在此基础上利用图像分割构建模板，再利用模板匹配进行初筛限定目标范围，分别求出模板与待识别显微图像中的动物的质心坐标和中心不变矩，分别比较模板的不变矩与图像中各动物的不变矩，差值小于设定的阈值即可初步判定为模板所示的动物。将场景中目标的质心与模板质心对齐，利用“缩放”和“旋转”等图像处理方法将模板与图像中的目标对齐以最终确定目标是否为模板所表示的动物，并计算出匹配到的目标动物的尺寸。

图1 基于显微图像的海洋浮游动物自动测量步骤Fig.1 Marine zooplankton automatic measurement steps based on microscopic images

1.2 图像预处理

浮游动物显微图像如图2所示，图像中存在大量噪声还夹杂一些碎片以及光源带来的较大干扰,即背景噪声，需要对图像进行预处理。

图2 浮游动物显微图像Fig.2 Zooplankton microscopic image

预处理主要是在创建模板之前对原始图像进行增强，目的是突出图像整体或者局部的特征，减少噪声，提高图像质量[6]。图像增强方法大体可以分为空域和频域增强两类。常用的空域增强法有中值滤波和均值滤波，而频率域图像增强法则是将经过预处理后的图像进行傅里叶变换，将其转换到频率域后再利用低通滤波器进行去噪处理[7]。

由图2可见，显微图像中的噪声主要为椒盐噪声，中值滤波对滤除椒盐类噪声具有较好的效果，实验中对比均值滤波、中值滤波和高斯低通滤波发现(图3)，中值滤波器在去除噪声的同时还能尽量保留浮游动物图像的边缘细节，因此本文采用一个9×9的模板的中值滤波器进行图像增强。最后,采用直方图均匀化增强动物目标与背景对比度，以利于阈值分割。

图3 中值滤波(a),高斯低通滤波(b)以及均值滤波(c)去噪效果Fig.3 The denoising result of median filter(a), Gaussian low-pass filter(b) and the average filter(c)

1.3 图像分割与模板创建

图像分割在整个过程中居于承前启后的位置，既可以检验图像预处理的效果，也是后续创建模板的前提。针对海洋浮游动物显微图像对比度较好的特点，本文综合采用基于灰度直方图确定阈值、Canny边缘检测算子和基于数学形态学闭运算等多种图像分割算法，得到的边缘可以提供很好的形状特征，其中的区域即可用来构建模板。

1.3.1 基于阈值的浮游动物目标分割

图2所示的海洋浮游动物显微图像的灰度直方图如图4所示，由图可见灰度值30～120区间为浮游动物体的像素，灰度值大于120的山峰对应背景像素。这里取120灰度值为阈值T，便可以将浮游动物与背景明显地分开，实现目标分割。

图4 浮游动物图像灰度直方图Fig.4 The gray histogram of zooplankton image

为了提高阈值分割精度，本文采用了OTSU[8]方法在30～120目标动物灰度值区间上自适应地再进行第二次阈值分割，得到比较满意的自适应阈值区间。OTSU算法根据灰度直方图的谷底选择100为阈值进行分割，第二次灰度值阈值分割后得到的区域如图5a所示。

图5a已经大致可以看出动物目标的外形，但是初步分割结果存在图像噪声和大量的杂质碎片。为此本文设计了一套基于特征直方图的阈值分割方案，利用面积、周长特征将背景噪声剔除。考虑到噪声碎片的面积小于最小的动物体面积，可以通过设置关于目标面积特征的阈值，排除所有面积过大或者面积过小的非动物目标，剩余的区域即为所要处理的目标区域(图5b)。

图5 初步阈值分割结果(a)和最终阈值分割结果(b)Fig.5 Preliminary threshold segmentation results(a) and the final threshold segmentation result(b)

1.3.2 基于Canny算子的浮游动物边缘检测

为精确定位海洋浮游动物体的轮廓边缘，需要对其进行边缘检测，本文采用了基于Canny算子的边缘检测算法，Canny边缘检测算子有低误判率，高定位精度，抗噪声能力强等优点[9]。Canny算子边缘检测具体实现步骤如下：

(1)高斯平滑滤波。

(2)非极大值抑制。通过抑制非极大值而保留局部梯度幅值最大的点，可以得到相对更加精确的边缘图像。

(3)双阈值和边缘连接。非极大值抑制处理后，图像边缘像素中可能还混杂着由噪声而引起的伪边缘。Canny边缘检测算法通过双阈值区分边缘与伪边缘。针对浮游动物图像的特点，绘制出梯度幅度图的梯度幅值直方图，然后通过设定非边缘点像素数与梯度幅度图中所有像素数目的比例来确定高阈值，取高阈值的0.4倍作为低阈值，程序自适应地设定双阈值，分别得到高阈值与低阈值图像。当阈值较高时，去噪效果好，但同时也会造成有效边缘信息的损失；反之则去噪效果不明显，但会保留边缘细节。在高阈值图像的基础上，通过低阈值图像的补充来连接浮游生物图像的边缘。图6为双阈值和边缘连接后得到的边缘效果图。

图6 Canny边缘检测结果Fig.6 Result of Canny edge detection operation

1.3.3 基于数学形态学的浮游动物图像处理

经过前面的处理后获得区域内仍存在孤立的小噪声点和目标区域内的细小空洞，在其周边存在凹坑和毛刺。主要采用腐蚀、膨胀及其组合方法即开运算与闭运算等数学形态学算子去除干扰。数学形态学是一种常用的数字图像处理和模式识别的理论和方法，主要内容是一套运算、概念和算法，用来描述图像的基本特征，它不同于常用的频域或时域的方法，而是建立在积分几何和随机理论的基础之上[10]。

本文首先对浮游动物图像轮廓进行膨胀操作，有效地填补图像中的空洞，并起到连接边缘检测后不连续区域的作用，再腐蚀图像边沿，消除一些边界点，使边界向内收缩。将膨胀操作得到的轮廓减去腐蚀操作得到的轮廓，在有效地平滑边界、消除目标间的粘连的同时保持其面积等特征信息。通过一系列的图像处理算法，最终得到图7a，在此基础上将区域裁剪出来，最终构成一个18×20的生物模板(图7b)。

图7 数学形态学处理后消除空洞的图像(a)和基于图像分割创建的模板(b)Fig.7 The image after mathematical morphology processing(a) and a template based on image segmentation(b)

1.4 不变矩与模板匹配

1.4.1 求不变矩特征

浮游动物的形态特征千姿百态，前面的方法能从生物图像中提取区域和轮廓构建模板，为了使模板能够匹配位姿不同的目标，本文针对海洋浮游动物形态特点，采用了一些对于平移、旋转和尺度变换的不变量作为特征，这样模板可以适应浮游动物目标的方位旋转或尺度变化，减少误差。

不变矩特征是形状识别中常用的特征，(m+n)阶的中心距具有平移变换不变性，定义为：

μmn=∑x∑y(x-x0)m(y-y0)nf(x,y)

(1)

(2)

因为数字图像中像元的离散性，这些中心距函数近似具有旋转、平移和缩放不变性。

1.4.2 模板匹配过程

将上述对于物体不变矩的描述方法与计算机视觉中的模板匹配算法相结合，提出一种适应旋转、平移和缩放的模板匹配方法，利用该方法匹配识别浮游动物显微图像，具体步骤如下：

(1)首先用创建好的模板来进行粗筛选，为了减少匹配复杂度，提高模板搜索速度，本文采用基于图形金字塔方法进行匹配。使用归一化相关系数作为匹配准则，对于每个滤波窗口，经过灰度分布标准化预处理后，使用动物模板进行匹配，若相关系数超过给定阈值(本文取0.25)，则认为通过了匹配初筛选，但是图中动物存在旋转或缩放等情况，这一阶段得到的匹配结果会包含其他非模板目标动物，误差较大，需要应用中心不变矩约束进一步提高匹配精度。

(2)分别计算图像中各匹配到的目标与模板的中心矩，比较并取相差值，差值小于某一阈值的即初步认为是模板表示的动物。对模板缩放操作，使其尺寸与图像中与之不变矩基本相同的动物一致，利用坐标的旋转变换来确定图像中动物对应的旋转角度，采用的算法分别为：

(3)

(4)

式中：(x2,y2)为像素新坐标，(x1,y1)为原坐标，(x0,y0)为形心坐标，N为缩放系数，θ为旋转角度。将两幅已经对准的图像进行“异或”操作，以剩余像素的多少来证实该物体是否为模板所示动物，若是，则给出匹配物体在图像中的坐标和个数。

2 实验验证

本文采用Leica M205A立体显微镜对在南海海域采集到的浮游动物样本进行显微成像，采用Halcon图像处理库实现上述所有的图像处理算法，并将匹配结果可视化显示在界面上。Halcon是一套完善的机器视觉算法包，拥有成熟配套的集成开发环境和专业的图像处理算子，能够胜任本文的研究。

为了验证方法的有效性，本实验选取了100幅浮游动物显微图像(主要包括水蚤类和桡足类)作为测试图像，其中20幅为无目标动物的图像，80幅同时含有非目标动物和目标动物的图像。图8为本文用到的部分显微图像示例，可以看到动物的尺寸、方位角度各不相同，按照上文方法创建动物模板。

利用本文提出的匹配方法对测试图像进行检测，共成功检测出71幅显微图像中存在模板对应的动物，准确率为88.75%，并统计出动物的数量和尺寸数据，部分实验结果如图9所示。

图8 部分浮游动物测试图像示例Fig.8 Examples of partial zooplankton test image

图9 匹配和计数结果Fig.9 Results of template matching and counting

将图9中检测出的与模板对应动物的尺寸数据与实验室显微镜测量软件测量数据进行对比验证。获得的浮游动物尺寸数据均是以像素为单位，需要将其转化为实际长度单位。根据显微镜CCD的参数，计算出像素长度转为实际尺寸(mm)的比例如下式：

1像素=0.003 5 mm

(5)

根据上式计算出浮游动物形态特征数据如表1所示。

表1 浮游动物形态特征测量数据Tab.1 Morphological characteristics data of zooplankton

由表1结果可知，实验测得动物尺寸与实测数据相比误差为0.032～0.106 mm。利用本文的方法能够快速地在大量显微图像中匹配到特定的浮游动物，这样可以大大提高海洋动物研究者的鉴定效率，还可以快速统计显微图像中的浮游动物个数和尺寸信息。

3 小结

本文根据处理浮游动物显微图像的需求，提出了将模板匹配方法与不变中心矩描述法相结合的方法对浮游动物的显微图像进行识别。与基于特征的识别方法相比，本文算法比较稳定，计算量小，识别耗时短，识别效率高。在实用性方面，还需要进一步提高模板的质量，同时应考虑动物体重叠和缺损等情况，使方法更具有实用价值。

[1]CULVERHOUSEPF,REGUERAB,GONZALEZGS,etal.Doexpertsmakemistakes[J].MarineEcologyProgress,2003,247(247):17-25.

[2]LOKERE,BUFJMHD,BAYERMM,etal.Diatomclassificationinecologicalapplications[J].PatternRecognition,2004,37(6):1 283-1 285.

[3]HUQ,DAVISC.Accurateautomaticquantificationoftaxa-specificplanktonabundanceusingdualclassificationwithcorrection[J].MarineEcologyProgress,2006,306(8):51-61.

[4]TANGX,STEWARTWK,HUANGH,etal.Automaticplanktonimagerecognition[J].ArtificialIntelligenceReview,1998,12(1):177-199.

[5]YANGR,ZHANGR,SUNS.Automatedclassificationofzooplanktonbasedondigitalimageprocessing[J].ComputerSimulation,2006,23(5):167-170.

[6]RAFAELCG,RICHARDEW.Digitalimageprocessing[M].Thirdedition.Beijing:PublishingHouseofElectronicsIndustry,2011.

[7]KENNETHRC.Digitalimageprocessing[M].Beijing:PublishingHouseofElectronicsIndustry,2002.

[8]OTSUNA.Athresholdselectionmethodfromgray-levelhistograms[J].IEEETransactionsSystemManandCybemectics,1979,9(1):62-66.

[9]CANNYJ.Acomputationalapproachtoedgedetection[J].IEEETransactionsonPatternAnalysis&MachineIntelligence,1986,8(6):679-98.

[10]CARSTENS,MARKUSU,CHRISTIANW.Machinevisionalgorithmsandapplication[M].Beijing:TsinghuaUniversityPress,2008.