陈蓓蓓,陈嘉玲,郑万挺
眩晕症是一种因对空间定位障碍而产生的运动性或位置性错觉,患者常常感觉自身或周围事物不停旋转、摇摆[1]。眩晕症的病因多样,主要为前庭病变,其表现为感觉到阵发性的外物或本身的旋转、倾倒感、堕落感,多伴有明显的恶心、呕吐等植物神经症状[2]。患者出现眩晕症时伴有异常的眼震,这种眼球震颤运动可以作为眩晕症诊断要素[3]。医生可以通过眼球运动轨迹图来分析病因。
移动医疗平台将互联网和医疗创造性地结合在一起,克服部分线下医疗服务等待时间长、即时性差、交互性差的缺点。诊断技术和移动技术的结合,使人们跨越时间、空间的障碍,享受更优质的医疗服务[4]。
如果将眩晕症诊断和移动医疗平台相结合,为医患提供一个初步诊断的平台[5],就能提高医生诊断效率,快速排查哪些是需要尽快就诊的患者,对患者来说也可以减少大量时间和精力。
移动远程眩晕症诊断平台的设计主要包括5 大部分:服务器、数据库、供用户操作的平台网页、视频处理程序、自动诊断程序,如图1所示。这5 个部分,由服务器上发布的Web Service实现连接,互相通信。一方面,Web Service 将患者录入的注册信息通信给数据库,进行录入,又通过数据库查询用户信息使用户成功登录,患者可以查询自己的诊断结果,医生可以查询其患者的联系方式。另一方面,患者将视频上传到服务器的文件系统上,Web Service 通信Matlab对获取的视频进行处理,得到瞳孔中心运动轨迹图和自动诊断的结果,得到的轨迹图存放在服务器上,自动诊断的结果录入数据库[6]。
图1 系统的总体设计Fig.1 Overall design of the system
网页端需满足患者注册、登录、录制视频、上传视频、查询自动诊断的需求,以及医生查询其名下病人瞳孔中心运动轨迹图和联系方式的需求[7],如图2所示,患者录制的视频通过Web Service 上传到服务器文件系统中,Web Service 通信Matlab 程序将视频进行处理,得到轨迹图和自动诊断的结果,轨迹图保存在服务器的文件系统中,以病人id 进行命名,自动诊断结果通过Web Service 存入数据库;医生查询轨迹图,Web Service 将服务器中对应患者的轨迹图调出,同时显示患者的联系方式。
图2 网页端流程Fig.2 Flowchart of web page
如图3所示,经Web Service 通信,调用Matlab 程序将上传到服务器的视频先进行解码,读取其中一部分的帧[8],对得到的图片进行瞳孔中心定位,最后得到一张瞳孔中心轨迹图,轨迹图中分别展示了水平方向和垂直方向瞳孔运动变化,轨迹图存放在服务器的文件系统中,以患者的id命名。
图3 视频处理流程Fig.3 Flowchart of video processing
本文设计的自动诊断程序建立在眩晕症临床表现——眼震的基础上。眼震是眼球一种非自主的节律性运动,根据震颤方向主要分为3种:水平性、垂直性、旋转性。本文主要研究水平性眼震以及垂直性眼震。眩晕症的临床分类主要分为以下两种,其临床表现有所差异[9-10],如表1、表2所示。
表1 周围性眩晕临床表现Tab.1 Clinical manifestations of peripheral vertigo
表2 中枢性眩晕临床表现Tab.2 Clinical manifestations of central vertigo
本文设计的眩晕症辅助自动诊断主要针对中枢性眩晕,根据其眩晕室眼球震颤方向不定、幅度大的特点来设计程序。通过网页提供触发眩晕症的图片,收集病人凝视眩晕症触发图片后的瞳孔中心轨迹图,分析水平方向眼震和垂直方向眼震的幅度,因此需要设定合理的幅度阈值,才能筛选出符合条件的瞳孔运动轨迹。
按眼球震颤幅度可分为三类:第一类幅度小,幅度在1 mm 以内;第二类为中等,幅度在1~3 mm 之间;第三类幅度大,幅度在3 mm 以上。根据中枢性眩晕眼球震颤幅度大这一特点,可确定合理阈值。
像素与毫米的转换需要知道DPI(每英寸点数),根据公式:象素数/DPI=英寸数;英寸数×25.4=毫米数。已知图片分辨率的情况下,水平方向和垂直方向眼震幅度转换为像素的公式分别如式(1)和式(2)所示:
遍历横坐标数组X与纵坐标数组Y,将每个数组的最大值最小值相减得到的差除以2,得到两个近似水平和垂直振幅的值Xzf和Yzf,与上面得到的阈值作比较,只要水平方向或垂直方向上的一个幅值大于其对应的阈值,则推断疑似患有眩晕症,两个方向上的幅值都小于阈值,则推断为正常人,如图4所示。
图4 自动诊断流程Fig.4 Flowchart of automatic diagnosis
在SQL serve 2008 中进行数据库的设计,主要设计了两张表来存放医生、患者的信息,如表3、表4所示。
数据库中的医生信息是预先录入的,医生只需登录即可,患者用户需要先进行注册,再进行登录。数据库还需要将自动诊断的结果录入,这种网页端和Matlab程序与数据库的交互都需要Web Service服务,如图5所示。
本文的移动远程眩晕症诊断平台的设计是基于不同平台不同语言进行开发的,因此使用ASP.NET自带的Web Service服务进行不同程序间的通信。在Visual Studio 2019 中创建Web Service 服务,VS 会自动搭建好Web Service 的框架,如图6所示。因此只需要在这个框架中写需要的函数即可,例如图6[WebMethod]中的HelloWorld函数。
表3 医生信息表Tab.3 Information sheet of doctors
表4 病患信息表Tab.4 Information sheet of patients
图5 数据库流程图Fig.5 Flowchart of databases
图6 Web Service框架Fig.6 Web Servive framework
UploadFile 函数:利用这个方法可以上传患者拍摄的眼球运动视频到服务器,并将眼球运动视频根据患者id 进行命名。先定义一个视频储存的文件路径,将患者上传的视频以其id命名存入该路径下。
TransToImage 函数:利用这个方法调用Matlab程序的TransToIamge 函数,将上传到服务器的患者视频进行处理,事先定义一个图片文件路径,得到瞳孔中心运动轨迹图并保存在该文件夹下,自动诊断的结果从Matlab 程序中的Diagnose 函数中得出,通过ExecuteSql函数,赋值给Patient表中的flag。
DownloadImage 函数:利用这个方法,根据选中的患者从服务器的图片文件夹中调取与病人id 相对应的图片,医生可以访问其病人的瞳孔中心运动轨迹图,进行诊断。
ExecuteSql 函数:利用这个方法可以执行sql 语句,对数据库中的数据进行增加,无需返回数据。
GetDataTable 函数:利用这个方法,通过sql 语句对数据库中的数据进行查询,将查询到的数据返回,获得需要的数据。
整个平台核心的功能在于从上传的眼球转动视频中提取出瞳孔中心运动轨迹图,因此在Matlab 中所设计的程序应当完成视频解码、高斯模糊去除图片噪声、灰度投影法确定人眼图像、离散灰度梯度微分法确定二值化阈值、进行二值化、Canny 算子边缘检测、选取最大轮廓、椭圆曲线拟合得到每张图片瞳孔中心,最后用For 循环得到一张完整的瞳孔中心运动轨迹图,如图7所示。
图7 视频处理及图像处理总流程Fig.7 Video processing and image processing
由于处理的图像是由手机摄像头获取的,难免因为手机摄像头感光元件的曝光导致不良光颗粒影响瞳孔中心的定位,为了淡化毛发和瞳孔中反光点的影响,本文利用高斯模糊的特点来减小图像的噪声,减少图像中不需要的细节。
对图像进行高斯模糊相当于加了一个低通滤波器,它的工作原理是将图像与正态分布做卷积,是一种图像模糊滤波器,它利用正态分布对图像中的每个像素进行变换。N维正态分布的方程如式(3)所示:
其中,r为高斯模糊的半径,σ为标准偏差。
本文采用的灰度投影法的原理是通过对被检测的人脸图像分别进行水平和垂直的灰度投影,并得到这两个方向上的像素灰度累加曲线图,根据曲线图上的波峰、波谷和人脸图像的先验知识来确定人眼的位置,具体为从水平灰度投影确定人眼区域的纵坐标范围,垂直灰度投影确定人眼区域的横坐标范围。
垂直灰度投影公式如式(4)所示:
由于本文研究瞳孔中心的定位只需要一只眼睛即可,因此先根据人脸的对称性,取人脸的左半边,对得到的左半边图像做垂直灰度投影,从投影曲线中找出背景、头发与人脸的边界。
水平灰度投影公式如式(5)所示:
对已经过垂直投影提取感兴趣部位的图像进行水平投影,得到水平方向上的像素灰度累加曲线图,从曲线图中波峰波谷找到头发、额头、眉毛、眼睛、脸颊之间的界限。
在获得眼部图像后,需要精确定位瞳孔,常用的二值化方法不能很好地将瞳孔与图像的其他部分分离。根据瞳孔部分像素灰度值与人眼其他部分像素具有非常明显的跳跃这一特征,离散灰度梯度微分法能从人眼图像的灰度直方图中找出这个跳跃边界,如图8所示,将小于这个阈值的部分填充为黑色,大于或等于这个阈值的部分为白色。
图8 人眼图像灰度直方图Fig.8 Gray histogram of human eye image
首先,统计眼部图像所有像素的灰度值,z表示灰度值,p(z)表示该灰度值上统计的像素个数。对{(z,p(z))}进行微分,以s为步长,得到w(z),如式(6)所示:
将w(z)中最小的一项记为w(z′),并记下此时的z=z′。得到阈值value,如式(7)所示:
将得到的value向上取整,得到最后的阈值。
本文采用的Canny 算子通过滞后阈值能够有效检测强边缘和弱边缘的连接,滞后阈值包括高阈值和低阈值,其步骤一般包括对图像去噪,利用一阶微分算子来计算图像M点处梯度幅值A及其方向a,对图像中局部像素抑制其非最大梯度点,利用双阈值得到边界点并连接。
Px(i,j)和Py(i,j)是对M点的x方向和y方向的一阶微分偏导,如式(8)和式(9)所示:
由此可以得出图像中各点的梯度幅值A和方向a,如式(10)和式(11)所示:
得到图像各点的梯度幅值和方向后,抑制非最大梯度点,主要步骤包括:以3×3 的模板在8 个方向上对所有像素按照梯度方向a插值。对于每个点,将邻域中心元素N(i,j) 与沿梯度方向的第二个梯度幅值插值的结果进行比较,如果N(i,j)的值小于梯度方向上的两个插值结果,就要将N(i,j)对应的边缘标志位赋予0值[11]。
最后,通过设定双阈值对经过上述处理的图片进行检测,得到需要的边缘点,将其连接。
经Canny 算子边缘检测后,由于睫毛、瞳孔中微弱光点的影响,得到的图像不止包括瞳孔的边缘,还会有小的不规则轮廓。因此要从中选择需要的瞳孔轮廓,瞳孔轮廓在所有轮廓中具有最大周长。本文引用了Matlab 中的Regionprops 函数,即get the properties of region,是Matlab中用来度量图像区域属性的函数,常用来统计被标记区域的面积分布,显示区域总数。调用这个函数中的“Perimeter”属性,这个属性表示图像各个区域边界地区的周长,即可筛选出需要的瞳孔轮廓。
传统的Hough 变换法找寻瞳孔中心将图像中的点代入圆的公式,如式(12)所示,找到最大的圆,圆心即为瞳孔中心:
本文采用椭圆公式根据已知的图像点主动拟合椭圆,利用椭圆拟合曲线方法得到的结果更加贴合实际情况中因眼睑、睫毛等遮挡致使瞳孔呈现非正圆的形状。
已知椭圆公式:
从椭圆公式可变形得到椭圆中心横坐标Xc,纵坐标Yc,如式(14)、式(15)所示,长轴倾角θ,如式(16)所示,长半轴a、短半轴b,如式(17)、式(18)所示:
原始测得N(N≥5)组数据(xi,yi),i=1,2,3,…,N中,根据椭圆方程通式和最小二乘法原理,求目标函数,如式(19)所示:
求目标函数的最小值来确定A、B、C、D、E5 个参数,令F(A,B,C,D,E)对各个参数偏导均为0,得到如下方程组:
求解此线性方程组可解出A、B、C、D、E,代入上述各公式即可解得拟合的椭圆参数,画出椭圆拟合曲线,得到椭圆中心坐标,即瞳孔中心坐标(Xc,Yc)。
将每一次循环中得到的(Xc,Yc)都分别存入数组X和数组Y中。将得到的瞳孔中心横坐标的数组X和纵坐标的数组Y分别输出成图片,便得到我们需要的瞳孔中心轨迹图。
平台的登录界面,如图9所示,在角色中选择医生或患者,输入姓名和密码[12],还未注册的患者点击登录键下的“还未注册?”进行注册,注册页面如图10所示,输入信息后,按“添加”完成注册,按“放弃”清空所输入的信息,按“返回”回到登录界面[13]。
如图11所示,界面左边是触发眩晕的动态旋转图,右边是录制视频的窗口,点击“打开摄像头”,摄像头打开,但没有开始录制,患者调整好姿势,如图12所示,使自己的脸正对摄像头,尽量保持脸部干净,减少背景入镜。点击“开始录制”,上方秒表开始计时,患者大致录制40 s的视频即可,点击“停止”,点击“保存”,即可上传视频,保存成功后,点击查询,可以知道初步辅助诊断的情况[14]。
图9 登录界面Fig.9 Login interface
图10 患者注册界面Fig.10 Patient registration interface
成功保存后,会跳出“保存成功”的提示消息,点击查询,会跳出诊断结果,如诊断结果为正常的则回复“正常”,诊断结果异常的则回复“疑似患有眩晕症”,视频处理程序异常的会跳出回复“未知”。
医生登陆后,进入查询页面,如图13所示,医生选择需要查询的病人进行查询,以自己的经验对瞳孔中心运动轨迹图进行分析[15],需要进行进一步检测的患者,通过其留下的联系方式进行联系。
3.4.1 高斯模糊经高斯模糊后的图片减轻了光照、毛发对图片产生的干扰,如图14、图15所示。本文肖像图均获得知情同意。
图11 拍摄视频界面Fig.11 Video shooting interface
图12 打开摄像头Fig.12 Turn on the camera
图13 医生查询界面Fig.13 Doctor query interface
图14 原图Fig.14 Original image
图15 高斯模糊后图片Fig.15 Image after Gaussian blurring
3.4.2 粗略定位人眼图像先根据人脸的对称性得到半张人脸,再将得到的图片转化为灰度图片,进行垂直方向的灰度投影,得到曲线图如图16所示,根据这张曲线图的波谷来确定人脸的边界,如图17所示。
图16 垂直灰度投影曲线Fig.16 Vertical grayscale projection curve
图17 经处理的人脸图片Fig.17 Processed face image
对得到的图17进行水平方向灰度投影,得到曲线图,如图18所示,从这条曲线的波峰波谷中辨别出额头、眉毛、人眼、脸颊的位置,得到人眼图像[16],如图19所示。
图18 水平灰度投影曲线Fig.18 Horizontal grayscale projection curve
图19 人眼图像Fig.19 Human eye image
3.4.3 人眼图像二值化得到人眼图像灰度直方图,如图20所示。利用离散灰度梯度微分法,取合理步长s得到我们需要的阈值,这里s=5,得到二值化后的图像,如图21所示。
图20 人眼灰度直方图Fig.20 Gray histogram of human eye image
3.4.4 边缘检测利用Canny 算子进行边缘检测,得到图22。
图21 二值化人眼图像Fig.21 Binarization of human eye image
图22 边缘检测后瞳孔Fig.22 Pupil after edge detection
3.4.5 筛选最大轮廓和椭圆曲线拟合从图22中筛选出最大周长的边缘就是瞳孔的边缘。将边缘的坐标带入椭圆曲线拟合算法,得到拟合的椭圆及其中心,如图23、图24所示。
图23 筛选出瞳孔边缘Fig.23 Outlined pupillary edge
图24 椭圆曲线拟合Fig.24 Fitting of elliptical curves
3.4.6 瞳孔运动曲线利用循环得到完整的瞳孔中心运动轨迹图[17],如图25所示。
图25 瞳孔中心运动轨迹图Fig.25 Motion trails of pupil center
该平台的设计主要满足两类用户,一类是患者,对于患者,他们可以进行注册登录、上传眼球运动视频、查询自动诊断结果[18]。其中上传视频,将视频转化成瞳孔中心运动轨迹是关键,诊断也是基于算出瞳孔中心运动轨迹才能够实现的。为实现这个功能,首先需要将视频分解成帧,取其中一帧进行分析[19]。先对得到的图片进行粗略的人眼定位,高斯模糊,减轻光点和毛发的影响,然后利用人脸的先验知识和灰度投影曲线定位人眼部分,虽然存在精度不够高的问题,但是应用于本算法中,效果良好。再对得到的人眼图像进行瞳孔中心定位,首先取合理阈值进行二值化,本文采用了离散灰度梯度微分法,该方法能区分出瞳孔和人眼其他部位,较其他二值化的方法具有更强的针对性;其次利用Canny算子进行边缘检测,筛选轮廓,将轮廓的坐标带入椭圆曲线拟合算法,得出瞳孔中心坐标,利用循环处理所有帧,得到最后的瞳孔中心运动轨迹图[20]。诊断部分分为两块,一个是医生通过查询该病患的瞳孔中心运动轨迹图,以自己的医学经验进行诊断,联系患者,另一个是依靠获得的瞳孔中心运动轨迹的各个点的坐标,其距离中值的平均值作为眼震幅度,与阈值进行判断,得出初步辅助的自动诊断结果,患者可自己查询。