はじめに
スマホゲームのようなアプリの運営において「ユーザーの継続率」はKPIの中でも特に重要な指標として上げられます。特に近年のゲームアプリはリリース直後に大型のプロモーションを行って大量のユーザーを一気に獲得してしまうことも多く、そこでゲームを始めたユーザーにゲームを中長期的に継続してもらうことがアプリの命運を握ると言っても過言ではない状況も生まれています。そのため、ユーザーの継続率の分析、ならびにその向上は非常に価値のある取り組みと言えるのです。
継続に関する分析は簡単ではない
一方、ユーザーの継続を分析するのはそうそう簡単な問題ではありません。以下のブログでも述べられていますが、時系列データのランダムウォークや平均への回帰といった問題を考慮せずに短期的な数字ばかりを追ってしまうと、単なる数値のランダムな変動を施策の効果として勘違いしてしまい、意味のない施策を連発するか、場合によっては継続とは逆効果の施策を連発して中長期的に振り返ってみると継続率を悪化させてしまったといったような目も当てられないような結果に陥ってしまうこともあります。意味のある施策を積み上げて着実に継続率を向上させていくためには、施策の効果を「実験的に分析する」ことが非常に重要性だと言えます。
一方で、継続率に対する施策の影響を分析することは必ずしも簡単ではありません。特にゲームのサービスにおいては「全ユーザーに対する公平性」の問題が付きまとうため、A/Bテスト的に一部のユーザーにだけランダムに施策を実施して効果を見るといったことは基本的にはできず、いわゆる「観察データ」を用いた分析がメインとなります*1。また、ゲームのアップデートは一度にたくさんの機能追加や改修が入ることが多く分析の際に様々な施策効果が混じってしまったり、もともとゲームに対する熱量が高く継続傾向の強いユーザーが施策に反応する(=施策の影響で継続に繋がったのか、もともと継続傾向の強いユーザーが施策に反応したかが鶏卵になる)といったケースがどちらかというと一般的で、かつ継続率という数字自体がそんなに劇的に変わることが珍しい(=一見すると上がったか上がっていないのか分からないような変化度合いケースが多い)ため、施策がユーザーの継続に効いたのかどうかをちゃんと評価すること自体の難易度が高かったりします。
統計的因果推論 x 継続分析
そもそも観察データの分析がメインになることから、ゲームの分析で施策の効果をちゃんと測るには統計的因果推論のフレームワークを活用した分析を行うことが理想的であると言えます。統計的因果推論の手法は様々ありますが、本記事では継続という観点を絡めて傾向スコア (IPW推定量) を用いた生存分析について触れたいと思います*2。
なお、傾向スコアないしはIPW推定量に関しては星野本あたりが参考になります。
IPWsurvival パッケージを用いた生存分析
IPW推定量で重みづけした上で処置によるイベント発生の期間を分析するパッケージとして、R では IPWsurvival が存在します。CRAN からこちらのパッケージをインストールすることができるため、まずはそちらをインストールしてライブラリを読み込みます。また、今回は分析する対象として同パッケージの中に含まれている DIVAT を用いることにします。
library(IPWsurvival) data(DIVAT)
同パッケージに含まれる DIVAT データは DIVAT データバンクにある腎臓移植のレシピエントのサンプルです。同テーブルは6つのカラムから構成されています。
- age: 腎臓移植時のレシピエントの年齢
- hla: ドナーとレシピエントの間の HLA 不適合性の指標
- retransplant: 再移植か否か (0=初回の腎臓移植、1=2回目以降)
- ecd: 拡張基準ドナー (Expended Criteria Donor) か否か (1=ECD、0=その他)
- times: 各患者の追跡期間
- failures: イベント発生の有無 (0=右側打ち切り、1=イベント発生)。イベントは人工透析に戻るないしは移植腎が機能したままの死亡が観測された時
私は医学分野に関しては全くの門外漢ではあるのですが、何となく ECD の効果を推定したいが、ECD を受けるにあたって患者の年齢や HLA や再移植の有無が影響するが、それが生存期間に対しても影響を及ぼすため純粋に ECD の効果を見るために傾向スコアを用いて分析をするといったものでしょうか。
IPW 推定量を用いずにそのまま生存分析を行う
IPWsurvaival パッケージの中にある adjusted.KM という関数が IPW推定量を用いた生存分析 (Kaplan-Meier 推定法) を行う際に用いる関数で、adjusted.LR という関数が log-rank 検定を行う関数になります。
adjusted.KM 関数は weight のパラメータを設定しない場合特に重みづけを行わずに生存分析を行った結果を返すため、まずは特に IPW 推定量を用わなかった場合の結果を見てみます。
# IPW 推定量を用いずに生存分析を実施 res.km <- adjusted.KM(times=DIVAT$times, failures=DIVAT$failures, variable=DIVAT$ecd, weights=NULL) # グラフを描画 plot(NULL, xlim=c(0,13), ylim=c(0,1), ylab="Graft and patient survival", xlab="Time post-transplantation(years)") lines(res.km$times[res.km$variable==1], res.km$survival[res.km$variable==1], type="s", col=2, lty=2, lwd=2) lines(res.km$times[res.km$variable==0], res.km$survival[res.km$variable==0], type="s", col=1, lty=2, lwd=2)
IPW 推定量を用いて生存分析を行う
次に、IPW 推定量を用いて重みづけを行った上で生存分析を行った結果を見てみます。
# 傾向スコアの算出 (ロジスティック回帰で実施) Pr1 <- glm(ecd ~ age + hla + retransplant, data=DIVAT, family=binomial(link="logit"))$fitted.values # IPW 推定量を算出 W_ATE <- (DIVAT$ecd==1) * (1/Pr1) + (DIVAT$ecd==0) * (1)/(1-Pr1) # ATE: 平均処置効果 W_ATT <- (DIVAT$ecd==1) + (DIVAT$ecd==0) * (Pr1)/(1-Pr1) # ATT: 処置群における平均処置効果 # IPW 推定量で重みづけした上で生存分析 res.akm.ATE <- adjusted.KM(times=DIVAT$times, failures=DIVAT$failures, variable=DIVAT$ecd, weights=W_ATE) # ATE res.akm.ATT <- adjusted.KM(times=DIVAT$times, failures=DIVAT$failures, variable=DIVAT$ecd, weights=W_ATT) # ATT # グラフを描画 par(mfrow=c(1,2)) plot(NULL, xlim=c(0,13), ylim=c(0,1), ylab="Graft and patient survival", xlab="Time post-transplantation(years)", main="Weighted Sample (ATE Weights)") lines(res.akm.ATE$times[res.akm.ATE$variable==1], res.akm.ATE$survival[res.akm.ATE$variable==1], type="s", col=2, lwd=2) lines(res.akm.ATE$times[res.akm.ATE$variable==0], res.akm.ATE$survival[res.akm.ATE$variable==0], type="s", col=1, lwd=2) legend("bottom", legend=c("Adjusted estimator for SCD", "Adjusted estimator for ECD"), col=c(1,2), lty=c(1,1), lwd=2, cex=0.6) plot(NULL, xlim=c(0,13), ylim=c(0,1), ylab="Graft and patient survival", xlab="Time post-transplantation(years)", main="Weighted Sample (ATT Weights)") lines(res.akm.ATT$times[res.akm.ATT$variable==1], res.akm.ATT$survival[res.akm.ATT$variable==1], type="s", col=2, lwd=2) lines(res.akm.ATT$times[res.akm.ATT$variable==0], res.akm.ATT$survival[res.akm.ATT$variable==0], type="s", col=1, lwd=2) legend("bottom", legend=c("Adjusted estimator for SCD", "Adjusted estimator for ECD"), col=c(1,2), lty=c(1,1), lwd=2, cex=0.6) par(mfrow=c(1,1))
IPW 推定量で重みづけしなかった場合と比べて、解析の結果が変わっているのが分かるかと思います。
log-rank検定を実施
群ごとの生存曲線の間の差の優位性の検定である log-rank 検定は adjusted.LR 関数を用いることで簡単に行うことができます。
> adjusted.LR(DIVAT$times, DIVAT$failures, DIVAT$ecd, W_ATE) # ATE $`statistic` [1] 2.173458 $p.value [1] 0.02974587 > adjusted.LR(DIVAT$times, DIVAT$failures, DIVAT$ecd, W_ATT) # ATT $`statistic` [1] 2.589789 $p.value [1] 0.009603476
いずれも p < 0.05 となっており、両群の生存曲線の間には p = 0.05 の水準で有意な差があると言えそうです。
統計的因果推論の関連記事
*1:商品のバナーを変えるといった程度のことであれば実施はできる範囲に入りますが、例えばランダムで一部のユーザーだけログインボーナスが増える、イベントが開催される、割引が実施される......といったことはゲームにおいては基本実施が不可能な部類に入ります。
*2:傾向スコアを用いた生存分析に関しては Peter C. Autsin の The use of propensity score methods with survival or time-to-event outcomes: reporting measures of effect similar to those used in randomized experiments がとても参考になります https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4285179/pdf/sim0033-1242.pdf
