光电工程  2018, Vol. 45 Issue (8): 180101      DOI: 10.12086/oee.2018.180101     
斑马鱼行为学自动观测装置关键算法研究
莫思特 , 刘天琪 , 何凌     
四川大学电气信息学院,四川 成都 610065
摘要:斑马鱼行为学研究已受到越来越多的关注,其空间坐标计算是斑马鱼行为分析的基础。本文设计了一种基于双目立体视觉技术的斑马鱼行为学自动观测装置。结合p率阈值化和模式阈值化,提出了图像阈值算法。计算斑马鱼图像轮廓像素坐标的平均值,得到单条斑马鱼的水平二维坐标XY。当不同斑马鱼的Y坐标差异较大时,对两台摄像机各自得到的斑马鱼的Y坐标分别按大小排序,同一序号即为同一目标;当Y坐标大小较为接近时,应用最小距离法识别同一斑马鱼。根据折射定理和观测装置结构,推导出斑马鱼三维空间坐标计算公式。算法复杂度分析显示本文提出的算法具有较少的运行时间。通过预设多条模型鱼的位置进行测试,计算结果与预设位置接近,验证了本文提出算法的正确性。
关键词斑马鱼行为学    双目立体视觉    三维坐标    
School of Electrical Engineering and Information, Sichuan University, Chengdu, Sichuan 610065, China
Mo Site, Liu Tianqi, He Ling     
School of Electrical Engineering and Information, Sichuan University, Chengdu, Sichuan 610065, China
Abstract: Studies on the zebrafish behavior have attracted more and more attention. The algorithm of 3D coordinate calculation is the basis of zebrafish behavioral analysis. An automatic device for observing the behavior of zebrafish was designed based on binocular stereo vision technology. According to methods of p-rate threshold and pattern threshold, the image threshold algorithm was proposed. Horizontal X and Y coordinates of single zebrafish were figured out by calculating the average of pixel coordinates of image contour. If the difference of Y coordinate values among different zebrafishes is large, two groups of zebrafishes captured by two cameras will be set up separately according to the Y coordinate values, and the same serial number in two groups means the same target. If the difference of Y coordinate is small, the same target in different cameras was identified by the method of minimum distance. The formula of three-dimensional coordinates was conducted based on the refraction theorem and the structure of the observation device. Algorithm analysis shows that the running time of the algorithm proposed in the paper is saved. Validating testing to a number of fishes is designed and carried out to show that the calculated coordinates are close to the preset locations.
Keywords: zebrafish behavior    binocular stereo vision    three-dimensional coordinate    

1 引言

斑马鱼是一种亚热带淡水鱼类,易于饲养,生殖周期短,繁殖能力强,具有完整器官、系统的脊椎生物,在生理学方面和基因遗传学方面与人类具有很高的同源性[1]。作为一种理想的模式动物,斑马鱼被广泛应用于发育学[2]、遗传学[3]、行为学和分子生物学等研究领域[4-6]。最近,斑马鱼在神经行为学方面越来越多的受到人们的关注。实验研究发现斑马鱼幼体、成体对不同药物有应激反应。已发表应用斑马鱼进行的行为学研究包括社会行为[7]、嗅觉行为[8]、焦虑行为[9]、成瘾行为[10]等等。

对于斑马鱼的行为学研究,需要分析斑马鱼的行为参数,用于表征斑马鱼的行为参数主要有垂直穿梭次数、顶部停留时间、水平移动距离和水平穿梭次数的变化等[11-12],这些行为参数可以通过采集视频信号,计算分析斑马鱼在水中的实时三维立体坐标,再对斑马鱼的实时三维立体坐标进行分析计算得到。

采用视频信号分析三维图像的方法主要是双目立体视觉法,现有的双目立体视觉法需要进行立体匹配和三维重建算法,立体匹配算法复杂,不适合于高分辨率高帧率视频信号的实时计算[13-16]。此外现有三维重建算法都只考虑了空气中的三维成像,水中的目标由于有折射,没法采用现有三维重建算法计算。

本文在现有双目立体三维视觉的基础上,设计了适用于斑马鱼行为学研究需求的自动观测装置,研究了两台摄像机中的目标匹配算法和斑马鱼在水中的三维立体坐标重建算法,分析实时记录的斑马鱼三维坐标可以得到斑马鱼的行为特征参数。

2 斑马鱼行为学自动观测装置组成结构

本文研究的斑马鱼行为学自动观测装置组成结构如图 1所示[17-18]

图 1 斑马鱼行为学自动观测装置组成结构 Fig. 1 Structure of the automatic device for observing the behavior of zebrafish

实验中,鱼缸的内壁被涂成白色。摄像机采集的视频数据通过通信接口输入计算机。与一般的双目立体成像设备不同的是,本文不需计算每个像素的三维坐标,只需将单条斑马鱼看成一点,计算每一条斑马鱼的三维坐标。由于拍摄物在水中,计算三维坐标时需考虑水的折射因素。

3 斑马鱼三维坐标计算

计算三维坐标,需先找出所拍摄的斑马鱼,分别计算各斑马鱼在两台摄像机中的水平二维坐标,匹配两台摄像机中属于同一斑马鱼的水平二维坐标,根据同一斑马鱼在两台摄像机中的水平二维坐标计算斑马鱼在水中的深度,从而得到斑马鱼的三维坐标。

3.1 单条斑马鱼水平二维坐标的确定

计算单条斑马鱼水平二维坐标,需先将图像阈值化,然后分割出每条鱼,计算每条鱼边沿像素的XY方向坐标平均值,将其平均值作为斑马鱼的水平二维坐标。

3.1.1 图像阈值化

图像阈值化的目的是分离背景信息和目标信息,就是将斑马鱼的图像从鱼缸中分离出来。本文采用自适应阈值化的方法,结合p率阈值化和模式阈值化技术的优点计算阈值。图像的主要特点是:背景鱼缸为高亮度,占较大面积;斑马鱼为低亮度,相对背景来说面积很小。根据图像特点确定阈值计算方法如下:计算图像的直方图,对直方图平滑处理,从低亮度开始寻找直方图的第一个频度峰值点,设该峰值点为AA点与最高亮度之间,寻找直方图频度峰值点(如果有多个频度峰值点,选择接近A点的峰值点),设其为B,在A点和B点之间,寻找频度最低点(如果有多个频度最低点,选择接近A点的最低点),设其为CC点的亮度即为阈值。以MATLAB函数为例,图像阈值化算法流程图如图 2所示。

图 2 图像阈值化算法流程图 Fig. 2 Flow chart of image threshold algorithm
3.1.2 斑马鱼平面坐标点计算方法

将图像阈值化后,再提取单条斑马鱼的坐标。单条斑马鱼坐标计算方法如下:先计算图像轮廓,通过图像轮廓将单条斑马鱼分割,然后对分割后的图像轮廓像素计算xy坐标的平均值,用xy的平均值作为对应斑马鱼的坐标。斑马鱼平面坐标点算法流程图如图 3所示。

图 3 斑马鱼平面坐标算法流程图 Fig. 3 Flow chart of zebrafish plane coordinate algorithm
3.2 两个摄像机中同一条斑马鱼的辨别

由于两个摄像机是平行安装,在没有误差的情况下,同一斑马鱼在水平Y方向有相同的坐标。但是由于两台摄像机存在安装误差,且镜头和感光传感器不完全相同,没法绝对保证同一斑马鱼在不同摄像机中的水平Y方向的坐标相同。如果将多条斑马鱼按水平Y方向坐标的大小依次排序的话,同一条斑马鱼在两摄像机中有相同的排序,由此可以判断两摄像机中的同一目标。

如果有多条斑马鱼位于相近的Y坐标,不能用排序的方法来识别同一目标,需根据最小距离的方法进行目标匹配,匹配方法为:对于Y坐标相近的几个目标,以坐标为中心,分别向上、下、左、右扩散a个像素,取边长为2a个像素的正方形区域,采用两台摄像机拍摄的灰度图像,分别计算各相近Y坐标在正方形区域内的距离,距离计算:

$d = \sum\limits_{i = 1}^{2a} {\sum\limits_{j = 1}^{2a} {\left| {A(i, j) - B(i, j)} \right|} } , $ (1)

式中:A(i, j),B(i, j)分别表示相近Y坐标点的目标在所计算正方形区域内的两个摄像机的灰度。根据式(1)计算某一目标点与另一摄像机相近Y坐标各目标点的距离,选择距离最小的点为目标匹配点。

3.3 三维坐标计算方法

图 4所示建立三维笛卡尔坐标系:鱼缸水面正中央为坐标原点O(0, 0, 0);鱼缸正中间的垂线为Z轴,向上为正方向;鱼缸的左右方向为X轴,向右方向为正方向;鱼缸前后方向为Y轴,向前方向为正方向。

图 4 三维坐标示意图 Fig. 4 Diagram of three-dimensional coordinates

摄像机垂直朝下,以Z轴为对称轴对称安装在XOZ平面上,左边摄像机的镜头位于L点,设其坐标为:(-c, 0, h),过L做一条与平面XOY垂直的直线,交X轴于A点;右边摄像机的镜头位于R点,设其坐标为:(c, 0, h),过R做一条与平面XOY垂直的直线,交X轴于B点。c为摄像机在水面的投影到坐标原点的距离,h为摄像机镜头距水面的高度。

3.3.1 根据右边摄像机确定鱼所在的直线

设两个摄像机镜头焦距为f,感光传感器像素尺寸为w×w,横向像素总数为M,纵向像素总数为N。左边摄像机Y坐标和X坐标记为LyLx;右边摄像机Y坐标和X坐标记为RyRx

图 4中连接右边摄像机感光传感器的感光点(Rx, Ry)与镜头的中心点的直线,与Z=0平面有一个交点T,设交点T坐标为(Tx+c, Ty, 0),则TxTy的值为

$\left\{ {\begin{array}{*{20}{l}} {{T_x} = \left( {{R_x} - \frac{M}{2}} \right) \times w \times \frac{h}{f}} \\ {{T_y} = \left( {{R_y} - \frac{N}{2}} \right) \times w \times \frac{h}{f}} \end{array}} \right.。$ (2)

设水中的斑马鱼位于Q点,过T点做直线KT垂直于Z=0平面,设直线QT与直线KT的夹角为θ,根据折射原理,θ的计算公式如下:

$\theta = \arcsin \frac{{\sqrt {T_x^2 + T_y^2} }}{{1.33 \times \sqrt {T_x^2 + T_y^2 + {h^2}} }}。$ (3)

如果直线QT的一个方向向量在Z轴投影为1,在X轴的投影为Px,在y轴的投影为Py,则PxPy的计算公式如下:

${P_x} = \tan \theta \frac{{{T_x}}}{{\sqrt {T_x^2 + T_y^2} }}, $ (4)
${P_y} = \tan \theta \frac{{{T_y}}}{{\sqrt {T_x^2 + T_y^2} }}, $ (5)

直线QT的方程可表示为

$\left\{ \begin{array}{l} x = {T_x} + c + t{P_x} \hfill \\ y = {T_y} + t{P_y} \hfill \\ z = t \hfill \\ \end{array} \right.。$ (6)
3.3.2 根据左边摄像机确定鱼所在的平面

连接左边摄像机感光传感器的感光点(Lx, Ly)与镜头的中心点的直线与Z=0平面有一个交点S,设交点S坐标为(Sx-c, Sy, 0),则SxSy的计算公式为

$\left\{ {\begin{array}{*{20}{l}} {{S_x} = \left( {{L_x} - \frac{M}{2}} \right) \times w \times \frac{h}{f}} \\ {{S_y} = \left( {{L_y} - \frac{N}{2}} \right) \times w \times \frac{h}{f}} \end{array}} \right.。$ (7)

连接点ALS三点,可确定一个平面ALS,平面ALS的法向量计算方法如下:

$\left| {\begin{array}{*{20}{l}} {\mathit{\boldsymbol{i}}}&{\mathit{\boldsymbol{j}}}&{\mathit{\boldsymbol{k}}} \\ 0&0&1 \\ {{S_x}}&{{S_y}}&0 \end{array}} \right|。$ (8)

由式(8)算出,平面ALS的法向量为:(-Sy, Sx, 0)。

所以,平面ALS的方程可表示为

$ - {S_y} \times (x + c - {S_x}) + {S_x}(y - {S_y}) = 0。$ (9)

由式(9)可以算出,平面方程为

$ - {S_y} \times (x + c) + {S_x} \times y = 0。$ (10)
3.3.3 确定鱼的坐标

直线QT与平面ALS相交,其交点为Q点。根据式(6)和式(10),可求出Q点的坐标,记为(x1, y1, z1)。

连接点Q和点S形成直线QS,根据与式(6)的类似推导方法,可推导出QS的方程,同样根据与式(10)类似的推导方法可以推导平面BRT的方程。直线QS与平面BRT相交,其交点也为Q点,根据直线QS方程与平面BRT方程计算得到Q点坐标记为(x2, y2, z2)。

在不存在误差的情况下,计算得到的(x1, y1, z1)与(x2, y2, z2)相等。由于存在误差,实际计算取(x1, y1, z1)和(x2, y2, z2)的平均值。最后的坐标值计算方法如下:

$\left\{ {\begin{array}{*{20}{l}} {x = \frac{{{x_1} + {x_2}}}{2}} \\ {y = \frac{{{y_1} + {y_2}}}{2}} \\ {z = \frac{{{z_1} + {z_2}}}{2}} \end{array}} \right.。$ (11)
4 算法分析

1) 阈值计算

阈值计算采用了直方图来计算,并结合p率阈值化和模式阈值化算法,根据2.1.1节的计算方法,阈值计算算法运行时间为:T(n)=40nn=256。

2) 斑马鱼平面坐标点计算方法

斑马鱼的轮廓计算采用了常规计算方法,在轮廓计算完成后,通过计算斑马鱼轮廓的XY坐标平均值计算斑马鱼的平面坐标点。除了轮廓计算的时间,计算XY坐标算法运行时间为斑马鱼轮廓像素个数,由于轮廓像素很少,因此运算时间很短。

3) 同一斑马鱼的识别

如果采用排序的方法识别,算法运行时间为排序算法的运行时间。由于排序的目标有限,算法时间也很少。对于相近Y坐标需要采用最小距离的方法进行识别,算法运行时间为T(n)=O((2a)2),由于a很小,所以运算时间也很短。

5 实验验证

图 1架设测量装置,将鱼缸四周涂成白色,校正摄像机安装角度和安装位置,鱼缸中加入水,随机在不同的位置放入5条静止模型鱼,分别计算每条模型鱼的坐标,计算值和设定值如表 1所示。从表 1可以看出,最大误差小于12 mm,最大误差率12%,可以满足斑马鱼实验需求。

表 1 多条模型鱼静止坐标实验数据 Table 1 The static coordinate experimental data of model fishs
设置坐标/mm 计算坐标/mm 计算误差/mm
X Y Z X Y Z X Y Z
1号鱼 84 -59 -82 96 -56 -88 12 3 -6
2号鱼 -33 -29 -44 -25 -31 -48 8 -2 -4
3号鱼 -77 64 -111 -77 70 -112 0 6 -1
4号鱼 -29 -96 -28 -31 -96 -33 -2 0 -5
5号鱼 -121 -91 -103 -123 -94 -96 -2 -3 7

采用单条模型鱼,用钓鱼线牵引,从坐标(0, 0, -150)处分别向前、后、左、右、上五个方向,用5 s时间移动100 mm,模拟模型鱼移动行为,并通过计算5 s内坐标轨迹变化,分析模拟斑马鱼行为,分析结果如表 2所示,图中方向用在XOY平面上的投影上与X正向的夹角α和与Z轴夹角β表示,αβ计算方法如下:

表 2 单条模型鱼行为分析实验数据 Table 2 Experimental data of behavior analysis in single model fish
开始坐标/mm 结束坐标/mm 方向 距离/mm 速度/(mm·s-1) 行为
X0 Y0 Z0 X1 Y1 Z1 α β
第1次 -4 3 -140 100 0 -145 -2 -2 104 20.7 向右
第2次 9 -1 -155 -90 3 -155 177 0 98 19.7 向左
第3次 5 5 -145 8 101 -140 91 3 97 19.3 向前
第4次 5 -8 147 5 -106 -151 269 -3 98 19.7 向后
第5次 -1 -6 -159 2 -3 -62 51 87 103 20.7 向上
$\alpha {\rm{ = }}\left\{ {\begin{array}{*{20}{l}} {\begin{array}{*{20}{l}} {\arctan \frac{{{y_1} - {y_0}}}{{{x_1} - {x_0}}}}&{\begin{array}{*{20}{l}} {}&{} \end{array}{x_1} - {x_0} \geqslant 0} \end{array}} \\ {\begin{array}{*{20}{l}} {{\rm{180}} - \arctan \frac{{{y_1} - {y_0}}}{{{x_1} - {x_0}}}}&{{x_1} - {x_0} < 0} \end{array}} \end{array}} \right., $ (12)
$\beta {\rm{ = }}\arcsin \frac{{{z_1} - {z_0}}}{{\sqrt {{{({x_1} - {x_0})}^2} + {{({y_1} - {y_0})}^2} + {{({z_1} - {z_0})}^2}} }}。$ (13)

根据αβ可以判断模型鱼行为。从表 2可以看出,模型鱼行为与预设动作一致。

6 结论及展望

本文设计了斑马鱼行为学自动观测装置,提出了自动观测装置的关键算法,推导了斑马鱼三维坐标计算公式。提出的算法主要包括阈值计算、斑马鱼平面坐标点计算和同一斑马鱼的识别算法,分析该算法的复杂度可看出,该算法具有较少的运行时间。通过实验测试数据,验证了上述算法和三维坐标计算方法的正确性。通过计算相邻时间的坐标变化,进行了模型鱼移动距离、游行速度、转向、上浮、下沉、悬停等行为分析,并通过实验验证了行为分析计算的正确性。

参考文献
[1]
Zon L I, Peterson R T. In vivo drug discovery in the zebrafish[J]. Nature Reviews Drug Discovery, 2005, 4(1): 35-44. DOI:10.1038/nrd1606
[2]
Sun Z H, Jia S J, Meng A M. Zebrafish: swimming in life sciences[J]. Chinese Bulletin of Life Sciences, 2006, 18(5): 431-436.
孙智慧, 贾顺姬, 孟安明. 斑马鱼:在生命科学中畅游[J]. 生命科学, 2006, 18(5): 431-436.
[3]
Wang Y X, Zhong T, Song H Y. Advances in genetics of zebrafish development[J]. Section of Genetics Foreign Medical Sciences, 2004, 27(4): 220-223.
王跃样, 钟涛, 宋后燕. 斑马鱼发育遗传学研究进展[J]. 国外医学遗传学分册, 2004, 27(4): 220-223.
[4]
Linney E, Upchurch L, Donerly S. Zebrafish as a neurotoxicological model[J]. Neurotoxicology and Teratology, 2004, 26(6): 709-718. DOI:10.1016/j.ntt.2004.06.015
[5]
Shin J T, Fishman M C. From zebrafish to human: modular medical models[J]. Annual Review of Genomics and Human Genetics, 2002, 3: 311-340. DOI:10.1146/annurev.genom.3.031402.131506
[6]
Ma J J, He Z X, Shu L P, et al. Expression of HOXC4 gene during early development in zebrafish[J]. Journal of Sichuan University (Medical Science Edition), 2013, 44(3): 371-374.
马健娟, 何志旭, 舒莉萍, 等. HOXC4基因在斑马鱼早期胚胎发育过程中的表达[J]. 四川大学学报(医学版), 2013, 44(3): 371-374.
[7]
Engeszer R E, Ryan M J, Parichy D M. Learned social preference in zebrafish[J]. Current Biology, 2004, 14(10): 881-884. DOI:10.1016/j.cub.2004.04.042
[8]
Vitebsky A, Reyes R, Sanderson M J, et al. Isolation and characterization of the laure olfactory behavioral mutant in the zebrafish, Danio rerio[J]. Developmental Dynamics, 2005, 234(1): 229-242. DOI:10.1002/(ISSN)1097-0177
[9]
Gil Barcellosa L J, Rittera F, Kreutz L C, et al. Whole-body cortisol increases after direct and visual contact with a predator in zebrafish, Danio rerio[J]. Aquaculture, 2007, 272(1-4): 774-778. DOI:10.1016/j.aquaculture.2007.09.002
[10]
Kily L J M, Cowe Y C M, Hussain O, et al. Gene expression changes in a zebrafish model of drug dependency suggest conservation of neuro-adaptation pathways[J]. Journal of Experimental Biology, 2008, 211(10): 1623-1634. DOI:10.1242/jeb.014399
[11]
Fang F, Yu L Z. Zebrafish-A model organism that can be used to study biological medicine anxiolytics[J]. Pharmacology and Clinics of Chinese Materia Medica, 2011, 27(1): 104-106.
方芳, 余林中. 斑马鱼——一种可用于中药抗焦虑药研究的模式生物[J]. 中药药理与临床, 2011, 27(1): 104-106.
[12]
Wang Z, Pu Y Z, Chen Y J, et al. Influence of dichlorvos on zebrafish behavior[J]. Journal of International Pharmaceutical Research, 2013, 40(3): 327-330.
王卓, 蒲韵竹, 陈怡君, 等. 敌敌畏对斑马鱼运动行为的影响[J]. 国际药学研究杂志, 2013, 40(3): 327-330.
[13]
Huang W Y, Xu X M, Wu F Q, et al. Research of underwater binocular vision stereo positioning technology in nuclear condition[J]. Opto-Electronic Engineering, 2016, 43(12): 28-33.
黄文有, 徐向民, 吴凤岐, 等. 核环境水下双目视觉立体定位技术研究[J]. 光电工程, 2016, 43(12): 28-33. DOI:10.3969/j.issn.1003-501X.2016.12.005
[14]
Liu J, Fu W P, Wang W, et al. Image matching based on improved SIFT algorithm[J]. Chinese Journal of Scientific Instrument, 2013, 34(5): 1107-1112.
刘佳, 傅卫平, 王雯, 等. 基于改进SIFT算法的图像匹配[J]. 仪器仪表学报, 2013, 34(5): 1107-1112.
[15]
Xu Y X, Chen F. Real-time stereo visual localization based on multi-frame sequence motion estimation[J]. Opto-Electronic Engineering, 2016, 43(2): 89-94.
许允喜, 陈方. 基于多帧序列运动估计的实时立体视觉定位[J]. 光电工程, 2016, 43(2): 89-94. DOI:10.3969/j.issn.1003-501X.2016.02.015
[16]
Bai T Z, Hou X B. An improved image matching algorithm based on SIFT[J]. Transactions of Beijing Institute of Technology, 2013, 33(6): 622-627.
白廷柱, 侯喜报. 基于SIFT算子的图像匹配算法研究[J]. 北京理工大学学报, 2013, 33(6): 622-627.
[17]
Mo S T, Li B X. Design of analog front end circuit for ICX274 based on AD9923A[J]. Opto-Electronic Engineering, 2009, 36(10): 141-145.
莫思特, 李碧雄. 基于AD9923A的ICX274模拟前端电路设计[J]. 光电工程, 2009, 36(10): 141-145. DOI:10.3969/j.issn.1003-501X.2009.10.027
[18]
Mo S T, Wu Z H. Design of 2.0 mega pixels digital video camera and key technology[J]. Opto-Electronic Engineering, 2009, 36(5): 117-121.
莫思特, 吴志红. 200万像素数字摄像机设计及关键技术研究[J]. 光电工程, 2009, 36(5): 117-121. DOI:10.3969/j.issn.1003-501X.2009.05.022