一种基于眼底镜图像的视网膜血管迂曲度计算方法及其应用
【技术领域】
[0001] 本发明设及医学图像处理及应用领域，具体设及一种基于眼底镜图像的视网膜血 管迂曲度计算方法及其应用。
【背景技术】
[0002] 典型的视网膜图像包括视盘结构，血管和神经由此进入视网膜、从视盘发出的血 管分支、视网膜中央的黄斑。正常的视网膜血管走形直或者轻度弯曲，但是在高血流、血管 新生和血管阻塞的病理情况下会扩张和迂曲。非光滑的血管走形被称作迂曲，通常是血管 形态的首要改变。临床上，迂曲度是反映糖尿病性视网膜病变严重程度的一个重要指标。曲 率是形态学分析的一个重要属性，迂曲度通常用血管的曲率变化进行描述，异常的迂曲度 是临床许多疾病的重要标志。曲率的估计对于整体眼底血管形状的分析和迂曲的计算具有 重要意义。在临床眼科工作中，眼科医生常用定性指标来描述该种改变，如轻度迂曲、中度、 重度和极度迂曲。而该种弯曲程度可W计算出数值，对于许多慢性疾病如糖尿病、高血压、 屯、脑血管疾病的疾病的眼底血管弯曲改变定量测量具有重要的意义。
[0003] 当前，在迂曲度的计算方面，Capowski等人通过测量直的血管上相对增加长 度，用一个柔和的曲线表示血管，测量血管积分曲线在一定空间频率的相对长度增加程 度，类似的方法也得到了临床应用。但对于测量相对增加长度的方法而言，仅测量曲 线的延伸而不设及形态学和血液动力学结果，而那些通过积分曲线的方法需要在采样 的坐标系内通过任意平滑范式来平滑噪声（CapwskiJJ,KylstraJA,化ee血anSF.A numericindexbasedonspatialfrequencyforthetortuosityofretinal vesselsanditsapplicationtoplusdiseaseinretinopathyofprematurity[J]. Retina, 1995, 15:490-500.)。Hart等人通过血管曲线的7个积分自动估计迂曲度， 但该种方法在某些肉眼看上去就有迂曲改变的图像上没有能分辨出迂曲度的差异 (HartWE,GoldbaumM,CoteB,etal.Automatedmeasurementofretinalvascular to;rtuosity[J].IntJMedInformatics, 1999, 53 (2-3): 239 - 252.)。Dou曲erty和 Varro等人通过血管中屯、线的二阶微分来计算迂曲度值ou曲ertyG,^hnsonMJ,Wiers MD.Measurementofretinalvasculartortuosityanditsapplicationtoretinal pathologies[J].MedBiolEngComput. 2010 ;48(1): 87-95.)。Grisan等人则提出 了一种 称为曲率相关测量的替代方法，使用曲率的变化来计算（GrisanE，ForacchiaM，Ruggeri A.Anovelmethodfortheautomaticgradingofretinalvesseltortuosity[J]. IE邸hansMedImaging.2008;27(3):310-9.)。这些研究都需要手工分割血管和放 置测量点。Sukkaew等人通过两个连续点么间的两个梯度差的和除臥采样间隔来计算 整个血管的迂曲度系数（SukkaewL，UyyanonvaraB，MaldianovSS，etal.Automatic tortuosity-basedretinopathyofprematurityscreeningsystem.lEICETransInf andSyst，2008，E91-D(12):2868 - 2874.)。另外，还有两种迂曲度测量矩阵用于对S 维结构的血管迂曲度定量分析，并已经被证实在二维和S维血管图像上有效（Johnson MJ,DoughertyG.民obustmeasuresofthree-dimensionalvasculartortuositybased ontheminimumcurvatureofapproximatingpolynomialsplinesfitstothe vesselmid-line.MedEngPhys,2007,29:677-690.Dou組ertyG,JohnsonMJ.Clinical validationofthree-dimensionalturtuositymetricsbasedontheminimum curvatureofapproximatingpolynomialsplines.MedEngPhys,2008,30:190-198.)。 Pal等人则提出一种基于链码的迂曲度计算方法，建立k-曲率的改进算法，主要用于提取 链码上的分叉点pi，k用作曲率计算的Pi的数目，一般根据实验W及不同分辨率设置为 4,6,8, 10, 12。K值越小，噪声越多；k值越大又会低估曲率值。传统的迂曲度计算方式T =L/C，L为曲线长度，C为曲线端点的直线距离。但该种方法不能判断曲线内部的凹凸程 度，因此重新定义迂曲度计算公式为：
[0004]
[0005] 式中nif表示血管中屯、线链码中扭转点的数量，L表示弧的长度，该个算法通过 对血管上每个像素的曲率求和W及考虑了每个血管亚支的扭转点的数量。除W弧的长度 后可W作为迂曲度密度参数，来比较不同长度血管的迂曲度，因此能够避免血管树的分割 结果影口向（PalS,BhowmicP.Estimationofdiscretecurvaturebasedonchaincode pairinganddigitalstrai曲tness[C].I邸ETransICIP, 2009. 1097 - 1100.)。由于正 常视网膜血管在视网膜上走行过程中有一定的生理弧度变化，而上述基于曲率的迂曲度表 征方法没有能考虑血管的正常弧度，对计算的结果临床意义产生一定的影响。
【发明内容】
[0006] 本发明的目的在于提供一种基于眼底镜图像的视网膜血管迂曲度计算方法。
[0007] 本发明结合视网膜血管在图像上可能存在正常生理弧度该一解剖特点，没有采 用传统的曲率计算方法，本发明提出了一种非监督学习的分割方法，并对分割后的视网膜 血管骨架进行了拓扑学分层，发明了一种基于随机过程的方法来对迂曲程度进行建模的方 法，并对该方法的临床结果进行了验证。
[0008] 本发明先获取视网膜血管的分割图像，再对视网膜血管进行拓扑分层，利用我们 发明的迂曲度计算模型计算拓扑层次上的血管迂曲度。
[0009] 本发明的技术方案是；经过绿色通道、灰度化等预处理操作后，利用UDWT方法对 获得的视网膜图像进行增强，然后视网膜灰度图像局部滴的纹理提取，利用FCM方法对视 网膜血管进行分割，对分割后的血管进行骨架化，并计算骨架的拓扑层次，并用我们开发出 迂曲度计算模型对血管骨架进行迂曲度计算。
[0010] 本发明提供了一种基于眼底镜图像的视网膜血管迂曲度计算方法，包括W下步 骤：
[0011] 步骤一：首先使用数字化眼底镜获取筛查人群的眼底图像；
[001引步骤二：利用非下采样的离散小波变换（Undecimateddiscretewavelet transform,UDWT)对图像进行增强；
[0013]步骤S;然后视网膜灰度图像局部滴的纹理提取；
[0014] 步骤四；然后利用模糊C聚类（FuzzyCmeanclustering,FCM)的方法对视网膜 血管进行分割；
[0015] 步骤五：对分割后的血管进行骨架化；
[0016] 步骤六：计算骨架的拓扑层次；
[0017] 步骤走：利用W下迂曲度计算模型对血管骨架每个层次进行迂曲度计算：
[0018] 根据血管曲线的每个像素点的X，Y坐标作为输入值，利用协方差矩阵的迹和曲线 的像素点个数对血管迂曲度进行建模；
[0019] 眼底血管待测曲线的每个像素点X，Y坐标为两个随机变量X和Y，根据协方差（记 为Cov狂，Y))定义；
[0020] Cov化Y) =E{技-E狂)][Y-E(Y) ]}
[00引]其中，E佩二化+X2+. . .Xn)/n;E(Y)=化+Y2+. . .+YJ/n
[0022]
称为rv〇(，Y)的协方差矩阵，该里X，Y分别对应于待测血管曲线上的每 个像素的横纵坐标值；
[002引 迂曲