R のベイジアン混合効果モデルから、パラメーターの異なるレベル間のコントラットを計算し、ベイズ ファクターを生成したいと考えています。私の結果 (Jud) はバイナリ (1 = はい/同期、0 = いいえ/非同期) であり、パラメーター SOAsF は 6 レベル (0、100、200、300、400、500) の因子です。
さまざまなチュートリアル/機能[#1]、[#2]、および[#3]に従って、3 つの異なる方法を使用したコードを次に示します。
library(emmeans)
library(brms)
library(modelbased)
brm_acc_1<-brm(Jud ~ SOAsF +(1|pxID),data =dat_long, family=bernoulli("logit"), prior = set_prior('normal(0,10)'), iter = 2000, chains=4, save_all_pars = TRUE)
summary(brm_acc_1)
brms::conditional_effects(brm_acc_1)
####1
groups <- emmeans(brm_acc_1, ~ SOAsF)
group_diff <- pairs(groups)
(groups_all <- rbind(groups, group_diff))
bayesfactor_parameters(groups_all, prior = brm_acc_1, direction = "two-sided", effects = c("fixed", "random", "all"))
####2
ppc <- pp_check(brm_acc_1, type = "stat_grouped", group = "SOAsF")
#contrast 200 - 300
contrast_300_200 <- ppc$data$value[ppc$data$group == "200"] - ppc$data$value[ppc$data$group == "300"]
quantile(contrast_300_200*100, probs = c(.5, .025, .975))
####3
h_1 <- hypothesis(brm_acc_1, "SOAsF200 < SOAsF300")
print(h1, digits = 4)
h2 <- hypothesis(brm_acc_1, "SOAsF200 > SOAsF300")
print(h2, digits = 4)
結果：
####1
# Bayes Factor (Savage-Dickey density ratio)
Parameter | BF
-----------------------
0, . | 8.08e-03
100, . | 0.61
200, . | 7.29e+03
300, . | 67.77
400, . | 21.81
500, . | 0.28
., 0 - 100 | 2.75
., 0 - 200 | 1.90e+05
., 0 - 300 | 410.42
., 0 - 400 | 570.11
., 0 - 500 | 1.03
., 100 - 200 | 0.5
., 100 - 300 | 0.05
., 100 - 400 | 0.02
., 100 - 500 | 7.13e-03
., 200 - 300 | 0.01
., 200 - 400 | 0.01
., 200 - 500 | 1.11
., 300 - 400 | 7.21e-03
., 300 - 500 | 0.1
., 400 - 500 | 0.04
* Evidence Against The Null: [0]
####2
50% 2.5% 97.5%
1.988631 -3.707585 7.680694
####3
Hypothesis Tests for class b:
Hypothesis Estimate Est.Error CI.Lower CI.Upper
1 (SOAsF200)-(SOAsF... < 0 0.0881 0.0903 -0.0612 0.2372
Evid.Ratio Post.Prob Star
1 0.1919 0.161
---
'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
'*': For one-sided hypotheses, the posterior probability exceeds 95%;
for two-sided hypotheses, the value tested against lies outside the 95%-CI.
Posterior probabilities of point hypotheses assume equal prior probabilities.
Hypothesis Tests for class b:
Hypothesis Estimate Est.Error CI.Lower CI.Upper
1 (SOAsF200)-(SOAsF... > 0 0.0881 0.0903 -0.0612 0.2372
Evid.Ratio Post.Prob Star
1 5.2112 0.839
---
'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
'*': For one-sided hypotheses, the posterior probability exceeds 95%;
for two-sided hypotheses, the value tested against lies outside the 95%-CI.
Posterior probabilities of point hypotheses assume equal prior probabilities.
200 対 300 の対比を例にとると、SOA 200 で提示されたサウンドは、SOA 300 で提示されたサウンドと比較して、同期している (はい) と同等に判断されますか?
方法 1 は帰無仮説 SOA 200 - SOA 300 = 0 with BF = 0.01; に対する証拠を提供しているようです。それで、それらは同期していると同等に判断されているように見えますか?
方法 2 は、1.988631% 95%CI [-3.707585, 7.680694] の証拠で、帰無仮説 SOA 200 = SOA 300 に対してほとんど証拠を提供しないようです。
方法 3 は、対立仮説 SOA 200 > SOA 300 または SOA 200 - SOA 300 < 0、BF = 5.2112 に対する証拠を提供しているようです。
#1 は両面であり、#3 は片面であるため、違いを見つけることができますか?
ただし、＃1の片側を実行することはできませんでした（方向=「左」または「右」）
bayesfactor_parameters(groups_all, prior = brm_acc_1, direction = ">", effects = c("fixed", "random", "all") )
Computation of Bayes factors: sampling priors, please wait...
Error in `$<-.data.frame`(`*tmp*`, "ind", value = 8L) :
replacement has 1 row, data has 0
または #3 両側 (仮説(brm_acc_1, " SOAsF200 - SOAsF300 = 0 ") )
Hypothesis Tests for class b:
Hypothesis Estimate Est.Error CI.Lower CI.Upper
1 (SOAsF200-SOAsF300) = 0 0.0881 0.0903 -0.0914 0.2682
Evid.Ratio Post.Prob Star
1 NA NA
---
'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
'*': For one-sided hypotheses, the posterior probability exceeds 95%;
for two-sided hypotheses, the value tested against lies outside the 95%-CI.
Posterior probabilities of point hypotheses assume equal prior probabilities.
私は立ち往生しています、どんな助けでも大歓迎です.Thanks.