μ's とAqours の人気の差

こんな記事があった。あるアニメショップでキャラの人気投票をしたら、ラブライブにおいてμ's のメンバーのほうが、Aqours のメンバーより総じて上位だったらしい。
というわけで、2グループの人気はどれくらいの差かを考える。
　
2グループ各9人、全部で18人のキャラの得票数がある。あるベース に各キャラの効果、グループ効果 があり、18 人の所属は であるとする。18人のハイパーパラメータは

投票確率 はディリクレ分布

得票数は多孔分布

でサンプリングされるとする。
　
結果としては 程度が多く、収束しなかった。また、n_eff が全然なかった。
また、 が何十万とかなって単純にμ's だと何倍人気になる、というのがわかりにくかったので、 の事後分布を各グループについて中央値を取って何倍人気に差があるか、にしている。すると2.5倍くらいμ's とAqours に人気の差があるようだった。

a <- apply(ex$alpha, 1, tapply, dat$group, median)
quantile(a[2,]/a[1,], c(0.025, 0.5, 0.975))

# 人気の差
    2.5%      50%    97.5% 
2.111500 2.469390 3.011622 

code <- "
data{
  #int<lower=1> N_vote;
  int<lower=1> Member;
  int<lower=1> Vote[Member];
  int<lower=0, upper=1> group[Member];
}
parameters{
  real<lower=-2, upper=2> base;
  real b[Member];
  real g;
  simplex[Member] p;
}
transformed parameters{
  vector<lower=0>[Member] alpha;
  for(i in 1:Member){
    #alpha[i] = inv_logit(g*group[i] + b[i] + base);
    alpha[i] = g*group[i] + b[i] + base;
  }
}
model{
  p ~ dirichlet(alpha);
  Vote ~ multinomial(p);
}
"

library(rstan)
rstan_options(auto_write = TRUE)
options(mc.cores = parallel::detectCores())
m <- stan_model(model_code=code)

standata <- list(Member=length(dat$name), Vote=dat$vote, group=c(dat$group)-1)
fit <- sampling(m, standata, warmup=5000, iter=20000)
ex <- extract(fit)

par(mar=c(4.5, 5.5, 2, 2), cex.lab=1.5)
plot(0, type="n", xlim=c(0, max(ex$p)), ylim=c(1, nrow(dat)), xlab="投票率", ylab="", yaxt="n")
abline(h=seq(nrow(dat)), lty=3)
for(i in 1:nrow(dat)) vioplot(ex$p[,i], at=i, add=TRUE, col=c("skyblue", "orange")[c(dat$group)[i]], horizontal=TRUE)
text(par()$usr[1]-0.027, seq(nrow(dat)), dat$name, xpd=TRUE, pos=4)

データ

name vote group
矢澤にこ 384 muse
西木野真姫 384 muse
南ことり 370 muse
東條希 336 muse
綾瀬絵里 283 muse
園田海未 263 muse
小泉花陽 251 muse
星空凛 226 muse
津島善子 190 aqours
黒澤ルビィ 171 aqours
国木田花丸 155 aqours
渡辺曜 141 aqours
高坂穂乃果 128 muse
松浦果南 110 aqours
黒澤ダイヤ 103 aqours
桜内梨子 93 aqours
小原鞠莉 80 aqours
高海千歌 74 aqours