PART I 缺失值的基本操作
- 基本处理
lactate <- c(0.2,3.6,4.2,NA,6.1,2.5)
is.na(lactate)
# [1] FALSE FALSE FALSE TRUE FALSE FALSE
which(is.na(lactate))
# [1] 4
x1 <- c(1,4,3,NA,7,8,NA,7,7)
x1 < 0
# [1] FALSE FALSE FALSE NA FALSE FALSE NA FALSE FALSE
NA 不可以和0进行比较
NA 也不可以和其它数值比较
函数is.na可以用于寻找和创建NA值
- 带有缺失值的计算
mean(lactate)
# [1] NA
sum(lactate)
# [1] NA
mean(lactate,na.rm = TRUE)
# [1] 3.32
sd(lactate,na.rm = TRUE)
# [1] 2.178761
- 带有缺失值的数据框(注意
<NA>和NA的区别)
ptid <- 1:6
sex <- c("m","f",NA,"f","m","m")
lactate <- c(0.2,3.3,4.5,NA,6.1,5.6)
data <- data.frame(ptid,sex,lactate)
data
# ptid sex lactate
# 1 1 m 0.2
# 2 2 f 3.3
# 3 3 <NA> 4.5
# 4 4 f NA
# 5 5 m 6.1
# 6 6 m 5.6
na.fail(data)
# Error in na.fail.default(data) : missing values in object
- 回归方程里的缺失值(注意omit和exclude的区别)
model.omit <- lm(ptid ~lactate,data = data,na.action = na.omit)
model.exclude <- lm(ptid ~lactate,data = data,na.action = na.exclude)
resid(model.omit)
# 1 2 3 5 6
# 0.5437376 -0.8962588 -0.8407735 -0.1001265 1.2934213
resid(model.exclude)
# 1 2 3 4 5 6
# 0.5437376 -0.8962588 -0.8407735 NA -0.1001265 1.2934213
fitted(model.omit)
# 1 2 3 5 6
# 0.4562624 2.8962588 3.8407735 5.1001265 4.7065787
fitted(model.omit)
# 1 2 3 5 6
# 0.4562624 2.8962588 3.8407735 5.1001265 4.7065787
na.omit(data)
# ptid sex lactate
# 1 1 m 0.2
# 2 2 f 3.3
# 5 5 m 6.1
# 6 6 m 5.6
complete.cases(data)
# [1] TRUE TRUE FALSE FALSE TRUE TRUE
complete.data <- data[complete.cases(data),]
complete.data
# ptid sex lactate
# 1 1 m 0.2
# 2 2 f 3.3
# 5 5 m 6.1
# 6 6 m 5.6
na.fail(complete.data)
# ptid sex lactate
# 1 1 m 0.2
# 2 2 f 3.3
# 5 5 m 6.1
# 6 6 m 5.6
- 缺失值的寻找、定位
is.na(data)
# ptid sex lactate
# [1,] FALSE FALSE FALSE
# [2,] FALSE FALSE FALSE
# [3,] FALSE TRUE FALSE
# [4,] FALSE FALSE TRUE
# [5,] FALSE FALSE FALSE
# [6,] FALSE FALSE FALSE
which(is.na(data)) #按列顺序
# [1] 9 16
unique(unlist(lapply (data,function(x) which (is.na(x)))))
# [1] 3 4
- 删除缺失超过10%的变量
missing.percent <- unlist(lapply(data, function(x) sum(is.na(x))))/nrow(data)
missing.percent
# ptid sex lactate
# 0.0000000 0.1666667 0.1666667
data.tenth <- data[,missing.percent <= 0.1]
data.tenth
# [1] 1 2 3 4 5 6
- 有意义的缺失值
(比如，乳酸值的缺失是由于患者的病情稳定)
data.discrete <- data.frame(ptid = c(1,2,3,4,5),
lactate = c("low",NA,"moderate",NA,"high"),
death = c("y","y",NA,"n","y"))
data.discrete
# ptid lactate death
# 1 1 low y
# 2 2 <NA> y
# 3 3 moderate <NA>
# 4 4 <NA> n
# 5 5 high y
na.fail(data.discrete)
# Error in na.fail.default(data.discrete) : missing values in object
以下代码无法重现，即使查看了章教授的原始文献: Missing values in big data research: some basic skills
# 对于分类变量
library(tree)
newdata.discrete <- na.tree.replace(data.discrete)
# newdata.discrete <- na.gam.replace(data.discrete)
na.fail(newdata.discrete)
# 经过这个处理，NA已经不是缺失值，而是新的一个类别
data
na.tree.replace(data)
# 对于连续变量
# install.packages("gam")
library(foreach)
library(gam)
newdata <- na.gam.replace(data)
na.fail(newdata)
image.png
经过
tree处理，NA已经不是缺失值，而是新的一个类别。注意NA前后的变化。
PART II 缺失值可视化
为何要探索缺失值
探索缺失数据有利于发现缺失形成的规律，从而找到插补方法
eg: 乳酸值的缺失可能是由于患者病情稳定
eg: CRP值缺失可以从WBC进行推断，从而进行插补
数据的缺失类型
image.png
模拟数据
set.seed(888)
age <- round(abs(rnorm(200,mean = 67,sd = 19)))
sex <- rbinom(200,1,0.45)
sex.miss.tag <- rbinom(200,1,0.3) # MCAR
sex[sex.miss.tag ==1] <- NA
sex[sex ==1] <- "male"
sex[sex ==0] <- "female"
lac <- round(abs(rnorm(200,mean = 3,sd = 4)),1)
lac.miss.tag <- rbinom(200,1,0.3)
lac[lac <= 3&lac.miss.tag == 1] <- NA #NMAR
wbc <- round(abs(rnorm(200,mean = 12,sd = 4)),1)
wbc.miss.tag <- rbinom(200,1,0.3)
wbc[wbc.miss.tag ==1] <- NA
crp <- round(abs(rnorm(200,mean = 50,sd = 100)),1)
crp.miss.tag <- rbinom(200,1,0.4)
crp[wbc <= 12& crp.miss.tag ==1] <- NA # MAR
data <- data.frame(age,sex,lac,wbc,crp)
data[10:21,]
image.png
library(mice)
md.pattern(data)
image.png
image.png
library(VIM)
matrixplot(data) #一种交互可视化图?(没有探索出章教授视频中提到的交互效果)
image.png
用图探索缺失值的特点
barMiss(data)
image.png
很奇怪，这里和视频里的图也有出入。视频里的x轴为lac。似乎应该有个交互界面。
table(complete.cases(data[lac <= 2&!is.na(lac),]))
# FALSE TRUE
# 33 11
table(complete.cases(data[is.na(lac),][c("age","sex","wbc","crp")]))
# FALSE TRUE
# 16 16
aggr(data,numbers = TRUE,prop = FALSE)
image.png
探索缺失变量之间的相互关系
marginplot(data[c("wbc","crp")],pch =c(20),col =c("green","red","blue"))
image.png
全部变量之间的缺失相互关系
marginmatrix(data)
image.png
spineMiss(data)
image.png
该图与视频中的图有出入。视频中的似乎是GUI交互界面，以lac为横轴。
scattmatrixMiss(data)
image.png
红色：在其它变量里有缺失值
淡蓝色：无缺失值的病例
星号：The red-cross symbols represent observations with missing values on any of the variables age, sex, lac, wbc and crp
Exploring missing data pattern by correlation matrix
一个变量的缺失与另一个变量的缺失的相关性
shadow <- as.data.frame(abs(is.na(data)))
miss.shadow <- shadow[,which(unlist(lapply(shadow,sum))!=0)]
round(cor(miss.shadow),3)
# sex lac wbc crp
# sex 1.000 -0.074 -0.048 0.023
# lac -0.074 1.000 0.022 -0.097
# wbc -0.048 0.022 1.000 -0.277
# crp 0.023 -0.097 -0.277 1.000
影子矩阵： 缺失值用1表示，非缺失值用0
相关性不强：一个变量的缺失与另一个变量的缺失不相关
某变量值的缺失会不会受其他变量数值的影响(仅针对数值型变量)？
round(cor(data[!names(data)%in%c("sex")],miss.shadow,use ="pairwise.complete.obs"),3)
# sex lac wbc crp
# age -0.167 -0.001 0.060 -0.015
# lac -0.071 NA -0.096 -0.040
# wbc -0.017 0.026 NA -0.326
# crp 0.006 0.109 0.005 NA
crp与wbc呈负相关(-0.326), 说明crp的缺失在低水平的wbc更容易出现
names(data)%in%c("sex") 逻辑值，在sex变量的位置出现true
data[!names(data)%in%c("sex")] ：排除sex变量，因其为二分类变量
参考资料
文中代码及部分截图来自章仲恒教授的丁香园课程：缺失数据的基本处理方法 及 缺失数据的可视化
网友评论