brmsパッケージでカーブフィッティング
計量言語学のMenzerath-Altmann法則をbrmsで。
(a)
s<-1:5 m<-c(3.1,2.53,2.29,2.12,2.09) dat <- data.frame(s, m) my.prior <- prior(normal(0,3), nlpar = a) + prior(normal(0,1), nlpar = b) + prior(normal(0,1), nlpar = c) fit <- brm(bf(m ~ (a*s^-b)*exp(1)^(-c*s), a + b + c ~ 1, nl = TRUE), data = dat, prior = my.prior, save_all_pars=T) summary(fit)
これ,データポイントが少ないから,事前分布の設定がすげえ問題だよね。
適当に,
- a ~ Normal(0,3)
- b ~ Normal(0,1)
- c ~ Normal(0,1)
とかにした。実際に-1とか1とかの値を取ることがほとんどないbやcはともかく,このaの設定は恣意的といわざるをえないけど。
ちなみに,より簡単な,
(a)
とベイズ因子やWAICで比べてみるとか。
my.prior <- prior(normal(0,3), nlpar = a) + prior(normal(0,1), nlpar = b) fit2 <- brm(bf(m ~ (a*s^-b), a + b ~ 1, nl = TRUE), data = dat, prior = my.prior, save_all_pars=T) summary(fit2)
bayes_factor(fit,fit2) waic(fit, fit2)
fitがfit2より高かったけど…微妙。