完全随机缺失 MCAR
所缺失的数据是完全随机的,缺失发生的概率与观察到的数据无关,这是一种比较理想的情况
随机缺失 MAR
数据的缺失并不是完全随机的,缺失数据发生的概率与所观察到的变量是有关的
非随机缺失 MNAR
缺失数据不仅依赖于其他变量,又依赖于变量本身,这种缺失即为不可忽略的缺失
假设数据缺失是完全随机缺失,但是如果数据缺失太多,也可能造成很大的问题。一般来说,数据的缺失量
小于数据总量的5%是可以接受的,那么如果某个特征或样本的缺失数据量超过5%,那么可以考虑是否需要留
下该特征,我们可以写一个简单的函数来检查某一特征或者某一样本的缺失量是否超过5%
使用 mice 中的 airquality 来进行介绍
library(mice)
## Loading required package: lattice
##
## Attaching package: 'mice'
## The following objects are masked from 'package:base':
##
## cbind, rbind
data <- airquality
pMiss <- function(x) {sum(is.na(x))/length(x)*100}
apply(data,2,pMiss)
## Ozone Solar.R Wind Temp Month Day
## 24.183007 4.575163 0.000000 0.000000 0.000000 0.000000
apply(data,1,pMiss)
## [1] 0.00000 0.00000 0.00000 0.00000 33.33333 16.66667 0.00000
## [8] 0.00000 0.00000 16.66667 16.66667 0.00000 0.00000 0.00000
## [15] 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
## [22] 0.00000 0.00000 0.00000 16.66667 16.66667 33.33333 0.00000
## [29] 0.00000 0.00000 0.00000 16.66667 16.66667 16.66667 16.66667
## [36] 16.66667 16.66667 0.00000 16.66667 0.00000 0.00000 16.66667
## [43] 16.66667 0.00000 16.66667 16.66667 0.00000 0.00000 0.00000
## [50] 0.00000 0.00000 16.66667 16.66667 16.66667 16.66667 16.66667
## [57] 16.66667 16.66667 16.66667 16.66667 16.66667 0.00000 0.00000
## [64] 0.00000 16.66667 0.00000 0.00000 0.00000 0.00000 0.00000
## [71] 0.00000 16.66667 0.00000 0.00000 16.66667 0.00000 0.00000
## [78] 0.00000 0.00000 0.00000 0.00000 0.00000 16.66667 16.66667
## [85] 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
## [92] 0.00000 0.00000 0.00000 0.00000 16.66667 16.66667 16.66667
## [99] 0.00000 0.00000 0.00000 16.66667 16.66667 0.00000 0.00000
## [106] 0.00000 16.66667 0.00000 0.00000 0.00000 0.00000 0.00000
## [113] 0.00000 0.00000 16.66667 0.00000 0.00000 0.00000 16.66667
## [120] 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
## [127] 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
## [134] 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
## [141] 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
## [148] 0.00000 0.00000 16.66667 0.00000 0.00000 0.00000
summary the dataset
summary(data)
## Ozone Solar.R Wind Temp
## Min. : 1.00 Min. : 7.0 Min. : 1.700 Min. :56.00
## 1st Qu.: 18.00 1st Qu.:115.8 1st Qu.: 7.400 1st Qu.:72.00
## Median : 31.50 Median :205.0 Median : 9.700 Median :79.00
## Mean : 42.13 Mean :185.9 Mean : 9.958 Mean :77.88
## 3rd Qu.: 63.25 3rd Qu.:258.8 3rd Qu.:11.500 3rd Qu.:85.00
## Max. :168.00 Max. :334.0 Max. :20.700 Max. :97.00
## NA's :37 NA's :7
## Month Day
## Min. :5.000 Min. : 1.0
## 1st Qu.:6.000 1st Qu.: 8.0
## Median :7.000 Median :16.0
## Mean :6.993 Mean :15.8
## 3rd Qu.:8.000 3rd Qu.:23.0
## Max. :9.000 Max. :31.0
##
Ozone变量缺失最多,有37个缺失
利用mice包查看数据缺失的模式
md.pattern(data)
## Wind Temp Month Day Solar.R Ozone
## 111 1 1 1 1 1 1 0
## 35 1 1 1 1 1 0 1
## 5 1 1 1 1 0 1 1
## 2 1 1 1 1 0 0 2
## 0 0 0 0 7 37 44
从上面的输出结果可以得到,有111个样本是完整的,35个样本缺失Ozone特征值,Solar.R变量有7个缺失值
Ozone有37个缺失值
## Loading required package: colorspace
## Loading required package: grid
## Loading required package: data.table
## Registered S3 methods overwritten by 'car':
## method from
## influence.merMod lme4
## cooks.distance.influence.merMod lme4
## dfbeta.influence.merMod lme4
## dfbetas.influence.merMod lme4
## VIM is ready to use.
## Since version 4.0.0 the GUI is in its own package VIMGUI.
##
## Please use the package to use the new (and old) GUI.
## Suggestions and bug-reports can be submitted at: https://github.com/alexkowa/VIM/issues
##
## Attaching package: 'VIM'
## The following object is masked from 'package:datasets':
##
## sleep
##
## Variables sorted by number of missings:
## Variable Count
## Ozone 0.24183007
## Solar.R 0.04575163
## Wind 0.00000000
## Temp 0.00000000
## Month 0.00000000
## Day 0.00000000
从图片中可以看出70%的样本没有缺失信息,23%的样本缺失Ozone特征
另一种有用的可视化方法
这个函数一次只能显示两个变量的缺失情况,左边红色boxplot代表了缺失Ozone数据项,但是Solar.R数据项
没有缺失的样本分布,而左边的蓝色boxplot则显示了剩余的数据点的分布。类似的,下面的两个图形也是这
么进行描述的
如果我们的数据服从完全随机分布的话,那么红色 boxplot 和蓝色 boxplot应该很相似
tempData <- mice(data,
m = 5, # 产生查补数据集的个数,默认是5个
maxit = 50,
method = "pmm", # 计算插值的方式,这里是用的预测均值
seed = 2019,
printFlag = F)
summary(tempData)
## Class: mids
## Number of multiple imputations: 5
## Imputation methods:
## Ozone Solar.R Wind Temp Month Day
## "pmm" "pmm" "" "" "" ""
## PredictorMatrix:
## Ozone Solar.R Wind Temp Month Day
## Ozone 0 1 1 1 1 1
## Solar.R 1 0 1 1 1 1
## Wind 1 1 0 1 1 1
## Temp 1 1 1 0 1 1
## Month 1 1 1 1 0 1
## Day 1 1 1 1 1 0
查看Ozone的插补情况
tempData$imp$Ozone
## 1 2 3 4 5
## 5 19 8 19 19 14
## 10 30 21 30 23 23
## 25 19 14 6 8 18
## 26 4 37 37 19 13
## 27 24 18 32 23 24
## 32 52 45 45 35 32
## 33 18 16 23 16 7
## 34 1 37 32 18 32
## 35 49 63 20 89 16
## 36 89 35 59 110 49
## 37 30 16 12 23 30
## 39 66 135 61 77 85
## 42 78 78 85 76 78
## 43 76 76 78 73 85
## 45 45 28 18 32 23
## 46 32 46 28 63 52
## 52 39 28 71 35 20
## 53 61 40 168 78 59
## 54 64 29 78 23 47
## 55 108 47 64 48 59
## 56 35 31 46 28 16
## 57 37 29 63 20 29
## 58 65 12 23 44 12
## 59 16 7 39 65 20
## 60 30 21 9 9 23
## 61 82 35 85 73 59
## 65 39 31 23 45 36
## 72 47 28 59 63 45
## 75 89 39 52 37 20
## 83 71 23 35 59 16
## 84 59 63 44 35 29
## 102 135 135 73 77 96
## 103 28 59 36 20 44
## 107 16 11 10 24 7
## 115 20 41 11 30 14
## 119 97 50 79 78 110
## 150 30 11 9 14 9
结果展示了每一个观测样本Ozone变量的插补值
查看每一个变量的插值方法
tempData$method
## Ozone Solar.R Wind Temp Month Day
## "pmm" "pmm" "" "" "" ""
现在,我们就可以对原始数据进行最后操作
completedata <- complete(tempData,1) # 1 表示使用5个插值数据集的第一个
首先,画一个 Ozone和其他变量间的散点图
xyplot(tempData,Ozone~Wind+Temp+Solar.R,pch = 18,cex = 1)
上图中,品红色的点(imputed)与蓝色的点(observed)几乎匹配,预示者插补方法很可靠
另一个有用的图是密度曲线图
densityplot(tempData)
密度图中,每一个插补数据集用品红色表示,实际观测值用蓝色表示
还有一个图用于展示预测的插值在变量中的分布情况
stripplot(tempData,pch = 20, cex = 1.2)
假设在后续的步骤中,你想对数据进行一个线性回归分析,那么应该选取那个插补数据集呢?MICE 包 提供了很方便的方法,先用每一个插补的数据集进行模型建立,最后将几个模型综合评估,得到最优的
结果
with-pool
modelfit <- with(tempData,lm(Temp~Ozone+Solar.R+Wind))
summary(pool(modelfit))
## estimate std.error statistic df p.value
## (Intercept) 72.39012292 2.934926474 24.665055 93.36193 0.000000e+00
## Ozone 0.16548481 0.028554047 5.795494 24.27635 5.406154e-06
## Solar.R 0.01054416 0.008871597 1.188530 17.18176 2.507933e-01
## Wind -0.34603937 0.208260501 -1.661570 115.30217 9.931570e-02