第15章 處理缺失數(shù)據(jù)的高級方法
15.1 處理缺失值的步驟
一個完整的處理方法通常包含以下幾個步驟：
(1)識別缺失數(shù)據(jù)；
(2)檢查導(dǎo)致數(shù)據(jù)缺失的原因书劝；
(3)刪除包含缺失值的實(shí)例或用合理的數(shù)值代替（插補(bǔ)）缺失值温算。
數(shù)據(jù)缺失的分類
統(tǒng)計(jì)學(xué)家通常將缺失數(shù)據(jù)分為三類。它們都用概率術(shù)語進(jìn)行描述，但思想都非常直觀石咬。我
們將用sleep研究中對做夢時長的測量（有12個動物有缺失值）來依次闡述三種類型俏扩。
(1)完全隨機(jī)缺失 若某變量的缺失數(shù)據(jù)與其他任何觀測或未觀測變量都不相關(guān)糜工，則數(shù)據(jù)為完全隨機(jī)缺失（MCAR）。若12個動物的做夢時長值缺失不是由于系統(tǒng)原因录淡，那么可認(rèn)為數(shù)據(jù)是MCAR捌木。注意，如果每個有缺失值的變量都是MCAR嫉戚，那么可以將數(shù)據(jù)完整的實(shí)例看做是對更大數(shù)據(jù)集的一個簡單隨機(jī)抽樣刨裆。
(2)隨機(jī)缺失 若某變量上的缺失數(shù)據(jù)與其他觀測變量相關(guān)澈圈，與它自己的未觀測值不相關(guān)，則數(shù)據(jù)為隨機(jī)缺失（MAR）帆啃。例如瞬女，體重較小的動物更可能有做夢時長的缺失值（可能因?yàn)檩^小的動物較難觀察），“缺失”與動物的做夢時長無關(guān)努潘，那么該數(shù)據(jù)就可以認(rèn)為是MAR诽偷。此時，一旦你控制了體重變量慈俯，做夢時長數(shù)據(jù)的缺失與出現(xiàn)將是隨機(jī)的渤刃。
(3)非隨機(jī)缺失 若缺失數(shù)據(jù)不屬于MCAR或MAR，則數(shù)據(jù)為非隨機(jī)缺失（NMAR）贴膘。例如卖子，做夢時長越短的動物也更可能有做夢數(shù)據(jù)的缺失（可能由于難以測量時長較短的事件），那么數(shù)據(jù)可認(rèn)為是NMAR刑峡。大部分處理缺失數(shù)據(jù)的方法都假定數(shù)據(jù)是MCAR或MAR洋闽。此時，你可以忽略缺失數(shù)據(jù)的生成機(jī)制突梦，并且（在替換或刪除缺失數(shù)據(jù)后）可以直接對感興趣的關(guān)系進(jìn)行建模诫舅。當(dāng)數(shù)據(jù)是NMAR時，想對它進(jìn)行恰當(dāng)?shù)胤治霰容^困難宫患，你既要對感興趣的關(guān)系進(jìn)行建模刊懈，還要對缺失值的生成機(jī)制進(jìn)行建模。（目前分析NMAR數(shù)據(jù)的方法有模型選擇法和模式混合法娃闲。由于NMAR數(shù)據(jù)的分析十分復(fù)雜虚汛，超出了本書的范疇，我們將忽略對它的討論皇帮。）
圖15-1列出了一系列可用來處理不完整數(shù)據(jù)的方法卷哩，以及相應(yīng)的R包。
15.2 識別缺失值
R使用NA（不可得）代表缺失值属拾，NaN（不是一個數(shù)）代表不可能的值将谊。另外，符號Inf和-Inf分別代表正無窮和負(fù)無窮渐白。函數(shù)is.na()尊浓、is.nan()和is.infinite()可分別用來識別缺失值、不可能值和無窮值纯衍。每個返回結(jié)果都是TRUE或FALSE眠砾。
這些函數(shù)返回的對象與其自身參數(shù)的個數(shù)相同。若每個元素的類型檢驗(yàn)通過，則由TRUE替換褒颈，否則用FALSE替換柒巫。
qsz <- c(1, 2, 3, NA) # 定義向量qsz，其中有一個缺失值谷丸。
is.na(qsz) # 識別向量中的缺失值堡掏。
## [1] FALSE FALSE FALSE TRUE
函數(shù)complete.cases()可用來識別矩陣或數(shù)據(jù)框中沒有缺失值的行。若每行都包含完整的實(shí)例刨疼，則返回TRUE的邏輯向量泉唁；若每行有一個或多個缺失值，則返回FALSE揩慕。
data(sleep, package = "VIM") # 加載數(shù)據(jù)集亭畜。
sleep[complete.cases(sleep),] # 列出沒有缺失值的行。
## BodyWgt BrainWgt NonD Dream Sleep Span Gest Pred Exp Danger
## 2 1.000 6.60 6.3 2.0 8.3 4.5 42.0 3 1 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
## 7 0.023 0.30 15.8 3.9 19.7 19.0 35.0 1 1 1
## 8 160.000 169.00 5.2 1.0 6.2 30.4 392.0 4 5 4
## 9 3.300 25.60 10.9 3.6 14.5 28.0 63.0 1 2 1
## 10 52.160 440.00 8.3 1.4 9.7 50.0 230.0 1 1 1
## 11 0.425 6.40 11.0 1.5 12.5 7.0 112.0 5 4 4
## 12 465.000 423.00 3.2 0.7 3.9 30.0 281.0 5 5 5
## 15 0.075 1.20 6.3 2.1 8.4 3.5 42.0 1 1 1
## 16 3.000 25.00 8.6 0.0 8.6 50.0 28.0 2 2 2
## 17 0.785 3.50 6.6 4.1 10.7 6.0 42.0 2 2 2
## 18 0.200 5.00 9.5 1.2 10.7 10.4 120.0 2 2 2
## 22 27.660 115.00 3.3 0.5 3.8 20.0 148.0 5 5 5
## 23 0.120 1.00 11.0 3.4 14.4 3.9 16.0 3 1 2
## 25 85.000 325.00 4.7 1.5 6.2 41.0 310.0 1 3 1
## 27 0.101 4.00 10.4 3.4 13.8 9.0 28.0 5 1 3
## 28 1.040 5.50 7.4 0.8 8.2 7.6 68.0 5 3 4
## 29 521.000 655.00 2.1 0.8 2.9 46.0 336.0 5 5 5
## 32 0.005 0.14 7.7 1.4 9.1 2.6 21.5 5 2 4
## 33 0.010 0.25 17.9 2.0 19.9 24.0 50.0 1 1 1
## 34 62.000 1320.00 6.1 1.9 8.0 100.0 267.0 1 1 1
## 37 0.023 0.40 11.9 1.3 13.2 3.2 19.0 4 1 3
## 38 0.048 0.33 10.8 2.0 12.8 2.0 30.0 4 1 3
## 39 1.700 6.30 13.8 5.6 19.4 5.0 12.0 2 1 1
## 40 3.500 10.80 14.3 3.1 17.4 6.5 120.0 2 1 1
## 42 0.480 15.50 15.2 1.8 17.0 12.0 140.0 2 2 2
## 43 10.000 115.00 10.0 0.9 10.9 20.2 170.0 4 4 4
## 44 1.620 11.40 11.9 1.8 13.7 13.0 17.0 2 1 2
## 45 192.000 180.00 6.5 1.9 8.4 27.0 115.0 4 4 4
## 46 2.500 12.10 7.5 0.9 8.4 18.0 31.0 5 5 5
## 48 0.280 1.90 10.6 2.6 13.2 4.7 21.0 3 1 3
## 49 4.235 50.40 7.4 2.4 9.8 9.8 52.0 1 1 1
## 50 6.800 179.00 8.4 1.2 9.6 29.0 164.0 2 3 2
## 51 0.750 12.30 5.7 0.9 6.6 7.0 225.0 2 2 2
## 52 3.600 21.00 4.9 0.5 5.4 6.0 225.0 3 2 3
## 54 55.500 175.00 3.2 0.6 3.8 20.0 151.0 5 5 5
## 57 0.900 2.60 11.0 2.3 13.3 4.5 60.0 2 1 2
## 58 2.000 12.30 4.9 0.5 5.4 7.5 200.0 3 1 3
## 59 0.104 2.50 13.2 2.6 15.8 2.3 46.0 3 2 2
## 60 4.190 58.00 9.7 0.6 10.3 24.0 210.0 4 3 4
## 61 3.500 3.90 12.8 6.6 19.4 3.0 14.0 2 1 1
sleep[!complete.cases(sleep),] # 列出有1個或多個缺失值的行迎卤。
## 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分別等價于數(shù)值1和0拴鸵，可用sum()和mean()函數(shù)來獲取關(guān)于缺失數(shù)據(jù)的有用信息。
sum(is.na(sleep$Dream)) # 統(tǒng)計(jì)數(shù)據(jù)集sleep中Dream變量缺失值個數(shù)蜗搔。
## [1] 12
mean(is.na(sleep$Dream)) # 統(tǒng)計(jì)數(shù)據(jù)集sleep中Dream變量缺失值平均值劲藐。
## [1] 0.1935484
mean(!complete.cases(sleep)) # 統(tǒng)計(jì)數(shù)據(jù)集sleep中Dream變量多少比例含一個或多個缺失值。
## [1] 0.3225806
15.3 探索缺失值模式
15.3.1 列表顯示缺失值
mice包中的md.pattern()函數(shù)可生成一個以矩陣或數(shù)據(jù)框形式展示缺失值模式的表格.
library(mice) # 調(diào)用mice包樟凄。
data(sleep, package = "VIM") # 調(diào)用數(shù)據(jù)集sleep聘芜。
md.pattern(sleep) # 探索數(shù)據(jù)集sleep缺失值。
## BodyWgt BrainWgt Pred Exp Danger Sleep Span Gest Dream NonD
## 42 1 1 1 1 1 1 1 1 1 1 0
## 9 1 1 1 1 1 1 1 1 0 0 2
## 3 1 1 1 1 1 1 1 0 1 1 1
## 2 1 1 1 1 1 1 0 1 1 1 1
## 1 1 1 1 1 1 1 0 1 0 0 3
## 1 1 1 1 1 1 1 0 0 1 1 2
## 2 1 1 1 1 1 0 1 1 1 0 2
## 2 1 1 1 1 1 0 1 1 0 0 3
## 0 0 0 0 0 4 4 4 12 14 38
表中1和0顯示了缺失值模式缝龄，0表示變量的列中有缺失值汰现，1則表示沒有缺失值。
15.3.2 圖形探究缺失數(shù)據(jù)
aggr()函數(shù)不僅繪制每個變量的缺失值數(shù)叔壤，還繪制每個變量組合的缺失值數(shù)瞎饲。
library(VIM) # 調(diào)用VIM包。
aggr(sleep, prop = FALSE, numbers = TRUE) # 圖形展示缺失值百新。選項(xiàng)prop=TRUE用比例作為y軸企软，選項(xiàng)numbers = FALSE（默認(rèn)）刪去數(shù)值型標(biāo)簽庐扫。
matrixplot()函數(shù)可生成展示每個實(shí)例數(shù)據(jù)的圖形饭望。數(shù)值型數(shù)據(jù)被重新轉(zhuǎn)換到[0,1]區(qū)間，并用灰度來表示大行瓮ァ：淺色表示值小铅辞，深色表示值大。默認(rèn)缺失值為紅色萨醒。
matrixplot(sleep) # 展示每個實(shí)例圖形斟珊。
marginplot()函數(shù)可生成一幅散點(diǎn)圖，在圖形邊界展示兩個變量的缺失值信息富纸。
marginplot(sleep[c("Gest", "Dream")], pch = c(20), col = c("darkgray", "red", "blue")) # 缺失值散點(diǎn)圖囤踩。
圖形的主體是Gest和Dream（兩變量數(shù)據(jù)都完整）的散點(diǎn)圖旨椒。左邊界的箱線圖展示的是包含（深灰色）與不包含（紅色）Gest值的Dream變量分布。注意堵漱，在灰度圖上紅色是更深的陰影综慎。
四個紅色的點(diǎn)代表著缺失了Gest得分的Dream值。在底部邊界上勤庐，Gest和Dream間的關(guān)系反過來了示惊。可以看到愉镰，妊娠期和做夢時長呈負(fù)相關(guān)米罚，缺失妊娠期數(shù)據(jù)時動物的做夢時長一般更長。兩個變量均有缺失值的觀測個數(shù)在兩邊界交叉處（左下角）用藍(lán)色輸出丈探。
15.3.3 用相關(guān)性探索缺失值
你可用指示變量替代數(shù)據(jù)集中的數(shù)據(jù)（1表示缺失录择，0表示存在），這樣生成的矩陣有時稱作影子矩陣类嗤。求這些指示變量間和它們與初始（可觀測）變量間的相關(guān)性糊肠，有助于觀察哪些變量常一起缺失，以及分析變量“缺失”與其他變量間的關(guān)系遗锣。
x <- as.data.frame(abs(is.na(sleep))) # 將sleep數(shù)據(jù)集中缺失值賦值货裹，若sleep的元素缺失，則數(shù)據(jù)框x對應(yīng)的元素為1精偿，否則為0弧圆。
head(sleep) # 查看數(shù)據(jù)集sleep前6行。
## 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
head(x) # 查看缺失值賦值后數(shù)據(jù)集的前6行笔咽。
## BodyWgt BrainWgt NonD Dream Sleep Span Gest Pred Exp Danger
## 1 0 0 1 1 0 0 0 0 0 0
## 2 0 0 0 0 0 0 0 0 0 0
## 3 0 0 1 1 0 0 0 0 0 0
## 4 0 0 1 1 0 1 0 0 0 0
## 5 0 0 0 0 0 0 0 0 0 0
## 6 0 0 0 0 0 0 0 0 0 0
y <- x[which(apply(x,2,sum)>0)] # 提取含缺失值的變量搔预。
cor(y) # 指示變量的相關(guān)系數(shù)。
## NonD Dream Sleep Span Gest
## NonD 1.00000000 0.90711474 0.48626454 0.01519577 -0.14182716
## Dream 0.90711474 1.00000000 0.20370138 0.03752394 -0.12865350
## Sleep 0.48626454 0.20370138 1.00000000 -0.06896552 -0.06896552
## Span 0.01519577 0.03752394 -0.06896552 1.00000000 0.19827586
## Gest -0.14182716 -0.12865350 -0.06896552 0.19827586 1.00000000
cor(sleep, y, use = "pairwise.complete.obs") # 缺失值變量與其他觀測變量的關(guān)系叶组。矩陣中行為可觀測變量拯田，列為表示缺失的指示變量。
## NonD Dream Sleep Span Gest
## BodyWgt 0.22682614 0.22259108 0.001684992 -0.05831706 -0.05396818
## BrainWgt 0.17945923 0.16321105 0.007859438 -0.07921370 -0.07332961
## NonD NA NA NA -0.04314514 -0.04553485
## Dream -0.18895206 NA -0.188952059 0.11699247 0.22774685
## Sleep -0.08023157 -0.08023157 NA 0.09638044 0.03976464
## Span 0.08336361 0.05981377 0.005238852 NA -0.06527277
## Gest 0.20239201 0.05140232 0.159701523 -0.17495305 NA
## Pred 0.04758438 -0.06834378 0.202462711 0.02313860 -0.20101655
## Exp 0.24546836 0.12740768 0.260772984 -0.19291879 -0.19291879
## Danger 0.06528387 -0.06724755 0.208883617 -0.06666498 -0.20443928
15.4 理解缺失數(shù)據(jù)的來由和影響
識別缺失數(shù)據(jù)的數(shù)目甩十、分布和模式有兩個目的：(1)分析生成缺失數(shù)據(jù)的潛在機(jī)制船庇；(2)評價缺失數(shù)據(jù)對回答實(shí)質(zhì)性問題的影響。具體來講侣监，我們想弄清楚以下幾個問題鸭轮。
?缺失數(shù)據(jù)的比例多大？
?缺失數(shù)據(jù)是否集中在少數(shù)幾個變量上橄霉，抑或廣泛存在窃爷？
?缺失是隨機(jī)產(chǎn)生的嗎？
?缺失數(shù)據(jù)間的相關(guān)性或與可觀測數(shù)據(jù)間的相關(guān)性，是否可以表明產(chǎn)生缺失值的機(jī)制呢按厘？
回答這些問題將有助于判斷哪種統(tǒng)計(jì)方法最適合用來分析你的數(shù)據(jù)医吊。
假使已經(jīng)知道了缺失數(shù)據(jù)的來源和影響，那么讓我們看看如何轉(zhuǎn)換標(biāo)準(zhǔn)的統(tǒng)計(jì)方法來適應(yīng)缺失數(shù)據(jù)的分析逮京。我們將重點(diǎn)學(xué)習(xí)三種非常流行的方法：恢復(fù)數(shù)據(jù)的推理方法遮咖、涉及刪除缺失值的傳統(tǒng)方法、涉及模擬的現(xiàn)代方法造虏。
15.5 理性處理不完整數(shù)據(jù)
推理方法會根據(jù)變量間的數(shù)學(xué)或者邏輯關(guān)系來填補(bǔ)或恢復(fù)缺失值御吞。推理研究法常常需要創(chuàng)造性和想法，同時還需要許多數(shù)據(jù)處理技巧漓藕，而且數(shù)據(jù)的恢復(fù)可能是準(zhǔn)確的或者近似的陶珠。
15.6 完整實(shí)例分析（行刪除）
只有每個變量都包含了有效數(shù)據(jù)值的觀測才會保留下來做進(jìn)一步的分析。實(shí)際上享钞，這樣會導(dǎo)致包含一個或多個缺失值的任意一行都會被刪除揍诽，因此常稱作行刪除法（listwise）、個案刪除（case-wise）或剔除栗竖。
newdata <- mydata[complete.cases(mydata),]
newdata <- na.omit(mydata)
兩行代碼表示的意思都是：mydata中所有包含缺失數(shù)據(jù)的行都被刪除暑脆，然后結(jié)果才存儲到newdata中。
options(digits = 1) # 有效數(shù)字保留1位狐肢。
cor(na.omit(sleep)) # 刪除缺失值后計(jì)算sleep數(shù)據(jù)集相關(guān)系數(shù)添吗。
## BodyWgt BrainWgt NonD Dream Sleep Span Gest Pred Exp Danger
## BodyWgt 1.00 0.96 -0.4 -0.07 -0.3 0.47 0.71 0.10 0.4 0.26
## BrainWgt 0.96 1.00 -0.4 -0.07 -0.3 0.63 0.73 -0.02 0.3 0.15
## NonD -0.39 -0.39 1.0 0.52 1.0 -0.37 -0.61 -0.35 -0.6 -0.53
## Dream -0.07 -0.07 0.5 1.00 0.7 -0.27 -0.41 -0.40 -0.5 -0.57
## Sleep -0.34 -0.34 1.0 0.72 1.0 -0.38 -0.61 -0.40 -0.6 -0.60
## Span 0.47 0.63 -0.4 -0.27 -0.4 1.00 0.65 -0.17 0.3 0.01
## Gest 0.71 0.73 -0.6 -0.41 -0.6 0.65 1.00 0.09 0.6 0.31
## Pred 0.10 -0.02 -0.4 -0.40 -0.4 -0.17 0.09 1.00 0.6 0.93
## Exp 0.41 0.32 -0.6 -0.50 -0.6 0.32 0.57 0.63 1.0 0.79
## Danger 0.26 0.15 -0.5 -0.57 -0.6 0.01 0.31 0.93 0.8 1.00
summary(lm(Dream ~ Span + Gest, data = na.omit(sleep))) # 按行刪除缺失值后進(jìn)行的回歸分析。
##
## Call:
## lm(formula = Dream ~ Span + Gest, data = na.omit(sleep))
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.333 -0.915 -0.221 0.382 4.183
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.480122 0.298476 8.31 3.7e-10 ***
## Span -0.000472 0.013130 -0.04 0.971
## Gest -0.004394 0.002081 -2.11 0.041 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1 on 39 degrees of freedom
## Multiple R-squared: 0.167, Adjusted R-squared: 0.125
## F-statistic: 3.92 on 2 and 39 DF, p-value: 0.0282
15.7 多重插補(bǔ)
多重插補(bǔ)（MI）是一種基于重復(fù)模擬的處理缺失值的方法份名。在面對復(fù)雜的缺失值問題時碟联，MI是最常選用的方法，它將從一個包含缺失值的數(shù)據(jù)集中生成一組完整的數(shù)據(jù)集（通常是3到10個）僵腺。每個模擬數(shù)據(jù)集中鲤孵，缺失數(shù)據(jù)將用蒙特卡洛方法來填補(bǔ)。此時辰如，標(biāo)準(zhǔn)的統(tǒng)計(jì)方法便可應(yīng)用到每個模擬的數(shù)據(jù)集上普监，通過組合輸出結(jié)果給出估計(jì)的結(jié)果，以及引入缺失值時的置信區(qū)間琉兜。R中可利用Amelia凯正、mice和mi包來執(zhí)行這些操作。本節(jié)中呕童，我們將重點(diǎn)學(xué)習(xí)mice包（利用鏈?zhǔn)椒匠痰亩嘣逖a(bǔ)）提供的方法漆际。
基于mice包的分析通常符合以下分析過程：
library(mice)
imp <- mice(data, m)
fit <- with(imp, analysis)
pooled <- pool(fit)
summary(pooled)
? mydata是一個包含缺失值的矩陣或數(shù)據(jù)框淆珊。
? imp是一個包含m個插補(bǔ)數(shù)據(jù)集的列表對象夺饲，同時還含有完成插補(bǔ)過程的信息。默認(rèn)地，m為5往声。
? analysis是一個表達(dá)式對象擂找，用來設(shè)定應(yīng)用于m個插補(bǔ)數(shù)據(jù)集的統(tǒng)計(jì)分析方法。方法包括做線性回歸模型的lm()函數(shù)浩销、做廣義線性模型的glm()函數(shù)贯涎、做廣義可加模型的gam()，以及做負(fù)二項(xiàng)模型的nbrm()函數(shù)慢洋。表達(dá)式在函數(shù)的括號中塘雳，~的左邊是響應(yīng)變量，右邊是預(yù)測變量（用+符號分隔開）普筹。
? fit是一個包含m個單獨(dú)統(tǒng)計(jì)分析結(jié)果的列表對象败明。
? pooled是一個包含這m個統(tǒng)計(jì)分析平均結(jié)果的列表對象。
library(mice) # 調(diào)用mice包太防。
data(sleep, package="VIM") # 調(diào)用數(shù)據(jù)集sleep妻顶。
imp <- mice(sleep, seed=1234) # 對sleep缺失值進(jìn)行插補(bǔ)。
##
## 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
fit25 <- with(imp, lm(Dream ~ Span + Gest)) # 對插補(bǔ)后數(shù)據(jù)集線性回歸蜒车。
pooled <- pool(fit25) # 統(tǒng)計(jì)分析平均結(jié)果列表讳嘱。
summary(pooled) # 返回結(jié)果。
## term estimate std.error statistic df p.value
## 1 (Intercept) 2.597 0.249 10.4 52 2e-14
## 2 Span -0.004 0.012 -0.3 56 7e-01
## 3 Gest -0.004 0.001 -3.0 55 5e-03
你可以通過檢查分析過程所創(chuàng)建的對象來獲取更多的插補(bǔ)信息酿愧。
imp # 查看對象imp沥潭。
## Class: mids
## Number of multiple imputations: 5
## Imputation methods:
## BodyWgt BrainWgt NonD Dream Sleep Span Gest Pred
## "" "" "pmm" "pmm" "pmm" "pmm" "pmm" ""
## Exp Danger
## "" ""
## PredictorMatrix:
## BodyWgt BrainWgt NonD Dream Sleep Span Gest Pred Exp Danger
## BodyWgt 0 1 1 1 1 1 1 1 1 1
## BrainWgt 1 0 1 1 1 1 1 1 1 1
## 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
## Number of logged events: 5
## it im dep meth out
## 1 3 2 Span pmm Sleep
## 2 3 2 Gest pmm Sleep
## 3 4 2 Span pmm Sleep
## 4 4 2 Gest pmm Sleep
## 5 4 4 Span pmm Sleep
通過提取imp對象的子成分，可以觀測到實(shí)際的插補(bǔ)值嬉挡。
imp$imp$Dream # 提取子成分Dream叛氨。
## 1 2 3 4 5
## 1 0.0 0.5 0.5 0.5 0.3
## 3 0.5 1.4 1.5 1.5 1.3
## 4 3.6 4.1 3.1 4.1 2.7
## 14 0.3 1.0 0.5 0.0 0.0
## 24 3.6 0.8 1.4 1.4 0.9
## 26 2.4 0.5 3.9 3.4 1.2
## 30 2.6 0.8 2.4 2.2 3.1
## 31 0.6 1.3 1.2 1.8 2.1
## 47 1.3 1.8 1.8 1.8 3.9
## 53 0.5 0.5 0.6 0.5 0.3
## 55 2.6 3.6 2.4 1.8 0.5
## 62 1.5 3.4 3.9 3.4 2.2
利用 complete() 函數(shù)可以觀察m個插補(bǔ)數(shù)據(jù)集中的任意一個。格式為：complete(imp, action=#)棘伴，其中 # 指定m個完整數(shù)據(jù)集中的一個來展示寞埠。
dataset3 <- complete(imp, action=3)
dataset3
## BodyWgt BrainWgt NonD Dream Sleep Span Gest Pred Exp Danger
## 1 7e+03 6e+03 3 0.5 3 39 645 3 5 3
## 2 1e+00 7e+00 6 2.0 8 4 42 3 1 3
## 3 3e+00 4e+01 11 1.5 12 14 60 1 1 1
## 4 9e-01 6e+00 13 3.1 16 7 25 5 2 3
## 5 3e+03 5e+03 2 1.8 4 69 624 3 5 4
## 6 1e+01 2e+02 9 0.7 10 27 180 4 4 4
## 7 2e-02 3e-01 16 3.9 20 19 35 1 1 1
## 8 2e+02 2e+02 5 1.0 6 30 392 4 5 4
## 9 3e+00 3e+01 11 3.6 14 28 63 1 2 1
## 10 5e+01 4e+02 8 1.4 10 50 230 1 1 1
## 11 4e-01 6e+00 11 1.5 12 7 112 5 4 4
## 12 5e+02 4e+02 3 0.7 4 30 281 5 5 5
## 13 6e-01 2e+00 8 2.7 10 4 42 2 1 2
## 14 2e+02 4e+02 3 0.5 3 40 365 5 5 5
## 15 7e-02 1e+00 6 2.1 8 4 42 1 1 1
## 16 3e+00 2e+01 9 0.0 9 50 28 2 2 2
## 17 8e-01 4e+00 7 4.1 11 6 42 2 2 2
## 18 2e-01 5e+00 10 1.2 11 10 120 2 2 2
## 19 1e+00 2e+01 5 1.3 6 34 210 1 2 1
## 20 6e+01 8e+01 12 6.1 18 7 14 1 1 1
## 21 5e+02 7e+02 12 0.3 12 28 400 5 5 5
## 22 3e+01 1e+02 3 0.5 4 20 148 5 5 5
## 23 1e-01 1e+00 11 3.4 14 4 16 3 1 2
## 24 2e+02 4e+02 11 1.4 12 39 252 1 4 1
## 25 8e+01 3e+02 5 1.5 6 41 310 1 3 1
## 26 4e+01 1e+02 10 3.9 13 16 63 1 1 1
## 27 1e-01 4e+00 10 3.4 14 9 28 5 1 3
## 28 1e+00 6e+00 7 0.8 8 8 68 5 3 4
## 29 5e+02 7e+02 2 0.8 3 46 336 5 5 5
## 30 1e+02 2e+02 8 2.4 11 22 100 1 1 1
## 31 4e+01 6e+01 7 1.2 8 16 33 3 5 4
## 32 5e-03 1e-01 8 1.4 9 3 22 5 2 4
## 33 1e-02 2e-01 18 2.0 20 24 50 1 1 1
## 34 6e+01 1e+03 6 1.9 8 100 267 1 1 1
## 35 1e-01 3e+00 8 2.4 11 2 30 2 1 1
## 36 1e+00 8e+00 8 2.8 11 5 45 3 1 3
## 37 2e-02 4e-01 12 1.3 13 3 19 4 1 3
## 38 5e-02 3e-01 11 2.0 13 2 30 4 1 3
## 39 2e+00 6e+00 14 5.6 19 5 12 2 1 1
## 40 4e+00 1e+01 14 3.1 17 6 120 2 1 1
## 41 2e+02 5e+02 12 1.0 12 24 440 5 5 5
## 42 5e-01 2e+01 15 1.8 17 12 140 2 2 2
## 43 1e+01 1e+02 10 0.9 11 20 170 4 4 4
## 44 2e+00 1e+01 12 1.8 14 13 17 2 1 2
## 45 2e+02 2e+02 6 1.9 8 27 115 4 4 4
## 46 2e+00 1e+01 8 0.9 8 18 31 5 5 5
## 47 4e+00 4e+01 11 1.8 12 14 63 2 2 2
## 48 3e-01 2e+00 11 2.6 13 5 21 3 1 3
## 49 4e+00 5e+01 7 2.4 10 10 52 1 1 1
## 50 7e+00 2e+02 8 1.2 10 29 164 2 3 2
## 51 8e-01 1e+01 6 0.9 7 7 225 2 2 2
## 52 4e+00 2e+01 5 0.5 5 6 225 3 2 3
## 53 1e+01 1e+02 3 0.6 3 17 150 5 5 5
## 54 6e+01 2e+02 3 0.6 4 20 151 5 5 5
## 55 1e+00 1e+01 8 2.4 11 13 90 2 2 2
## 56 6e-02 1e+00 8 2.2 10 4 30 3 1 2
## 57 9e-01 3e+00 11 2.3 13 4 60 2 1 2
## 58 2e+00 1e+01 5 0.5 5 8 200 3 1 3
## 59 1e-01 2e+00 13 2.6 16 2 46 3 2 2
## 60 4e+00 6e+01 10 0.6 10 24 210 4 3 4
## 61 4e+00 4e+00 13 6.6 19 3 14 2 1 1
## 62 4e+00 2e+01 8 3.9 11 13 38 3 1 1
15.8 處理缺失值的其他方法
15.8.1 成對刪除
處理含缺失值的數(shù)據(jù)集時，成對刪除常作為行刪除的備選方法使用焊夸。對于成對刪除仁连，觀測只是當(dāng)它含缺失數(shù)據(jù)的變量涉及某個特定分析時才會被刪除。
15.8.2 簡單（非隨機(jī)）插補(bǔ)
所謂簡單插補(bǔ)阱穗，即用某個值（如均值饭冬、中位數(shù)或眾數(shù)）來替換變量中的缺失值。若使用均值替換揪阶，Dream變量中的缺失值可用1.97來替換昌抠，NonD中的缺失值可用8.67來替換（兩個值分別是Dream和NonD的均值）。注意這些替換是非隨機(jī)的鲁僚，這意味著不會引入隨機(jī)誤差（與多重插補(bǔ)不同）炊苫。
簡單插補(bǔ)的一個優(yōu)點(diǎn)是裁厅，解決“缺失值問題”時不會減少分析過程中可用的樣本量。雖然簡單插補(bǔ)用法很簡單侨艾，但是對于非MCAR的數(shù)據(jù)會產(chǎn)生有偏的結(jié)果执虹。若缺失數(shù)據(jù)的數(shù)目非常大，那么簡單插補(bǔ)很可能會低估標(biāo)準(zhǔn)差唠梨、曲解變量間的相關(guān)性袋励，并會生成不正確的統(tǒng)計(jì)檢驗(yàn)的p值。與成對刪除一樣当叭，我建議在解決缺失數(shù)據(jù)的問題時盡量避免使用該方法茬故。
15.9 小結(jié)
參考資料：
- 《R語言實(shí)戰(zhàn)》（中文版），人民郵電出版社蚁鳖，2013.
- 如何理解主成分分析法 （PCA）均牢，https://zhuanlan.zhihu.com/p/170398464
- 主成分分析法，https://blog.csdn.net/weixin_43914260/article/details/99585202
