应用预测建模第三章数据预处理习题3.2【退化分布、缺失值的缺失模式探索、处理缺失值】

本文链接：https://blog.csdn.net/totobey/article/details/103219066

《应用预测建模》Applied Predictive Modeling (2013) by Max Kuhn and Kjell Johnson，林荟等译

第三章数据预处理

3. 2 UC机器学习数据库里还有一组大豆数据。数据收集自683 个大豆样本，其目的是用来预测大豆疾病。35 个预测变量大部分是分类变量，包括环境条件的信息（如温度和降水量）和种植条件（如霉菌生长）。结果变量是19 个不同类的分类变量。

( a）检查分类预测变量的频数分布。其中是否存在某些变量服从退化分布（如本章之前讨论过的）？
( b )数据大约有18% 的缺失。缺失的部分是否集中于某些预测变量？缺失数据的模式是否与结果变量有关？
( c）提出一个处理缺失值的策略。可以删除或填补。

library(mlbench)
data(Soybean)
str(Soybean)
head(Soybean)

> str(Soybean)
'data.frame':	683 obs. of  36 variables:
 $ Class          : Factor w/ 19 levels "2-4-d-injury",..: 11 11 11 11 11 11 11 11 11 11 ...
 $ date           : Factor w/ 7 levels "0","1","2","3",..: 7 5 4 4 7 6 6 5 7 5 ...
 $ plant.stand    : Ord.factor w/ 2 levels "0"<"1": 1 1 1 1 1 1 1 1 1 1 ...
 $ precip         : Ord.factor w/ 3 levels "0"<"1"<"2": 3 3 3 3 3 3 3 3 3 3 ...
 $ temp           : Ord.factor w/ 3 levels "0"<"1"<"2": 2 2 2 2 2 2 2 2 2 2 ...
 $ hail           : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 2 1 1 ...
 $ crop.hist      : Factor w/ 4 levels "0","1","2","3": 2 3 2 2 3 4 3 2 4 3 ...
 $ area.dam       : Factor w/ 4 levels "0","1","2","3": 2 1 1 1 1 1 1 1 1 1 ...
 $ sever          : Factor w/ 3 levels "0","1","2": 2 3 3 3 2 2 2 2 2 3 ...
 $ seed.tmt       : Factor w/ 3 levels "0","1","2": 1 2 2 1 1 1 2 1 2 1 ...
 $ germ           : Ord.factor w/ 3 levels "0"<"1"<"2": 1 2 3 2 3 2 1 3 2 3 ...
 $ plant.growth   : Factor w/ 2 levels "0","1": 2 2 2 2 2 2 2 2 2 2 ...
 $ leaves         : Factor w/ 2 levels "0","1": 2 2 2 2 2 2 2 2 2 2 ...
 $ leaf.halo      : Factor w/ 3 levels "0","1","2": 1 1 1 1 1 1 1 1 1 1 ...
 $ leaf.marg      : Factor w/ 3 levels "0","1","2": 3 3 3 3 3 3 3 3 3 3 ...
 $ leaf.size      : Ord.factor w/ 3 levels "0"<"1"<"2": 3 3 3 3 3 3 3 3 3 3 ...
 $ leaf.shread    : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
 $ leaf.malf      : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
 $ leaf.mild      : Factor w/ 3 levels "0","1","2": 1 1 1 1 1 1 1 1 1 1 ...
 $ stem           : Factor w/ 2 levels "0","1": 2 2 2 2 2 2 2 2 2 2 ...
 $ lodging        : Factor w/ 2 levels "0","1": 2 1 1 1 1 1 2 1 1 1 ...
 $ stem.cankers   : Factor w/ 4 levels "0","1","2","3": 4 4 4 4 4 4 4 4 4 4 ...
 $ canker.lesion  : Factor w/ 4 levels "0","1","2","3": 2 2 1 1 2 1 2 2 2 2 ...
 $ fruiting.bodies: Factor w/ 2 levels "0","1": 2 2 2 2 2 2 2 2 2 2 ...
 $ ext.decay      : Factor w/ 3 levels "0","1","2": 2 2 2 2 2 2 2 2 2 2 ...
 $ mycelium       : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
 $ int.discolor   : Factor w/ 3 levels "0","1","2": 1 1 1 1 1 1 1 1 1 1 ...
 $ sclerotia      : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
 $ fruit.pods     : Factor w/ 4 levels "0","1","2","3": 1 1 1 1 1 1 1 1 1 1 ...
 $ fruit.spots    : Factor w/ 4 levels "0","1","2","4": 4 4 4 4 4 4 4 4 4 4 ...
 $ seed           : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
 $ mold.growth    : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
 $ seed.discolor  : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
 $ seed.size      : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
 $ shriveling     : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
 $ roots          : Factor w/ 3 levels "0","1","2": 1 1 1 1 1 1 1 1 1 1 ...
> head(Soybean)
                  Class date plant.stand precip temp hail crop.hist area.dam sever
1 diaporthe-stem-canker    6           0      2    1    0         1        1     1
2 diaporthe-stem-canker    4           0      2    1    0         2        0     2
3 diaporthe-stem-canker    3           0      2    1    0         1        0     2
4 diaporthe-stem-canker    3           0      2    1    0         1        0     2
5 diaporthe-stem-canker    6           0      2    1    0         2        0     1
6 diaporthe-stem-canker    5           0      2    1    0         3        0     1
  seed.tmt germ plant.growth leaves leaf.halo leaf.marg leaf.size leaf.shread
1        0    0            1      1         0         2         2           0
2        1    1            1      1         0         2         2           0
3        1    2            1      1         0         2         2           0
4        0    1            1      1         0         2         2           0
5        0    2            1      1         0         2         2           0
6        0    1            1      1         0         2         2           0
  leaf.malf leaf.mild stem lodging stem.cankers canker.lesion fruiting.bodies
1         0         0    1       1            3             1               1
2         0         0    1       0            3             1               1
3         0         0    1       0            3             0               1
4         0         0    1       0            3             0               1
5         0         0    1       0            3             1               1
6         0         0    1       0            3             0               1
  ext.decay mycelium int.discolor sclerotia fruit.pods fruit.spots seed mold.growth
1         1        0            0         0          0           4    0           0
2         1        0            0         0          0           4    0           0
3         1        0            0         0          0           4    0           0
4         1        0            0         0          0           4    0           0
5         1        0            0         0          0           4    0           0
6         1        0            0         0          0           4    0           0
  seed.discolor seed.size shriveling roots
1             0         0          0     0
2             0         0          0     0
3             0         0          0     0
4             0         0          0     0
5             0         0          0     0
6             0         0          0     0

( a）检查分类预测变量的频数分布。其中是否存在某些变量服从退化分布（如本章之前讨论过的）？

退化分布：具有唯一值的预测变量（即零方差预测变量，变量只有一个取值）或同时具有以下两个特征的预测变量：
相对于样本数量，它们具有很少的唯一值；最常见值的频率与次最常见值的频率之比很大。

可见，预测变量"leaf.mild" "mycelium" "sclerotia"服从退化分布。

建议从数据集中移除这些变量，这些变量会削弱一些模型，删除这些变量能显著提高模型的表现与稳定性（p31)

#nearZeroVar可诊断具有唯一值的预测变量（即零方差预测变量）或同时具有以下两个特征的预测变量：
#相对于样本数量，它们具有很少的唯一值；最常见值的频率与次最常见值的频率之比很大。
#nearZeroVar(x,freqCut = 95/5,uniqueCut = 10,saveMetrics = FALSE,names = FALSE,foreach = FALSE,allowParallel = TRUE)
#freqCut 最常见值与第二常见值之比的临界值，默认95/5
#uniqueCut 样本总数中不同值的百分比的临界值，默认10
near.zero.ind<-nearZeroVar(Soybean)
names(Soybean[near.zero.ind])# "leaf.mild" "mycelium"  "sclerotia"

( b )数据大约有18% 的缺失。缺失的部分是否集中于某些预测变量？缺失数据的模式是否与结果变量有关？

#计算每个变量的缺失值
na.var<-apply(Soybean,2,function(x) sum(is.na(x)))
#按缺失值数量排序
sort(na.var,decreasing = TRUE)

将每个变量的缺失值排序。可见hail 、sever、seed.tmt、lodging、 germ 等变量都含有较多缺失值。

> sort(na.var,decreasing = TRUE)
           hail           sever        seed.tmt         lodging            germ 
            121             121             121             121             112 
      leaf.mild fruiting.bodies     fruit.spots   seed.discolor      shriveling 
            108             106             106             106             106 
    leaf.shread            seed     mold.growth       seed.size       leaf.halo 
            100              92              92              92              84 
      leaf.marg       leaf.size       leaf.malf      fruit.pods          precip 
             84              84              84              84              38 
   stem.cankers   canker.lesion       ext.decay        mycelium    int.discolor 
             38              38              38              38              38 
      sclerotia     plant.stand           roots            temp       crop.hist 
             38              36              31              30              16 
   plant.growth            stem            date        area.dam           Class 
             16              16               1               1               0 
         leaves 
              0


#查看因变量的取值，总共有19个类别
unique(Soybean$Class)
#做频数交叉分布图
#新建一个na.flg向量，原始数据每一行是否有缺失值
na.flg<-apply(Soybean[,-1],1,function(x) sum(is.na(x))>0  )
table(Soybean$Class,na.flg)

可见，某些因变量的取值类别无缺失值，如alternarialeaf-spot、anthracnose 等。

某些因变量的取值类别的所有观测都有缺失值，如 cyst-nematode、 diaporthe-pod-&-stem-blight 。

                              FALSE TRUE
  2-4-d-injury                    0   16
  alternarialeaf-spot            91    0
  anthracnose                    44    0
  bacterial-blight               20    0
  bacterial-pustule              20    0
  brown-spot                     92    0
  brown-stem-rot                 44    0
  charcoal-rot                   20    0
  cyst-nematode                   0   14
  diaporthe-pod-&-stem-blight     0   15
  diaporthe-stem-canker          20    0
  downy-mildew                   20    0
  frog-eye-leaf-spot             91    0
  herbicide-injury                0    8
  phyllosticta-leaf-spot         20    0
  phytophthora-rot               20   68
  powdery-mildew                 20    0
  purple-seed-stain              20    0
  rhizoctonia-root-rot           20    0

#用相关性探索缺失值
#可用指示变量替代数据集中的数据（1表示缺失，0表示存在），这样生成的矩阵有时被称作影子矩阵。
#求这些指示变量之间和它们与初始（可观测）变量之间的相关性，有助于观察哪些变量常一起缺失，以及分析变量“缺失”与其他变量间的关系。
#若数据集的元素缺失，则数据框x对应的元素为1，否则为0
x <- as.data.frame(abs(is.na(Soybean)))
#可提取含（但不全部是）缺失值的变量
y <- x[which(apply(x,2,sum)>0)]
#列出这些指示变量间的相关系数：
correlations<-cor(y)#求相关性
library(corrplot)#可视化
#可通过图观测预测变量之间的缺失值是否有相关性
corrplot(correlations,order='hclust', tl.cex = .35)

以上是探究预测变量之间缺失值的相关性，但题目主要是想探究预测变量缺失值是否与响应变量有关：

#看题，主要是为了观察自变量的缺失与因变量有无关系
#因为响应变量为因子型（且为无序），因此使用卡方检验（否则可使用相关系数）
library(vcd)
mytable <- xtabs(~na.flg+Soybean$Class)
chisq.test(mytable)

从全局角度出发，p值接近于0，说明预测变量的缺失与响应变量有关:

> chisq.test(mytable)

	Pearson's Chi-squared test

data:  mytable
X-squared = 576.98, df = 18, p-value < 2.2e-16

预测变量分开看:

#预测变量分开看
myChisqTest <- function(x) {
  mytable <- xtabs(~y[,x]+Soybean$Class)
  chisq.test(mytable)$p.value
}
result<-sapply(names(y), myChisqTest)
summary(result)

可见，p值均小于0.05，说明所有预测变量的缺失都与响应变量有关。

> result
           date     plant.stand          precip            temp            hail 
   1.198859e-03   4.268778e-118   2.292469e-133   2.292469e-133   6.272918e-111 
      crop.hist        area.dam           sever        seed.tmt            germ 
  2.292469e-133    1.198859e-03   6.272918e-111   6.272918e-111   5.674080e-104 
   plant.growth       leaf.halo       leaf.marg       leaf.size     leaf.shread 
  2.292469e-133    5.567557e-93    5.567557e-93    5.567557e-93    1.733025e-98 
      leaf.malf       leaf.mild            stem         lodging    stem.cankers 
   5.567557e-93   1.296453e-100   2.292469e-133   6.272918e-111   2.292469e-133 
  canker.lesion fruiting.bodies       ext.decay        mycelium    int.discolor 
  2.292469e-133   1.995997e-108   2.292469e-133   2.292469e-133   2.292469e-133 
      sclerotia      fruit.pods     fruit.spots            seed     mold.growth 
  2.292469e-133   4.553154e-103   1.995997e-108   2.542290e-105   2.542290e-105 
  seed.discolor       seed.size      shriveling           roots 
  1.995997e-108   2.542290e-105   1.995997e-108   2.292469e-133 
> summary(result)
     Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
0.000e+00 0.000e+00 0.000e+00 7.052e-05 0.000e+00 1.199e-03

( c）提出一个处理缺失值的策略。可以删除或填补。

可用caret包中的preProcess()函数，其中有三种缺失值插补方法：K近邻、装袋树、中位数。

#处理缺失值
library( caret )

#预测变量为因子，需要转换为数值型
Soybean_num<-apply(Soybean[,-1],2,as.numeric)
as.data.frame(Soybean_num)
head(Soybean_num)
summary(Soybean_num)

#中位数插补
trans<-preProcess(Soybean_num, method=c("medianImpute"), na.remove=FALSE )
transformed<-predict(trans,Soybean_num)
summary(transformed)

#装袋法插补
trans<-preProcess(Soybean_num, method=c("bagImpute"), na.remove=FALSE )
transformed<-predict(trans,Soybean_num)
summary(transformed)