例題

以下A国とB国の市民の、「収入」、「性別」に関するデータがあるとする。はじめに、A国のデータのみを取り出し、収入の平均と分散、性別の度数分布表を計算する。次に、収入のヒストグラム収入と性別の箱ひげ図を描く。

0.1 データの読み込みと準備

# データを読み込む
dat <- data.frame(
  income = c(150, 213, 221, 284, 314, 240, 609, 612, 703, 1660, 
             487, 493, 494, 497, 496, 502, 502, 502, 516, 517),
  nation = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
             2, 2, 2, 2, 2, 2, 2, 2, 2, 2),
  gender = c(1, 2, 2, 1, 1, 1, 2, 1, 1, 2, 
             2, 2, 2, 1, 2, 1, 1, 1, 2, 2)
)

# 因子型の変数を定義する
dat$nation <- factor(dat$nation, 
                     levels = c(1, 2), labels = c("A", "B")) # 1 = A国, 2 = B国
dat$gender <- factor(dat$gender, 
                     levels = c(1, 2), 
                     labels = c("male", "female")) #  1 = 男性、2 = 女性

データの抽出

# A国のデータのみ取り出す
dat_A <- subset(dat, subset = dat$nation=="A")

なお、上述のスクリプトではsubset関数を用いているが、以下のようにデータフレームから選択条件を工夫することで、直接取り出すこともできる。

dat_A <- dat[dat$nation=="A", ]   # nationがA国である行のみ取り出すという意味

1 記述統計

1.1 方法1: summary関数を用いる。

summary関数を用いると簡便に記述統計を出すことができる。

summary(dat_A)　# 平均、中央値

##      income       nation    gender 
##  Min.   : 150.0   A:10   male  :6  
##  1st Qu.: 225.8   B: 0   female:4  
##  Median : 299.0                    
##  Mean   : 500.6                    
##  3rd Qu.: 611.2                    
##  Max.   :1660.0

1.2 方法2: 個別に計算する

# 収入
mean(dat_A$income)　　　# 平均
var(dat_A$income)   　　# 分散
sd(dat_A$income)    　　# 標準偏差
max(dat_A$income)    　 # 最大値
min(dat_A$income)       # 最小値
median(dat_A$income)    # 中央値

# 性別
tab_A_gender <- table(dat_A$gender)     # 度数分布
tab_A_gender
addmargins(tab_A_gender)                # 周辺分布付き
prop.table(tab_A_gender)                # 割合表示


#　特定の変数のみから直接テーブルをマニュアルで作りたい場合の例
data.frame(
  mean = mean(dat_A$income),
  sd   = mean(dat_A$income)
)

1.3 方法3: describe関数を用いる

psychパッケージのdescribe関数を用いることもできる。ただしこちらの関数では、質的変数（ここではgender）についてはうまく扱えず、そのまま量的変数として扱われる。

library(psych)

describe(dat_A)             # 因子型もそのまま扱われる。

##         vars  n  mean     sd median trimmed    mad min  max range skew kurtosis
## income     1 10 500.6 452.50    299  399.50 174.21 150 1660  1510 1.56     1.40
## nation*    2 10   1.0   0.00      1    1.00   0.00   1    1     0  NaN      NaN
## gender*    3 10   1.4   0.52      1    1.38   0.00   1    2     1 0.35    -2.05
##             se
## income  143.09
## nation*   0.00
## gender*   0.16

2 記述統計のプロット

2.1 ヒストグラム

# 収入のヒストグラム
hist(dat_A$income, 
     main = "A国の収入",   #　題名 
     xlab = "収入（円）",  #  x軸名
     ylab = "n")　　　　　 #  y軸名

より凝ったヒストグラムの例

hist(dat_A$income,
     main = "A国の収入",
     xlab = "収入（円）",
     ylab = "n",
     right = FALSE,           #　右側を「未満」まで含むようにする
     breaks = c(0,1500,3000), # 「0以上1500未満」、「1500以上3000未満」に分割
     ylim = c(0,10),　　　　　# 　yの範囲（０～１０）
     label = TRUE,　　　　　　# 　棒の上に数字を入れる　　　
     col = "BLUE"             # 　色
)

2.2 箱ひげ図

# 箱ひげ図（y = 収入；ｘ＝性別）
boxplot(dat_A$income ~ dat_A$gender, 
        ylab = "収入", 
        xlab = "性別", 
        names=c("男性", "女性"),　# x軸のラベル
        main = "A国における男女の収入の箱ひげ図")

2.3 ggplotを用いたプロット

上記ではRのデフォルトのプロット関数を使ったが、ggplotを使うとより美しいプロットが得られる。

# ライブラリを読み込む
library(ggplot2)

ヒストグラム

ggplot(data = dat_A,
       aes(x = income)) +  # まずggplotでデータを参照し、次にグラフの種類を+でつないで指定する。
  geom_histogram(bins = 3)+ # この数を変えるといくつの棒に区切るかを変更できる。
　labs(title =  "A国の収入の分布", x = "収入", y = "度数")

バイオリンプロット

ボックスプロットはgeom_boxplot()でできるが、ggplotの場合にはバイオリンプロットも可能(点の分布密度を凹凸で表現。さらに実際の点をgeom_point()で重ね書きすることができる。

ggplot(data = dat_A,                # ggplotでデータを参照
       aes(y = income, x = gender), 
       fill = gender) +
  geom_violin() + 
  labs(title = "A国における男女の収入の箱ひげ図", x = "性別", y = "収入") + 
  scale_x_discrete(labels = c("male"="男性", "female"="女性")) +
  geom_point(position = position_jitter(width = 0.1, seed = 1))

記述統計の計算とグラフ化

佐藤圭一

2023-09-12

例題