R中數據缺失值的處理--基於mice包

簡介

缺失值是一個實際數據處理中常見的問題。其缺失機制大致可以分為完全隨機缺失（MCAR），非隨機缺失（MNAR）

• 缺失的例子說明：假設一個數據集有3個變數 X1 , X2 , Y ,假設 X1 , X2 是完全變數, Y存在缺失值，那麼

1. 當 Y 以0.5概率缺失，為MCAR

2. 當 X1<0 或者當 Y<0 , Y 發生缺失，為MNAR

對於完全隨機缺失，直接剔除一般不會帶來偏誤，而對於非完全隨機缺失，由於缺失值出現的位置可能和其他變數有關聯，帶來的後果則相對複雜得多。

常見處理方法

刪除法

1. 直接刪除有缺失值的樣本

2. 刪除存在大面積缺失值的變數。或完全變數分析：如果研究的問題只涉及到全部變數中的一部分變數，這部分變數是完整的，那麼可以只分析這幾個完整變數之間的關係

填補法

• 單變數填補

1. 簡單隨機填補：對於每一個缺失值，從已有的該變數數據中隨機抽樣作為填補值，填補進缺失位置。僅僅考慮到了缺失變數本身，而並沒有考慮到相關變數的信息。因此，信息量的利用少。

2. 均值/中位數/分位數填補：用存在缺失值的變數的已有值的均值/中位數/分位數，作為填補值。這種方法顯然會導致方差偏小。

3. 回歸填補：將缺失變數作為因變數，相關變數（其他變數）作為自變數，進行回歸擬合，用預測值作為填補值。用於作為自變數的變數最好是具有完全數據（無缺失）。

4. 熱平台和冷平台：熱平台法又稱匹配插補法，思路是在完全數據樣本中，找到一個和具有缺失值的樣本相似的完全數據樣本，用完全數據樣本值作為填充值，其過程有點類似於K階近鄰的思想。冷平台法又稱條件均值插補法，思路是先將總體分層（聚類），採用樣本所在層（類）的完全數據的均值來替代缺失值。

• 多變數填補

回歸插補法：對缺失變數和完全數據變數擬合多元回歸模型來預測缺失值。是多重填補法的一種應用。多重填補法（Multiple Imputation Missing Data）的具體技術方法眾多且相對複雜，限於篇幅此處不一一展開。

基於mice包的缺失值填補方法

mice即是基於多重填補法構造的。基本思想是對於一個具有缺失值的變數，用其他變數的數據對這個變數進行擬合，再用擬合的預測值對這個變數的缺失值進行填補。

數據準備（可跳過）

set.seed(2016)ndata <- airqualityndata[sample(nrow(data),7),3] <- NAndata[sample(nrow(data),7),4] <- NAndata <- data[-c(5,6)]n

採用R自帶的airquality數據，其第一和第二列數據已經有很多缺失值，現在我們再人為地在第三個第四列中加入隨機的7個缺失值。由於最後兩列是月份和日期，不適合作為自變數，所以在缺失值填充中先剔除掉。

一個簡單可用的版本

以下代碼起到一個初始化的作用，使用時只需要改data為你自己需要填補缺失值的數據集名稱

library(mice)ninit = mice(data, maxit=0)nmeth = init$methodnpredM = init$predictorMatrixn

如果不想要使全部變數都參與擬合，比如排除掉Ozone變數，加入以下代碼：

predM[, c("Ozone")]=0n

當然，即使你排除了擬合變數，並不意味著簡單的將其排除，他只是不作為擬合變數，但仍然會作為被擬合變數，進行缺失值填充。

如果你想要跳過某個變數,如Temp，不對其填補，加入以下代碼：

meth[c("Temp")]=""n

在這種情況下，雖然這個變數不會被作為被擬合變數進行填充，但仍然會作為擬合變數用於擬合其他變數的缺失值。

對於每個變數的擬合，可以指定所用的擬合方法：

meth[c("Ozone")]="norm"nmeth[c("Solar.R")]="logreg"nmeth[c("Wind")]="polyreg"n

=後面的雙引號內的即為方法的名字，=前面的雙引號中的為所需指定的變數。norm代表貝葉斯線性回歸，logreg代表logit回歸擬合，polyreg代表多項式擬合。

設定完成後執行填充：

imputed = mice(data, method=meth, predictorMatrix=predM, m=5)n

## n## iter imp variablen## 1 1 Ozone Solar.R Wind Tempn## 1 2 Ozone Solar.R Wind Tempn## 1 3 Ozone Solar.R Wind Tempn## 1 4 Ozone Solar.R Wind Tempn## 1 5 Ozone Solar.R Wind Tempn...n

這裡我們選擇讓所有變數都進入擬合，每個變數的方法用默認方法（即不運行前面的設定代碼），注意變數擬合時採用方法需要根據不同變數的情況選擇，否則可能導致效果不佳或者演算法無法計算。這需要事先對數據情況有一定了解。另外作為一個例子，本例僅僅作為代碼演示，並不代表正確的參數選擇方法。

輸出填充結果到imputed：

imputed <- complete(imputed)n

檢查是否存在缺失值

sapply(imputed, function(x) sum(is.na(x)))n

## Ozone Solar.R Wind Temp n## 0 0 0 0n

到這裡就已經可以解決大部分缺失值問題了，如果想繼續深入了解mice包的用法，請看下一節。

進階版本

對分類變數的填充效果一般不是非常好，不建議對分類變數（categorical variables ）做填充，本例中只對連續變數進行填充。

假設數據是MCAR的，一般數據量的5%的樣本存在缺失值是相對安全的比例。

計算數據的缺失率

miss <- function(x){sum(is.na(x))/length(x)*100}napply(data,2,miss)n

## Ozone Solar.R Wind Temp n## 24.183007 4.575163 4.575163 4.575163n

apply(data,1,miss)n

## [1] 25 0 0 0 50 25 0 0 25 25 25 0 0 0 0 0 0 25 0 0 25 25 0n## [24] 0 25 25 50 25 0 0 0 25 50 25 25 25 25 0 25 0 25 25 25 0 25 25n## [47] 0 0 0 0 0 25 25 25 25 25 25 25 50 25 25 0 0 0 25 0 0 0 0n## [70] 0 0 50 0 0 25 0 0 0 0 0 0 0 25 25 0 0 0 0 0 0 25 0n## [93] 0 0 0 25 25 25 0 0 0 25 25 0 0 0 25 0 0 0 0 0 0 0 25n## [116] 0 0 0 25 0 0 0 0 0 0 0 0 25 25 0 0 0 0 0 0 0 25 0n## [139] 0 0 0 0 0 0 0 0 0 0 0 25 0 0 0n

從變數角度看，可以看到Qzone的缺失率高達25%，因此我們可能需要剔除這個變數。

對於樣本，由於只有4個變數，缺失一個就有高達25%的缺失率，權衡一下樣本量，我們可以選擇剔除缺失高達50%和以上的樣本。

使用mice判斷缺失情況

對於上述缺失率等缺失情況的計算，在mice包中還有更簡便的函數：

md.pattern(data)n

## Solar.R Wind Temp Ozone n## 100 1 1 1 1 0n## 32 1 1 1 0 1n## 5 0 1 1 1 1n## 6 1 0 1 1 1n## 5 1 1 0 1 1n## 2 0 1 1 0 2n## 1 1 0 1 0 2n## 2 1 1 0 0 2n## 7 7 7 37 58n

最左邊一列是樣本數，右邊的0-1矩陣的1代表列對應的變數沒有缺失，0代表有缺失。最右邊一列是對應模式累計缺失值個數，最下面一列是變數累積存在的缺失模式個數。

假如 第三行的數據410111，表示有4個樣本滿足這行的缺失模式，1011表示這行的缺失模式是缺失了Solar.R變數，最右邊的1表示這種缺失模式缺失了1個變數數據。

可以用VIM包獲得缺失值的可視化表示

library(VIM)naggr_plot <- aggr(data, col=c(navyblue,red), numbers=TRUE, labels=names(data), cex.axis=.7, gap=3, ylab=c("Histogram of missing data","Pattern"))n

marginplot(data[c(1,2)])n

左邊的紅色箱線圖表示有Ozone缺失的樣本的Solar.R的分布，藍色的箱線圖表示的是剩下的數據點的分布。下方的2個箱線圖的含義類似。

在MCAR的假設下，藍色和紅色的箱線圖應該非常接近。

填充缺失值

這次我們詳細介紹mice函數的用法

tempData <- mice(data,m=5,maxit=50,meth=pmm)n

summary(tempData)n

## Multiply imputed data setn## Call:n## mice(data = data, m = 5, method = "pmm", maxit = 50)n## Number of multiple imputations: 5n## Missing cells per column:n## Ozone Solar.R Wind Temp n## 37 7 7 7 n## Imputation methods:n## Ozone Solar.R Wind Temp n## "pmm" "pmm" "pmm" "pmm" n## VisitSequence:n## Ozone Solar.R Wind Temp n## 1 2 3 4 n## PredictorMatrix:n## Ozone Solar.R Wind Tempn## Ozone 0 1 1 1n## Solar.R 1 0 1 1n## Wind 1 1 0 1n## Temp 1 1 1 0n## Random generator seed value: NAn

mice()函數用於生成填充矩陣tempData，其中的參數有：

• data，需要填充的數據集
• m，多重填補法的填補矩陣數。默認為5
• method，填補用的方法，pmm代表預測均值匹配(predictive mean matching),用 methods(mice) 可以看到有哪些可用的方法
• maxit，迭代次數，默認50次

如果想要查看填充的是那些值可以用以下代碼：

tempData$imp$Tempn

## 1 2 3 4 5n## 1 83 85 84 79 83n## 9 72 77 72 57 80n## 33 68 93 82 68 84n## 41 87 78 81 82 78n## 59 82 82 68 75 58n## 129 82 79 71 76 76n## 137 69 75 56 73 81n

最左邊的一列表示的是被填充樣本序號，之後的5列是多重填補法生成的5個填補矩陣對這個變數產生的填充數。

調用meth值可以查看每個變數用的是什麼方法（前面我們指定了pmm方法）

tempData$methn

## Ozone Solar.R Wind Temp n## "pmm" "pmm" "pmm" "pmm"n

當然也可以如前一節一樣對每個變數指定各自的擬合方法。

最後生成完全數據集：

completedData <- complete(tempData,1)n

1表示用tempData$imp中的第一個矩陣來填充，如果想要用其他矩陣則可以改成2，3，4等等。

Reference

Imputing Missing Data with R; MICE package

Handling missing data with MICE package; a simple approach