実習8.5 Posterior Prediction

「実習8.1,4 "Hot Hand" は存在するか？」からの続き。


実習8.5「JAGS, 離散型分析, Metropolis Algorithm , Gibbsサンプリング」のいずれかの事後確率を使って、2つのコインを投げて、結果が y₁ = 1, y₂ = 0 となる確率を求める。

ここでは、ランダムウォークの結果の事後確率から、予測値を算出する（posterior rediction）。

左は、スクリプト BernTwoJags.R の実行結果。このスクリプトの末尾に、以下を追加。

chainLength = length( theta1Sample )
yPred = matrix( NA , nrow=2 , ncol=chainLength ) 
for ( stepIdx in 1:chainLength ) {
  pHead1 = theta1Sample[stepIdx]
  yPred[1,stepIdx] = sample( x=c(0,1), prob=c(1-pHead1,pHead1), size=1 )
  pHead2 = theta2Sample[stepIdx]
  yPred[2,stepIdx] = sample( x=c(0,1), prob=c(1-pHead2,pHead2), size=1 )
}
pY1eq1andY2eq0 = sum( yPred[1,]==1 & yPred[2,]==0 ) / chainLength


以下、この追加分の要点を説明。

行列 yPred にはサンプリング値を使った「予測値」、1 行目にθ₁ 、2 行目に θ₂ の予測値（0 or 1）が格納。 よって、求める確率は「y₁ = 1, y₂ = 0」は

　sum( yPred[1,]==1 & yPred[2,]==0 ) / chainLength

今回は 0.3797324 。

yPred への予測値の格納について。

JAGSサンプリング後、theta1Sample, theta2Sample には、左のようにθ₁, θ₂ のサンプリング値が格納。

> par(mfrow=c(2,1))
> hist(theta1Sample, breaks = seq(0,1,by=0.1),xlab="θ",ylab="")
> hist(theta2Sample, breaks = seq(0,1,by=0.1),xlab="θ",ylab="")

このサンプリング値（ 1 である確率）を使って、θ₁ と θ₂ の「予測値」を生成して yPred に格納。

yPred[1,stepIdx] = sample( x=c(0,1), prob=c(1-pHead1,pHead1), size=1 )
yPred[2,stepIdx] = sample( x=c(0,1), prob=c(1-pHead2,pHead2), size=1 )