R语言学习系列16异常值处理.docx

1、R语言学习系列16异常值处理16. 异常值处理目录：1.用箱线图检测异常值2.使用局部异常因子法（ LOF 法）检测异常值三 . 用聚类方法检测异常值四 . 检测时间序列数据中的异常值五. 基于稳健马氏距离检测异常值正文：异常值， 是指测量数据中的随机错误或偏差， 包括错误值或偏离 均值的孤立点值。 在数据处理中， 异常值会极大的影响回归或分类的 效果。为了避免异常值造成的损失， 需要在数据预处理阶段进行异常值 检测。另外，某些情况下，异常值检测也可能是研究的目的，例如， 数据造假的发现、电脑入侵的检测等。一、用箱线图检测异常值在一条数轴上，以数据的上下四分位数（Q1-Q3）为界画一个矩形盒子

2、（中间 50%的数据落在盒） ；在数据的中位数位置画一条线段为中位线；用标记数据的均值；默认延长线不超过盒长的 1.5倍,之外的点认为是异常值（用O标记）。QH 肢小值1 1L !盒须 11 彳 ! !h y I 5 x IQR QI Q2 QI I.5XIQH盒形图的主要应用就是，剔除数据的异常值、判断数据的偏态和尾重。R语言实现，使用函数boxplot.stats（） ，基本格式为：stats, n, conf, out=boxplot.stats（x, coef=1.5, do.co nf=TRUE, do.out=TRUE）其中，x为数值向量（NA、NaN值将被忽略）；coef为盒须的

3、长度 为几倍的IQR （盒长），默认为1.5； do.conf和do.out设置是否输出 conf 禾口 out返回值：stats返回5个元素的向量值，包括盒须最小值、盒最小 值、中位数、盒最大值、盒须最大值；n返回非缺失值的个数；conf 返回中位数的95%置信区间；out返回异常值。单变量异常值检测：set.seed(2016)x-rnorm(100) #生成100个服从N(0,1)的随机数summary(x) #x的汇总信息Min. 1st Qu. Media n Mean 3rd Qu. Max.-2.7910 -0.7173 -0.2662 -0.1131 0.5917 2.1940

4、boxplot.stats(x) #用箱线图检测x中的异常值$stats1 -2.5153136 -0.7326879 -0.2662071 0.5929206 2.1942200 $n1 100$conf1 -0.47565320 -0.05676092$out1 -2.791471boxplot(x) #绘制箱线图多变量异常值检测：x-rnorm(100)yv-rnorm(100) dfv-data.frame(x,y) #用x,y生成两列的数据框head(df)x y10.41452353 0.48522682-0.47471847 0.696768830.06599349 0.1855

5、1394-0.50247778 0.70073355-0.82599859 0.311681060.16698928 0.7604624#寻找x为异常值的坐标位置a-which(x %in% boxplot.stats(x)$out)a1 78 81 92#寻找y为异常值的坐标位置b-which(y %in% boxplot.stats(y)$out)b1 27 37intersect(a,b) #寻找变量x,y都为异常值的坐标位置in teger(O)plot(df) #绘制x, y的散点图p2-union(a,b) #寻找变量x或y为异常值的坐标位置1 78 81 92 27 37poi

6、nts(dfp2,col=red,pch=x,cex=2) # 标记异常值二、使用局部异常因子法（LOF法）检测异常值局部异常因子法（LOF法），是一种基于概率密度函数识别异常 值的算法。LOF算法只对数值型数据有效。算法原理：将一个点的局部密度与其周围的点的密度相比较， 若前者明显的比后者小（LOF值大于1），贝S该点相对于周围的点来说就处于一个相对比较稀疏的区域，这就表明该点是一个异常值R语言实现：使用DMwR或dprep包中的函数lofactor() ，基本格式为：lofactor(data, k)其中，data为数值型数据集；k为用于计算局部异常因子的邻居数量library(DMwR)

7、#只选数值型的前4列iris2-iris,1:4head(iris2)Sepal.Length Sepal.Width Petal.Length Petal.Widthplot(de nsity(out.scores)#LOF值排前5的数据作为异常值，out-order(out.scores,decreasi ng=TRUE)1:5out1 42 107 23 16 99iris2out, #异常值数据Sepal.Length Sepal.Width Petal.Length Petal.Width424.52.31.30.31074.92.54.51.7234.63.61.00.2165.7

8、4.41.50.4995.12.53.01.1density.default(x = out.scores)对鸢尾花数据进行主成分分析，并利用产生的前两个主成分绘制成双标图来显示异常值：*-n row(iris2) # 样本数n1 150labels-1:n #用数字1-n标注labels-out-. #非异常值用.标注biplot(prcomp(iris2),cex=0.8,xlabs=labels)PC1说明：函数prcompO对数据集iris2做主成份分析，biplot() 取主成份分析结果的前两列数据即前两个主成份绘制双标图。上图中， x轴和y轴分别代表第一、二主成份，箭头指向了原始变

9、量名，其中 5个异常值分别用对应的行号标注。也可以通过函数pairs() 绘制散点图矩阵来显示异常值，其中异常值用红色的+标注：pchs-rep(., n) pchsout=+ cols-rep(black ,n) colsout-red pairs(iris2,pch=pchs,col=cols)注：另外，Rlof包中函数lof()可实现相同的功能，并且支持并行 计算和选择不同距离。三、用聚类方法检测异常值通过把数据聚成类，将那些不属于任何一类的数据作为异常值。 比如，使用基于密度的聚类 DBSCAN，如果对象在稠密区域紧密相连，则被分组到一类；那些不会被分到任何一类的对象就是异常值。也可以

10、用k-means算法来检测异常值：将数据分成 k组，通过把 它们分配到最近的聚类中心。然后，计算每个对象到聚类中心的距离 (或相似性)，并选择最大的距离作为异常值。36.850000 3.073684 5.742105 2.071053kmea ns.result$cluster #输出聚类结果1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 230 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 1 1 3 1 1 1 1 159 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1

11、1 1 1 1 3 1 1 1 1 1 1 1 1 188 1 1 1 1 1 1 1 1 1 1 1 1 1 3 1 3 3 3 3 1 3 3 3 3 3 3 1 1 3117 3 3 3 1 3 1 3 1 3 3 1 1 3 3 3 3 3 1 3 3 3 3 1 3 3 3 1 3 3146 3 1 3 3 1#centers返回每个样本对应的聚类中心样本cen ters - kmea ns.result$ce nterskmea ns.result$cluster,#计算每个样本到其聚类中心的距离dista nces-sqrt(rowSums(iris2-ce nters)A2)#

12、找到距离最大的5个样本，认为是异常值out-order(dista nces,decreasi ng=TRUE)1:5 out #异常值的样本号1 99 58 94 61 119 iris2out, #异常值Sepal.Length Sepal.Width Petal.Length Petal.Width995.12.53.01.1584.92.43.31.0945.02.33.31.0615.02.03.51.01197.72.66.92.3#绘制聚类结果plot(iris2,c(Sepal.Le ngth,Sepal.Width),pch=o,col=kmea ns.result$clus

13、ter,cex=0.3)#聚类中心用*标记poin ts(kmea ns.result$ce nters,c(Sepal.Le ngth, Sepal.Width), col=1:3, pch=8, cex=1.5)#异常值用+标记poi nts(iris2out,c(Sepal.Le ngth, Sepal.Width), pch四、检测时间序列数据中的异常值对时间序列数据进行异常值检测，先用函数Stl()进行稳健回归 分解，再识别异常值。函数stl()，基于局部加权回归散点平滑法(LOESS),对时间序 列数据做稳健回归分解，分解为季节性、趋势性、不规则性三部分。3, 9.009778 1

14、23.6572 -0.6670047sts-f$time.series#用红色x标记异常值points(time(sts)out, 0.8*sts,remainderout, pch=x, col=red)par(op)五、基于稳健马氏距离检测异常值检验异常值的基本思路是观察各样本点到样本中心的距离， 若某些样本点的距离太大，就可以判断是异常值。若使用欧氏距离，则具有明显的缺点：将样本不同属性（即各指 标变量）之间的差别等同看待。而马氏距离则不受量纲的影响，并且 在多元条件下，还考虑到了变量之间的相关性。对均值为卩，协方差矩阵为 艺的多变量向量，其马氏距离为（x-卩）艺匕-卩）但是传统的马氏距

15、离检测方法是不稳定的，因为个别异常值会把 均值向量和协方差矩阵向自己方向吸引，这就导致马氏距离起不了检测异常值的所用。解决方法是利用迭代思想构造一个稳健的均值和协方差矩阵估计量，然后计算稳健马氏距离，这样异常值就能正确地被识别出来。用mvoutlier包实现,library(mvoutlier)set.seed(2016)x-cb in d(rnorm(80),rnorm(80)yv-cb in d(rnorm(10,5,1), rnorm(10,5,1) #噪声数据zv-rbi nd(x,y)res1-u ni.plot(z) #维数据的异常值检验#返回outliers 标记各样本是否为异常

16、值，md返回数据的稳健马氏距离 which(res1$outliers=TRUE) #返回异常值的样本号1 81 82 83 84 85 86 87 88 89 90 res2-aq.plot(z) #基于稳健马氏距离的多元异常值检验which(res2$outliers=TRUE) #返回异常值的样本号1 81 82 83 84 85 86 87 88 89 90上图为在一维空间中观察样本数据。9寸CNlO说明：图1-1为原始数据；图1-2的X轴为各样本的稳健马氏距 离排序，Y轴为距离的经验分布，红色曲线为卡方分布，蓝色垂线表 示阀值，在阀值右侧的样本判断为异常值；图 2-1和2-2均是用不

17、同颜色来表示异常值，只是阀值略有不同。若数据的维数过高，则上述距离不再有很大意义(例如基因数据 有几千个变量，数据之间变得稀疏)。此时可以融合主成份降维的思 路来进行异常值检验。mvoutlier包中提供了函数pcout()来对高维数 据进行异常值检验。data(swiss) #使用 swiss 数据集res3-pcout(swiss)#返回wfinal01标记是否为异常值，0表示是 which(res3$wfi nal01=0) #返回异常值的样本号Delem ont Fran ches-Mnt Porre ntruy Broye2 3 6 7Gla ne Gruyere Sarine Veveyse8 9 10 11La Vallee Con they En trem ont Here ns19 31 32 33Martigwy Mon they St Maurice Sierre34 35 36 37Sio n V. De Gen eve38 45注：对于分类数据，一个快速稳定的异常检测的策略是 AVF(Attribute Value Freque ncy)算法。主要参考文献：R语言-异常值处理1-3，银河统计学，博客园.cn blogs./cloudtj/category/780800.html