ドラマの視聴率を予測する（２）〜stanを用いて、キャストの視聴率を算出する。
こんにちは、進撃のアナリストです。
前回の記事から早4ヶ月、今回もドラマの視聴率データを用いて、キャストの視聴率を算出したいと思います。
今回は前回とは異なりstanを用います。前回はキャスト視聴率は単なる按分平均を用いて算出していました。
今回はMCMCによるベイズ推定により、キャスト自体が持つ視聴率を直接モデルに組み込んで算出します。
参考 ATNDリンク BUGS/stan勉強会 #1
ちなみにこれが私のstanデビューです。
先に概念的な話として、今回のモデルや結果について解説してから、後半で技術的な話に触れます。
概念的な話
モデルのイメージは以下の通りです。
ドラマの初回視聴率＝B1×年度にとる効果＋B2×時間帯による効果＋B3×前回ドラマ最終話（ある場合は前シリーズ）最終回視聴率＋CP木村拓哉×出演ランク＋CP中居正広×出演ランク＋・・・（俳優890人分続く）・・・
（出演ランクは出演順を表す。主演なら5、2番目は4、3番目は3、2番目は2、1番目は1となります。）
要は単なる重回帰モデルです。これまでに使った時間帯等の効果によるモデルに加え、各キャストが出ているか出ていないかのフラグに対してキャストの効果となる係数、CP○○○○を付加したものです。
前回に記載した通り、重回帰やその他のロバストであるとされている回帰モデル手法等を試したのですが全く効果ありませんでした。
CPはキャストパワーの略です。ある特定のキャストが出演している際に、キャストパワー分の視聴率が付加されるイメージです。
今回はこのキャストパワーを俳優別に知ることを目的としています。
キャストパワーは一部のキャストのみ非常に高い値を持つ、ほとんどのキャストは0に近いと考えました（対数正規分布に従う）。
その結果、各キャストが持つキャストパワーBEST30は以下のようになりました。
結果
| 順位 | 俳優 | CP(中央値) | CP標準偏差 |
| 1 | 木村拓哉 | 2.16 | 0.23 |
| 2 | 黒木瞳 | 1.47 | 0.23 |
| 3 | 米倉涼子 | 1.44 | 0.19 |
| 4 | 仲間由紀恵 | 1.36 | 0.19 |
| 5 | 堂本剛 | 1.32 | 0.22 |
| 6 | 小雪 | 1.29 | 0.2 |
| 7 | 鈴木京香 | 1.29 | 0.21 |
| 8 | 豊川悦司 | 1.25 | 0.22 |
| 9 | 大塚寧々 | 1.25 | 0.24 |
| 10 | 竹野内豊 | 1.17 | 0.16 |
| 11 | 高橋克典 | 1.17 | 0.17 |
| 12 | 藤原紀香 | 1.16 | 0.2 |
| 13 | 中居正広 | 1.12 | 0.2 |
| 14 | 夏川結衣 | 1.11 | 0.2 |
| 15 | 坂口憲二 | 1.11 | 0.19 |
| 16 | 志田未来 | 1.1 | 0.19 |
| 17 | 香取慎吾 | 1.08 | 0.19 |
| 18 | 佐藤浩市 | 1.07 | 0.2 |
| 19 | 福山雅治 | 1.07 | 0.21 |
| 20 | ともさかりえ | 1.01 | 0.19 |
| 21 | 木村佳乃 | 1.01 | 0.2 |
| 22 | 森山未來 | 1 | 0.2 |
| 23 | 手越祐也 | 0.98 | 0.2 |
| 24 | 大沢たかお | 0.98 | 0.2 |
| 25 | 山田孝之 | 0.98 | 0.17 |
| 26 | 西島秀俊 | 0.97 | 0.19 |
| 27 | 黒木メイサ | 0.96 | 0.18 |
| 28 | 酒井若菜 | 0.93 | 0.19 |
| 29 | 天野ひろゆき | 0.9 | 0.22 |
| 30 | 田中美佐子 | 0.9 | 0.2 |
前回に比べると少しましになったという印象がありますが、ここに掲載していないものも含めて１を切るあたりから妥当性が怪しくなってきます。
とりあえずドラマに出演している回数が少ない俳優さんは排除したほうが良いかなというのが私の感想です。
しかしながら、その他の係数が想定していた数値と異なっていました。
以下はB1、B2、B3の中央値です。
> median(extract(fitmerged)$beta1)
[1] -1.425487
> median(extract(fitmerged)$beta2)
[1] -1.741203
> median(extract(fitmerged)$beta3)
[1] 2.380221
全て、正の値で１付近で収束することを期待していたのですが、時間帯と年度の効果が思いっきり負の値に振り切れています。
それを補うように前回ドラマの視聴率の係数が２を越えています。。。
時間帯や年度の効果をキャスト効果と同様に１を基準に考えたいのでこの辺りは完全に失敗です。
技術的な話
モデルを表すstanファイルは以下の通りです。
data {
int N;
int C;
vector[N] rating;
vector[N] trend;
vector[N] timezone;
vector[N] ratinglastdrama;
matrix[N, C] castmat;
}
parameters {
real beta1;
real beta2;
real beta3;
vector<lower=0>[C] mu;
}
transformed parameters {
vector[C] castbeta;
castbeta <- log(mu);
}
model {
vector[N] q;
beta1 ~ normal(1,1);
beta2 ~ normal(1,1);
beta3 ~ normal(1,1);
for(j in 1:C){
mu[j] ~ normal(-2.5,0.5);
}
for (i in 1:N) {
q[i] <- beta1 * timezone[i] + beta2 * trend[i] + beta3 * ratinglastdrama[i] + dot_product(castmat[i],castbeta);
rating[i] ~ normal(q[i], 2);
}
}
最初はBUGSを使っていたのですが、遅すぎて話にならずstanに変更しました。
実行時間の目安は以下の通りです。
| 言語 | 試行回数 | チェイン数 | 時間 |
| BUGS | 1500 | 3 | 36h |
| stan | 1500 | 4 | 30m |
圧倒的に早くなり、また収束のスピードもstanの方が早く感じました。
環境は以下の通りです。
用いた環境
課題
今回のはまだまだモデルとして未成熟で、stanを動かすテストです。次回までに以下のような課題を少しずつクリアーしたいきたいと考えています。
- 変数の取り入れ方・・・時間帯や年度の効果に関して今は代表値を用いて一つの変数にしていますが、フラグを利用した変数に変換して各それぞれの時間帯や年度のフラグに基づく効果を純粋に抽出してもよいかと考えています。色々工夫が必要なんだろうなと思います。
- 状態空間モデル化・・・俳優が旬であるかどうかをモデルに組み込みんだり、月9などに主演した俳優が高い評価を得てさらに次のドラマに抜擢されるような仕組みを導入したいと思っています。
stanを動かすまでにかなり苦労したのですが、これであれば現実的な計算時間でモデルをいじることができるので今からこれを使ってモデルの妥当性を上げていく作業を行いたいと思います。
