第6章欠損値の扱い

欠損値 (missing value) とは観測ができなかった、あるいはその時点において観察対象が存在していなかった、などの理由から変数が欠損していることを示しています。

Rでは主にNAとして表示されます。ここでは、Rにおける欠損値の扱い方について説明します。

ここでは欠損値の対処法（例、多重代入法）は触れません。

library(tidyverse)
library(readxl)

6.1 欠損値を読み込む

Rでは読み込むデータにおいて空欄である場合に欠損値として認識し、データセット内ではNAと表示されます。例えば、dyplyrに含まれているstarwarsというデータセットを見てみると、身長や重量、誕生年にNAがいくつかあるのが分かります。

starwars <- starwars %>% 
  select(1:10)
summary(starwars)

##      name               height           mass          hair_color       
##  Length:87          Min.   : 66.0   Min.   :  15.00   Length:87         
##  Class :character   1st Qu.:167.0   1st Qu.:  55.60   Class :character  
##  Mode  :character   Median :180.0   Median :  79.00   Mode  :character  
##                     Mean   :174.4   Mean   :  97.31                     
##                     3rd Qu.:191.0   3rd Qu.:  84.50                     
##                     Max.   :264.0   Max.   :1358.00                     
##                     NA's   :6       NA's   :28                          
##   skin_color         eye_color           birth_year         sex           
##  Length:87          Length:87          Min.   :  8.00   Length:87         
##  Class :character   Class :character   1st Qu.: 35.00   Class :character  
##  Mode  :character   Mode  :character   Median : 52.00   Mode  :character  
##                                        Mean   : 87.57                     
##                                        3rd Qu.: 72.00                     
##                                        Max.   :896.00                     
##                                        NA's   :44                         
##     gender           homeworld        
##  Length:87          Length:87         
##  Class :character   Class :character  
##  Mode  :character   Mode  :character  
##                                       
##                                       
##                                       
##

問題は、もとのデータセットにおいて欠損値が空欄になっているとは限らないという点です。例えば、Polity IVというデータセットを開いてみます。

polity <- read_excel("data/p4v2017.xls") %>% 
  select(ccode, country, year, democ, autoc)
summary(polity)

##      ccode         country               year          democ          
##  Min.   :  2.0   Length:17395       Min.   :1800   Min.   :-88.00000  
##  1st Qu.:200.0   Class :character   1st Qu.:1894   1st Qu.:  0.00000  
##  Median :366.0   Mode  :character   Median :1959   Median :  1.00000  
##  Mean   :404.2                      Mean   :1940   Mean   : -0.08928  
##  3rd Qu.:630.0                      3rd Qu.:1991   3rd Qu.:  7.00000  
##  Max.   :950.0                      Max.   :2017   Max.   : 10.00000  
##      autoc         
##  Min.   :-88.0000  
##  1st Qu.:  0.0000  
##  Median :  4.0000  
##  Mean   :  0.3837  
##  3rd Qu.:  7.0000  
##  Max.   : 10.0000

一見すると欠損値はないように見えますが、よく見るとdemocとautocの最小値が-88となっています。コードブックをよく見ると分かるのですが、実は-66、-77、-88はそれぞれ被介入、無政府、移行期を意味しています。それらの違いを無視することは問題があるものの、これらを欠損値に置き換えることがしばしばあります。

polity$democ[polity$democ %in% c(-66, -77, -88)] <- NA
polity$autoc[polity$autoc %in% c(-66, -77, -88)] <- NA
summary(polity)

##      ccode         country               year          democ       
##  Min.   :  2.0   Length:17395       Min.   :1800   Min.   : 0.000  
##  1st Qu.:200.0   Class :character   1st Qu.:1894   1st Qu.: 0.000  
##  Median :366.0   Mode  :character   Median :1959   Median : 1.000  
##  Mean   :404.2                      Mean   :1940   Mean   : 3.564  
##  3rd Qu.:630.0                      3rd Qu.:1991   3rd Qu.: 7.000  
##  Max.   :950.0                      Max.   :2017   Max.   :10.000  
##                                                    NA's   :776     
##      autoc       
##  Min.   : 0.000  
##  1st Qu.: 0.000  
##  Median : 4.000  
##  Mean   : 4.059  
##  3rd Qu.: 7.000  
##  Max.   :10.000  
##  NA's   :776

無事、欠損値が正しく認識されて、最小値が0になりました。

なお、欠損値が-や.のような文字列になっている場合、欠損値に入れ替えた後にas.numeric()などで数値データに変換するのを忘れないようにしましょう。

tidyverseのdplyrで欠損値を扱うときの方法を説明します。まず。元データにおける欠損値が空欄ではない場合は、データを読み込む際にオプションnaで指定します。

polity <- read_excel("data/p4v2017.xls", na = c("-66", "-77", "-88")) %>% 
  select(ccode, country, year, democ, autoc)
summary(polity)

##      ccode         country               year          democ       
##  Min.   :  2.0   Length:17395       Min.   :1800   Min.   : 0.000  
##  1st Qu.:200.0   Class :character   1st Qu.:1894   1st Qu.: 0.000  
##  Median :366.0   Mode  :character   Median :1959   Median : 1.000  
##  Mean   :404.2                      Mean   :1940   Mean   : 3.564  
##  3rd Qu.:630.0                      3rd Qu.:1991   3rd Qu.: 7.000  
##  Max.   :950.0                      Max.   :2017   Max.   :10.000  
##                                                    NA's   :776     
##      autoc       
##  Min.   : 0.000  
##  1st Qu.: 0.000  
##  Median : 4.000  
##  Mean   : 4.059  
##  3rd Qu.: 7.000  
##  Max.   :10.000  
##  NA's   :776

文字列として指定する点に注意してください。

あるいは、面倒くさいですがna_if()という関数を使うことも可能です。

polity <- read_excel("data/p4v2017.xls", na = c("-66", "-77", "-88")) %>% 
  select(ccode, country, year, democ, autoc) %>% 
  mutate(democ = na_if(democ, -66),
         democ = na_if(democ, -77),
         democ = na_if(democ, -88),
         autoc = na_if(autoc, -66),
         autoc = na_if(autoc, -77),
         autoc = na_if(autoc, -88))
summary(polity)

##      ccode         country               year          democ       
##  Min.   :  2.0   Length:17395       Min.   :1800   Min.   : 0.000  
##  1st Qu.:200.0   Class :character   1st Qu.:1894   1st Qu.: 0.000  
##  Median :366.0   Mode  :character   Median :1959   Median : 1.000  
##  Mean   :404.2                      Mean   :1940   Mean   : 3.564  
##  3rd Qu.:630.0                      3rd Qu.:1991   3rd Qu.: 7.000  
##  Max.   :950.0                      Max.   :2017   Max.   :10.000  
##                                                    NA's   :776     
##      autoc       
##  Min.   : 0.000  
##  1st Qu.: 0.000  
##  Median : 4.000  
##  Mean   : 4.059  
##  3rd Qu.: 7.000  
##  Max.   :10.000  
##  NA's   :776

6.2 欠損値を含む観察を除外する

後述するように基本的には欠損値がデータセットに含まれていても問題ないのですが、しばしば欠損値を含む観察を取り除く必要があります。

特定のベクトルのどの要素が欠損しているかを判定する関数はis.na()です。これは名前の通り、欠損している場合にTRUEが返ってくるので、次のようにすることでdemocとautocが欠損しているサンプルを除外できます。

no_na <- polity[!(is.na(polity$democ) | is.na(polity$autoc)),]
summary(no_na)

##      ccode         country               year          democ       
##  Min.   :  2.0   Length:16619       Min.   :1800   Min.   : 0.000  
##  1st Qu.:200.0   Class :character   1st Qu.:1892   1st Qu.: 0.000  
##  Median :365.0   Mode  :character   Median :1959   Median : 1.000  
##  Mean   :403.2                      Mean   :1940   Mean   : 3.564  
##  3rd Qu.:630.0                      3rd Qu.:1991   3rd Qu.: 7.000  
##  Max.   :950.0                      Max.   :2017   Max.   :10.000  
##      autoc       
##  Min.   : 0.000  
##  1st Qu.: 0.000  
##  Median : 4.000  
##  Mean   : 4.059  
##  3rd Qu.: 7.000  
##  Max.   :10.000

同様のことは次のようにもできます。

no_na <- polity[!is.na(polity$democ) & !is.na(polity$autoc),]
summary(no_na)

##      ccode         country               year          democ       
##  Min.   :  2.0   Length:16619       Min.   :1800   Min.   : 0.000  
##  1st Qu.:200.0   Class :character   1st Qu.:1892   1st Qu.: 0.000  
##  Median :365.0   Mode  :character   Median :1959   Median : 1.000  
##  Mean   :403.2                      Mean   :1940   Mean   : 3.564  
##  3rd Qu.:630.0                      3rd Qu.:1991   3rd Qu.: 7.000  
##  Max.   :950.0                      Max.   :2017   Max.   :10.000  
##      autoc       
##  Min.   : 0.000  
##  1st Qu.: 0.000  
##  Median : 4.000  
##  Mean   : 4.059  
##  3rd Qu.: 7.000  
##  Max.   :10.000

tidyverseで欠損値を除く場合はdrop_na()を使います。

polity %>% 
  drop_na() %>% 
  summary()

##      ccode         country               year          democ       
##  Min.   :  2.0   Length:16619       Min.   :1800   Min.   : 0.000  
##  1st Qu.:200.0   Class :character   1st Qu.:1892   1st Qu.: 0.000  
##  Median :365.0   Mode  :character   Median :1959   Median : 1.000  
##  Mean   :403.2                      Mean   :1940   Mean   : 3.564  
##  3rd Qu.:630.0                      3rd Qu.:1991   3rd Qu.: 7.000  
##  Max.   :950.0                      Max.   :2017   Max.   :10.000  
##      autoc       
##  Min.   : 0.000  
##  1st Qu.: 0.000  
##  Median : 4.000  
##  Mean   : 4.059  
##  3rd Qu.: 7.000  
##  Max.   :10.000

drop_na()の中で変数を指定することで、その変数が欠損している観察だけを除外することもできます。

polity %>% 
  drop_na(democ) %>% 
  summary()

##      ccode         country               year          democ       
##  Min.   :  2.0   Length:16619       Min.   :1800   Min.   : 0.000  
##  1st Qu.:200.0   Class :character   1st Qu.:1892   1st Qu.: 0.000  
##  Median :365.0   Mode  :character   Median :1959   Median : 1.000  
##  Mean   :403.2                      Mean   :1940   Mean   : 3.564  
##  3rd Qu.:630.0                      3rd Qu.:1991   3rd Qu.: 7.000  
##  Max.   :950.0                      Max.   :2017   Max.   :10.000  
##      autoc       
##  Min.   : 0.000  
##  1st Qu.: 0.000  
##  Median : 4.000  
##  Mean   : 4.059  
##  3rd Qu.: 7.000  
##  Max.   :10.000

6.3 欠損値を置き換える

is.na()を使うと欠損値を他の数値に書き換えることも可能です。例えば、欠損値を全て-1に置き換えるには次のようにします。

polity[is.na(polity$democ),]$democ <- -1
polity[is.na(polity$autoc),]$autoc <- -1
summary(polity)

##      ccode         country               year          democ       
##  Min.   :  2.0   Length:17395       Min.   :1800   Min.   :-1.000  
##  1st Qu.:200.0   Class :character   1st Qu.:1894   1st Qu.: 0.000  
##  Median :366.0   Mode  :character   Median :1959   Median : 1.000  
##  Mean   :404.2                      Mean   :1940   Mean   : 3.361  
##  3rd Qu.:630.0                      3rd Qu.:1991   3rd Qu.: 7.000  
##  Max.   :950.0                      Max.   :2017   Max.   :10.000  
##      autoc       
##  Min.   :-1.000  
##  1st Qu.: 0.000  
##  Median : 4.000  
##  Mean   : 3.834  
##  3rd Qu.: 7.000  
##  Max.   :10.000

tidyverseで欠損値を他の値に置き換えるときはreplace_na()を使います。

polity <- polity %>% 
  replace_na(list(democ = -1,
                  autoc = -1))
summary(polity)

##      ccode         country               year          democ       
##  Min.   :  2.0   Length:17395       Min.   :1800   Min.   :-1.000  
##  1st Qu.:200.0   Class :character   1st Qu.:1894   1st Qu.: 0.000  
##  Median :366.0   Mode  :character   Median :1959   Median : 1.000  
##  Mean   :404.2                      Mean   :1940   Mean   : 3.361  
##  3rd Qu.:630.0                      3rd Qu.:1991   3rd Qu.: 7.000  
##  Max.   :950.0                      Max.   :2017   Max.   :10.000  
##      autoc       
##  Min.   :-1.000  
##  1st Qu.: 0.000  
##  Median : 4.000  
##  Mean   : 3.834  
##  3rd Qu.: 7.000  
##  Max.   :10.000

複数の変数を一括で変更する場合はlist()で選択します。

6.4 欠損値がある場合の分析*

基本的に、欠損値がある場合でも分析をする上で問題はありません。

polity <- read_excel("data/p4v2017.xls") %>% 
  select(ccode, country, year, democ, autoc)
polity$democ[polity$democ %in% c(-66, -77, -88)] <- NA
polity$autoc[polity$autoc %in% c(-66, -77, -88)] <- NA

6.4.1 一変数の記述統計

例えば、平均値などの記述統計を求める際にはna.rm = TRUEというオプションをつけます。

mean(polity$democ)

## [1] NA

mean(polity$democ, na.rm = TRUE)

## [1] 3.564174

6.4.2 二変数の記述統計

相関係数や共分散など二変数以上の場合は若干ややこしくなります。 useオプションで以下の対処法のうち一つを決定します。

everything：欠損値を含む場合はNAを返す（デフォぱルト）。
all.obs：欠損値を含む場合はエラーとなる。
complete.obs：全ての変数が欠損していないものだけを使って計算する。
pairwise.complete.obs：二変数が欠損していないサンプルを使って計算する。

いまいちよくわからないのでデモデータを作って確認してみます。 irisの上から10のサンプルで品種以外のデータを使います。

demo <- iris[1:10,1:4]
demo[1,1] <- NA
demo[2,2] <- NA
demo[3,2] <- NA
demo

6.4.2.1 everythingの場合

欠損値を含まない変数Petal.LengthとPetal.WIdthのみ計算されています。

cor(demo, use = "everything")

##              Sepal.Length Sepal.Width Petal.Length Petal.Width
## Sepal.Length            1          NA           NA          NA
## Sepal.Width            NA           1           NA          NA
## Petal.Length           NA          NA    1.0000000   0.5216405
## Petal.Width            NA          NA    0.5216405   1.0000000

6.4.2.2 all.obsの場合

欠損値が含まれているのでエラーが起こります。

cor(demo, use = "all.obs")

## Error in cor(demo, use = "all.obs"): missing observations in cov/cor

6.4.2.3 complete.obsの場合

これはどの変数も欠損していないサンプル（つまり、4番目から10番目のサンプル）の相関係数になります。

cor(demo, use = "complete.obs")

##              Sepal.Length Sepal.Width Petal.Length Petal.Width
## Sepal.Length    1.0000000   0.8545604    0.7624091   0.4247114
## Sepal.Width     0.8545604   1.0000000    0.5684454   0.7269985
## Petal.Length    0.7624091   0.5684454    1.0000000   0.5385368
## Petal.Width     0.4247114   0.7269985    0.5385368   1.0000000

cor(demo[4:10,])

##              Sepal.Length Sepal.Width Petal.Length Petal.Width
## Sepal.Length    1.0000000   0.8545604    0.7624091   0.4247114
## Sepal.Width     0.8545604   1.0000000    0.5684454   0.7269985
## Petal.Length    0.7624091   0.5684454    1.0000000   0.5385368
## Petal.Width     0.4247114   0.7269985    0.5385368   1.0000000

6.4.2.4 pairwise.completeの場合

先ほどと一致しているところと一致していないところがあります。

cor(demo, use = "pairwise.complete")

##              Sepal.Length Sepal.Width Petal.Length Petal.Width
## Sepal.Length    1.0000000   0.8545604    0.6853616   0.4225771
## Sepal.Width     0.8545604   1.0000000    0.4845437   0.6915641
## Petal.Length    0.6853616   0.4845437    1.0000000   0.5216405
## Petal.Width     0.4225771   0.6915641    0.5216405   1.0000000

例えば、Sepal.LengthとSepal.Widthの場合、どちらも欠損していないのは4番目から10番目のサンプルになります。なので、complete.obsのときと同じ結果になります。

cor(demo[4:10,]$Sepal.Length, demo[4:10,]$Sepal.Width)

## [1] 0.8545604

しかし、Petal.LengthとPetal.Widthの場合、欠損しているサンプルはないので、1番目から10番目のサンプル全てを使った値になります。

cor(demo$Petal.Length, demo$Petal.Width)

## [1] 0.5216405

everythingのときと同じであることも分かります。

6.4.3 回帰分析

回帰分析の際には欠損値を含むサンプルを除外した上で計算します。

summary(lm(democ ~ autoc, data = polity))

## 
## Call:
## lm(formula = democ ~ autoc, data = polity)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.5316 -1.5316  0.2009  2.1118  2.5576 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  7.442414   0.023232   320.4   <2e-16 ***
## autoc       -0.955418   0.004321  -221.1   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.964 on 16617 degrees of freedom
##   (776 observations deleted due to missingness)
## Multiple R-squared:  0.7463, Adjusted R-squared:  0.7463 
## F-statistic: 4.888e+04 on 1 and 16617 DF,  p-value: < 2.2e-16

776 observations deleted due to missingnessと表示されています。

Rで計量政治学入門

第6章欠損値の扱い

6.1 欠損値を読み込む

6.2 欠損値を含む観察を除外する

6.3 欠損値を置き換える

6.4 欠損値がある場合の分析*

6.4.1 一変数の記述統計

6.4.2 二変数の記述統計

6.4.2.1 everythingの場合

6.4.2.2 all.obsの場合

6.4.2.3 complete.obsの場合

6.4.2.4 pairwise.completeの場合

6.4.3 回帰分析

6.5 欠損値のクラス*

第6章 欠損値の扱い

6.1 欠損値を読み込む

6.2 欠損値を含む観察を除外する

6.3 欠損値を置き換える

6.4 欠損値がある場合の分析*

6.4.1 一変数の記述統計

6.4.2 二変数の記述統計

6.4.2.1 everythingの場合

6.4.2.2 all.obsの場合

6.4.2.3 complete.obsの場合

6.4.2.4 pairwise.completeの場合

6.4.3 回帰分析

6.5 欠損値のクラス*

第6章欠損値の扱い