====== 【社会情報処理】第5回 母平均の仮説検定 ======

===== 中心極限定理 =====
母平均 $\mu$、母分散 $\sigma^{2}$ の母集団から $n$ 個の標本 $X_{1}, X_{2},\cdots, X_{n}$ を無作為抽出する。
どのような母集団((コーシー分布のように中心極限定理が成り立たない分布もあるが、ほとんどの場合は中心極限定理が成り立つと考えてよい。))でも、標本のサイズ $n$ が大きければ標本平均 $\bar{X}$ の標本分布は（近似的に）$N(\mu,\sigma^{2}/n)$ になる。
これを__**中心極限定理**__という。
母集団が正規分布の場合は $n$ に関わらず厳密だが、それ以外の場合は $n$ が大きいほど近似は良くなる。

===== 仮説検定とは =====

__**仮説検定**__とは、母集団について仮定された命題（仮説）を標本に基づき統計的に検証することである。
検証したい仮説（__**対立仮説**__）とそれに反する仮説（__**帰無仮説**__)を立て、もし帰無仮説が正しいと仮定するならばまず起こらない結果が標本から得られることで **<fc #ff0000>帰無仮説を棄却して対立仮説を採択すること</fc>**を目的とする。

この、まず起こらない結果の範囲を__**棄却域**__といい、__**有意水準**__ $\alpha$ という小さな確率で決められる。
$\alpha$ は恣意的に与える値で、0.1（10%）, 0.05（5%）, 0.01（1%）がよく使われる。
帰無仮説が正しいと仮定すると、標本から得られた検定統計量の値が棄却域の範囲に入る確率が $\alpha$ になる。
$\alpha$ を小さく取るほど棄却域の範囲に入る可能性は低くなるが、もしそんな小さい確率でも入ったのであれば、帰無仮説が正しいという仮定が間違っていたと考えていいだろうということで、帰無仮説を棄却するという判断をする。

仮説検定は確率的に判断するので、この判断が間違いがである可能性はゼロではなく、次の二通りの間違いが起こりうる。
  * __**第一種の過誤**__：帰無仮説が正しいのに帰無仮説を棄却してしまう
  * __**第二種の過誤**__：帰無仮説は正しくないのに帰無仮説を棄却できない
 
$\alpha$ が小さいほど第一種の過誤が起こりにくいが、逆に第二種の過誤が起こりやすくなってしまうので、小さいほどよいということにならない。
第一種の過誤と第二種の過誤の両方を起こりにくくしたい場合は、$\alpha$ を小さくし、かつ、標本サイズ $n$ を大きくすることが必要となる。



===== 母平均の仮説検定 =====

ここでは母平均の仮説検定を行う。

母平均 $\mu$、母分散 $\sigma^{2}$ の母集団から $n$ 個の標本 $X_{1},X_{2},\cdots,X_{n}$ を抽出したとする。

母平均の仮説検定に用いる検定統計量は標本平均 $\bar{X}$ である。
\[
\bar{X}=\frac{1}{n}\sum_{i=1}^{n}X_{i}
\]

----
==== 仮説検定の手順 ====

仮説検定の手順は以下の通り。

=== 1. 仮説の設定 ===

  * 右片側検定の場合
    * 帰無仮説 $H_{0}:\mu=\mu_{0}$
    * 対立仮説 $H_{1}:\mu>\mu_{0}$
  * 左片側検定の場合
    * 帰無仮説 $H_{0}:\mu=\mu_{0}$
    * 対立仮説 $H_{1}:\mu<\mu_{0}$
  * 両側検定の場合
    * 帰無仮説 $H_{0}:\mu=\mu_{0}$
    * 対立仮説 $H_{1}:\mu\neq\mu_{0}$

ここで $\mu_{0}$ は定数である。

=== 2. 検定統計量と標本分布の決定 ===
母平均の仮説検定に用いる検定統計量は $\bar{X}$ である。
棄却域を求めるためには $\bar{X}$ の標本分布を知らなくてはならない。
中心極限定理によれば $\bar{X}$ の標本分布は（近似的に）$N(\mu,\sigma^{2}/n)$ である。

母分散が既知の場合は正規化して標本分布 $N(0,1)$ から棄却域を求められるが、母分散が未知の場合は、母分散 $\sigma^{2}$ の代わりに不偏分散 $s^{2}$ を使わざるを得ないため、標本分布として自由度 $n-1$ のスチューデントのt分布 $t(n-1)$ を使うことになる（小標本かつ正規母集団の場合）。

  * 母分散が既知の場合
    * 検定統計量 $Z=\frac{\bar{X}-\mu_{0}}{\sigma/\sqrt{n}}$
    * 標本分布 $\mathrm{N}(0,1)$

  * 母分散が未知の場合（小標本かつ正規母集団）
    * 検定統計量 $t=\frac{\bar{X}-\mu_{0}}{s/\sqrt{n}}$
    * 標本分布 $t(n-1)$ 

  * 母分散が未知の場合（大標本）
    * 検定統計量 $Z=\frac{\bar{X}-\mu_{0}}{s/\sqrt{n}}$
    * 標本分布 $\mathrm{N}(0,1)$ （大標本では $s\simeq\sigma$）


標本から得られる $t$ の値は __**t 値**__ と呼ばれることがある。

=== 3. 有意水準と棄却域の決定 ===

有意水準を $\alpha$ とする。
棄却域 $R$ はそれぞれ次のようになる。

  * 母分散が既知の場合
    * 右片側検定の場合 $R=\{Z;Z>Z_{\alpha}\}$
    * 左片側検定の場合 $R=\{Z;Z< -Z_{\alpha}\}$
    * 両片側検定の場合 $R=\{Z;|Z|>Z_{\alpha/2}\}=\{Z;Z< -Z_{\alpha/2},Z_{\alpha/2}<Z\}$

{{ socialip_2021:z-test.png?nolink | σ^2 が既知の場合}}

  * 母分散が未知の場合（小標本かつ正規母集団）
    * 右片側検定の場合 $R=\{t;t>t_{\alpha}\}$
    * 左片側検定の場合 $R=\{t;t< -t_{\alpha}\}$
    * 両片側検定の場合 $R=\{t;|t|>t_{\alpha/2}\}=\{t;t< -t_{\alpha/2},t_{\alpha/2}<t\}$

  * 母分散が未知の場合（大標本）
    * 右片側検定の場合 $R=\{Z;Z>Z_{\alpha}\}$
    * 左片側検定の場合 $R=\{Z;Z< -Z_{\alpha}\}$
    * 両片側検定の場合 $R=\{Z;|Z|>Z_{\alpha/2}\}=\{Z;Z< -Z_{\alpha/2},Z_{\alpha/2}<Z\}$

$Z_{\alpha}$ は $\mathrm{N}(0,1)$の上側確率 $\alpha$ のパーセント点、
$t_{\alpha}$ は $t(n-1)$ の上側確率 $\alpha$ のパーセント点である。 

=== 4. 検定の実施 ===

標本から求めた検定統計量が棄却域 $R$ に入れば帰無仮説を棄却し、対立仮説を採択する。
棄却域に入らなければ帰無仮説は棄却されないので、対立仮説は採択されない。


検定統計量に $Z$ を用いる場合を __**Z 検定**__、$t$ を用いる場合を __**t 検定**__という。


===== 小テスト =====

[[KMS>|Moodle Server（非公式）]]で第5回の小テストを受験しなさい。


===== p 値 =====
上の手順では、有意水準 $\alpha$ からパーセント点を求めて棄却域を決定したが、
$\alpha$ がどの程度まで小さくても帰無仮説が棄却されるかを知るためには少々手間である。

そこで、標本から得られた検定統計量の値から（仮に帰無仮説が正しいとしたときの）標本分布の確率を求めることが一般に行われる。
この確率を__**p 値**__という。
標本から得られた検定統計量の値がパーセント点となる有意水準ともいえる。

帰無仮説を棄却するかどうかを、棄却域の代わりに p 値を使って判断することがよく行われる。

下図は Z 検定の場合の p 値である（図中の $Z$ は標本から得られた検定統計量の値）。
{{ socialip_2021:p-value.png?nolink |p 値}}

p 値より大きい有意水準 $\alpha$ で帰無仮説は棄却される。

例えば、p 値が $p=0.012$ なら帰無仮説は有意水準 $\alpha=0.05$（5%）では棄却されるが、$\alpha=0.01$（1%）では棄却されない。
\[
0.01 < 0.012 < 0.05
\]

p 値の定義は片側検定か両側検定かによって異なるので注意が必要である。
===== 母平均の仮説検定の例 =====

母集団は R のデータセット airquality の Wind とする。
ヒストグラムを描くと正規分布に近いので、これを正規母集団とみなす。
<code rsplus>
> hist(airquality$Wind)
</code>

{{ socialip_2021:wind_histogram.png?nolink |airqualityのWind}}

<code rsplus>
> mean(airquality$Wind)
[1] 9.957516
> sd(airquality$Wind) * sqrt(1 - 1 / length(airquality$Wind))
[1] 3.511469
</code>

このように母平均 $\mu=9.957516\mathrm{[mph]}$ と分かるが、分かっていないものとして仮説検定を行う。


母集団から $n=10$ 個の標本を抽出し、$\mu_{0}=12$ の仮説検定を行う。
次のプログラム test_mean1.R と test_mean2.R はそれぞれ Z 検定と t 検定で p値を計算するプログラムである。
ただし、標本抽出は ''set.seed(1)'' で固定しているので、毎回ランダムにしたければ ''set.seed(1)'' を削除すること。

=== Z検定のプログラム ===
<file rsplus test_mean1.R[enable_line_numbers="true"]>
# 母平均の仮説検定（Z検定）

# 母集団データ
pdata <- na.omit(airquality$Wind)

# 母標準偏差 σ
sigma <- sd(pdata) * sqrt(1 - 1 / length(pdata))

# 標本データ
n <- 10
set.seed(1)
sdata <- sample(pdata, n)

# 標本平均
xbar <- mean(sdata)

print(sprintf("母集団 : σ = %g", sigma), quote = FALSE)
print(sprintf("標本 : n = %d, xbar = %g", n, xbar), quote = FALSE)

# 帰無仮説
mu0 <- 12
print(sprintf("帰無仮説 : μ = %g", mu0), quote = FALSE)

# 検定統計量
z.value <- (xbar - mu0) / (sigma / sqrt(n))
# 右片側検定
z.right.p <- pnorm(z.value, lower.tail = FALSE)
# 左片側検定
z.left.p <- pnorm(z.value)
# 両側検定
z.twosided.p <- pnorm(abs(z.value), lower.tail = FALSE) * 2

print(sprintf("Z = %g, 右片側 p = %g, 左片側 p = %g, 両側 p = %g",
              z.value, z.right.p, z.left.p, z.twosided.p),
      quote = FALSE)
</file>

=== t検定のプログラム ===
<file rsplus test_mean2.R[enable_line_numbers="true"]>
# 母平均の仮説検定（t検定）

# 母集団データ
pdata <- na.omit(airquality$Wind)

# 標本データ
n <- 10
set.seed(1)
sdata <- sample(pdata, n)

# 標本平均
xbar <- mean(sdata)

# 不偏標準偏差
s <- sd(sdata)

print(sprintf("標本 : n = %d, xbar = %g, s = %g",
              n, xbar, s), quote = FALSE)

# 帰無仮説
mu0 <- 12
print(sprintf("帰無仮説 : μ = %g", mu0), quote = FALSE)

# 右片側検定
t.right <- t.test(sdata, alternative = "greater", mu = mu0)
# 左片側検定
t.left <- t.test(sdata, alternative = "less", mu = mu0)
# 両側検定
t.twosided <- t.test(sdata, alternative = "two.sided", mu = mu0)

print(sprintf("t = %g, 右片側 p = %g, 左片側 p = %g, 両側 p = %g",
              t.right$statistic, t.right$p.value, t.left$p.value, t.twosided$p.value),
      quote = FALSE)
</file>

test_mean1.R の実行結果
<code rsplus>
[1] 母集団 : σ = 3.51147
[1] 標本 : n = 10, xbar = 10.25
[1] 帰無仮説 : μ = 12
[1] Z = -1.57597, 右片側 p = 0.942484, 左片側 p = 0.0575159, 両側 p = 0.115032
</code>

test_mean2.R の実行結果
<code rsplus>
[1] 標本 : n = 10, xbar = 10.25, s = 3.04202
[1] 帰無仮説 : μ = 12
[1] t = -1.81918, 右片側 p = 0.948882, 左片側 p = 0.051118, 両側 p = 0.102236
</code>


----
==== Z検定の例 ====

airquality の Wind を母集団とする。
この母集団から $n=10$ 個の標本を抽出したところ、標本平均 $\bar{X}=10.25\mathrm{[mph]}$ だった。
この標本を使って**<fc #ff0000>母平均は 12[mph] より小さい</fc>**といえるかを **<fc #ff0000>有意水準 10%</fc>** で検定する。
ただし、**<fc #ff0000>母分散は既知</fc>**とする（母標準偏差 $\sigma=3.511469\mathrm{[mph]}$）。

=== 1. 仮説の設定 ===

  * 帰無仮説 $H_{0}:\mu=12\mathrm{[mph]}$
  * 対立仮説 $H_{1}:\mu<12\mathrm{[mph]}$

これは左片側検定である。

=== 2. 検定統計量と標本分布の決定 ===
母分散が既知なので、検定統計量は $Z$ を用いる。
$Z$ の標本分布は $N(0,1)$ である。

=== 3. 有意水準の決定 ===
有意水準 10% なので $\alpha=0.1$ である。

=== 4. 検定の実施 ===
標本から $Z=-1.57597$ で $p=0.0575159<0.1$ だった。
よって有意水準 10% で帰無仮説は棄却されるので、対立仮説が採択される。
つまり、有意水準 10% で母平均は 12[mph] より小さいといえる。

----
==== t検定の例 ====

airquality の Wind を母集団とする。
この母集団から $n=10$ 個の標本を抽出したところ、標本平均 $\bar{X}=10.25\mathrm{[mph]}$、標本標準偏差 $s=3.04202\mathrm{[mph]}$ だった。
この標本を使って**<fc #ff0000>母平均は 12[mph] ではない</fc>**といえるかを **<fc #ff0000>有意水準 5%</fc>** で検定する。
ただし、**<fc #ff0000>母分散は未知</fc>**とする。

=== 1. 仮説の設定 ===

  * 帰無仮説 $H_{0}:\mu=12\mathrm{[mph]}$
  * 対立仮説 $H_{1}:\mu\neq12\mathrm{[mph]}$

これは両側検定である。

=== 2. 検定統計量と標本分布の決定 ===
母分散が未知なので、検定統計量は $t$ を用いる。
$t$ の標本分布は $t(9)$ である。

=== 3. 有意水準の決定 ===
有意水準 5% なので $\alpha=0.05$ である。

=== 4. 検定の実施 ===
標本から $t=-1.81918$ で $p=0.102236>0.05$ だった。
よって有意水準 5% で帰無仮説は棄却されない。
つまり、有意水準 5% で母平均は 12[mph] ではないとはいえない。

===== t検定の関数 =====

t 検定については関数 ''t.test()'' で p 値を求めるのが便利である。

<code rsplus
> t.test(データ, alternative, mu)
</code>
引数 alternative は検定方式を文字列で指定する。
''"greater"'' は右片側検定、''"less"'' は左片側検定、''"two.sided"'' は両側検定になる。

引数 mu は帰無仮説と対立仮説に用いる $\mu_{0}$ の値である（デフォルトは 0）。

戻り値に含まれる主な要素は以下の通りである。
^  ラベル    ^  説明                                                                   ^
| statistic  | t 値。                                                                  |
| parameter  | 自由度。df（degrees of freedom）とも書く。                              |
| p.value    | p 値。                                                                  |
| conf.int   | 区間推定の信頼区間。デフォルトの信頼係数（conf.level）は 0.95 である。  |
| estimate   | 標本平均（1標本の場合）。                                               |


先ほどの [[#t検定のプログラム]] では ''t.test()'' を使っている。
例えば右片側検定の結果 t.right は以下の通り。
<code rsplus>
> t.right

	One Sample t-test

data:  sdata
t = -1.8192, df = 9, p-value = 0.9489
alternative hypothesis: true mean is greater than 12
95 percent confidence interval:
 8.486598      Inf
sample estimates:
mean of x 
    10.25 
</code>

  * 「''alternative hypothesis: true mean is greater than 12''」は「対立仮説 : $\mu>12$」という意味である。
  * 「''95 percent confidence interval''」は区間推定における 95%信頼区間の下限値と上限値である。
===== 課題  =====
[[KMS>|Moodle Server（非公式）]] で第5回の課題を行いなさい。

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

----
==== 解答例 ====

山形県の市町村データの場合（ファイル名は yamagata.csv にしている）

<file rsplus report.R[enable_line_numbers="true"]>
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年度の人口増減率を標本とする
sdata <- df[df$調査年 == 2015,][["#A05101_人口増減率【％】"]]
n <- length(sdata)

# 標本平均
xbar <- mean(sdata)

# 不偏標準偏差
s <- sd(sdata)

print(sprintf("標本 : n = %d, xbar = %g, s = %g",
              n, xbar, s), quote = FALSE)

# 帰無仮説
mu0 <- 0
print(sprintf("帰無仮説 : μ = %g", mu0), quote = FALSE)

# t検定（母分散が未知）
# 右片側検定
t.right <- t.test(sdata, alternative = "greater", mu = mu0)
# 左片側検定
t.left <- t.test(sdata, alternative = "less", mu = mu0)
# 両側検定
t.twosided <- t.test(sdata, alternative = "two.sided", mu = mu0)

print(sprintf("t = %g, 右片側 p = %g, 左片側 p = %g, 両側 p = %g",
              t.right$statistic, t.right$p.value, t.left$p.value, t.twosided$p.value),
      quote = FALSE)
</file>

実行結果
<code rsplus>
[1] 標本 : n = 35, xbar = -6.568, s = 3.55321
[1] 帰無仮説 : μ = 0
[1] t = -10.9357, 右片側 p = 1, 左片側 p = 5.62276e-13, 両側 p = 1.12455e-12
</code>

  * 「''5.62276e-13''」 は $5.62276\times 10^{-13}$ のことである。

これより解答は以下のようになる。

=== 1. 仮説の設定 ===

  * 帰無仮説 H0:μ=0
  * 対立仮説 H1:μ<0

これは左片側検定である。

=== 2. 検定統計量と標本分布の決定 ===
母分散が未知なので、検定統計量は t を用いる。
t の標本分布は t(34) である。

=== 3. 有意水準の決定 ===
有意水準 10% なので α=0.1 である。

=== 4. 検定の実施 ===
標本から t = -10.9357 で p = 5.62276e-13 < 0.1 だった。
よって有意水準 10% で帰無仮説は棄却されるので、対立仮説が採択される。
つまり、有意水準 10% で母平均は 0 より小さいといえる。