一.什么是缺失值,NA与NULL的区别
(1)NA表示数据集中的该数据遗失、不存在。在针对具有NA的数据集进行函数操作的时候,该NA不会被直接剔除。如x<-c(1,2,3,NA,4),取mean(x),则结果为NA,如果想去除NA的影响,需要显式告知mean方法,如mean(x,na.rm=T);NA是没有自己的mode的,在vector中,它会“追随”其他数据的类型,比如刚刚的x,mode(x)为numeric,mode(x[4])亦然
(2)NULL表示未知的状态。它不会在计算之中,如x<-c(1,2,3,NULL,4),取mean(x),结果为2.5。NULL是不算数的,length(c(NULL))为0,而length(c(NA))为1。可见NA“占着”位置,它存在着,而NULL没有“占着”位置,或者说,“不知道”有没有真正的数据。
二.处理缺失值的步骤:
一个完整的处理方法通常包含以下几个步骤:
(1) 识别缺失数据;
(2) 检查导致数据缺失的原因;
(3) 删除包含缺失值的实例或用合理的数值代替(插补)缺失值。
缺失数据的分类
统计学家通常将缺失数据分为三类。尽管它们都用概率术语进行描述,但思想都非常直观。我们将用sleep研究中对做梦时长的测量(12种动物有缺失值)来依次阐述三种类型。
library("VIM") data(sleep) dim(sleep) [1] 62 10 head(sleep) BodyWgt BrainWgt NonD Dream Sleep Span Gest Pred Exp Danger 1 6654.000 5712.0 NA NA 3.3 38.6 645 3 5 3 2 1.000 6.6 6.3 2.0 8.3 4.5 42 3 1 3 3 3.385 44.5 NA NA 12.5 14.0 60 1 1 1 4 0.920 5.7 NA NA 16.5 NA 25 5 2 3 5 2547.000 4603.0 2.1 1.8 3.9 69.0 624 3 5 4 6 10.550 179.5 9.1 0.7 9.8 27.0 180 4 4 4
① 完全随机缺失 若某变量的缺失数据与其他任何观测或未观测变量都不相关,则数据为完全随机缺失(MCAR)。若12种动物的做梦时长值缺失不是出于系统原因,那么可以认为数据是MCAR。注意,如果每个有缺失值的变量都是MCAR,那么可以将数据完整的实例看作对更大数据集的一个简单随机抽样。
② 随机缺失 若某变量上的缺失数据与其他观测变量相关,与它自己的未观测值不相关,则数据为随机缺失(MAR)。例如,如果体重较小的动物更可能有做梦时长的缺失值(可能因为较小的动物更难观察),而且该“缺失”与动物的做梦时长无关,那么就可以认为该数据是MAR。此时,一旦控制了体重变量,做梦时长数据的缺失与出现将是随机的。
③ 非随机缺失 若缺失数据不属于MCAR和MAR,则数据为非随机缺失(NMAR) 。例如,做梦时长越短的动物更可能有做梦数据的缺失(可能由于难以测量时长较短的事件),那么可认为数据是NMAR。大部分处理缺失数据的方法都假定数据是MCAR或MAR。此时,你可以忽略缺失数据的生成机制,并且(在替换或删除缺失数据后)可以直接对感兴趣的关系进行建模。
当数据是NMAR时,想对它进行恰当的分析比较困难,你既要对感兴趣的关系进行建模,又要对缺失值的生成机制进行建模。(目前分析NMAR数据的方法有模型选择法和模式混合
法。)
(1) 识别缺失值NA
在R语言中缺失值通常以NA表示,判断是否缺失值的函数是is.na。
另一个常用到的函数是complete.cases,它对数据框进行分析,判断某一观测样本是否完整。下面我们读取VIM包中的sleep数据作为例子,它的样本数为62,变量数为10,由complete.cases函数计算可知完整的样本个数为42。
sum(complete.cases(sleep)) [1] 42 sleep[complete.cases(sleep),] #列出没有缺失值的行 sleep[!complete.cases(sleep),] #列出有一个或多个缺失值的行 BodyWgt BrainWgt NonD Dream Sleep Span Gest Pred Exp Danger 1 6654.000 5712.0 NA NA 3.3 38.6 645 3 5 3 3 3.385 44.5 NA NA 12.5 14.0 60 1 1 1 4 0.920 5.7 NA NA 16.5 NA 25 5 2 3 13 0.550 2.4 7.6 2.7 10.3 NA NA 2 1 2 14 187.100 419.0 NA NA 3.1 40.0 365 5 5 5 19 1.410 17.5 4.8 1.3 6.1 34.0 NA 1 2 1 20 60.000 81.0 12.0 6.1 18.1 7.0 NA 1 1 1 21 529.000 680.0 NA 0.3 NA 28.0 400 5 5 5 24 207.000 406.0 NA NA 12.0 39.3 252 1 4 1 26 36.330 119.5 NA NA 13.0 16.2 63 1 1 1 30 100.000 157.0 NA NA 10.8 22.4 100 1 1 1 31 35.000 56.0 NA NA NA 16.3 33 3 5 4 35 0.122 3.0 8.2 2.4 10.6 NA 30 2 1 1 36 1.350 8.1 8.4 2.8 11.2 NA 45 3 1 3 41 250.000 490.0 NA 1.0 NA 23.6 440 5 5 5 47 4.288 39.2 NA NA 12.5 13.7 63 2 2 2 53 14.830 98.2 NA NA 2.6 17.0 150 5 5 5 55 1.400 12.5 NA NA 11.0 12.7 90 2 2 2 56 0.060 1.0 8.1 2.2 10.3 3.5 NA 3 1 2 62 4.050 17.0 NA NA NA 13.0 38 3 1 1
由于逻辑值TRUE和FALSE分别等价于数值1和0,可用sum()和mean()函数来获取关于缺失数据的有用信息。如:
> sum(is.na(sleep$Dream)) [1] 12 > mean(is.na(sleep$Dream)) [1] 0.19 > mean(!complete.cases(sleep)) [1] 0.32
结果表明变量Dream有12个缺失值, 19%的实例在此变量上有缺失值。另外,数据集中32%的实例包含一个或多个缺失值。
对于识别缺失值,有两点需要牢记。第一, complete.cases()函数仅将NA和NaN识别为缺失值,无穷值(Inf和-Inf)被当作有效值。第二,必须使用与本章中类似的缺失值函数来识别R数据对象中的缺失值。像myvar == NA这样的逻辑比较无法实现。
三、识别缺失数据的模式
存在缺失数据情况下,需进一步判断缺失数据的模式是否随机。
①在数据量不大的情况下可以用列表显示缺失值,在R中是利用mice包中的md.pattern函数来实现的。md.pattern()函数可生成一个以矩阵或数据框形式展示缺失值模式的表格。
library(mice) md.pattern(sleep) BodyWgt BrainWgt Pred Exp Danger Sleep Span Gest Dream NonD 42 1 1 1 1 1 1 1 1 1 1 0 2 1 1 1 1 1 1 0 1 1 1 1 3 1 1 1 1 1 1 1 0 1 1 1 9 1 1 1 1 1 1 1 1 0 0 2 2 1 1 1 1 1 0 1 1 1 0 2 1 1 1 1 1 1 1 0 0 1 1 2 2 1 1 1 1 1 0 1 1 0 0 3 1 1 1 1 1 1 1 0 1 0 0 3 0 0 0 0 0 4 4 4 12 14 38
表中的1和0显示了缺失值模式: 0表示变量的列中有缺失值, 1则表示没有缺失值。第一行表述了“无缺失值”的模式(所有元素都为1) 。第二行表述了“除了Span之外无缺失值”的模式。第一列表示各缺失值模式的实例个数,最后一列表示各模式中有缺失值的变量的个数。此处可以看到,有42个实例没有缺失值,仅2个实例缺失了Span。9个实例同时缺失了NonD和Dream的值。数据集包含了总共(42×0)+(2×1)+…+(1×3)=38个缺失值。最后一行给出了每个变量中缺失值的数目。
②在数据量大的情况下可以用图形显示缺失值,如aggr()、matrixplot()和scattMiss()。这里只说明aggr()函数,也是我觉得比较好理解的。aggr()函数不仅绘制每个变量的缺失值数,还绘制每个变量组合的缺失值数。
library("VIM") aggr(sleep, prop=FALSE, numbers=TRUE)
marginplot()函数可生成一幅散点图,在图形边界展示两个变量的缺失值信息。以做梦时长与哺乳动物妊娠期时长的关系为例,来看下列代码:
marginplot(sleep[c("Gest","Dream")], pch=c(20), col=c("darkgray", "red", "blue"))
图形的主体是Gest和Dream(两变量数据都完整)的散点图。左边界的箱线图展示的是包含(深灰色)与不包含(红色) Gest值的Dream变量分布。注意,在灰度图上红色是更深的阴影。四个红色的点代表缺失了Gest得分的Dream值。在底部边界上, Gest和Dream间的关系反过来了。可以看到,妊娠期和做梦时长呈负相关,缺失妊娠期数据时动物的做梦时长一般更长。两个变量均有缺失值的观测个数在两边界交叉处用蓝色输出(左下角的0)。
四、理解缺失值的由来
如果缺失数据集中在几个相对不太重要的变量上,那么你可以删除这些变量,然后再进行正常的数据分析。如果有一小部分数据(如小于10%)随机分布在整个数据集中(MCAR) ,那么你可以分析数据完整的实例,这样仍可以得到可靠且有效的结果。如果可以假定数据是MCAR或者MAR,那么你可以应用多重插补法来获得有效的结论。如果数据是NMAR,你则需要借助专门的方法,收集新数据,或者加入一个相对更容易、更有收益的行业。
五、处理缺失数据
对于缺失数据通常有三种应付手段:
(1)当缺失数据较少时直接删除相应样本删除缺失数据样本,其前提是缺失数据的比例非常少,而且缺失数据是随机出现的,这样删除缺失数据后对分析结果影响不大。这种方法已过时,一般不建议使用。
(2)对缺失数据进行插补
用变量均值或中位数来代替缺失值,其优点在于不会减少样本信息,处理简单。但是缺点在于当缺失数据不是随机出现时会产成偏误。多重插补法(Multiple imputation):多重插补是通过变量间关系来预测缺失数据,利用蒙特卡罗方法生成多个完整数据集,再对这些数据集分别进行分析,最后对这些分析结果进行汇总处理。可以用mice包实现。这是我们关注的重点。
(3)使用对缺失数据不敏感的分析方法,例如决策树。
多重插补(MI)是一种基于重复模拟的处理缺失值的方法。在面对复杂的缺失值问题时,MI是最常选用的方法,它将从一个包含缺失值的数据集中生成一组完整的数据集(通常是3到10个)。每个模拟数据集中,缺失数据将用蒙特卡洛方法来填补。此时,标准的统计方法便可应用到每个模拟的数据集上,通过组合输出结果给出估计的结果,以及引入缺失值时的置信区间。R中可利用Amelia、mice和mi包来执行这些操作。下面重点介绍mice包(利用链式方程的多元插补)提供的方法。
基于mice包的分析通常符合以下分析过程:
library(mice) imp <- mice(data, m) fit <- with(imp, analysis) pooled <- pool(fit) summary(pooled)
其中,
data是一个包含缺失值的矩阵或数据框。
imp是一个包含m个插补数据集的列表对象,同时还含有完成插补过程的信息。 默认m为5。
analysis是一个表达式对象,用来设定应用于m个插补数据集的统计分析方法。方法包
括做线性回归模型的lm()函数、做广义线性模型的glm()函数、做广义可加模型的
gam(),以及做负二项模型的nbrm()函数。 表达式在函数的括号中, ~的左边是响应变量,
右边是预测变量(用+符号分隔开)。
fit是一个包含m个单独统计分析结果的列表对象。
pooled是一个包含这m个统计分析平均结果的列表对象。
函数mice()首先从一个包含缺失数据的数据框开始,然后返回一个包含多个(默认为5个)完整数据集的对象。每个完整数据集都是通过对原始数据框中的缺失数据进行插补而生成的。由于插补有随机的成分,因此每个完整数据集都略有不同。然后, with()函数可依次对每个完整数据集应用统计模型(如线性模型或广义线性模型),最后, pool()函数将这些单独的分析结果整合为一组结果。最终模型的标准误和p值都将准确地反映出由于缺失值和多重插补而产生的不确定性。
> library(mice) > data(sleep, package="VIM") > imp <- mice(sleep, seed=1234) iter imp variable 1 1 NonD Dream Sleep Span Gest 1 2 NonD Dream Sleep Span Gest 1 3 NonD Dream Sleep Span Gest 1 4 NonD Dream Sleep Span Gest 1 5 NonD Dream Sleep Span Gest 2 1 NonD Dream Sleep Span Gest 2 2 NonD Dream Sleep Span Gest 2 3 NonD Dream Sleep Span Gest 2 4 NonD Dream Sleep Span Gest 2 5 NonD Dream Sleep Span Gest 3 1 NonD Dream Sleep Span Gest 3 2 NonD Dream Sleep Span Gest 3 3 NonD Dream Sleep Span Gest 3 4 NonD Dream Sleep Span Gest 3 5 NonD Dream Sleep Span Gest 4 1 NonD Dream Sleep Span Gest 4 2 NonD Dream Sleep Span Gest 4 3 NonD Dream Sleep Span Gest 4 4 NonD Dream Sleep Span Gest 4 5 NonD Dream Sleep Span Gest 5 1 NonD Dream Sleep Span Gest 5 2 NonD Dream Sleep Span Gest 5 3 NonD Dream Sleep Span Gest 5 4 NonD Dream Sleep Span Gest 5 5 NonD Dream Sleep Span Gest > fit <- with(imp, lm(Dream ~ Span + Gest)) > pooled <- pool(fit) > summary(pooled) est se t df Pr(>|t|) lo 95 hi 95 nmis fmi lambda (Intercept) 2.546199168 0.254689696 9.997260 52.12563 1.021405e-13 2.035156222 3.0572421151 NA 0.08710301 0.05273554 Span -0.004548904 0.012039106 -0.377844 51.94538 7.070861e-01 -0.028707741 0.0196099340 4 0.08860195 0.05417409 Gest -0.003916211 0.001468788 -2.666287 55.55683 1.002562e-02 -0.006859066 -0.0009733567 4 0.05442170 0.02098354 > imp Multiply imputed data set Call: mice(data = sleep, seed = 1234) Number of multiple imputations: 5 Missing cells per column: BodyWgt BrainWgt NonD Dream Sleep Span Gest Pred Exp 0 0 14 12 4 4 4 0 0 Danger 0 Imputation methods: BodyWgt BrainWgt NonD Dream Sleep Span Gest Pred Exp "" "" "pmm" "pmm" "pmm" "pmm" "pmm" "" "" Danger "" VisitSequence: NonD Dream Sleep Span Gest 3 4 5 6 7 PredictorMatrix: BodyWgt BrainWgt NonD Dream Sleep Span Gest Pred Exp Danger BodyWgt 0 0 0 0 0 0 0 0 0 0 BrainWgt 0 0 0 0 0 0 0 0 0 0 NonD 1 1 0 1 1 1 1 1 1 1 Dream 1 1 1 0 1 1 1 1 1 1 Sleep 1 1 1 1 0 1 1 1 1 1 Span 1 1 1 1 1 0 1 1 1 1 Gest 1 1 1 1 1 1 0 1 1 1 Pred 0 0 0 0 0 0 0 0 0 0 Exp 0 0 0 0 0 0 0 0 0 0 Danger 0 0 0 0 0 0 0 0 0 0 Random generator seed value: 1234
从输出结果可以看到,五个数据集同时被创建,预测均值(pmm)匹配法被用来处理每个含缺失数据的变量。 BodyWgt、 BrainWgt、 Pred、 Exp和Danger没有进行插补(" "),因为它们并没有缺失数据。 VisitSequence从左至右展示了插补的变量,从NonD开始,以Gest结束。最后,预测变量矩阵(PredictorMatrix)展示了进行插补过程的含有缺失数据的变量,它们利用了数据集中其他变量的信息。(在矩阵中,行代表插补变量,列代表为插补提供信息的变量, 1和0分别表示使用和未使用。)
imp$imp$Dream #提取imp对象的子成分,观测实际的插补值 1 2 3 4 5 1 1.0 0.5 0.5 0.5 0.3 3 2.6 2.1 1.5 1.8 1.3 4 3.4 3.1 3.4 1.2 3.4 14 0.3 0.5 0.5 0.3 1.2 24 1.8 1.3 3.6 0.9 5.6 26 2.3 3.1 2.0 2.6 2.1 30 1.2 0.3 3.4 2.6 2.3 31 3.4 0.5 0.6 1.0 0.5 47 0.5 1.5 1.5 2.2 3.4 53 0.3 0.5 0.5 0.5 0.6 55 0.5 0.9 2.6 2.7 2.4 62 1.0 2.1 0.5 3.9 3.6
展示了在Dream变量上有缺失值的12种动物的5次插补值。检查该矩阵可以帮助你判断插补值是否合理。若睡眠时长出现了负值,插补将会停止(否则结果将会很糟糕)。
dataset3 <- complete(imp, action=3) #展示了多重插补过程中创建的第三个完整数据集 dead(dataset3) BodyWgt BrainWgt NonD Dream Sleep Span Gest Pred Exp Danger 1 6654.000 5712.00 2.1 0.5 3.3 38.6 645.0 3 5 3 2 1.000 6.60 6.3 2.0 8.3 4.5 42.0 3 1 3 3 3.385 44.50 11.0 1.5 12.5 14.0 60.0 1 1 1 4 0.920 5.70 13.2 3.4 16.5 2.0 25.0 5 2 3 5 2547.000 4603.00 2.1 1.8 3.9 69.0 624.0 3 5 4 6 10.550 179.50 9.1 0.7 9.8 27.0 180.0 4 4 4
六、处理缺失值的其他方法
R还支持其他一些处理缺失值的方法。虽然它们不如之前的方法应用广泛,但在一些专业领域非常有用。