Rによる欠損値処理

今回は機械学習などでモデルを作る際に必須となる欠損値の扱いをtidyverseを用いて実装します。
欠損値の生成メカニズムなどの理論的な話は、参考サイトや参考文献をご覧ください。

#データセットの確認
今回はR組み込みデータセットのairqualityデータセットを使用します。

#使用パッケージの読み込み
library(tidyverse)
library(skimr)
library(modelr)

#データの確認
data("airquality")
df <- airquality
df %>% head()

#出力
  Ozone Solar.R Wind Temp Month Day
1    41     190  7.4   67     5   1
2    36     118  8.0   72     5   2
3    12     149 12.6   74     5   3
4    18     313 11.5   62     5   4
5    NA      NA 14.3   56     5   5
6    28      NA 14.9   66     5   6

欠損値が確認できますね。

↓skimrパッケージは欠損データの確認に便利です

df %>% skimr::skim()

#出力
Skim summary statistics
 n obs: 153 
 n variables: 6 

-- Variable type:integer -------------------------------------------------------
 variable missing complete   n   mean    sd p0    p25   p50    p75 p100
      Day       0      153 153  15.8   8.86  1   8     16    23      31
    Month       0      153 153   6.99  1.42  5   6      7     8       9
    Ozone      37      116 153  42.13 32.99  1  18     31.5  63.25  168
  Solar.R       7      146 153 185.93 90.06  7 115.75 205   258.75  334
     Temp       0      153 153  77.88  9.47 56  72     79    85      97
     hist
 ▇▇▇▇▆▇▇▇
 ▇▇▁▇▁▇▁▇
 ▇▆▃▃▂▁▁▁
 ▃▃▃▃▅▇▇▃
 ▂▂▃▆▇▇▃▃

-- Variable type:numeric -------------------------------------------------------
 variable missing complete   n  mean   sd  p0 p25 p50  p75 p100     hist
     Wind       0      153 153  9.96 3.52 1.7 7.4 9.7 11.5 20.7 ▁▃▇▇▅▅▁▁

この結果から、Ozone、Solar.Rに欠損値が含まれていることがわかります。
以下、簡単のためOzoneの欠損値に絞って記載します。

#リストワイズ法(lisewise deletion)
欠損値が一つでもあるデータを除きます。

df_lw <- 
  na.omit(df)

df_lw %>%  head()

#出力
  Ozone Solar.R Wind Temp Month Day
1    41     190  7.4   67     5   1
2    36     118  8.0   72     5   2
3    12     149 12.6   74     5   3
4    18     313 11.5   62     5   4
7    23     299  8.6   65     5   7
8    19      99 13.8   59     5   8

#平均値代入法(mean imputation)
欠損値を平均値で置き換えます。

df_mi <-
  df %>% 
  mutate(Ozone_new=if_else(condition = is.na(Ozone),
                           true = mean(Ozone,na.rm = T),
                           false = Ozone %>% as.double()))

df_mi %>%  head()

#出力
  Ozone Solar.R Wind Temp Month Day Ozone_new
1    41     190  7.4   67     5   1  41.00000
2    36     118  8.0   72     5   2  36.00000
3    12     149 12.6   74     5   3  12.00000
4    18     313 11.5   62     5   4  18.00000
5    NA      NA 14.3   56     5   5  42.12931
6    28      NA 14.9   66     5   6  28.00000

#回帰代入法(regression imputation)
回帰モデルを用いて予測した値で置き換えます。

#回帰モデル作成(今回はWindを用いて1次回帰しています。)
model_lm=
  lm(formula = Ozone~Wind,data = df)

df_ri <-
df %>% 
  add_predictions(model = model_lm,var = "Ozone_pred") %>% 
  mutate(Ozone_new=if_else(condition = is.na(Ozone),
                           true = Ozone_pred,
                           false = Ozone %>% as.double())) %>% 
  select(-Ozone_pred)

df_ri %>%  head()

#出力
  Ozone Solar.R Wind Temp Month Day Ozone_new
1    41     190  7.4   67     5   1   41.0000
2    36     118  8.0   72     5   2   36.0000
3    12     149 12.6   74     5   3   12.0000
4    18     313 11.5   62     5   4   18.0000
5    NA      NA 14.3   56     5   5   17.4947
6    28      NA 14.9   66     5   6   28.0000

回帰モデルでなくても、KNNやRandomForestなども使えます。

###参考
データ分析プロセス (シリーズ Useful R 2)
Rで実践！欠損データ分析入門
データ前処理〜欠損値への対応〜
データ前処理〜データの要約〜
skimrパッケージ