賭ケグルイを見ている。
黒髪ロングのプリーツスカート黒タイツのはやみんなので視聴意欲がやばい。
ここで、賭ケグルイの1話で、投票じゃんけんという変則じゃんけんで勝負する。
クラス 人がグー
、チョキ
、パー
のいずれかの手を投票する。
蛇喰さんと芽亜里はその投票箱から重複なしで3枚カードを引く。
各々の持ち手はその3枚となり、あとはじゃんけんのルールに従って勝敗が決まる。出した手は捨てられる。
ということである。ここで、両者120万円(120枚のチップ)から始めて、出した手の順に(他はアニメ演出上、蛇喰さんが把握できた手について)
1回戦：2枚賭け 蛇喰勝(G??) 芽亜里(C??) アニメ絵では残りの2枚はG であることが分かっているが蛇喰さんには見えない
2回戦：50枚賭け 蛇喰負(GCC) 芽亜里(GG?)
3回戦：2枚賭け 蛇喰勝(P??) 芽亜里(G??) アニメ絵では残りの2枚はC と(たぶん) G であることが分かっているが蛇喰さんには見えない。
4回戦：50枚賭け 蛇喰負(PG?) 芽亜里(PP?)
その後、じゃんけんぽん、が3回は行われて、蛇喰さんがすべてのチップを奪われるが、即金で
N回戦：1000万賭け 蛇喰勝(C??) 芽亜里(PPG)
というふうになって最終的に蛇喰さんが勝利する。
ここでこの投票じゃんけんイカサマのネタバレだが、芽亜里はクラスメイト30人中21人の投票を操り、主人公の鈴井くんの提示した手(N回戦ではチョキ) 以外の手(つまりグー
とパー
)のどちらかに投票する。こうすると、クラスの
の割合が投票を操作されているとき、
の投票確率は
となる。
ですべて
、
のときグーとパーのみ出るようになって各々
である。
このとき、 が出る確率はかなり低いから、
を出しておけば少なくとも負けることはないだろう、というのが芽亜里の作戦である、たぶん。
ここで蛇喰さんは、2枚賭けのときは芽亜里が適当に手を出しているのに50枚賭けという大勝負のときには策を講じた手を出していること、自分の真後ろにたっていた主人公の手を手鏡で観察していたことから、N回戦のときに芽亜里にイカサマの存在を伝え、なおかつそれが、「10人から20人が芽亜里の指示通りに投票する」ということを見抜いた。
実際にこれを見抜けるかやってみる。
いま、 がある確率で出るとき(等確率
とする)、自分の手札の揃え方は10通りある。この10通りが出る確率は多項分布であり、以下の通りである。
x <- expand.grid(0:3, 0:3, 0:3) x <- x[rowSums(x) == 3,] colnames(x) <- c("G", "C", "P") prob <- c(1, 1, 1)/3 dm <- apply(x, 1, dmultinom, prob=prob) cbind(x, dm)
G C P dm 4 3 0 0 0.03703704 7 2 1 0 0.11111111 10 1 2 0 0.11111111 13 0 3 0 0.03703704 19 2 0 1 0.11111111 22 1 1 1 0.22222222 25 0 2 1 0.11111111 34 1 0 2 0.11111111 37 0 1 2 0.11111111 49 0 0 3 0.03703704
ここで、蛇喰さんがイカサマを見抜けたのは、「の手札が揃う勝負が極端に少ない」からという仮説にしておく。というのも、投票が偏らなければ投票箱は各々
になるはずだから、とする。たぶん。
ここで、上の各手の出方のうち、以外に絞れば、余事象的に「手札がすべて同じ(3枚)」と「2枚と1枚(と0枚)」は計算できる。
そしてこれはか、それ以外にすれば、二項分布で計算できる。
二項分布から95%信頼区間を計算してもよいが、蛇喰さんは「10人から、万全を期すなら20人」といっているので、適当にrstan を使って事後分布で評価してみる。
例えば、1回の勝負で が来たときの、投票箱の手の分布はこのような感じになる。
ここで、 の出る確率は足して1 なので、
の3次元平面になるが、平面なので2次元平面に写像した図にしている。三角形頂点のGCP の点に近いほど、その手の出る確率が1 に近くなる。
1回の勝負で2枚と1枚、例えば のときはこうなる。
1回の勝負で2枚と1枚、例えば のときはこうなる。
library(rstan) code <- " data{ int trial; int gcp[3, trial]; } parameters{ simplex[3] a; } model{ for(i in 1:trial) gcp[, i] ~ multinomial(a); } " m <- stan_model(model_code=code) trial <- 3 standata <- list(trial=trial, gcp=replicate(trial, c(1, 1, 1))) fit <- sampling(m, standata, iter=30000) M <- cbind(c(0,1,0.5), c(0,0,sqrt(3)/2)) # 2次元に写像する行列 triangle <- diag(1, 3) %*% M xy <- ex$a %*% M kd <- kde2d(xy[,1], xy[,2], c(bandwidth.nrd(xy[,1]), bandwidth.nrd(xy[,2])), n=500) cols <- matlab::jet.colors(100) i <- which(kd$z == max(kd$z), arr.ind=TRUE) out <- c(kd$x[ i[1] ], kd$y[ i[2] ], 1) %*% solve(cbind(M, 1)) par(mar=c(3, 3, 2, 2), cex.axis=1.5) image(kd$z, col=cols) polygon(rbind(triangle[c(1:3,1),],c(-5,0),c(-5,5),c(5,5),c(5,-5),c(-5,-5),c(-5,0),triangle[1,]), col="white", border=NA) polygon(triangle) points(triangle, pch=16, cex=3, xpd=TRUE, col="hotpink") text(triangle, c("G", "C", "P"), col="yellow", font=2, xpd=TRUE) points(out %*% M, pch="★", col="hotpink", cex=1.5, xpd=TRUE) legend("topright", legend=mapply(sprintf, "%s: %.4f", c("G", "C", "P"), out), bty="n", cex=2)
ここで大事なのが、確率的な問題なので、たとえ が出る確率が低くても(青の部分)、引く手が
になることはありえる、ということである。ましてや、等確率
はそこそこありうる確率点である。
しかし、この分布も、勝負の回数が増えることで、等確率 はおかしい、という状況が出てくる。
さてここで、買収されている(奴隷)の割合を とすると、
が出る確率は
、それ以外(2枚同じか、3枚とも同じ) 手が出る確率は
となる。勝負が進んで蛇喰さんが手の内訳を確認できた、すなわち、蛇喰さん自身の手(これは勝負の回数分)と、芽亜里の手(初手で勝負が決まれば不明だが、2回目以降じゃんけんが進めばわかる)を足した回数、ということにする。つまり、アニメの4回戦が終わったときには、手を6回確認できていることとする。
芽亜里の奴隷の割合 が一様な事前分布を仮定したときの、N回手を確認したあとの事後分布は以下のようになる。
回が増えると分布は極端に、つまりみんな奴隷に傾く。
% 信用区間をプロットした。12回くらい手のうちがわかると、10/30 人が95% CI、20/30 人が50% CI、つまり中央値になるので、ここらへんのことを言っているのだろうか。
とはいっても、クラスメイト全員が芽亜里の奴隷であることはあるのはあるだろうが、ちょっと現実的ではない、かもしれない。このとき、奴隷割合 の事前分布をちょっといじってみる。
いま、一様分布を事前分布として仮定すると、も
も等確率で生じてしまう。30人のクラスのうち、さすがに芽亜里が奴隷にできるのは半分もいかないだろう、という仮定をおけば、違う事前分布を設定できる。ここで、
は確率なので、適当にベータ分布からサンプリングしてみよう。ベータ分布のパラメータはふたつで決まる。ここで、ふたつのパラメータをどちらも1 にすれば、一様分布になる。
今回は適当に3と10を選んだ。おそらく あたりが最頻値になる。
これで同じようにすると、50回手のうちを確認してもたかだか 程度までしか事後確率は変動しない。
% CI をプロットしてみると、10/30 人の奴隷、20/30 人の奴隷を推定値として採用するのは難しいような気がする。というのも、50枚賭けを2回負けた時点で、持ちチップが20枚くらいしかないので、1枚ずつ賭けて負けても50回手のうちを確認することができないからである。
とすれば、 より多い割合で奴隷がいる、と事前分布を仮定するのが無難っぽい。というのも、転校してきて最初の勝負を自信満々にふっかけられたら何か裏があると思っておいたほうがいいっぽい。
統計ガチ勢の人がいたらぜひとも数式的にごにょごにょ解いてください。
code <- " data{ int n; } parameters{ real<lower=0, upper=1> p[n]; } model{ p ~ beta(3, 10); for(i in 1:n) i ~ binomial(i, (2+p[i])/3); } " m <- stan_model(model_code=code) standata <- list(n=50) fit <- sampling(m, standata, iter=6000, warmup=1000) ex <- extract(fit) med <- colMeans(ex$a) cint <- c(99.9, 99, 95, 90, 80, 75, 50)/100 # 信用区間 ci <- apply(ex$p, 2, quantile, (1-cint)/2) # 累積確率分布 cm <- apply(ex$p, 2, sort) par(mar=c(5, 5, 2, 2), cex.axis=1.5, cex.lab=1.5) matplot(seq(0, 1, length=nrow(cm)), cm, type="l", col=matlab::jet.colors(ncol(cm)), lwd=3, lty=1, xlab="index", ylab="cumulative probability") abline(0, 1, lty=3) # 信用区間幅のプロット meari <- 10:20 cols <- matlab::jet.colors(length(cint)) par(mar=c(5, 5, 2, 2), cex.axis=1.5, cex.lab=1.5) matplot(t(ci), type="l", lwd=3, col=cols, lty=1, ylim=c(0, 1), xlab="カード内訳を確認した回数", ylab="クラス内の寝返り割合") abline(h=meari/N, lty=3) abline(h=n0/N) legend("topleft", legend=sprintf("%.1f%s", cint*100, "%"), title="CI", lty=1, lwd=3, col=cols, bg="white") # 分布のベタ書き dxy <- mapply(function(z){ dx <- density(unlist(z)) i0 <- 0 < dx$x & dx$x < 1 return(list(x=dx$x[i0], y=dx$y[i0])) }, apply(ex$p, 2, list), SIMPLIFY=FALSE) par(mar=c(5, 5, 2, 2), cex.axis=1.5, cex.lab=1.5) yl <- range(0, dxy) cols <- jet.colors(length(dxy)) plot(0, type="n", ylim=yl, xlim=c(0, 1), xlab="probability", ylab="density") for(i in seq(dxy)){ lines(dxy[[i]]$x, dxy[[i]]$y, col=cols[i], lwd=2) }
