53-R语言中缺失值处理方法
2020-03-21 本文已影响0人
wonphen
1、概览
缺失值被认为是预测建模的首要障碍，尽管一些机器学习算法声称能够从根本上解决这个问题，但是谁又能知道究竟在“黑盒子”里能解决得多好。
缺失值填补方法的选择，在很大程度上影响了模型的预测能力。一般处理方法是直接删除相关行，但这样并不好，因为会造成信息丢失。
2、Hmice包
Hmice是一个多用途的软件包，可用于数据分析、高级图形、缺失值处理、高级表格制作、模型拟合和诊断(线性回归、 Logit模型和cox回归)等。 该软件包包含的功能范围广泛，它提供了两个强大的函数，用于处理缺失值。分别为 impute ()和 aregImpute ()。
impute()函数使用用户定义的统计方法(中间值，最大值，平均值等)来估算缺失值。 默认是使用中位数。另一方面，aregImpute()允许使用加性回归、自举和预测平均匹配进行填补（additive regression, bootstrapping, and predictive mean matching）。
bootstrapping对替代原始数据的样本拟合了一个柔性可加模型(非参数回归方法) ，并利用非缺失值(自变量)对缺失值(因变量)进行了预测。然后，使用预测均值匹配(缺省值)来估算缺失值。
> library(pacman)
> p_load(Hmisc)
> data(iris)
>
> # 制造10%的缺失值
> iris.mis <- missForest::prodNA(iris, noNA = 0.1)
> sum(is.na(iris.mis))/(dim(iris)[1] * dim(iris)[2])
## [1] 0.1
使用平均值填充：
> iris.mis$imputed.1 <- with(iris.mis, impute(Sepal.Length, mean))
使用随机值填充：
> iris.mis$imputed.2 <- with(iris.mis, impute(Sepal.Length, "random"))
> head(iris.mis,1)
## Sepal.Length Sepal.Width Petal.Length Petal.Width Species imputed.1 imputed.2
## 1 NA 3.5 1.4 0.2 setosa 5.836364 6.0
> summary(iris.mis$Sepal.Length)
## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## 4.300 5.100 5.800 5.836 6.400 7.900 18
同样，还可以使用min,max,median来估算缺失值。
aregImpute ()自动识别变量类型并相应地处理它们：
> impute.arg <- aregImpute(~Sepal.Length + Sepal.Width + Petal.Length + Petal.Width +
+ Species, data = iris.mis, n.impute = 6)
> impute.arg
##
## Multiple Imputation using Bootstrap and PMM
##
## aregImpute(formula = ~Sepal.Length + Sepal.Width + Petal.Length +
## Petal.Width + Species, data = iris.mis, n.impute = 6)
##
## n: 150 p: 5 Imputations: 6 nk: 3
##
## Number of NAs:
## Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## 14 14 19 13 15
##
## type d.f.
## Sepal.Length s 2
## Sepal.Width s 2
## Petal.Length s 2
## Petal.Width s 2
## Species c 2
##
## Transformation of Target Variables Forced to be Linear
##
## R-squares for Predicting Non-Missing Values for Each Variable
## Using Last Imputations of Predictors
## Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## 0.897 0.754 0.980 0.956 0.991
输出显示预测缺失值的 R2值， 数值越高，预测的数值越好。还可以使用以下命令查看估算值：
> impute.arg$imputed$Sepal.Length
## [,1] [,2] [,3] [,4] [,5] [,6]
## 12 4.9 5.0 5.1 4.9 5.1 4.9
## 24 4.7 5.2 5.0 5.3 5.0 5.2
## 34 5.1 4.8 5.7 5.4 5.2 5.6
## 37 4.6 5.0 5.3 4.5 5.1 5.0
## 39 4.4 5.0 4.8 5.0 4.4 4.4
## 56 6.4 5.8 5.9 6.1 5.6 5.7
## 57 6.4 6.3 6.5 6.3 6.4 6.7
## 80 5.8 5.4 5.5 5.2 5.1 5.1
## 90 5.5 5.8 5.6 5.8 5.8 5.8
## 92 6.7 6.3 6.5 6.3 6.0 6.6
## 98 5.7 6.4 6.3 5.6 6.0 6.0
## 115 7.0 6.7 6.5 5.4 5.7 6.5
## 148 6.7 6.3 7.0 6.3 6.4 6.1
## 150 5.7 6.3 6.3 6.8 5.8 6.9
3、mice包
> p_load(mice)
> iris.mis2 <- missForest::prodNA(iris, noNA = 0.1)
> iris.mis2 <- iris.mis2[, -5]
>
> # 查看缺失值分布
> mis.tab <- md.pattern(iris.mis2, rotate.names = T, plot = T)
缺失值分布
> mis.tab
## Sepal.Length Petal.Width Sepal.Width Petal.Length
## 100 1 1 1 1 0
## 12 1 1 1 0 1
## 14 1 1 0 1 1
## 1 1 1 0 0 2
## 10 1 0 1 1 1
## 2 1 0 1 0 2
## 1 1 0 0 1 2
## 6 0 1 1 1 1
## 3 0 1 1 0 2
## 1 0 0 1 1 2
## 10 14 16 18 58
画个好看一点的图：
> mice.plot <- VIM::aggr(iris.mis2, col = c("navyblue", "yellow"),
+ numbers = TRUE, sortVars = TRUE,
+ labels = names(iris.mis), cex.axis = 0.7, gap = 3,
+ ylab = c("Missing data", "Pattern"))
缺失值分布
##
## Variables sorted by number of missings:
## Variable Count
## Petal.Length 0.12000000
## Sepal.Width 0.10666667
## Petal.Width 0.09333333
## Sepal.Length 0.06666667
估算缺失值：
> impute.data <- mice(iris.mis2,
+ # 生成的数据集数量
+ m=5,
+ # 估算缺失值的迭代次数
+ maxit = 50,
+ # 估算方法pmm,logreg,polyreg,polr
+ method = 'pmm',
+ # 随机数种子
+ seed = 500)
> summary(impute.data)
## Class: mids
## Number of multiple imputations: 5
## Imputation methods:
## Sepal.Length Sepal.Width Petal.Length Petal.Width
## "pmm" "pmm" "pmm" "pmm"
## PredictorMatrix:
## Sepal.Length Sepal.Width Petal.Length Petal.Width
## Sepal.Length 0 1 1 1
## Sepal.Width 1 0 1 1
## Petal.Length 1 1 0 1
## Petal.Width 1 1 1 0
pmm：预测均值匹配(PMM)-用于数值变量
logreg： (Logit模型)-二元变量
polyreg(Bayesian polytomous regression)：因子变量(>=2个水平)
polr：Proportional odds model(ordered, >= 2 levels)
查看估算的缺失值
> impute.data$imp$Sepal.Length
## 1 2 3 4 5
## 2 4.7 4.8 4.4 4.4 4.8
## 28 4.8 5.0 5.1 4.9 5.1
## 37 4.7 5.0 5.0 5.4 5.1
## 75 6.0 5.7 5.5 5.5 5.9
## 76 6.0 6.4 6.7 5.8 5.9
## 104 5.4 6.9 6.0 6.4 6.3
## 131 6.7 6.7 6.5 6.5 6.8
## 134 5.4 5.6 6.2 6.4 6.3
## 142 5.8 5.4 6.1 6.1 5.8
## 147 6.3 5.5 5.7 5.9 5.6
由于生成有5个输入数据集，您可以使用 complete ()函数选择任何数据集：
> complete.data <- complete(impute.data, 2)
> dim(iris)
## [1] 150 5
> dim(complete.data)
## [1] 150 4
还可以使用5个数据集构建模型，最后将结果合并：
> fit <- with(data = impute.data, expr = lm(Sepal.Width ~ Sepal.Length + Petal.Width))
> combine <- pool(fit)
> summary(combine)
## term estimate std.error statistic df p.value
## 1 (Intercept) 1.9718135 0.33748673 5.842640 101.9862 6.184242e-08
## 2 Sepal.Length 0.2821694 0.06911890 4.082377 109.4387 8.504625e-05
## 3 Petal.Width -0.4655174 0.07232868 -6.436138 133.8076 2.017213e-09
> # 选取其中一个估算数据集回归
> fit2 <- with(data = complete(impute.data, 3), expr = lm(Sepal.Width ~ Sepal.Length +
+ Petal.Width))
> summary(fit2)
##
## Call:
## lm(formula = Sepal.Width ~ Sepal.Length + Petal.Width)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.01692 -0.21227 0.01469 0.25813 1.00384
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.83945 0.32187 5.715 5.90e-08 ***
## Sepal.Length 0.30724 0.06633 4.632 7.91e-06 ***
## Petal.Width -0.48635 0.07116 -6.835 2.04e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3888 on 147 degrees of freedom
## Multiple R-squared: 0.2505, Adjusted R-squared: 0.2403
## F-statistic: 24.57 on 2 and 147 DF, p-value: 6.228e-10
> # 使用原始数据回归
> fit3 <- with(data = iris, expr = lm(Sepal.Width ~ Sepal.Length + Petal.Width))
> summary(fit3)
##
## Call:
## lm(formula = Sepal.Width ~ Sepal.Length + Petal.Width)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.99563 -0.24690 -0.00503 0.23354 1.01131
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.92632 0.32094 6.002 1.45e-08 ***
## Sepal.Length 0.28929 0.06605 4.380 2.24e-05 ***
## Petal.Width -0.46641 0.07175 -6.501 1.17e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3841 on 147 degrees of freedom
## Multiple R-squared: 0.234, Adjusted R-squared: 0.2236
## F-statistic: 22.46 on 2 and 147 DF, p-value: 3.091e-09
对比一下：
> df <- data.frame(fit.pool = combine$pooled$estimate, fit.2 = coef(fit2), fit = coef(fit3))
>
> df
## fit.pool fit.2 fit
## (Intercept) 1.9718135 1.8394542 1.9263208
## Sepal.Length 0.2821694 0.3072358 0.2892867
## Petal.Width -0.4655174 -0.4863496 -0.4664143
使用生成的6个数据集合并后的回归系数与原始数据的回归系数还是非常接近的。
