高校数学の質問スレ（医者・東大卒専用） Part438 (991ﾚｽ)
上下前次1-新
抽出解除 必死ﾁｪｯｶｰ(本家) (べ) 自ID ﾚｽ栞 あぼーん

---

838: 05/17(土)02:21 ID:zAzyVzie(1) AAS
>>837
# --- 必要パッケージ ---
library(rjags)
library(coda)
library(HDInterval)

# --- データ定義 ---
data_list <- list(
r0 = 19, n0 = 202, # プラセボ
r1 = 5, n1 = 201, # 旧薬（試験1）
r2 = 9, n2 = 203, # 新薬
r3 = 5, n3 = 204 # 旧薬（試験2）
)

# --- 階層モデル定義 ---
model_hier <- "
model {
r0 ~ dbin(p0, n0)
r1 ~ dbin(p1, n1)
r2 ~ dbin(p2, n2)
r3 ~ dbin(p3, n3)

p0 ~ dbeta(1, 1)
p2 ~ dbeta(1, 1)

mu_old ~ dbeta(1, 1)
tau ~ dgamma(0.001, 0.001) # 弱情報事前分布
p1 ~ dbeta(mu_old * tau, (1 - mu_old) * tau)
p3 ~ dbeta(mu_old * tau, (1 - mu_old) * tau)

p2_worse_than_p0 <- step(p2 - p0)

rd_p0_p1 <- p0 - p1
rd_p1_p2 <- p1 - p2
rd_p0_p3 <- p0 - p3

}
"

jags_model <- jags.model(textConnection(model_hier),
data = data_list, n.chains = 2, quiet=TRUE)
update(jags_model, 3000, progress.bar="none")
jags_samples <- coda.samples(jags_model,
c("p0","p1","p2","p3",
"p2_worse_than_p0", "rd_p0_p1","rd_p1_p2", "rd_p0_p3"),
n.iter=10000, progress.bar="none")
gelman.plot(jags_samples)
plot(jags_samples)

js <- as.data.frame(as.matrix(jags_samples))
mean(js$p2_worse_than_p0)
hdi(js$rd_p0_p1) # 旧薬(試験1) vs プラセボ
hdi(js$rd_p1_p2) # 旧薬(試験2) vs 新薬
hdi(js$rd_p0_p3) # 旧薬(試験2) vs プラセボ　仮想
hdi(js$p0-js$p2) # 新薬　　　　 vs プラセボ　仮想

source("plotpost.R")
layout(matrix(c(1,2,3,3), 2, 2, byrow=TRUE))
plotpost(js$p2, col='lightcoral',xlab="新薬",cex.lab=1.5,main="")
plotpost(js$p0, col='lightgreen', xlab="プラセボ",cex.lab=1.5,main="")
plotpost(js$p0 - js$p2, compVal = 0, col=c('lightcoral', 'lightgreen'),
xlab="プラセボ - 新薬", cex.main=2,main="二項分布階層モデル",cex.lab=1.5, breaks="scott")

---

上下前次1-新書関写板覧索設栞歴

---

ｽﾚ情報 赤ﾚｽ抽出 画像ﾚｽ抽出 歴の未読ｽﾚ

---

ぬこの手 ぬこTOP 1.718s*