第12講 - 複数のグループ間の違いを検証する
(Press ? for help, n and p for next and previous slide)
村田 昇
分散分析における仮説の例
- ある小売店について“売上高は月によって差があるか”
- ある銘柄の株価について“収益率は曜日によって差があるか”
2種類のデータの分散が等しいか否かを検定する
\begin{equation} H_0:\sigma_1^2=\sigma_2^2\quad\text{vs}\quad H_1:\sigma_1^2\neq\sigma_2^2 \end{equation}
\(F\) 検定
\begin{equation} \text{検定統計量 : } F=\frac{s_1^2}{s_2^2} \left( = \frac {((m{-}1)s_1^2/\sigma_1^2)/(m{-}1)} {((n{-}1)s_2^2/\sigma_2^2)/(n{-}1)} \right) \end{equation}
\(H_0\) の下では以下が成立
\begin{equation} P(F < F_{\alpha/2}(m{-}1,n{-}1) \text{ または }F > F_{1{-}\alpha/2}(m{-}1,n{-}1))=\alpha \end{equation}
第一種過誤の上限が \(\alpha\) の棄却域
\begin{equation} R_{\alpha}= \left(-\infty,F_{\alpha/2}(m{-}1,n{-}1)\right) \cup\left(F_{1{-}\alpha/2}(m{-}1,n{-}1),\infty\right) \end{equation}
以下の場合，帰無仮説を棄却
\begin{equation} F < F_{\alpha/2}(m{-}1,n{-}1) \quad\text{または}\quad F>F_{1{-}\alpha/2}(m{-}1,n{-}1) \end{equation}
各グループ \(i\) ごとに \(n_i\) 個のデータを観測
\begin{equation} Y_{i1},Y_{i2},\dotsc,Y_{in_i},\quad (i=1,2,\dotsc,p) \end{equation}
小売店の売上高の問題
- \(A_1,A_2,\dotsc,A_p\) : 月
- \(Y_{i1},Y_{i2},\dotsc,Y_{in_i}\) : \(i\) 月の各日における売上高
観測データのモデル
\begin{equation} Y_{ij}=\mu_i+\varepsilon_{ij}\quad(i=1,\dotsc,p;\;j=1,\dotsc,n_i). \end{equation}
検定問題
各水準 \(A_1,A_2,\dotsc,A_p\) の平均 \(\mu_1,\mu_2,\dotsc,\mu_{p}\) に差があるか否かを検定する
\begin{equation} H_0:\mu_1=\cdots=\mu_p\quad\text{vs}\quad H_1:\text{ある\(i,j\)に対して}\mu_i\neq\mu_j. \end{equation}
各種平均 (\(n=\sum_{i=1}^pn_i\) は全サンプル数)
\begin{align} &\text{(全データの標本平均)} &\bar{Y}_{\cdot\cdot}&=\frac{1}{n}\sum_{i=1}^p\sum_{j=1}^{n_i}Y_{ij},\\ &\text{(水準\(A_i\)の標本平均)} &\bar{Y}_{i\cdot}&=\frac{1}{n_i}\sum_{j=1}^{n_i}Y_{ij}\quad (i=1,\dotsc,p). \end{align}
各種変動
\begin{align} &\text{(全変動)} &SS_T&=\sum_{i=1}^p\sum_{j=1}^{n_i}(Y_{ij}-\bar{Y}_{\cdot\cdot})^2,\\ &\text{(級内変動)} &SS_W&=\sum_{i=1}^p\sum_{j=1}^{n_i}(Y_{ij}-\bar{Y}_{i\cdot})^2,\\ &\text{(級間変動)} &SS_B&=\sum_{i=1}^p\sum_{j=1}^{n_i}(\bar{Y}_{i\cdot}-\bar{Y}_{\cdot\cdot})^2 =\sum_{i=1}^pn_i(\bar{Y}_{i\cdot}-\bar{Y}_{\cdot\cdot})^2. \end{align}
変動の分解
\begin{equation} \underbrace{\text{(全変動) }SS_T}_{\text{自由度 }n{-}1} =\underbrace{\text{(級内変動) }SS_W}_{\text{自由度 }n{-}p} +\underbrace{\text{(級間変動) }SS_B}_{\text{自由度 }p{-}1} \end{equation}
変動の分解
\begin{equation} \underbrace{\text{(全変動) }SS_T}_{\text{自由度 }n{-}1} =\underbrace{\text{(級内変動) }SS_W}_{\text{自由度 }n{-}p} +\underbrace{\text{(級間変動) }SS_B}_{\text{自由度 }p{-}1} \end{equation}
検定統計量
\begin{equation} F=\frac{SS_B/(p{-}1)}{SS_W/(n{-}p)} \end{equation}
\(H_0\) の下で以下が成立
\begin{equation} P(F>F_{1{-}\alpha}(p{-}1,n{-}p))=\alpha \end{equation}
第一種過誤の上限が \(\alpha\) となる棄却域
\begin{equation} R_{\alpha}= \left(F_{1{-}\alpha}(p{-}1,n{-}p),\infty\right) \end{equation}
以下の場合，帰無仮説を棄却
\begin{equation} F>F_{1{-}\alpha}(p{-}1,n{-}p) \end{equation}
\(p\) 値を計算 (右片側検定の場合の計算方法)
\begin{equation} \text{(\(p\) 値)}=\int_{F}^\infty f(x)dx \end{equation}
| 自由度 | 平方和 | 平均平方和 | \(F\) 値 | \(p\) 値 | ||
|---|---|---|---|---|---|---|
| 級間 | \(p{-}1\) | \(SS_B\) | \(SS_B/(p{-}1)\) | \(F\) | \(\mathstrut\int_{F}^\infty f(x)dx\) | |
| 級内 | \(n{-}p\) | \(SS_W\) | \(SS_W/(n{-}p)\) | |||
| 全変動 | \(n{-}1\) | \(SS_T\) |
平均値の別表現
\begin{equation} \mu_i=\mu+\alpha_i,\quad \mu=\frac{1}{n}\sum_{i=1}^pn_i\mu_i,\quad \sum_{i=1}^pn_i\alpha_i=0 \end{equation}
帰無仮説 \(H_0\) は以下と同等
\begin{equation} H_{0}:\alpha_1=\cdots=\alpha_p=0 \end{equation}
基本書式
aov(formula, data)
#' formula: 式，一元配置の場合は (観測値 ~ 因子)
#' data: データフレーム
#' 分析の結果を参照する関数がいくつかあるので，多くの場合適当なオブジェクトに代入
分散分析表の表示
aov(formula, data) |> summary() # 分析表形式での表示
aov(formula, data) |> anova() # データフレーム(data.frame形式)
aov(formula, data) |> broom::tidy() # データフレーム(tibble形式)
aov(formula, data) |> model.table(type = "means") # 平均値
aov(formula, data) |> model.table(type = "effects") # 効果(既定値)
基本書式
oneway.test(formula, data, subset, na.action, var.equal = FALSE)
#' formula: 式
#' data: データフレーム
#' subset: subsetの指定
#' na.action: 欠損値の扱い
#' var.equal: 等分散を仮定する場合は TRUE，標準では Welch の近似が用いられる
tokyo_weather.csv) の
気温の項目について以下の問に答えよ．
曜日ごとの気温の平均と分散を求めよ．
#' 週日を因子(factor)とするには例えば以下のような項目を加えればよい
tw_data <- read_csv("data/tokyo_weather.csv") |>
mutate(day_of_week = ordered(day_of_week, # 曜日を順序付き因子にする
levels = wday(1:7, label = TRUE)))
#' 関数 wday() の出力は使用言語に依存するので，強制的に英語にするには locale を指定して
#' wday(..., locale = "en_US") とすれば良い
#' 省略形を用いない場合は abbr を指定する．以下はドイツ語の例
#' wday(1:7, label = TRUE, abbr = FALSE, locale = "de_DE")
#' 利用可能なlocale(言語)を調べるには stringi::stri_locale_list() を用いる
複数の薬の治験の問題
複数の薬の効能を複数の被験者に投与する実験
- “薬の種類”と“被験者番号”の2種類の因子
- “薬の種類”という因子での薬の効能の差を検証したい
- “被験者番号”という因子で効能に差があることは許容したい (薬の効き目には個人差があると考えられるため)
観測データのモデル
\begin{equation} Y_{ij}=\alpha_i+\beta_j+\varepsilon_{ij}\quad (i=1,\dotsc,a;\;j=1,\dotsc,b). \end{equation}
検定問題
因子 \(A\) の各水準の効果に差があるか否かを検定する
(因子 \(B\) の効果は除いて検定したい)
対立仮説 : 因子 \(A\) の効果に差のあるものがある
\begin{equation} H_0:\alpha_1=\cdots=\alpha_a\quad\text{vs}\quad H_1:\text{ある\(i_1,i_2\)に対して}\alpha_{i_1}\neq\alpha_{i_2}. \end{equation}
各種平均 (\(n=ab\))
\begin{align} &\text{(全データの標本平均)} &\bar{Y}_{\cdot\cdot}&=\frac{1}{n}\sum_{i=1}^a\sum_{j=1}^{b}Y_{ij},\\ &\text{(水準\(A_i\)の標本平均)} &\bar{Y}_{i\cdot}&=\frac{1}{b}\sum_{j=1}^bY_{ij}\quad (i=1,\dotsc,a),\\ &\text{(水準\(B_j\)の標本平均)} &\bar{Y}_{\cdot j}&=\frac{1}{a}\sum_{i=1}^aY_{ij}\quad (j=1,\dotsc,b). \end{align}
因子効果の推定量としての標本平均
\begin{align} \bar{Y}_{i\cdot} &\to\alpha_i+\frac{1}{b}\sum_{j=1}^b\beta_j,\\ \bar{Y}_{\cdot j} &\to\frac{1}{a}\sum_{i=1}^a\alpha_i+\beta_j,\\ \bar{Y}_{\cdot\cdot} &\to\frac{1}{ab}\sum_{i=1}^a\sum_{j=1}^b(\alpha_i+\beta_j)=\frac{1}{a}\sum_{i=1}^a\alpha_i+\frac{1}{b}\sum_{j=1}^b\beta_j \end{align}\begin{equation} Y_{ij}-\bar{Y}_{i\cdot} -\bar{Y}_{\cdot j}+\bar{Y}_{\cdot\cdot} \to\varepsilon_{ij} \end{equation}
各種変動
\begin{align} &\text{(行間変動)} &SS_A&=b\sum_{i=1}^a(\bar{Y}_{i\cdot}-\bar{Y}_{\cdot\cdot})^2,\\ &\text{(列間変動)} &SS_B&=a\sum_{j=1}^b(\bar{Y}_{\cdot j}-\bar{Y}_{\cdot\cdot})^2,\\ &\text{(誤差変動)} &SS_E&=\sum_{i=1}^a\sum_{j=1}^b(Y_{ij}-\bar{Y}_{i\cdot} -\bar{Y}_{\cdot j}+\bar{Y}_{\cdot\cdot})^2. \end{align}
変動の分解
\begin{align} \underbrace{\text{(全変動) }SS_T}_{\text{自由度 }ab-1} =&\underbrace{\text{(行間変動) }SS_A}_{\text{自由度 }a-1}\\ &+\underbrace{\text{(列間変動) }SS_B}_{\text{自由度 }b-1} +\underbrace{\text{(誤差変動) }SS_E}_{\text{自由度 }(a-1)(b-1)} \end{align}
変動の分解
\begin{align} \underbrace{\text{(全変動) }SS_T}_{\text{自由度 }ab-1} =&\underbrace{\text{(行間変動) }SS_A}_{\text{自由度 }a-1}\\ &+\underbrace{\text{(列間変動) }SS_B}_{\text{自由度 }b-1} +\underbrace{\text{(誤差変動) }SS_E}_{\text{自由度 }(a-1)(b-1)} \end{align}
検定統計量
\begin{equation} F_A=\frac{SS_A/(a{-}1)}{SS_E/(a{-}1)(b{-}1)} \end{equation}
\(H_0\) の下で以下が成立
\begin{equation} P(F_A>F_{1{-}\alpha}(a{-}1,(a{-}1)(b{-}1)))=\alpha \end{equation}
第一種過誤の上限が \(\alpha\) となる棄却域
\begin{equation} R_{\alpha}= (F_{1{-}\alpha}(a{-}1,(a{-}1)(b{-}1)),\infty) \end{equation}
以下の場合，帰無仮説を棄却
\begin{equation} F_A>F_{1{-}\alpha}(a{-}1,(a{-}1)(b{-}1)) \end{equation}
\(p\) 値を計算 (右片側検定の場合の計算方法)
\begin{equation} \text{(\(p\) 値)}=\int_{F_A}^\infty f(x)dx \end{equation}
| 自由度 | 平方和 | 平均平方和 | \(F\text{値}\) | \(p\text{値}\) | |
|---|---|---|---|---|---|
| \(\text{因子}A\) | \(a{-}1\) | \(SS_A\) | \(\displaystyle\frac{SS_A}{a{-}1}\) | \(F_A\) | \(\displaystyle\int_{F_A}^\infty f(x)dx\) |
| \(\text{因子}B\) | \(b{-}1\) | \(SS_B\) | \(\displaystyle\frac{SS_B}{b{-}1}\) | \(F_B\) | \(\displaystyle\int_{F_B}^\infty f(x)dx\) |
| \(\text{誤差}\) | \((a{-}1)(b{-}1)\) | \(SS_E\) | \(\displaystyle\frac{SS_E}{(a{-}1)(b{-}1)}\) |
相対効果による定式化
モデルの書き換え
\begin{align} \mu^*&=\bar{\alpha}+\bar{\beta},& \alpha^*_i&=\alpha_i-\bar{\alpha},& \beta^*_j&=\beta_j-\bar{\beta}\\ &&\bar{\alpha}&=\frac{1}{a}\sum_{i=1}^a\alpha_i, &\bar{\beta}&=\frac{1}{b}\sum_{j=1}^b\beta_j \end{align}\begin{equation} Y_{ij}=\mu^*+\alpha^*_i+\beta^*_j+\varepsilon_{ij}\quad (i=1,\dotsc,a;\;j=1,\dotsc,b) \end{equation}\begin{equation} \text{ただし }\sum_{i=1}^a\alpha^*_i=\sum_{j=1}^b\beta^*_j=0 \end{equation}
帰無仮説 \(H_0\) は以下と同等
\begin{equation} H_{0}:\alpha^*_1=\cdots=\alpha^*_a=0 \end{equation}
基本書式
aov(formula, data)
#' formula: 式，二元配置の場合は (観測値 ~ 因子 + 因子)
#' data: データフレーム
package::datarium に含まれている
jobsatisfaction データについて以下の問に答えよ．
性別 (
gender) と 学歴 (education_level) の違いによる 仕事の満足度 (score) を収集したデータ
#' パッケージのインストールと読み込みは以下のように行うことができる
install.packages("datarium") # インストール
library("datarium") # 読み込み