【譯文】R語言中的缺失值處理
作者 Selva Prabhakaran
譯者 錢亦欣
在處理一些真實資料時，樣本中往往會包含缺失值（Missing values)。我們需要對缺失值進行適宜的處理，才能建立更為有效的模型，使得後續預測分析能有更小的偏差。本文將羅列不同的缺失值處理方法，並進行具體應用。
資料準備和缺失模式設定
本文使用mlbench包中的BostonHousing資料集作為示例來演示不同的缺失值處理方法。由於原始的資料集並不包含缺失值，我們需要隨機刪除一些資料。通過這種方法，我們不僅可以評估由資料缺失帶來的精度損失，也可以比較不同處理方式的效果好壞。讓我們先載入這個資料集，並隨機刪除一些資料。
# 載入資料集
data ("BostonHousing", package="mlbench")
original <- BostonHousing # backup original data
# 引入缺失值
set.seed(100)
BostonHousing[sample(1:nrow(BostonHousing), 40), "rad"] <- NA
BostonHousing[sample(1:nrow(BostonHousing), 40), "ptratio"]
#> crim zn indus chas nox rm age dis rad tax ptratio b lstat medv
#> 1 0.00632 18 2.31 0 0.538 6.575 65.2 4.0900 1 296 15.3 396.90 4.98 24.0
#> 2 0.02731 0 7.07 0 0.469 6.421 78.9 4.9671 2 242 17.8 396.90 9.14 21.6
#> 3 0.02729 0 7.07 0 0.469 7.185 61.1 4.9671 2 242 17.8 392.83 4.03 34.7
#> 4 0.03237 0 2.18 0 0.458 6.998 45.8 6.0622 3 222 18.7 394.63 2.94 33.4
#> 5 0.06905 0 2.18 0 0.458 7.147 54.2 6.0622 3 222 18.7 396.90 5.33 36.2
#> 6 0.02985 0 2.18 0 0.458 6.430 58.7 6.0622 3 222
缺失值已經生成好了，儘管我們已經知道哪些位置的資料缺失，但還是用mice包中的md.pattern函式快速檢查下。
# 缺失值的模式
library(mice)
md.pattern(BostonHousing) # 返回資料的缺失值的模式
#> crim zn indus chas nox rm age dis tax b lstat medv rad ptratio
#> 431 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0
#> 35 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1
#> 35 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1
#> 5 1 1 1 1 1 1 1 1 1 1 1 1 0 0 2
#> 0 0 0 0 0 0 0 0 0 0 0 0 40 40 80
缺失值處理方法
目前共有四種方法來處理缺失值：
1. 刪除觀測（記錄）
如果你的資料集擁有大量觀測，足以用來建立模型，那你可以把包含缺失值的觀測刪去（或者在建模時選擇不納入這些觀測，如設定na.action=na.omit）。在刪去相應觀測後，請確保：
- 你有足夠的樣本點可以用來建模。
- 沒有引入偏差（譯者注：即認為這些缺失值是隨機產生的，刪除對應觀測後，樣本總體還是一個隨機樣本而非選擇樣本）。
# 例子
lm(medv ~ ptratio + rad, data=BostonHousing, na.action=na.omit)
2.刪除變數（欄位）
如果某個變數包含大量的缺失值，我們可以直接刪除這個變數來保留更多的觀測，除非這個變數對於模型而言特別重要。應用這個方法需要我們在變數的重要性和觀測的數量之間做權衡。
3.用均值、中位數或眾數插值
把缺失值用相應變數的均值、中位數或眾數替換是一種比較粗糙的處理方法。其可行性也要取決於具體情境，如果變數的數值本身波動比較小或者對相應變數的影響較小，使用這種粗略的插值法才可以得到使人滿意的結果。
library(Hmisc)
impute(BostonHousing$ptratio, mean) # 均值替代
impute(BostonHousing$ptratio, median) # 中位數替代
impute(BostonHousing$ptratio, 20) # 用特殊值替代（20）
# 也可以手動插值
BostonHousing$ptratio[is.na(BostonHousing$ptratio)] <- mean(BostonHousing$ptratio, na.rm = T)
讓我們看看均值插值的效果
library(DMwR)
actuals <- original$ptratio[is.na(BostonHousing$ptratio)]
predicteds <- rep(mean(BostonHousing$ptratio, na.rm=T), length(actuals))
regr.eval(actuals, predicteds)
#> mae mse rmse mape
#> 1.62324034 4.19306071 2.04769644 0.09545664
4.用預測值插值
用預測值插值是一種比較前沿的方法，我們有很多模型可以實現這個過程，比如KNN插值，rpart還有mice。
4.1. KNN插值
DMwR包中的knnImputation函式會使用k近鄰方法來填補缺失值。具體流程如下：對於每個需要插值的觀測，先基於歐氏距離找到k個和它最近的觀測。再將這k個近鄰的資料利用距離逆加權得到插補值，最後用該值替代缺失值。
這種方式的優勢在於你只要呼叫一次函式就能把所有缺失值插補好。該函式會把整個資料框作為引數，你不需要做其他設定。但在使用時請不要把響應變數也一併輸入，因為在你對測試集做處理時，你無法用未知的響應變數來插值。
library(DMwR)
knnOutput <- knnImputation(BostonHousing[, !names(BostonHousing) %in% "medv"]) # 使用KNN插值.
anyNA(knnOutput)
#> FALSE
檢驗該方法的精度
actuals <- original$ptratio[is.na(BostonHousing$ptratio)]
predicteds <- knnOutput[is.na(BostonHousing$ptratio), "ptratio"]
regr.eval(actuals, predicteds)
#> mae mse rmse mape
#> 1.00188715 1.97910183 1.40680554 0.05859526
與均值插值相比，mape的值降低了39個百分點。總體還不錯。
4.2 rpart
利用knn插值的侷限在於它對於因子類變數的插補效果可能不盡如人意。這種情況下rpart和mice就提供了更靈活的解決方案。rpart的優勢是你只需要一個未缺失值就可以插補整個樣本。
插值思路是利用rpart（決策樹）替代knn來預測缺失值。對於因子類變數而言，我們在呼叫rpart函式式可以把method設為class（譯者注：即用分類樹)，數值型變數就設定method=anova（迴歸樹）。當然，我們也要避免把響應變數傳入函式。
library(rpart)
class_mod <- rpart(rad ~ . - medv, data=BostonHousing[!is.na(BostonHousing$rad), ], method="class", na.action=na.omit) # 因為rad是因子
anova_mod <- rpart(ptratio ~ . - medv, data=BostonHousing[!is.na(BostonHousing$ptratio), ], method="anova", na.action=na.omit) # ptratio是數值變數
rad_pred <- predict(class_mod, BostonHousing[is.na(BostonHousing$rad), ])
ptratio_pred <- predict(anova_mod, BostonHousing[is.na(BostonHousing$ptratio), ])
ptratio的插補精度
actuals <- original$ptratio[is.na(BostonHousing$ptratio)]
predicteds <- ptratio_pred
regr.eval(actuals, predicteds)
#> mae mse rmse mape
#> 0.71061673 0.99693845 0.99846805 0.04099908
與knn相比，mape值又額外下降了30%，可喜可賀。
rad的插補精度
actuals <- original$rad[is.na(BostonHousing$rad)]
predicteds <- as.numeric(colnames(rad_pred)[apply(rad_pred, 1, which.max)])
mean(actuals != predicteds) # 計算誤分類比率
#> 0.25
僅有25%的缺失值被誤分類，這個結果也不壞。
4.3 mice
mice是鏈式方程多元插值的簡寫（Multivariate Imputation by Chained Equations）。R中有個同名包提供了多種先進的缺失值處理方法。它使用一種頗不常見的方法來進行兩步插值：先利用mice函式建模再用complete函式生成完整資料。mice(df)操作會返回df的多個完整副本，每個副本都對缺失的資料插補了不同的值。complete()函式則會返回這些資料集中的一個（預設）或多個。讓我們看看如何對rad和ptratio兩個變數插值：
library(mice)
miceMod <- mice(BostonHousing[, !names(BostonHousing) %in% "medv"], method="rf") # 基於隨機森林模型進行mice插值
miceOutput <- complete(miceMod) # 生成完整資料
anyNA(miceOutput)
#> FALSE
計算ptratio的插值精度：
actuals <- original$ptratio[is.na(BostonHousing$ptratio)]
predicteds <- miceOutput[is.na(BostonHousing$ptratio), "ptratio"]
regr.eval(actuals, predicteds)
#> mae mse rmse mape
#> 0.36500000 0.78100000 0.88374204 0.02121326
mape值與rpart相比又提升了48個百分點，亦可賽艇。
再看看rad的插值效果：
actuals <- original$rad[is.na(BostonHousing$rad)]
predicteds <- miceOutput[is.na(BostonHousing$rad), "rad"]
mean(actuals != predicteds) # compute misclass error.
#> 0.15
誤分類比率降低到了15%，也就是說40個缺失觀測裡插補錯誤的只有6個。相較於rpart的錯誤率（25%），這是一個了不起的提升。
如果你想了解的更深入，這裡是mice包的手冊和DataScience+上另一篇關於mice包的文章。
儘管通過本文你已經對各類處理方法有了初步瞭解，可這些還不足以幫助你判斷每種方法的優劣。但當你下次處理缺失值的時候，逐一測試這些方法是值得一試的。
如果你有任何問題請在下方留言或訪問我的LinkedIn。
原文刊載於datascience+網站