【社会情報処理】第4回 基本統計量

基本統計量とは分布の特徴を表すのに使われる基本的な統計量である。

• 分布の中心を表す統計量：平均、中央値、最頻値
• 分布の広がりを表す統計量：分散、標準偏差
• その他の統計量：最大値、最小値、四分位数

以下の数式では $X_{1},X_{2},\cdots,X_{n}$ の $n$ 個の観測値のデータがあるとする。

R の例として用いる標本データは、データセット airquality の Ozone から $n=21$ 個の標本を抽出したものとする。

> n <- 21
> set.seed(12); (data <- sample(na.omit(airquality$Ozone), n))
 [1]  35  73   9  91  27  89   7  47  12   8  30 168  13  32  18  48  18  20 122  45  52

set.seed() は乱数の種（seed）を設定するためのもので、毎回同じ標本抽出を再現するために使用している。

それぞれの観測値がある範囲（階級）で何回観測されたかをカウントした数を度数といい、それぞれの階級で上限値と下限値を足して2で割ったものを階級値という。 度数を全データ数で割ったものを相対度数という。

階級ごとの度数のデータを級別データといい、級別データを表にしたものが度数分布表である。

度数分布表の例（階級幅は 20）

横軸に階級、縦軸に度数をとって描いた柱状グラフをヒストグラムという。

R でヒストグラムを描くには関数 hist() を使う。

> hg <- hist(data)

戻り値は変数 hg に代入しているが、代入しなくてもヒストグラムは描ける。

デフォルトでは観測値がその階級の下限値より大きい、かつ上限値以下の場合に度数としてカウントする。 \[ 下限値<観測値\le 上限値 \]

日本では、観測値がその階級の下限値以上、かつ上限値未満の場合に度数としてカウントするのが一般的である。 \[ 下限値\le 観測値<上限値 \] これを行うには、引数に right = FALSE を入れる（breaks をベクトルで与える場合は、include.lowest = TRUE が必要かもしれない）。

> hg <- hist(data, right = FALSE)

階級数はデフォルトではスタージェスの公式（$k$ は階級数、$n$ はデータ数） \[ k\simeq \log_{2}n+1 \] を使って決められるが、階級の度数次第で微調整されるようである。 スタージェスの公式は関数 nclass.Sturges() で計算できる。

> nclass.Sturges(data)

hist() の戻り値には以下の要素が含まれる。

最大値

> max(data)

最小値

> min(data)

平均は分布の重心を表す統計量である。

\[ \bar{X}=\frac{1}{n}\sum_{i=1}^{n}X_{i} \]

関数 mean() で計算できる。

> mean(data)

合計の関数 sum() を使うと以下のようになる。

> sum(data) / length(data)

中央値は分布をちょうど二等分する点の値である。 $X_{i}$ が小さい順（もしくは大きい順）に並んでいるとすると、定義は以下のようになる。 \[ \mathrm{Me}= \left\{ \begin{array}{ll} X_{(n+1)/2} & （n が奇数）\\ (X_{n/2}+X_{n/2+1})/2 & （n が偶数） \end{array} \right. \]

例えば、$n=5$ だと \[ X_{1},X_{2},X_{3},X_{4},X_{5} \] の中で中央値は $\mathrm{Me}=X_{3}$ になる。

$n=4$ だと \[ X_{1},X_{2},X_{3},X_{4} \] の中で中央値は $\mathrm{Me}=(X_{2}+X_{3})/2$ になる。

関数 median() で求められる。

> median(data)

最頻値は級別データで度数が最大となるところの階級値である。 ヒストグラムの最も山が高いところでもある。

hist() の戻り値 hg を使って以下のようにして求められる。

> hg$mids[which.max(hg$counts)]

関数 which.max() はデータの中で最大値の要素の番号を返す。

Moodle Server（非公式）で第4回の小テストを受験しなさい。

四分位数 は分布を四等分したときの三つの分割点の値である。 $X_{i}$ を小さい順に並べて四等分し、その三つの分割点を小さい順に 第一四分位数 $Q_{1}$、 第二四分位数 $Q_{2}$、 第三四分位数 $Q_{3}$ という。 $Q_{1}$ は 25%パーセンタイル、$Q_{2}$ は 50%パーセンタイル、$Q_{3}$ は 75%パーセンタイルにあたる¹⁾ 。

データを小さい順に並べたときの分布を四分位数を使って「箱」と「ひげ」の長さで表したグラフを箱ひげ図という。 データがどの付近に集中して分布しているのかを見やすくするためのものである。

四分位数の定義には複数あるが、ここでは二つの定義を挙げる。

John Wilder Tukey が提案した定義である。

1. データの中央値を $Q_{2}$ としてデータを半分に分ける。
   • $n$ が奇数なら $Q_{2}$ を両半分に含める。
   • $n$ が偶数ならちょうど半分に分けられる。
2. 小さい方の半分の中央値を $Q_{1}$ とする。
3. 大きい方の半分の中央値を $Q_{3}$ とする。

これに最小値 $Q_{0}$ と最大値 $Q_{4}$ を加えて五数要約という。

五数要約は fivenum() で計算できる（もしくは summary() や quantile()）。

> fivenum(data)

このテューキーの五数要約による箱ひげ図は boxplot() で簡単に描くことができる。

> boxplot(data)

デフォルトでは縦方向の箱ひげ図だが、引数に horizontal = TRUE を入れると横方向の箱ひげ図になる。

データのプロットと並べると以下のようになる。 データのプロットは見やすいように jitter() を使って y 方向にランダムに少しずらしてある。

> layout(1 : 2)
> plot(x = data, y = jitter(numeric(length(data))), ylim = c(-0.2, 0.2), ylab = "", yaxt = "n", xlab = "オゾン濃度[ppb]", main = "データ")
> boxplot(data, horizontal = TRUE, xlab = "オゾン濃度[ppb]", main = "テューキーの方法")
> layout(1)

IQR（Inter Quantile Range） は四分位範囲（箱の長さ）である。 この中にはデータの 50% が含まれている。 \[ \mathrm{IQR}=Q_{3}-Q_{1} \]

ひげの下端は $Q_{1}$ から $1.5\times\mathrm{IQR}$ 以内の範囲にあるデータの最小値、 ひげの上端は $Q_{3}$ から $1.5\times\mathrm{IQR}$ 以内の範囲にあるデータの最大値である。 ひげの下端より小さいデータや上端より大きいデータは外れ値として○で描かれる。

外れ値がなければ、ひげの下端は $Q_{0}$、上端は $Q_{4}$ に等しい。

1. データの中央値を $Q_{2}$ としてデータを半分に分ける。
   • $n$ が奇数なら $Q_{2}$ を両半分に含めない。
   • $n$ が偶数ならちょうど半分に分けられる。
2. 小さい方の半分の中央値を $Q_{1}$ とする。
3. 大きい方の半分の中央値を $Q_{3}$ とする。

R ではこの方法で四分位数を求める関数は組み込まれていないが、以下のようにして計算できる。

> q <- quantile(data, type = ifelse(length(data) %% 2 == 0, 5, 6))

中学校指導要領で定義されている箱ひげ図は、 箱の上端を $Q_{3}$、中の横線を $Q_{2}$、下端を $Q_{1}$ とし、 ひげの上端を最大値、下端を最小値にする。 この定義だと外れ値はない。

データのプロットと箱ひげ図を並べると以下のようになる。 boxplot() だとテューキーの方法になってしまうので、代わりに $Q_{1},Q_{2},Q_{3}$ などを指定してbxp() で箱ひげ図を描いている。

> layout(1 : 2)
> plot(x = data, y = jitter(numeric(length(data))), ylim = c(-0.2, 0.2), ylab = "", yaxt = "n", xlab = "オゾン濃度[ppb]", main = "データ")
> bp.japan <- list(stats = matrix(q, 5, 1), n = length(data))
> bxp(bp.japan, horizontal = TRUE, xlab = "オゾン濃度[ppb]", main = "中学校指導要領の方法")
> layout(1)

boxplot() は複数の種類のデータをまとめて箱ひげ図として描画できる。

次のプログラム boxplot.R は第2回に e-Stat からダウンロードした山形県の市町村データ（ファイル名は yamagata.csv にしている）において、2015年度の第1次、第2次、第3次就業者比率をまとめて箱ひげ図で表示するプログラムである。

boxplot.R

library(readr)
 
# データの読み込み
df <- read_csv("yamagata.csv", locale = locale(encoding = "SHIFT-JIS"))
 
df$`/項目` <- NULL
df$調査年 <- as.numeric(sapply(df$調査年, gsub, pattern="年度", replacement = ""))
df <- df[order(df$調査年),]
 
# 2015年度の就業者比率
data <- df[df$調査年 == 2015,
           c("#F01201_第1次産業就業者比率【％】",
             "#F01202_第2次産業就業者比率【％】",
             "#F01203_第3次産業就業者比率【％】")]
 
# ラベルの書き換え
colnames(data) <- c("第1次産業", "第2次産業", "第3次産業")
 
# 箱ひげ図
boxplot(data, horizontal = TRUE, main = "山形県市町村の就業者比率", xlab = "就業者比率【％】")

分散は分布が平均からどれくらい広がっているかを表す統計量である。

\[ \sigma^{2}=\frac{1}{n}\sum_{i=1}^{n}(X_{i}-\bar{X})^{2}=s^{2}\left(1-\frac{1}{n}\right) \]

R には分散を直接計算する関数は組み込まれていないので以下で計算する。

> var(data) * (1 - 1 / length(data)) 

標本から計算し、母分散の推定に使われる。

\[ s^{2}=\frac{1}{n-1}\sum_{i=1}^{n}(X_{i}-\bar{X})^{2} \]

> var(data)

標準偏差は分散の平方根をとったものである。

\[ \sigma=\sqrt{\sigma^{2}}=s\sqrt{1-\frac{1}{n}} \]

> sd(data) * sqrt(1 - 1 / length(data))

\[ s=\sqrt{s^{2}} \]

> sd(data)

もしくは

> sqrt(var(data))

Moodle Server（非公式） で第4回の課題を行いなさい。

締め切り：2021年7月15日（木）20時