鱼类行为是指鱼类进行的各种运动，是鱼类对环境变化的外在反应[1]。鱼类游泳行为数据是研究鱼类生活习性的重要信息来源之一[2]。鱼类游泳轨迹在一定程度上反映了鱼类对所处的环境是否适应，速度数据是鱼类在特定环境下游泳能力的重要指标[3]。传统的鱼类游泳行为监测主要以人为观察、手工记录为主，劳动强度大、耗时长，且效率较低[4]。机器视觉技术主要利用计算机模仿人的视觉功能，从目标物体的图像序列中提取特征信息，对其进行处理并加以分析理解，用于实际的检测、跟踪和识别[5]。早期的目标跟踪算法以运动模型为主，经典的算法有卡尔曼滤波[6]等，后期人们意识到主要特征运用于物体跟踪的重要性[7]，如基于颜色特征的均值漂移算法[8]等。而有些运动目标形状时刻发生着变化，研究人员发现，通过训练所得的分类器对目标跟踪会产生较好的效果，如利用支持向量机[9-10]训练出分类器进行目标跟踪等。针对鱼类游泳跟踪算法，许多学者进行了大量的研究。徐盼麟等[11]提出一种基于单摄像机的鱼类三维观测方法，即将防水镜面安装在鱼缸上方，达到等效于另一台摄像机拍摄的功能。该方法虽然实现了鱼类三维坐标、速度等信息的自动输出，但由于防水镜得到的视角偏差较大，重建的跟踪轨迹误差较大。黄一凡等[12]提出一种最大重合面积跟踪算法，研究斑马鱼在不同水质的应激行为。该方法只适用于单个鱼类目标或多个鱼类目标之间距离较大的条件，当鱼类游泳过快或者两帧之间时间间隔太大导致重合面积极小或为零时，无法继续跟踪。颜鹏东等[13]将鱼体图像连通域面积最大的轮廓作为个体鱼区域，把区域的形心点作为个体鱼的位置，对每一帧图像执行相同的操作，得到鱼的整个运动轨迹。但对多帧图像进行预处理及形态学处理获取质点，使计算量增大，无法保证鱼类跟踪的实时性。基于核相关滤波器的跟踪算法较少运用到鱼类跟踪，该算法的优势在于高速的跟踪速度，使用多通道梯度方向直方图特征训练出核相关滤波器，提高了运动物体跟踪的实时性。
本研究中，对核相关滤波跟踪算法[14]进行了优化升级，实现了对青鳉Oryzias latipes游泳行为的跟踪研究，可为研究青鳉在特定环境下的生活习性提供数据支持。
1 材料与方法
1.1 材料
青鳉是国际医学用鱼，其体型较小，对水质、环境变化特别敏感，且活动性大，有利于进行视频跟踪分析[15-16]。为实现对不同大小的青鳉游泳领地进行分析研究，本试验中选用3尾青鳉为研究对象，鱼体体长分别为1.8、2.9、3.6 cm，鱼龄为5～6个月。试验前将待观测青鳉在温度为26.1 ℃、pH值为7.5的相同水体中暂养7 d，持续供氧，以使待观测青鳉适应试验环境。
1.2 方法
1.2.1 试验设计 试验环境由硬件部分和软件部分组成。硬件部分如图1所示，主要包括鱼缸(38 cm×26 cm×27 cm)、加热棒(20 W)、环形灯、摄像头、图像采集卡、pH与温度传感器和PC机。在鱼缸的正上方安装一个摄像头，通过图像采集卡将其连接到PC机，采集青鳉游泳视频图像。其中摄像头为奥速X200，采集视频帧分辨率为640×480像素，帧率为30帧/s，输出视频帧为RGB三通道彩色图像[17]。
图1 青鳉游泳跟踪硬件平台
Fig.1 Experimental device diagram for tracking medaka swimming
软件部分采用OpenCV3.2.0和MFC实现[18]。其整体流程如图2所示，通过高斯滤波图像处理，使用增强的核相关滤波跟踪算法对青鳉的游泳行为进行跟踪，最后利用所得到青鳉游泳参数数据结合有关理论对青鳉的运动行为进行分析。
1.2.2 高斯滤波 在青鳉游泳跟踪过程中，外界光照、鱼缸中的杂质、气泡等因素均会引入大量噪声，因此,需在青鳉游泳跟踪前先对视频图像进行滤波处理，以提高青鳉游泳跟踪效果。本系统采用高斯滤波对视频图像滤波处理。
高斯滤波是一种线性平滑滤波，其原理是用一个模板扫描图像中的每一个像素，用模板领域内像素加权的平均灰度值代替模板中心像素点的值[19]。在图像预处理中，一般采用二维高斯函数作为平滑滤波器[20]，其计算公式为
图2 青鳉游泳跟踪软件平台
Fig.2 Swimming tracking software platform of the medaka
(1)
其中：x、y分别为图像中像素的行列值；σ为标准差。
二维高斯函数具有旋转对称性，即在各个方向上的平滑程度是相同的。高斯滤波器宽度由参数σ表征，σ越大高斯滤波器的频带越宽，平滑程度就越好(图3)。
1.2.3 增强核相关滤波跟踪算法 核相关滤波跟踪算法使用多通道的梯度方向直方图特征[21]进行训练，具有高速的跟踪速度。在线目标跟踪评估基准[22]中整体性能表现优异，平均跟踪准确率达到73.2%(图4)。
核相关滤波跟踪算法把跟踪问题抽象为一个线性回归模型求解问题，设代表目标图像的输入为z，权重系数为w，输出为f(z)=wTz，其核心就是要找到f(z)的值，利用训练样本(xi,yi)，最小化平方误差目标函数,找到系数也就是滤波器w，即
(2)
其中：λ为正则化系数,用于控制回归分析中的过适应问题；f(xi)为输出值；yi为回归期望值。解得：
w=(XHX+λI)-1XHy。
(3)
其中：XH为复共轭转置矩阵;I为单位矩阵。
该算法引入循环矩阵，循环矩阵能够被离散傅里叶矩阵对角化，使得矩阵求逆转换为特征值求逆，将公式(3)转换到频域进行计算，应用离散傅里叶变换提高运算速度，然后再解逆变换回空域，从而得到响应最大的解ymax。
图3 二维高斯分布图
Fig.3 Gaussian 2D distribution
图4 50种图像序列集试验结果[14]
Fig.4 Experimental results of 50 image sequences[14]
循环矩阵X的构造是由基样本x的循环偏移量组成，即核相关滤波算法的训练样本均由目标样本的循环偏移量得到，这样就增加了训练样本的数量，提高了训练的分类器效果。
(4)
将式(4)代入式(3)得到：
(5)
对式(5)两边同时傅里叶变换得：
(6)
解得：
(7)
可由式(7)求傅里叶逆变换得到线性回归系数w。
由于大多数情况下解决的是非线性问题，于是引入核函数[23]来将问题扩展到非线性空间。把低维空间的计算映射到高维的核空间，上述问题的回归系数w和对偶空间α的线性组合为
(8)
回归问题就转化为
f(z)=
(9)
其中，k(z,xi)=ψ(xi)ψ(zi)为核函数，常用的有高斯核、多项式核和线性核，核相关滤波算法用的是高斯核，因此，将求回归系数w转化为求α。可得非线性问题的解为
α=(K+λI)-1y。
(10)
结合之前w的推导，可得:
(11)
其中，Kxx为核矩阵K=C(Kxx)的第一行。
这样对于输入z，同样构造矩阵K=C(Kxx),Kxx就是x和z的核相关，据此可以对所有部分进行回归检测：
(12)
然后将响应值
傅里叶逆变换到空域，找到最大值对应的部分就是跟踪的结果。
增强的核相关滤波跟踪算法原理与传统核相关滤波算法跟踪基本原理相同，增强的核相关滤波跟踪算法比传统核相关滤波跟踪算法在第一帧多训练出m-1个核相关滤波器(m为大于1的整数)，然后与新输入帧进行相关性计算，得到的响应值就是预测的跟踪结果，响应值最大处是跟踪对象的可能性最高。增强的核相关滤波算算法流程如图5所示。
图5 增强的核相关滤波算法
Fig.5 Enhanced correlation filter tracking algorithm
1.2.4 核密度分析及热力图绘制 核密度分析用于计算要素在其周围邻域中的密度，通过考察规则区域中点密度的空间变化来研究点的分布特征[24]。其计算公式为
(13)
其中：f(s)为空间位置s处的核密度计算函数；h为带宽；n为与位置s的距离小于或等于h的要素点数；k为核函数；s-ci为估计点s到样本ci处的距离。然后，可以根据核密度分析方法，绘制出青鳉游泳轨迹热力图[25-26]。
2 结果与分析
利用上述中搭建的青鳉游泳跟踪硬件平台进行试验，将摄像头固定在鱼缸的正上方，以获得较好的视野。使用 OpenCV3.2.0和MFC软件开发的青鳉游泳跟踪应用程序，实现对青鳉游泳跟踪，程序界面实时显示青鳉游泳的当前帧数、位置坐标、距离、游泳速度、pH和温度值，应用程序主界面截图如图6所示。
图6 青鳉游泳跟踪主界面
Fig.6 Main interface of the swimming tracking system in medaka
2.1 青鳉游泳跟踪
在水体环境温度为26.1 ℃、pH为7.5的条件下，选用5～6月龄长度约为1.8、2.9、3.6 cm的3尾青鳉作为观测目标。结果表明，即使青鳉在游泳过程中出现少许视频帧遮挡时，所用算法也能够实现实时稳定跟踪青鳉游泳行为，如图7(a)所示。
对得到的位置数据进行分析，可得到青鳉游泳轨迹曲线。从图7(b)可见，利用位置数据画出的青鳉游泳轨迹与利用增强的核相关滤波跟踪算法实时获取的青鳉游泳轨迹基本吻合，说明所用算法能有效地提取青鳉的游泳位置信息，实现对青鳉的跟踪。测试2万帧视频序列中，丢失806帧，丢失率为4.03%，因此，可有效地对青鳉游泳进行跟踪。
提取100帧青鳉视频图像序列跟踪速度数据进行常规统计分析(表1)，绘制出青鳉游泳速度曲线如图8所示。从图8可见：青鳉游泳速度随着帧数的变化而变化，试验鱼1的游泳速度开始时接近匀速，到第18帧左右时持续增加，第60帧左右时平稳下降；试验鱼2的游泳速度开始时呈上升趋势，到第25帧左右时开始持续下降，第50帧左右时又开始持续增加，到第96帧左右时趋于平稳；试验鱼3的游泳速度先平稳上升，到第30帧左右时开始下降，到第63帧左右时趋于平稳。从青鳉游泳速度曲线可以清晰地看出青鳉游泳速度状态。
图7 青鳉游泳轨迹
Fig.7 Swimming track of the medaka
图8 青鳉游泳速度曲线
Fig.8 Swimming speed curve of the medaka
表1 青鳉100帧游泳速度数据特点
Tab.1 Swimming speed data of the medaka at 100 frames
青鳉体长/cmbody length平均速度/(cm·s-1)average speed极差/(cm·s-1)range方差σ2/(cm·s-1)2square variance标准差σ/(cm·s-1)standard deviation1.81.4971.0300.3870.6222.91.4251.4300.1140.3373.60.7160.0620.8320.249
2.2 传统与增强核相关滤波算法的跟踪
利用青鳉游泳跟踪硬件平台进行试验，分别使用传统的核相关滤波跟踪算法和增强的核相关滤波跟踪算法实现对青鳉游泳行为的跟踪，在青鳉游泳视频图像序列第1032帧，两种跟踪算法的跟踪效果如图9所示。传统的核相关滤波跟踪算法跟踪目标丢失，而增强的核相关滤波跟踪算法未出现跟踪目标丢失的状况(图9)。测试5万帧视频图像，增强的核相关滤波算法比传统核相关滤波算法跟踪率提升4.45%。
2.3 青鳉游泳轨迹分析
首先利用青鳉游泳跟踪硬件平台，分别获取3种规格青鳉各1万帧游泳视频图像序列，然后利用增强的核相关滤波跟踪算法提取青鳉的位置数据，继而利用核密度分析方法，绘制出3种规格青鳉游泳轨迹热力图和青鳉游泳领地等值线分布图，最后分别对3种规格青鳉游泳领地进行分析，结果如图10所示。
3种规格青鳉游泳轨迹热力图如图10-A～C所示，其中图例表示单位面积青鳉点位数量，即红色区域表示单位面积青鳉游泳位置数据密集区域的核密度值较大。不同规格青鳉游泳位置相似，集中在鱼缸右侧和下侧，其中右侧靠上部分和下侧靠左部分居多。青鳉游泳活跃度在鱼缸右上侧和左下侧居高(图10-D)，体长1.8 cm的青鳉活动范围最大，体长3.6 cm的青鳉活动范围最小。由于在青鳉游泳视频采集时，获取的是刚放入试验鱼缸的青鳉游泳视频图像序列，在受到外界环境刺激青鳉会有短期的逃避行为，使得其在鱼缸边缘聚集，致使青鳉在鱼缸边缘的活跃度居高。
3 讨论
本研究中基于核相关滤波跟踪算法，在试验环境下实现了对青鳉游泳轨迹实时跟踪。利用传统的核相关滤波跟踪算法可以较好地实现青鳉游泳跟踪，但在青鳉游泳过快和长时间跟踪时，会出现目标跟踪漂移和丢失的情况，从而使跟踪不能持续进行。而增强的核相关滤波跟踪算法实现青鳉游泳跟踪，改善了青鳉游泳过快和长时间跟踪目标出现漂移和丢失的情况，但是会出现轻微的跟踪延迟现象，基本满足了信息实时采集的要求。
本研究中利用增强的核相关滤波跟踪算法进行青鳉游泳跟踪过程中，在目标严重多帧遮挡的情况下，会出现青鳉跟踪漂移和丢失的问题，使后续的跟踪不能持续进行。在进行多尾青鳉跟踪的过程中，因为鱼体规格差异较大，互相遮挡对试验结果也影响较大。针对目标遮挡问题，邢运龙等[27]提出融合卡尔曼滤波器形成预测跟踪校准的跟踪机制，以克服跟踪目标的完全遮挡问题。熊昌镇等[28]提出了一种自适应特征融合的核相关滤波目标跟踪算法，该算法结合物体颜色特征和方向梯度直方图特征，在目标严重遮挡情况下，仍然能够稳定地跟踪目标。在后续的研究中，考虑利用核相关滤波跟踪算法和青鳉检测相结合，以改善青鳉游泳严重多帧遮挡的问题，即在青鳉不存在遮挡时，使用核相关滤波跟踪算法进行跟踪，当存在严重多帧遮挡时，使用青鳉检测分类器进行检测，找到鱼体后，再利用核相关滤波跟踪算法进行跟踪。
4 结语
本研究中通过设计青鳉游泳硬件平台及试验环境，对获取的视频图像进行高斯滤波等预处理，再利用增强的核相关滤波跟踪算法进行青鳉游泳跟踪，对传统的核相关滤波算法进行了优化升级，提高了青鳉跟踪效果。基于核相关滤波算法的鱼类跟踪系统能够有效提取青鳉的游泳轨迹，并实时获取青鳉游泳位置、速度数据、游泳环境温度值和pH值，同时对不同大小的青鳉在鱼缸中的游泳领地进行了空间分析，为研究青鳉游泳行为提供有力支撑，对探究青鳉在特定环境下的生活习性提供参考。
注：A为传统的核相关滤波算法实现青鳉游泳跟踪过程；B为增强的核相关滤波算法实现青鳉游泳跟踪过程
Note: A, The tracking process based on the kernelized correlation filter tracking algorithm; B, The tracking process based on the enhanced kernelized correlation filter tracking algorithm
图9 青鳉游泳跟踪过程
Fig.9 Swimming tracking process in the medaka
注：A、B、C分别为体长1.8、2.9、3.6 cm青鳉游泳轨迹热力图；D为青鳉游泳领地等值线分布图
Note: A,Swimming trajectory thermodynamic diagram of the medaka (body length 1.8 cm); B, Swimming trajectory thermodynamic diagram of the medaka (body length 2.9 cm); C,Swimming trajectory thermodynamic diagram of the medaka (body length 3.6 cm); D, Territory distribution contour line diagrams in the medaka
图10 青鳉游泳轨迹热力图及领地等值线分布图
Fig.10 Swimming trajectory thermodynamic diagrams and the territory distribution contour line diagrams of the medaka
[1] 冯春雷,李志国,黄洪亮,等.鱼类行为研究在捕捞中的应用[J].大连水产学院学报,2009,24(2):166-170.
[2] 于晓明,张秀梅.鱼类游泳能力测定方法的研究进展[J].南方水产科学,2011,7(4):76-84.
[3] 范良忠,刘鹰,余心杰,等.基于计算机视觉技术的运动鱼检测算法[J].农业工程学报,2011,27(7):226-230.
[4] 张胜茂,张衡,唐峰华,等.计算机视觉技术在监测鱼类游泳行为中的研究进展[J].大连海洋大学学报,2017,32(4):493-500.
[5] 郭静,罗华,张涛.机器视觉与应用[J].电子科技,2014,27(7):185-188.
[6] Kalman R E.A new approach to linear filtering and prediction problems[J].Journal of Basic Engineering,1960,82(1):35-45.
[7] Wang Naiyan,Shi Jianping,Yeung D Y,et al.Understanding and diagnosing visual tracking systems[C]//Proceedings of 2015 IEEE International Conference on Computer Vision.Santiago,Chile:IEEE,2015:3101-3109.
[8] Comaniciu D,Meer P.Meanshift:a robust approach toward feature space analysis[J].IEEE Transactions on Pattern Analysis and Machine Intelligence,2002,24(5):603-619.
[9] Johnson A E,Hebert M.Using spin images for efficient object recognition in cluttered 3D scenes[J].IEEE Transactions on Pattern Analysis and Machine Intelligence,1999,21(5):433-449.
[10] Avidan S.Support vector tracking[J].IEEE Transactions on Pattern Analysis and Machine Intelligence,2004,26(8):1064-1072.
[11] 徐盼麟,韩军,童剑锋.基于单摄像机视频的鱼类三维自动跟踪方法初探[J].水产学报,2012,36(4):623-628.
[12] 黄一凡,陈欣,袁飞.鱼类应激行为作用下的水质视频监测分析系统[J].厦门大学学报:自然科学版,2017,56(4):584-589.
[13] 颜鹏东,谭均军,高柱,等.基于视频跟踪的竖缝式鱼道内鱼类运动行为分析[J].水生生物学报,2018,42(2):250-254.
[14] Henriques J F,Caseiro R,Martins P,et al.High-speed tracking with kernelized correlation filters[J].IEEE Transactions on Pattern Analysis and Machine Intelligence,2015,37(3):583-596.
[15] 李治国,路娜,郭晓芝,等.基于青鳉鱼行为的在线生物监测系统及综合预警平台建设研究[J].科学技术与工程,2017,17(7):288-292,298.
[16] 张融,郑宏远,李录,等.一种计量青鳉鱼胸鳍和尾鳍摆动频率和幅值的计算机视觉算法[J].生态毒理学报,2015,10(4):154-161.
[17] 张胜茂,赵申,张衡,等.基于视频分析的鱼类运动目标提取[J].渔业信息与战略,2017,32(1):44-50.
[18] 伊俊超,刘直芳.基于OpenCV的运动目标检测与跟踪[J].计算机工程与设计,2011,32(8):2817-2820.
[19] Nixon M S,Aguado A S.Feature extraction and image processing[M].Amsterdam:Academic Press,2008.
[20] Haddad R A,Akansu A N.A class of fast Gaussian binomial filters for speech and image processing[J].IEEE Transactions on Signal Processing,1991,39(3):723-727.
[21] Dalal N,Triggs B.Histograms of oriented gradients for human detection[C]//Proceedings of 2005 IEEE Computer Society Conference on Computer Vision and Pattern Recognition.San Diego,CA,USA:IEEE,2005.
[22] Wu Yi,Lim J,Yang M H.Online object tracking:a benchmark[C]//Proceedings of 2013 IEEE Conference on Computer Vision and Pattern Recognition.Portland,OR,USA:IEEE,2013.
[23] 汪廷华,陈峻婷.核函数的选择研究综述[J].计算机工程与设计,2012,33(3):1181-1186.
[24] 禹文豪,艾廷华,杨敏,等.利用核密度与空间自相关进行城市设施兴趣点分布热点探测[J].武汉大学学报:信息科学版,2016,41(2):221-227.
[25] 董浩洋,张东戈,万贻平,等.战场态势热力图构建方法研究[J].指挥控制与仿真,2017,39(5):1-8.
[26] 邓超,宋金伟,孙瑞志,等.基于热力图的卷烟市场数据可视分析系统[J].烟草科技,2016,49(12):91-97.
[27] 邢运龙,李艾华,崔智高,等.改进核相关滤波的运动目标跟踪算法[J].红外与激光工程,2016,45(S1):214-221.
[28] 熊昌镇,赵璐璐,郭芬红.自适应特征融合的核相关滤波跟踪算法[J].计算机辅助设计与图形学学报,2017,29(6):1068-1074.