xiangze's sparse blog

機械学習、ベイズ統計、コンピュータビジョンと関連する数学について

東方プロジェクト人気投票の統計解析　ベイズ統計モデリング編と本

時間が空いてしまいましたが以前ブログに書いた　
東方プロジェクト人気投票の統計解析 記述統計編 - xiangze's sparse blog

東方人気投票の統計モデルの事後分布が計算できたので記します。

また同人誌を作りました。

キャラクターのモデルと事後分布

和モデル

キャラクターを

  • 整数作品の登場キャラ M
  • 整数作品のボスキャラ bossとそのレベル Lv
  • 非整数作品 Sub
  • 秘封倶楽部、その他

に分類しそれぞれ対応する影響力のパラメーター係数sigma b,s,と各キャラクターがそれぞれの属性に当たるかを示す変数との積の和によって構成されるモデルを考えます。i番目のキャラクターのt回目の投票での正規化した投票数 nVote_{i,t}がディリクレ分布Dirを用いて

nVote_{i,t} ~ Dir(\sum_{l=1}^{TM} M_{j(i,t-l),t-l} \sigma_{j,l} + \sum_{l=1}^{TM} boss_{j(i,t-1),t-l}Lv_i b_{j(i,t-l),l} +\sum_{l=1}^{TM} Sub_{j(i,t-l)} s_{j,l}+\epsilon_i )

と表されるとします。ディリクレ分布Dirは各要素の合計が1になるようなベクトルに値を持つ確率分布関数であり正規化した投票結果に対して用いることができます。カッコ内の値が大きければ大きいほど高い得票比率になりますが1を超えることはありません。また0を下回ることもありません。j(i,t-l)という表記はi番目のキャラクターのt-l回目の投票の登場作品を表す関数です。
コードは
https://github.com/xiangze/toubou_vote_analysis/blob/main/analyse_char.py

https://github.com/xiangze/toubou_vote_analysis/blob/main/model/charpower_template_hyper.stan

にありj(i,t-l)という表記はiとtに対応するfor文の中で条件に一致したところで代入している部分に相当します。

非整数作品でのみ登場するキャラと再登場キャラは区別していませんが、同一の影響を受けると仮定しました。また以前の投票結果(人気)に応じて再登場の機会が変わると考えられますがその効果も明示的に入れていません。

このようにデータをグループ分けするモデルは階層モデルと呼ばれます。今回のキャラとそれが登場した整数作品のような階層的な関係が実際にある場合に自然な過程となります。データの分布の観点からは階層モデル、特に個体(individual)ごとのバラツキを考慮した混合モデルはデータのばらつきが大きくなってしまう過分散と呼ばれる問題に対処できるとされます(「個体差」の統計モデリング http://hdl.handle.net/2115/26401 )。

結果として計算された事後分布は以下のようになります。

キャラクター和モデルの係数事後分布


実行結果は
https://github.com/xiangze/toubou_vote_analysis/blob/main/stan_analyse.ipynb
でも見れます。

ベクトルのパラメーターは各要素を重ね合わせて描いています。chainによって分布が大きく変わらないことも確認できます。

個別のキャラクターに起因する項(indivisual)の影響が非常に大きくまたわかりづらいですが紅魔郷が強いことなどがそこに現れています。

個体の効果が上位のキャラクターですが初期の投票回での咲夜の強さの寄与が大きく、実質3強状態になっています。それ以下とは大きな差があり、また幽々子、紫間で大きな差があります。その次のグループでは小傘がトップにいます。全キャラクターの結果に興味がある方は
https://github.com/xiangze/toubou_vote_analysis/blob/main/book/img/char_indivisual_mean.png

をご覧ください。
平均値が大きいものが標準偏差が大きい傾向、近年の作品では投票回数が少ないため標準偏差が大きいのもは想定通りです。これらの平均、標準偏差の値は箱ひげ図やバイオリンプロットと呼ばれる手法で重ね書きしたほうがわかりやすいのですが今回は断念しました*1。

キャラクター和モデルの個体項平均値上位
キャラクター和モデルの個体項平均値上位の標準偏差

和モデル(整数、非整数作品固有の寄与を入れた場合)

微妙に違うモデルとして整数、非整数作品固有の寄与の項を追加したバージョン
nVote_{i,t} ~ Dir(\sum_{l=1}^{TM} M_{j(i,t-l),t-l} \sigma_{j,l} + \sum_{l=1}^{TM} boss_{j(i,t-1),t-l}Lv_i b_{j(i,t-l),l} +\sum_{l=1}^{TM} Sub_{j(i,t-l)} s_{j,l} + titl_{j(i)} titlpow_j + noninttitl_{j(i)} noninttitlpow_j + \epsilon_i )

https://github.com/xiangze/toubou_vote_analysis/blob/main/model/charpower_template_withtitlesubtitle.stan

を試してみますと以下のように異なる結果が得られます。

キャラクター和モデル(整数、非整数作品固有の寄与を入れた場合)の係数事後分布

最初の和モデルと同様にindivisualの影響が強く、それと比較するとタイトルによる影響はあまり見えませんでした。boss,非整数作品による寄与はばらつきが少なく一様になっています。作品タイトル固有の効果titlepow,noninttitlepowは整数作品ボスとしての登場の効果と同程度で、非整数や自機としての登場寄りはその効果は弱いです。

個体の効果が上位のキャラクターとその値も和モデルと大きな違いはありません。

キャラクター和モデル(整数、非整数作品固有の寄与を入れた場合)平均値上位
キャラクター和モデル(整数、非整数作品固有の寄与を入れた場合)標準偏差

積和モデル

和モデルではindividualの影響が強く、タイトルによる影響はあまり見えませんでした。そこで個体の影響stem:[\epsilon_i]に登場作品など属性の影響がかけ合わさるようなモデルを考えました。

nVote_{i,t} ~ Dir( (\sum_{l=1}^{TM} M_{j(i,t-l),t-l} \sigma_{j,l} + \sum_{l=1}^{TM} boss_{j(i,t-1),t-l}Lv_i b_{j(i,t-l),l} +\sum_{l=1}^{TM} Sub_{j(i,t-l)} s_{j,l})\epsilon_i )

https://github.com/xiangze/toubou_vote_analysis/blob/main/model/charpower_template_prodsum.stan

キャラクター積和モデルの係数事後分布

しかしながらこのモデルは収束せず事後分布のtraceはバラバラになってしまいました。すなわち妥当ではないモデルと言えるのではないでしょうか。しかし後述するように音楽に関してはそうはなりませんでした。

音楽のモデルと結果

キャラと同様に和モデル、積和モデルを考えます。違いは旧作、秘封倶楽部、音楽CDリリース作品が大きな割合を占めていることでそこでは登場順序はあまり重要ではないという仮定をしています。また非整数作品での再録曲の情報は用いませんでしたが実際には人気に影響があると考えられます。

和モデル

キャラクターの場合と同様に

  • 整数作品のボスキャラ、道中曲とそのレベル inttitlepow,平均 mu_t,分散 sigma_t
  • 非整数作品 noninttitlepow,平均 mu_i, 分散 sigma_i
  • 書籍についているCD bookpow
  • 秘封倶楽部 hifuupow
  • 音楽CDでのリリース(秘封含む) CDpow
  • 旧作 oldpow
  • 西方 seihoupow
  • オリジナル orgpow
  • その他 otherpow

の各寄与を表す係数の個別の項ε(indivisual)を加えたものをディリクレ分布の引数としています。

そのコードとモデルは以下のようになります。

https://github.com/xiangze/toubou_vote_analysis/blob/main/analyse_music.py
https://github.com/xiangze/toubou_vote_analysis/blob/main/model/music_template_sum.stan

またキャラクター解析の時と異なりある属性に含まれるかどうかというflag(is*という変数)をdataframeとして持っていますがこの方がデータサイズは大きくなるものの計算が非常に高速化することがわかりました(キャラクターの場合は約2日、楽曲は30分程度)。

事後分布の計算結果は以下のようになります。
https://github.com/xiangze/toubou_vote_analysis/blob/main/stan_analyse_music.ipynb

楽曲和モデルの係数事後分布

各係数(pow)と整数、非整数作品の平均、分散(mu_t,mu_i,sigma_t,sigma_i)の分布がなめらかになっていてより一致性が高いように見えます。キャラクターの場合と同様個体(indivisual)の項が最も大きくまたばらついており次に整数作品そして非整数作品の寄与係数が大きく、またばらついています。それと比べるとその他（書籍、秘封倶楽部、CD、旧作、西方、オリジナル、その他）の寄与は小さいです。またそれらの事後分布は0に近い方が大きな指数分布に近い形をしています。

individualの平均値上位はU.N.オーエンやセプテットが突出しており、その次に墨染の桜、ネクロファンタジアと続き紅魔郷、妖々夢のラストに近い曲であり妥当、王道というと言えるのではないでしょうか。しかしスモーキングドラゴンはかなり以外でもとの人気投票のデータを見てもそれほど得票は大きくありません。データ処理かプログラムへの渡し方にバグがあるのではないかという懸念がありますが一方でそれ以外の楽曲は妥当なので不思議さが残ります。標準偏差も大まかには平均と同じで新作はばらつきしかしこれを見てもスモーキングドラゴンの大きさの説明は付きません。全楽曲の結果に興味がある方は
https://github.com/xiangze/toubou_vote_analysis/blob/main/book/img/indivisual_music_sum_mean.png
をご覧ください。

楽曲和モデルの個体項平均値上位
楽曲和モデルの個体項平均値上位の標準偏差

積和モデル

個別楽曲の項εiを他の項の総和にかけ合わせたモデルです。このモデルはキャラクターの場合とは異なり収束しました。整数タイトル、非整数タイトルの係数のばらつきは和モデルより小さく、タイトルによってはその力が秘封や書籍、CDよりも小さいです。非整数タイトルのほうが平均して寄与が大きいことがわかります。
また旧作は平均的寄与は小さいのですが、以下に示すようにその中には作品の寄与に依らない突出した個性を持つ楽曲があります。

https://github.com/xiangze/toubou_vote_analysis/blob/main/model/music_template_sumprod.stan

楽曲積和モデルの係数事後分布

individualの値上位を見るとレインカーネーションを筆頭に旧作のものがいくつかあり実力派という感じの並びです。U.N.オーエンやセプテットもあり作品効果だけではない力を感じます。
全楽曲の結果に興味がある方は
https://github.com/xiangze/toubou_vote_analysis/blob/main/book/img/indivisual_music_sumprod_mean.png

をご覧ください。

楽曲積和モデルの個体項平均値上位
楽曲積和モデルの個体項平均値上位の標準偏差

事前分布

パラメーターに加えられる仮定である事前分布は計算の収束に大きな影響を与えます。stanでは事前分布を指定しない場合は一様分布が用いられますが、無限の区間幅を持つため現実的ではありません。現実的な範囲内に収めることが必要でありまたそうすることにより事後分布が収束するスピードが早まります。
正規分布があるいはパラメータが正の値を持つという仮定を置くのであれば指数分布、整数のパラメーターに対しては対してはポアッソン分布などが用いられます。「StanとRでベイズ統計モデリング」にも事前分布の選択について書かれており、逆ガンマ分布あるいは半コーシー分布という裾の広い分布関数を使用すると良いとされています。

https://stats.stackexchange.com/questions/237847/what-are-the-properties-of-a-half-cauchy-distribution

発展的話題、展望

上記のモデルに組み入れられなかった特徴として

  • 新作の影響の持続性のパラメーター化
  • 楽曲の再録
  • 作品、キャラ、楽曲関係と因果関係

などがあります。今、推し活が盛んです。東方でも何が入り口になるのかを知ることはコミュニティの恒常的な発展にとって重要です。アンケートデータからは投票者の入れ替わりが激しいことがわかりましたが*2、作品から特定のキャラクターを推すようになるのか、音楽から入っていくのかなどの説明が可能になるかもしれません。このような因果関係の推論は統計学でも一分野を築いています。また非整数作品での再登場、再録はその時点までの人気で決定されていると考えられそこまで説明できれば有益です。

他のデータの利用

pixivの画像やニコニコ動画、静画ではタグに対応したコンテンツが豊富にあり、時間的解像度も高いです。カップリングなどの組み合わせがどこまでも詳しく解析できます。またpixiv百科事典やニコニコ大百科には著名な東方アレンジの情報が豊富にあり解析しがいがあるかもしれません。

予測

機械学習は予測のための手法であるのに対し、統計分析は説明のための手法であると言われることがあります。しかしながら今回の投票データのような時系列に関するデータがあれば、そこから未来を予測したいという願望が出てきます。記述統計のところで触れましたが徐々に新作登場時の人気の勢いがなくなり上位が固定化しています。そうすると中、下位のキャラクターのテコ入れがあると活性化すると予測できますが、そんな中2023年東方獣王園が発表され、獣キャラの活躍が期待されています。

別の観点でのモデリング

上ではキャラクター、音楽を主体としたモデリングをしましたが、投票者の行動に基づいたモデリングも考えられます。すなわち個々の票に対して作品選択、キャラクターの選択を確率変数として推しへの票もモデリングに取り入れられそうですが、その場合票数分だけ計算の繰り返しが必要になり時間がかかりそうです。

ためになる本

pythonについての記述はありませんが英語版と著者のgithubサイト(https://github.com/MatsuuraKentaro/Bayesian_Statistical_Modeling_with_Stan_R_and_Python )にはpython,stanのコードがあります。

同人誌

5/7(日)東京ビッグサイト

booth、技術書展にもあります。
overfitting.booth.pm

techbookfest.org

その他

static.chunichi.co.jp

競艇は他の公営競技と違い発走前にエンジンを調整することができその効果も推定できます。

*1:arvizのさサイト https://python.arviz.org/en/stable/examples/index.htmlに豊富な実例があります。

*2:東方カンファレンスで考察がされています。 https://twitter.com/blue_comment/status/1644603082469429249