2017年12月6日 更新

ベイズ統計学を用いた平均値の差の推測

753 view お気に入り 0
本記事では,最近注目を集めているベイズ統計について,Rでの実装を行います.具体的な内容としては,平均値の差の推定を行います.今回は,Rstan をというライブラリを用いるのですが,環境構築に関してはRstanの公式Githubが日本語訳されているので,こちらを参考して頂けると良いかと思います.https://github.com/stan-dev/rstan/wiki/RStan-Getting-Started-(Japanese)

そして今回は,以下の例題を考えます.
[例題] とある学校の X組と Y組におけるテストの平均点を比較する. X組に属する20人の点数は [49 , 66 , 69 , 55 , 54 , 72 , 51 , 76 , 40 , 62 , 66 , 51 , 59 , 68 , 66 , 57 , 53 , 66 , 58 , 57] 平均59.8点で, Y組に属する20人の点数は [41 , 55 , 21 , 49 , 53 , 50 , 52 , 67 , 54 , 69 , 57 , 48 , 31 , 52 , 56 , 50 , 46 , 38 , 62 , 59] 平均50.5点であった.このとき X組と Y組の平均点および点数差の事後分布を求めよ.
まず前提条件として,標本X および Y は,平均がそれぞれ μx と μy (母平均),かつ分散 σ² (母分散) が等しい正規分布に従う母集団から無作為抽出されたと考えます.それぞれを式で表すと,
 (4429)

となります.また,母数の集まりは,
 (4430)

であります.次に尤度は,「母平均および母分散が与えられたとき (を知っていたとすると),今回の標本が得られる確率」であり,次の式で表されます.
 (4431)

となります.最後に,事前分布は,
 (4432)

で表されます.以上より,事後分布を,
 (4433)

の式で表すことができます.この式を利用して, MCMC法を用いることで,正規化定数の解析的計算は必要なくなり,母数 θ の事後分布に従う乱数の生成が可能となります.MCMC法によって得られた事後分布を以下に示します.この事後分布から様々な統計量が得られます.まず,X組およびY組の平均値の事後分布 (a) をみると,期待値 (X組:59.72,Y組:50.47) と分散 (X組:2.33,Y組:2.32) となりました.今回の場合,90,000回サンプリング操作を繰り返しています.その中で,x < y を満たす回数を数え,サンプリングの回数 (今回の場合は 90,000) で割れば,「X組の真の平均値 < Y組の真の平均値」となる確率を求めることができ,0.997 という値が算出されました.続いて,真の平均値の差 が 95 % の確率で含まれる区間は,平均値の差の事後確率分布 (b) で 95 % のデータが含まれる範囲を計算すると算出でき,これを 95 % ベイズ確信区間と呼び,本例題では,2.77 > μx - μy > 15.69 が得られました.
 (4426)

# データを受け取る
data { 
  int<lower=0>  N_x;    
  int<lower=0>  N_y;    
  real<lower=0> x[N_x]; 
  real<lower=0> y[N_y]; 
}

# 平均および分散を定義
parameters {
  real mu_x;
  real mu_y;
  real<lower=0> sigma;
}
 
#事前分布を定義
model {
  mu_x ~ normal(0, 100);
  mu_y ~ normal(0, 100);
  sigma ~ cauchy(0, 5);
  
  x ~ normal(mu_x, sigma);
  y ~ normal(mu_y, sigma);
}
 
# 統計量を定義.今回は平均値の差
generated quantities {
  real diff;
  diff = mu_x - mu_y;
}
bayesian_inference.stan
16 件

関連する記事 こんな記事も人気です♪

この記事のキーワード

この記事のキュレーター

井上 大輝 井上 大輝