最強のM-1漫才師は誰だ
はじめに
この記事はStan Advent Calendar2017,12月11日のエントリー記事です。
Stanというかベイズ統計を勉強するには,犬4匹本もいいんですが,コワい本もまたいいですね。
コワい本の中に,7人の科学者というお題がありまして。7人がある実験でスコアをとったのだけど,どうも最初の二人は未熟者で変なスコアを出している。残りの5人の熟練者から考えるとスコアは10.00点ぐらいになるはずなんだけど・・・という話があります。これをモデリングするにあたって,真の値は一つなんだけど,評定者ごとに「測定誤差」を持っている,というモデルを立ててます。すなわち,科学者$$k$$によって得られた計測値を$$X_{k}$$とすると,$$X_{k} \sim N(\mu,\sigma_k)$$というわけです。
これを応用して,「評定者のくせを考慮したスコア」というのを算出することを考えました。 データはそう,M-1グランプリです。
データの出典はこちらです。 そして今年の分はビデオを見ながら手入力しました。
結論;ブラマヨが一番すごいんじゃないか
この結論に至る手続きを,以下解説いたします。
手続き
データを取り込みます。審査員は毎回変わりますから,欠損値が多いデータセットです。 元のファイルも置いておきますね。UTF-8のcsvファイルをこちらからダウンロードしてください。
## 年代 演者 上沼恵美子 松本人志 ## Min. : 1 笑い飯 : 9 Min. :81.00 Min. :70.00 ## 1st Qu.: 4 麒麟 : 5 1st Qu.:89.00 1st Qu.:85.00 ## Median : 7 ハライチ : 4 Median :90.50 Median :89.00 ## Mean : 7 フットボールアワー: 4 Mean :91.17 Mean :88.12 ## 3rd Qu.:10 千鳥 : 4 3rd Qu.:95.00 3rd Qu.:93.00 ## Max. :13 POISON GIRL BAND : 3 Max. :98.00 Max. :97.00 ## (Other) :90 NA's :73 NA's :46 ## 博多大吉 春風亭小朝 中川礼二 渡辺正行 ## Min. :50.00 Min. :65.00 Min. :87.00 Min. :75.00 ## 1st Qu.:70.00 1st Qu.:80.00 1st Qu.:89.00 1st Qu.:87.00 ## Median :84.00 Median :88.00 Median :90.00 Median :89.00 ## Mean :79.13 Mean :84.28 Mean :90.71 Mean :88.16 ## 3rd Qu.:91.00 3rd Qu.:90.00 3rd Qu.:92.25 3rd Qu.:90.50 ## Max. :97.00 Max. :95.00 Max. :95.00 Max. :95.00 ## NA's :72 NA's :90 NA's :91 NA's :64 ## オール巨人 増田英彦 岩尾望 吉田敬 ## Min. :79.00 Min. :85 Min. :85.00 Min. :83.00 ## 1st Qu.:86.00 1st Qu.:87 1st Qu.:89.00 1st Qu.:85.00 ## Median :88.50 Median :89 Median :90.00 Median :89.00 ## Mean :88.46 Mean :89 Mean :90.67 Mean :87.56 ## 3rd Qu.:91.00 3rd Qu.:91 3rd Qu.:92.00 3rd Qu.:90.00 ## Max. :96.00 Max. :93 Max. :96.00 Max. :93.00 ## NA's :73 NA's :110 NA's :110 NA's :110 ## 徳井義実 富澤たけし 石田明 佐藤哲夫 ## Min. :88.00 Min. :89.00 Min. :83 Min. :88.00 ## 1st Qu.:89.00 1st Qu.:91.00 1st Qu.:87 1st Qu.:89.00 ## Median :89.00 Median :92.00 Median :88 Median :90.00 ## Mean :90.56 Mean :91.89 Mean :89 Mean :90.56 ## 3rd Qu.:91.00 3rd Qu.:93.00 3rd Qu.:92 3rd Qu.:92.00 ## Max. :96.00 Max. :94.00 Max. :94 Max. :93.00 ## NA's :110 NA's :110 NA's :110 NA's :110 ## 哲夫 島田紳助 南原清隆 大竹まこと ## Min. :87.00 Min. : 50.00 Min. :76.00 Min. :75.00 ## 1st Qu.:88.00 1st Qu.: 80.25 1st Qu.:84.75 1st Qu.:82.00 ## Median :90.00 Median : 86.00 Median :88.00 Median :85.00 ## Mean :89.78 Mean : 84.77 Mean :87.78 Mean :85.18 ## 3rd Qu.:91.00 3rd Qu.: 90.00 3rd Qu.:90.25 3rd Qu.:89.00 ## Max. :93.00 Max. :100.00 Max. :98.00 Max. :97.00 ## NA's :110 NA's :37 NA's :83 NA's :47 ## 宮迫博之 中田カウス 東国原英夫 ラサール石井 ## Min. :88 Min. :79.00 Min. :85.00 Min. :68.00 ## 1st Qu.:90 1st Qu.:86.00 1st Qu.:86.00 1st Qu.:82.00 ## Median :91 Median :90.00 Median :88.00 Median :86.00 ## Mean :92 Mean :89.98 Mean :87.67 Mean :86.84 ## 3rd Qu.:93 3rd Qu.:95.00 3rd Qu.:89.00 3rd Qu.:92.00 ## Max. :98 Max. :98.00 Max. :92.00 Max. :96.00 ## NA's :110 NA's :38 NA's :110 NA's :64 ## 島田洋七 西川きよし 立川談志 鴻上尚史 ## Min. :75.00 Min. :75.00 Min. :50 Min. :73.00 ## 1st Qu.:84.00 1st Qu.:80.00 1st Qu.:70 1st Qu.:82.25 ## Median :90.00 Median :86.00 Median :70 Median :83.50 ## Mean :88.09 Mean :85.68 Mean :70 Mean :81.80 ## 3rd Qu.:92.00 3rd Qu.:90.00 3rd Qu.:70 3rd Qu.:84.00 ## Max. :98.00 Max. :97.00 Max. :80 Max. :85.00 ## NA's :74 NA's :100 NA's :110 NA's :109 ## 青島幸男 ## Min. :75.00 ## 1st Qu.:76.25 ## Median :80.00 ## Mean :81.50 ## 3rd Qu.:85.00 ## Max. :90.00 ## NA's :109
データは次のようにして,縦長に整形しました。
m1 %>% tidyr::gather(審査員,val,-年代,-演者) %>% na.omit ->m1.long m1.long$審査員 <- factor(m1.long$審査員) m1.long$演者 <- factor(m1.long$演者) r.name <- levels(m1.long$審査員) p.name <- levels(m1.long$演者) datastan <- list(L=nrow(m1.long), N=max(as.numeric(m1.long$演者)), M=max(as.numeric(m1.long$審査員)), idX=as.numeric(m1.long$演者), idY=as.numeric(m1.long$審査員), X=m1.long$val)
モデルは先ほどと同様で,漫才師(演者)$$k$$の真のお笑いの実力$$\theta_k$$を持っていたとして,審査員(評定者)$$j$$がつけた得点$$X_{jk}$$は$$X_{jk}\sim N(\theta_k,\sigma_j)$$としています。
data{ int<lower=1> L; //data Length int<lower=1> N; //number of players int<lower=1> M; //number of rators int idX[L]; //player ID index int idY[L]; //rator ID index real X[L]; // scores } parameters{ real<lower=0> theta[N]; real<lower=0> sig[M]; } model{ for(l in 1:L){ X[l] ~ normal(theta[idX[l]],sig[idY[l]]); } theta ~ normal(50,100); sig ~ cauchy(0,5); } generated quantities{ real PredScore[N,M]; for(n in 1:N){ for(m in 1:M){ PredScore[n,m] = normal_rng(theta[n],sig[m]); } } }
これで推定した結果をグラフにします。まずは各演者の「漫才力$$\theta_k$$」のプロット。 点がEAP推定値,横幅は50%と95%のHDIです。
かわいそうなDonDokoDon。っていうか皆さん覚えてます?ぐっさんと平畠のコンビなんですが。
ただ,このモデルの中に含まれていることは,「評定者のスコアは100点満点で絶対的な感覚を持って採点している」ということ。だんだんスコアがインフレしているかも・・・というのは仮定していません。初期のメンバーはスコアが低くなりがちですが,ひょっとしたら審査員もみんな慣れてきて,スコアが上がってきたとかってことがあるかもです。
次に審査員の$$\sigma_j$$の結果。これは歪んでいるのでMED推定値にしています。