梁少军 ,张世荣 ,孙澜琼
(1.陆军工程大学 军械士官学校,湖北 武汉 430075;2.武汉大学 电气与自动化学院,湖北 武汉 430072)
随着信息技术的飞速发展,数据产生与获取方式的渠道增多,数据样式呈现多样化,如音视频数据、图像数据、文本数据等.众多研究领域的数据维度也越来越高,如人脸识别[1]、信息检索[2]、气候变化、恒星光谱、基因分布[3]等.高维数据会导致数据计算的高度复杂化,给后期数据处理工作带来巨大挑战,这通常被称为“维数灾难”[4].高维数据包含大量的冗余信息,体现为数据的稀疏性,数据降维是应对高维稀疏数据的有效手段.常用的线性降维方法有主成分分析(principal components analysis,PCA)算法[5]、线性判别(linear discriminant analysis,LDA)算法[6]、多维尺度变换(multi-dimensional scaling,MDS)算法[7]等.在实际应用中大多数数据各维度间呈现非线性相关特征,因此出现了一些非线性降维方法,如拉普拉斯映射(Laplacian eigenmaps,LEIGS)算法[8]、局部线性嵌入(locally linear embedding,LLE)算法[9]、局部切空间排列(local target space alignment,LTSA)算法、等距映射 (isometric mapping,ISOMAP)算法[10]等.其中,ISOMAP算法是一种使用广泛的典型非线性降维算法,该算法基于局部邻域距离计算各数据点之间的全局测地距离,在尽量保持高维数据测地距离与低维数据空间距离对等关系的基础上实现降维.但ISOMAP算法本身缺乏噪声应对能力,噪声会破坏算法的拓扑稳定性影响降维效果[11];局部邻域距离的计算没有方向性,当高维流形存在较大曲率时将导致“短路现象”[12],降维后数据不能很好保持高维拓扑结构.
为解决以上问题,提出了一些改进方法,例如,监督的ISOMAP (supervised isometric mapping,S–ISOMAP)算法[13]、多流形判别ISOMAP (multi-manifold dimensional isometric mapping,MMD–ISOMAP)算法[14]等.但此类方法无法适用于无标签数据集.基于密度缩放因子的ISOMAP(density scaling factor based isometric mapping,D–ISOMAP)算 法[15]虽然降低了对噪声的敏感性,但仍无法应对较大曲率引起的“短路现象”.本文针对现有算法存在的上述问题,提出了一种基于最优密度方向的流形降维方法(optimal density direction based ISOMAP,ODD–ISOMAP),该算法用向量表征高维数据与近邻关系,通过寻找高维数据最优密度方向构建距离缩放矩阵,旨在增强算法应对强的噪声抗干扰的能力并避免高维流形较大曲率导致的“短路现象”.文章将对比ODD–ISOMAP算法与ISOMAP,HLLE,LTSA,LEIGS,LLE和PCA算法的降维性能,验证了算法的有效性.
设原始 高维数据为X=[x1x2···xm]T∈Rm×n,降维后的数据矩阵为Y=[y1y2···ym]T∈Rm×d;其中:m为采样数,n为原始变量维度,d为降维后变量维度.经典的ISOMAP算法包括以下几个核心环节:1) 构建无向邻域距离矩阵;2) 获取测地距离矩阵;3) 获取低维嵌入表示矩阵.
无向邻域距离的计算方式有两种:一种基于距离,将xi与其邻域半径r范围内的数据间距离作为邻域距离;另一种基于数量,将xi与其k个最近邻集合k(xi)数据间的距离作为邻域距离.本文选择后者构建邻域距离矩阵,将xi与X中数据间的无向邻域距离定义为
则X数据间的无向邻域距离矩阵D可表示为
再使用Dijkstra或Floyd最短路径算法即可得到任意两个数据间的最短路径距离,即测地距离dG(xi,j).
获取低维嵌入表示可通过求解以下最优化问题来实现:
使用MDS算法求解该问题,即得到高维数据X的低维嵌入表示矩阵Y.
ODD–ISOMAP算法的基本思想是首先确定每个高维数据xi的最优密度方向,然后获取xi的k个近邻在最优密度方向上的投影.再以投影角度、方向和相对大小修订xi与k个近邻的局部邻域距离,即得到有向邻域距离,从而引导数据沿着流形方向计算测地距离,降低短路风险,增强算法抵制噪声干扰的能力.
高维流形上的任一数据xi都有一个指向局部最高数据密度且处于流形表面的方向,将该方向定义为最优密度方向,并记为.为了确保最优密度方向沿流形方向,同时具有较强的抗噪声干扰能力,本文引入自然邻居(natural neighbor,NN)概念.
先按照经典ISOMAP算法计算高维数据X中任一采样 数据xi的k近 邻集 合k(xi).在k(xi)基 础上xi的自然邻居集合定义为上式中ψ(xi)表示xi的自然邻居集合.由该式可知,数据xj属于ψ(xi)表示xi与xj都在对方的k近邻集合中.
图1展示了三维流形中xi自然邻居选取过程及自然邻居的抗干扰能力,此例中近邻数k=4.k(xi)包含4个数 据xj,xp,xt与xn,xi与xn由于流 形紧密 折叠出现了短路现象.如图1(a)所示,在噪声干扰前xj与xp被选为xi的自然邻居,xn和xt虽然属于xi的最近邻但xt局部范围内有更相近的数据,xn在流形的另一端其近邻数据大概率也不包含xi,故xn和xt均未入选xi的自然邻居.可见自然邻居的选取标准比较严格,这一特征能确保xi筛选出的自然邻居沿局部最大密度方向,由于数据xi的局部最大密度方向一般沿着局部流形的方向分布,故可以基于自然邻居得到的最优密度方 向也会 沿着局 部流形方向.另一方 面,如 图1(b)所示,当部分数据受到噪声干扰发生位移时,自然邻居的筛选结果不容易改变.这是由于自然邻居选择过程需要检查数据双方相互的近邻关系,这一近邻关系是基于近邻位置顺序而非直接距离的,因此具有一定稳定性.所以基于自然邻居得到的最优密度方向也将是稳定的.
图1 自然邻居的抗干扰能力Fig.1 Resistance of noise of natural neighbors
在示意图2中,差向量如红色箭头所示.
图2 差向量Fig.2 Difference vector
获取差向量后,可进一步获得xi对应的最优密度方向向量
差向量 的模标识数据间距离的远近.式(6)中,将差向量的模取平方是为了确保xi的自然邻居集合中距离其越近的数据对最优密度方向的影响越大.在图3所示示例中,红色箭头表示了xi未归一化的最优密度方向,结合前文分析可知xi的最优密度方向将指向xi的局部流形方向并且具有稳定性,下文将利用最优密度方向对xi与其最近邻的位置关系进行缩放,并设计缩放因子的大小与噪声对数据的影响程度相关,进而降低噪声影响,还原数据与其近邻间的原始分布情况.
图3 未归一化的最优密度方向Fig.3 Unnormalized optimal density direction
式(8)中:abs(·)表示取绝对值,|·|表示取向量模,ε为系数调节因子.
图4 K 近邻在最优密度方向上映射Fig.4 Projection ofK -nearest neighbors on optimal density direction
cosγj的符号代表了k(→xi)中向量与的方向关系,这里进一步定义符号函数sg以区分同向与反向数据
其中sgn(·)表示符号函数运算.
则xi相对于xj的密度缩放因子可按照下式计算:
其中:exp(·)表示指数运算,引入密度缩放因子ρ是为了方便控制沿最优密度方向的数据间距离的缩放程度,ρ越大缩放程度越剧烈.
遍历X中所有数据后,将得到所有数据相对于其k个最近邻的密度缩放因子矩阵DM m×n.X中某一数据xi相对于 其k最近邻xj的密度缩放因子与xj相对于xi的密度缩放因子并不相同,故DM为非对称矩阵.D M需要按下式进行对称处理:
使用D M ′对高维数据X数据间的无向邻域距离矩阵D进行缩放,可得有向邻域距离矩阵O M
使用密度缩放因子对高维数据的局部邻域距离进行缩放,会使得与最优密度方向同向的近邻点之间的距离缩小,且投影向量与最优密度方向夹角余弦绝对值越小或与当前高维数据点距离越近收缩程度越大;相反的,与最优密度方向反向的近邻点之间的距离会被放大,且投影向量与最优密度方向夹角余弦绝对值越大或与当前高维数据点距离越远则放大程度越大.图5示例中展示了数据缩放后的距离与图1所示原始距离的对比.
图5 缩放后距离Fig.5 Scaled distance
步骤1输入待处理原始高维数据矩阵X,设定最近邻集合的个数k,设定密度缩放因子ρ及系数调节因子ε(一般取10−4<ε≤ 10−1),设定密度缩放因子矩阵DM m×n为全1矩阵.从X中任选一高维数据xi作为待处理数据.
步骤2按 照式 (4) 计 算xi的自然邻居集 合ψ(xi).若ψ(xi)为空集,则xi为离群点,若实际情况允许可将xi剔除,否则可将k(xi)视为ψ(xi).
步骤3基于ψ(xi)按照式(5)–(6)获取xi最优密度方向向量
步骤4从中任选一向量,计算的差向量,按照式(7)–(8)分别计算差向量与最优密度方向向量夹角的余弦cosγj及在上的相对投 影.重复此步骤,直到获取中所有向量与夹角的余弦及所有向量在上的相对投影集合.
步骤5按照式(9)–(10)计 算xi相对于k(xi)中每一个近邻数据的密度缩放因子.
步骤6将xi相对于k(xi)的密度缩放因子存入矩阵DM对应位置.从X中任选一没有处理过的高维数据,重复步骤2–5,直到X中所有数据被处理完毕.
步骤7按照式(11)对DM对称处理,在获取有向邻域距离矩阵D基础上按照式(12)计算距离缩放后的有向邻域距离矩阵OD.
步骤8基于O D使用Dijkstra或Floyd最短路径算法获取测地距离矩阵,最后使用多维尺度变换算法MDS得到原始高维数据矩阵X降维后数据矩阵Y.
需要说明的是,文中最近邻及自然邻居的求取仍然使用欧式距离.而随着数据集维度的增加该距离计算方式将呈现出一定的弊端.
ODD–ISOMAP算法步骤2中需要计算原空间数据的自然邻居,应用了KD树优化后的自然邻居选择时间复杂度为O(mlogm)[16].步骤3计算所有数据的最优密度方向的复杂度为O(m).步骤4与步骤5中计算空间所有数据与其k个最近邻的差向量、相对最优密度方向夹角余弦和投影以及计算所有数据相对其k个近邻的密度缩放因子的复杂度均为O(km),故此步骤总体时间复杂度为O(4km).步骤6中获取DM矩阵的复杂度为O(1).步骤7与步骤8是经典的ISOMAP算法步骤,该步准确复杂度为O(m3+dm2logm)[17].
故ODD–ISOMAP算法准 确时间 复杂度 为O(m·logm+m+4km+1+1 +dm2logm+m3),考虑到k ≪ m及d ≪ m,则ODD–ISOMAP总体时间复杂度与经典ISOMAP算法的时间复杂度大致相同,为O(m3).
另外,经典ISOMAP算法中获取无向邻域矩阵D(也称为 构建邻域图)需要筛 选数据 的k个最近 邻,ODD–ISOMAP算法中自然邻居的选取只需要在此基础上进一步筛选出自然邻居即可,额外付出的计算单位为 O(m).故相对经典 ISOMAP 算法而言,ODD–ISOMAP 算法中 需要额 外付出 的总体 计算单元为O(2m+4km+1+1).相比O(m3)的总体时间复杂度,ODD–ISOMAP算法付出了较少的额外算力消耗取得了较好的噪声抗干扰能力.
为检验文章 所提算法的算法降 维效 果,本节在Swiss roll,Gauss,Iris,Seeds,Wine,Vertebral column-2c,QCM sensor alcohol共7种 数据集上对ODD–ISOMAP算法进行了测试,并将降维效果与其他6种常见降维算 法HLLE,LTSA,LEIGS,LLE,PCA,ISOMAP进行了对比,实验环境为MATLAB 2019a.
7种 数据 集中 前2种为人工 合成数据集,后5种 来自UCI machine learning实测数据库,去除异常值后各数据集基本参数如表1所示.
为直观展示各算法降维效果,将所有数据降维为2维,并使用 统一参 数设置 的经典K-mediods算法对第2至7种多类别数据集降维后数据进行聚类分析.由于数据集的类别数量已知,因此若降维算法能够保持高维数据的拓扑结构,经各算法降维后的2维数据聚类结果应该与真实数据类标相同,否则说明算法降维效果较 差.基于以 上规律,引入聚 类正确 率指标(clustering accuracy,CA)定量描述聚类结果与数据真实类别相似程度,CA表示为
上式中:yi表示高维数据xi映射后对应的低维数据,[xi,yi]表示高低维数据对的数据量,Cp表示类标,t为数据簇类别数量,m为数据簇数据的总量.CA表示所有类别中映射前后属于同一个类标的高低维数据对数据量总和与数据簇总数据量的比值.CA的值介于0到1之间,取值越大表示聚类结果与高维数据真实类别相似度越高,即算法降维效果越好.
表1 测试数据Table 1 Test datasets
Swiss roll数据为3维人工合成数据,常用来进行流形降维测试,本节为了模拟较大曲率及噪声对数据的影响,对1000组Swiss roll数据的维度3进行了适度压缩,并添加了60%的随机噪声,最终生成的Swiss roll数据如图6(a)所示.使用本文提出的ODD–ISOMAP算法及其他6种降维算法将该数据降为2维,其中各算法近邻数k取最优,各算法降维效果如图6(b)–6(h)所示.
分析图6 可 知,在模拟强噪声和大曲率情况下,HLLE算法降维效果出现了重叠,LTSA,LLE,PCA算法降维效果未能保持 3 维数据流形,LEIGS 与 ISOMAP算法 的降维效果受到了一定影 响,只 有ODD–ISOMAP保持了较好的降维效果.
Gauss数据为3维5簇人工合成数据,各数据簇边缘贴近但不重叠,如图7(a)所示.图7(b)–7(h)展示了各算法在该数据上的降维表现,从效果来 看 ODD–ISOMAP算法能够有效促使各数据簇向簇中心靠拢,进而减少降维后簇间数据重叠度,而其他算法降维后都不可避免的出现了数据重叠.
图6 各算法在Swiss roll数据上降维表现Fig.6 Dimensionality reduction performance of each algorithm on Swiss roll dataset
图7 各算法在Gauss数据上降维表现Fig.7 Dimensionality reduction performance of each algorithm on Gauss dataset
为进一步研究各算法 对噪声的敏感 度差别,对Gauss数据集逐步添加3%至15%的高斯噪声,统计各算法在各噪声幅度下的CA值,如表2所示.图8以曲线形式展示了Gauss数据各算法CA值随噪声幅度变化规律,分析表2与图8可知,随着高斯噪声幅度的增加各算法对应的CA值呈下降趋势,说明噪声确实对降维效果有负面影响,但不同算法受影响程度不同.在不同幅度噪声影响 下ODD–ISOMAP算法的CA值 均高于经典ISOMAP算法.表明在最优密度方向引导下,噪声干扰被密度缩放因子有效修正,在一定程度上保持了数据间的真实位置关系.另外,在噪声幅度不大于15%情况下,只 有LEIGS与ODD–ISOMAP算法的CA值保持在0.8以上,并且总体来看后者优于前者;其他算法均受到了较大影响,其中LLE算法效果最差.从噪声角度来看,只有在添加3%噪声时LEIGS算法降维效果略优于ODD–ISOMAP算法,其他情况下后者效果均优于其他6种算法,说明尽管噪声影响了数据的原始分布,但是多个数据间近邻关系存在一定的稳定性,而利用了此规律的ODD–ISOMAP算法对噪声有较好的抵抗能力.
表2 Gauss数据值随CA噪声变化Table 2 CA values of algorithms vary with noises on Gauss dataset
图8 Gauss数据集各算法CA值随噪声变化曲线Fig.8 CA curves of algorithms vary with noises on Gauss dataset
图9展示了ODD–ISOMAP算法与ISOMAP算法在添加3%高斯噪声的Gauss数据上降维后聚类效果.图中不同形状表示原始高维数据真实类标,不同颜色表示降维后数据聚类类标,红色圈中数据为分类错误数据.从图9 中可看出,相比 ISOMAP 算法,ODD–ISOMAP算法能够有效放大各数据簇间距离,并使各簇数据向簇心收缩,从而减少分类错误,较好保持了高维数据原始样式.
本节在优化近邻数值k、系数调节因子ε和密度缩放因子ρ前提下,对比了 ODD–ISOMAP 算法与 ISOMAP,HLLE,LTSA,LEIGS,LLE,PCA算法在5类 实测数据上的降维表现.表3统计了各算法在实测数据集上的CA值,图10 至图14分别展示了各数据集上CA值最高的前2种算法降维效果.各图中不同形状表示原始高维数据的真实类标,不同颜色表示降维后数据聚类类标,红色圈中数据为错误分类数据.
从图10来看,对Iris数据降维效果最好的为ODD–ISOMAP算法和ISOMAP算法.对比图10(a)与图10(b)的聚类簇1,可以发现ODD–ISOMAP算法能够使真实簇更加聚集.对比两图聚类簇2–3可知,两种算法降维后聚类效果均出现了重叠,但ODD–ISOMAP算法能够使两个不同簇沿着各自簇内数据密集方向收拢,从而减少了错分数据,提高了分类正确率.图11至图13展示了Seeds,Wine和Vertebral column2c数据集的结果,获得了与图10同样的结论.这是由 于ODD–ISOMAP算法能够根据高维数据的原始分布情况沿最优密度方向合理缩放了数据间距离,使得降维后各数据集同一簇内数据更加紧密,不同簇间距离更加分散,簇边缘数据向簇内靠拢,较好保持了高维数据真实簇分类特征,从而有效改善了聚类效果.图14 中虽然ODD–ISOMAP,ISOMAP 和 PCA 算法都能在 QCM Sensor Alcohol数据上取得极高正确率(CA=1),但是前者通过数据缩放功能使得分类效果更加显著.从表3数据角 度来看,在Wine数据集 上ODD–ISOMAP算法和HLLE算法均能得到最高 的CA值,其 他数据集上ODD–ISOMAP算法的CA值为最高,说明本文所提算法降维效果最好.
为进 一步研究各算法对 噪声的敏感度差 别,在Seeds数据集上逐步添加5%至25%的高斯噪声,统计各算法在各噪声幅度下的CA值,如表4所示.图15以曲线形式展示了Seeds 数据集上各算法CA值随噪声幅度变化规律.从表4据图15可以看出,除了在25%强噪声干扰下 ODD–ISOMAP 算法降 维效果 略低于LTSA算法外;在其他情况下ODD–ISOMAP算法都能保持更好的降维效果.
图9 Gauss数据3%噪声下降维后聚类效果Fig.9 Clustering effect after dimensionality reduction of Gauss dataset with 3% noise
图10 Iris数据二维映射Fig.10 Two-dimensional maps of Iris dataset
图11 Seeds数据二维映射Fig.11 Two-dimensional maps of Seeds dataset
图12 Wine数据二维映射Fig.12 Two-dimensional maps of Wine dataset
图13 Vertebral column 2c数据二维映射Fig.13 Two-dimensional maps of Vertebral column 2c dataset
图14 QCM Sensor Alcohol数据二维映射Fig.14 Two-dimensional maps of QCM Sensor Alcohol dataset
表3 各算法在多个实测数据上的CA值Table 3 CA values of algorithms on difference datasets
表4 各算法在不同噪声幅度Seeds数据上的CA值Table 4 CA values of algorithms vary with noises on seeds dataset
图15 Seeds数据集各算法CA值随噪声变化曲线Fig.15 CA curves of algorithms vary with noises on Seeds dataset
由于数据的多个近邻间的关系能够反映数据的本质分布情况,在一定噪声的影响下,数据与其自然邻居间的关系存在一定的稳定性,并且由于自然邻居的获取过程取决于近邻数据间的相互位置关系,因此数据的自然邻居一般沿着高维流形方向分布.本文利用了此种规律引入了最优密度方向概念,提出了具有较好噪声抗干扰能力的ODD–ISOMAP算法.该算法从高维数据点的k个最近邻中筛选出自然邻居集合,进一步得到了各个数据的最优密度方向,此最优密度方向指向流形方向且不易受噪声干扰.之后通过计算高维数据的k个最近邻在最优密度方向上投影的角度、方向和长度,获取各近邻数据相对高维数据距离的密度缩放因子.此密度缩放因子的大小与噪声对数据的影响程度密切相关,因此使用该密度缩放因子对数据间的局部邻域距离进行缩放,能够在一定程度上还原数据与其近邻间的原始分布情况,从而降低噪声对降维过程 的影响.在实验 环节对比了 ODD–ISOMAP算法与ISOMAP,HLLE,LTSA,LEIGS,LLE和PCA算法在2类人工合成数据集和5类实测数据集上的降维表现,并测试了各算法在不同噪声压力下的降维性能,验证了本文所提算法的实用性和有效性.