第9章 データの要約

データを読み込んで最初にすることは、データの概形を掴むことだと思います。 分析に入る前に、データがどのような特徴を持っているのか、データを作成した場合はミスをしていないかを確認するべきでしょう。

データから代表的な値を計算することでデータを要約することができます。 このような方法を記述統計 (descriptive statistics) と呼びます。

もう一つのアプローチはグラフなどを描いて視覚的に把握することですが、これは後述します。

library(tidyverse)
library(skimr)
library(summarytools)
## Registered S3 method overwritten by 'pryr':
##   method      from
##   print.bytes Rcpp
## For best results, restart R session and update pander using devtools:: or remotes::install_github('rapporter/pander')
## 
## Attaching package: 'summarytools'
## The following object is masked from 'package:tibble':
## 
##     view
tri <- read_csv("data/triangle.csv")
## Parsed with column specification:
## cols(
##   statea = col_double(),
##   stateb = col_double(),
##   year = col_double(),
##   dependa = col_double(),
##   dependb = col_double(),
##   demauta = col_double(),
##   demautb = col_double(),
##   allies = col_double(),
##   dispute1 = col_double(),
##   logdstab = col_double(),
##   lcaprat2 = col_double(),
##   smigoabi = col_double(),
##   opena = col_double(),
##   openb = col_double(),
##   minrpwrs = col_double(),
##   noncontg = col_double(),
##   smldmat = col_double(),
##   smldep = col_double(),
##   dyadid = col_double()
## )

9.1 データを見る

データを把握する最もナイーブな方法は生のデータを見ることです。 例えば、.csv.xlsファイルであれば適当な表計算ソフトで開くことで見ることができます。 また、View()という関数でRに取り込んだデータを見ることもできます。

View(tri)

しかし、大抵の場合、観察数が多すぎてデータをそのまま見てもなんだかよく分からないことが多いです。

  • データが大きすぎると開くのに時間がかかることもあります。

まずは、head()によって、データの冒頭だけ見るのがよいでしょう。 これによって、データの形式が想定通りか暫定的にチェックできます。

head(tri)
  • tail()を使うとデータの最後を見ることがきでます。

9.2 記述統計

記述統計とは変数をいくつかの指標によって表現することを言います。

9.2.1 平均

平均 (mean) は次のように定義されます。

\[ \bar{x} = \sum_{i=1}^n \frac{x_i}{n} \]

Rでは、まずデータフレームから変数を抜き出し、mean()に入れます。

mean(tri$smldep)
## [1] 0.002328343
  • ちなみに、欠損値がある場合はNAとなるので、na.rm = TRUEというオプションを付けます。
  • smldepは二国間の相互依存度を示しています。

tidyverse風に書くなら、次のようになります。

tri %>% 
  pull(smldep) %>% 
  mean()
## [1] 0.002328343
  • pull()でベクトルを抜き出すことができます。

9.2.2 分散と標準偏差

分散 (variance) は次のように定義され、変数の散らばりを意味します。

\[ \sigma^2 = \sum_{i=1}^n \frac{(x - \bar{x})^2}{n} \]

Rではvar()によって求めます。

  • 厳密に言えば、不偏分散と呼ばれるもので、\(n\)の代わりに\(n-1\)で割られています。
var(tri$smldep)
## [1] 5.220904e-05

分散の平方根\(\sigma\)標準偏差 (standard error) と呼びsd()で計算します。

sd(tri$smldep)
## [1] 0.007225582

実際に二乗すると分散と同じ値になります。

sd(tri$smldep)^2
## [1] 5.220904e-05

9.2.3 分位点

\(q\)%分位点 (quantile) とは、下から数えて\(q\)%に当たる点を意味します。 例えば、50%分位点は、中央値 (median) とも呼ばれ、上から見ても下から見てもちょうど真ん中に位置します。

中央値はmedian()で求めることができ、一般的に分位点はquantile()で求めることができます。

median(tri$smldep)
## [1] 0.0002466768
quantile(tri$smldep, 0.5)
##          50% 
## 0.0002466768
  • どっちも同じであることが分かります。

quantile()は一度に複数の分位点を求めることができます。

quantile(tri$smldep, c(0.25, 0.5, 0.75))
##          25%          50%          75% 
## 9.789278e-06 2.466768e-04 1.486773e-03
  • 25%と75%分位点は第1四分位点、第3四分位点とも呼ばれ、中央値と共に報告されることが多いです。

9.2.4 記述統計表

通常は、特定の変数だけの記述統計を報告するということはなく、分析で使用する変数全てについて報告します。 summary()にデータフレームを入れると、全体の記述統計を出してくれます。

summary(tri)
##      statea        stateb           year         dependa         
##  Min.   :  2   Min.   : 20.0   Min.   :1885   Min.   :0.0000000  
##  1st Qu.:100   1st Qu.:300.0   1st Qu.:1935   1st Qu.:0.0000233  
##  Median :212   Median :438.0   Median :1967   Median :0.0006987  
##  Mean   :248   Mean   :473.5   Mean   :1957   Mean   :0.0068491  
##  3rd Qu.:365   3rd Qu.:710.0   3rd Qu.:1979   3rd Qu.:0.0038495  
##  Max.   :900   Max.   :950.0   Max.   :1991   Max.   :0.8484982  
##                                                                  
##     dependb             demauta           demautb            allies      
##  Min.   :0.0000000   Min.   :-10.000   Min.   :-10.000   Min.   :0.0000  
##  1st Qu.:0.0000317   1st Qu.: -7.000   1st Qu.: -7.000   1st Qu.:0.0000  
##  Median :0.0019763   Median :  6.000   Median : -5.000   Median :0.0000  
##  Mean   :0.0175759   Mean   :  2.538   Mean   : -1.165   Mean   :0.2155  
##  3rd Qu.:0.0153196   3rd Qu.: 10.000   3rd Qu.:  8.000   3rd Qu.:0.0000  
##  Max.   :0.6529297   Max.   : 10.000   Max.   : 10.000   Max.   :1.0000  
##                                                                          
##     dispute1        logdstab        lcaprat2           smigoabi     
##  Min.   :0.000   Min.   :1.872   Min.   :0.000103   Min.   :  0.00  
##  1st Qu.:0.000   1st Qu.:6.740   1st Qu.:1.486062   1st Qu.: 13.00  
##  Median :0.000   Median :8.042   Median :2.830705   Median : 22.00  
##  Mean   :0.047   Mean   :7.626   Mean   :2.978976   Mean   : 25.33  
##  3rd Qu.:0.000   3rd Qu.:8.630   3rd Qu.:4.246247   3rd Qu.: 34.00  
##  Max.   :1.000   Max.   :9.407   Max.   :9.494078   Max.   :130.00  
##                                                                     
##      opena            openb           minrpwrs         noncontg     
##  Min.   :0.0008   Min.   :0.0008   Min.   :0.0000   Min.   :0.0000  
##  1st Qu.:0.1261   1st Qu.:0.1149   1st Qu.:0.0000   1st Qu.:0.0000  
##  Median :0.2346   Median :0.2439   Median :0.0000   Median :1.0000  
##  Mean   :0.3083   Mean   :0.3237   Mean   :0.2356   Mean   :0.6307  
##  3rd Qu.:0.4263   3rd Qu.:0.4550   3rd Qu.:0.0000   3rd Qu.:1.0000  
##  Max.   :1.9985   Max.   :1.9985   Max.   :1.0000   Max.   :1.0000  
##  NA's   :1204     NA's   :2731                                      
##     smldmat            smldep              dyadid      
##  Min.   :-10.000   Min.   :0.000e+00   Min.   :  2020  
##  1st Qu.: -8.000   1st Qu.:9.790e-06   1st Qu.:100365  
##  Median : -7.000   Median :2.467e-04   Median :212220  
##  Mean   : -3.311   Mean   :2.328e-03   Mean   :248516  
##  3rd Qu.:  1.000   3rd Qu.:1.487e-03   3rd Qu.:365439  
##  Max.   : 10.000   Max.   :1.719e-01   Max.   :900910  
## 

skimrskim()のほうがカッコよく見えます。

skim(tri)
表 9.1: Data summary
Name tri
Number of rows 39996
Number of columns 19
_______________________
Column type frequency:
numeric 19
________________________
Group variables None

Variable type: numeric

skim_variable n_missing complete_rate mean sd p0 p25 p50 p75 p100 hist
statea 0 1.00 248.04 195.04 2.00 100.00 212.00 365.00 900.00 ▇▇▃▂▁
stateb 0 1.00 473.48 223.77 20.00 300.00 438.00 710.00 950.00 ▃▇▃▆▂
year 0 1.00 1957.43 27.60 1885.00 1935.00 1967.00 1979.00 1991.00 ▂▂▂▅▇
dependa 0 1.00 0.01 0.02 0.00 0.00 0.00 0.00 0.85 ▇▁▁▁▁
dependb 0 1.00 0.02 0.04 0.00 0.00 0.00 0.02 0.65 ▇▁▁▁▁
demauta 0 1.00 2.54 7.70 -10.00 -7.00 6.00 10.00 10.00 ▅▁▁▁▇
demautb 0 1.00 -1.16 7.52 -10.00 -7.00 -5.00 8.00 10.00 ▇▂▂▁▅
allies 0 1.00 0.22 0.41 0.00 0.00 0.00 0.00 1.00 ▇▁▁▁▂
dispute1 0 1.00 0.05 0.21 0.00 0.00 0.00 0.00 1.00 ▇▁▁▁▁
logdstab 0 1.00 7.63 1.20 1.87 6.74 8.04 8.63 9.41 ▁▁▂▅▇
lcaprat2 0 1.00 2.98 1.86 0.00 1.49 2.83 4.25 9.49 ▇▇▅▂▁
smigoabi 0 1.00 25.33 17.57 0.00 13.00 22.00 34.00 130.00 ▇▅▁▁▁
opena 1204 0.97 0.31 0.25 0.00 0.13 0.23 0.43 2.00 ▇▂▁▁▁
openb 2731 0.93 0.32 0.28 0.00 0.11 0.24 0.46 2.00 ▇▃▁▁▁
minrpwrs 0 1.00 0.24 0.42 0.00 0.00 0.00 0.00 1.00 ▇▁▁▁▂
noncontg 0 1.00 0.63 0.48 0.00 0.00 1.00 1.00 1.00 ▅▁▁▁▇
smldmat 0 1.00 -3.31 6.68 -10.00 -8.00 -7.00 1.00 10.00 ▇▂▁▁▂
smldep 0 1.00 0.00 0.01 0.00 0.00 0.00 0.00 0.17 ▇▁▁▁▁
dyadid 0 1.00 248516.26 195149.01 2020.00 100365.00 212220.00 365439.00 900910.00 ▇▇▃▂▁

ただし、これをそのまま論文などに載せるのはカッコよくありません。 そこで、summarytoolsというパッケージを使います。

  • データフレームにすることにこだわらないのであれば、stargazerというパッケージでもいいと思います。

descr()という関数でmarkdownのテーブル形式で表示します。

descr(tri, transpose = TRUE)
## Warning: `funs()` is deprecated as of dplyr 0.8.0.
## Please use a list of either functions or lambdas: 
## 
##   # Simple named list: 
##   list(mean = mean, median = median)
## 
##   # Auto named with `tibble::lst()`: 
##   tibble::lst(mean, median)
## 
##   # Using lambdas
##   list(~ mean(., trim = .2), ~ median(., na.rm = TRUE))
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_warnings()` to see where this warning was generated.

その結果をさらにtb()に入れることで、データフレームにできるので、.csvファイルなどとして書き出すこともできます。

tri %>% 
  descr() %>% 
  tb()

また、dfSummary()view()でファンシーに表示してくれます。

tri %>% 
  dfSummary() %>% 
  summarytools::view()
  • summarytools::view()としているのは、tidyrにもview()という関数があるので、summarytoolsのもであることを明示するためです。

9.3 データのカウント

例えば、戦争が起こった・起こらなかった、のようにカテゴリカル変数の場合、上記の記述統計を用いるより、数えてしまう方が良いことがあります。 table()で数えることができます。

table(tri$dispute1)
## 
##     0     1 
## 38116  1880

table()に複数の便数を入れることで、クロス表を作ることもできます。

table(tri$dispute1, tri$noncontg)
##    
##         0     1
##   0 13492 24624
##   1  1279   601

tidyverseにはcount()という関数があり、数えることができます。

tri %>% 
  count(dispute1)
tri %>% 
  count(dispute1, noncontg)
  • count()の結果をクロス表にするにはpivot_wider()を使います。

summarytools()であれば、freq()ctable()を使います。

freq(tri$dispute1)
## Frequencies  
## tri$dispute1  
## Type: Numeric  
## 
##                Freq   % Valid   % Valid Cum.   % Total   % Total Cum.
## ----------- ------- --------- -------------- --------- --------------
##           0   38116     95.30          95.30     95.30          95.30
##           1    1880      4.70         100.00      4.70         100.00
##        <NA>       0                               0.00         100.00
##       Total   39996    100.00         100.00    100.00         100.00
ctable(tri$dispute1, tri$noncontg)
## Cross-Tabulation, Row Proportions  
## dispute1 * noncontg  
## Data Frame: tri  
## 
## ---------- ---------- --------------- --------------- ----------------
##              noncontg               0               1            Total
##   dispute1                                                            
##          0              13492 (35.4%)   24624 (64.6%)   38116 (100.0%)
##          1               1279 (68.0%)     601 (32.0%)    1880 (100.0%)
##      Total              14771 (36.9%)   25225 (63.1%)   39996 (100.0%)
## ---------- ---------- --------------- --------------- ----------------

9.4 グループごとの要約

特定のグループごとに記述統計を求めたい場合があるとします。 このような場合、tidyversegroup_by()summarise()を組み合わせます。

例えば、次のようにして各年のsmldepの平均値を求めることができます。

tri %>% 
  group_by(year) %>% 
  summarise(mean_smldep = mean(smldep))
## `summarise()` ungrouping output (override with `.groups` argument)

group_by()summarytoolsでも使えます。 紛争が起こっているかどうかでグループ分けをします。

tri %>% 
  group_by(dispute1) %>% 
  descr() %>% 
  tb()

9.5 ユニークなデータ*