アクチュアリーの技術をマーケティング (離脱 予測) に応用してみる @ill-identified 26.10.2019 1

自己紹介 • Twitter: @ill_identified → • ブログ: http://ill-identified.hatenablog.com/ • LinkedIn: https://www.linkedin.com/in/satoshi-katagiri/ • Twitter: https://twitter.com/ill_Identified • github: https://github.com/Gedevan-Aleksizde 2

自己紹介 • 経済学 (修士) → エンジニア (?) → 機械学習エンジニア (?) • 現在の勤務先: Web 広告の会社 • 創業 20 周年なのでロゴが変わった →

アジェンダ イントロダクション 先行研究 今回の提案 実験 まとめ 4

イントロダクション

事の発端 • 今年 6 月の人工知能学会の発表 (森下, 2019) を見て唐突に 思いつく • 発表者の所属と今日の会場が同じなのは全くの偶然 5

今回話すこと • 森下 (2019) の離脱予測の問題 • Abema TV ユーザーの離脱要因を見るためのコホート分析 • アクチュアリーや人口統計学で使われる方法を応用できな いか • その方法の名はリー‧カーター (Lee-Carter) 法 6

先行研究

先行研究 (森下, 2019) の問題意識 • マーケティングには顧客の生涯価値 (LTV) の計算が重要 • 計算には顧客がいつ離脱するかの予測が必要 • Abema TV での離脱予測したい • ユーザ登録がないサービス • ユーザの利用開始タイミングがばらばら • ユーザの離脱タイミング不明 7

先行研究での解法 1/2 • 訪問間隔をモデリング • 冪分布にしたがう場合もあると言う先例 • 下図のようにユニークユーザー数の相対頻度に対数線形回 帰モデルを当てはめる Figure 1: 訪問間隔が冪分布に従う例 (森下, 2019 図 1 より) 8

先行研究での解法 2/2 • 期間‧訪問間隔別に UU 数集計 • 推定した冪分布による理論値との残差を計算しプロット • 残差を特異値分解で要因分解 •「ユーザー由来の継続要因」 「施策由来の復帰要因」 Figure 2: 残差のヒートマップ (森下, 2019 図 2 より) 9

今回の提案

先行研究の気になる点 1. 離脱判定の定義には根拠がない 2. 大域的に冪分布であることが正しいという前提 • (1) は今回の代案でも未解決 • (2) 冪分布が真という前提で, 残差分析が意味をなす. • 以前の発表, Causal Impact の話 でも書いた • 参加者中の順位で健康度を判定する健康診断に意味はない • 代案: Lee-Carter 法を利用してモデリングする 10

Lee-Carter 法とは • Lee and Carter (1992) によって提案された • アクチュアリー業務や人口統計学で利用される • 生年‧年齢コホートごとの死亡率の確率モデル • 生命保険や年金のリスク計算に応用 (井川, 2014) • 将来人口推計へ応用 (Wiśniowski et al., 2015) • 縁起が悪い? • 出生率予測モデルにも使われる (Lee, 1993) から安心 (?) 11

Lee-Carter モデルの定義 • 生年 t, 年齢 s のグループ (コホート) に対する死亡率曲線 µ(t, s) を近似. • 死亡率の対数が以下のような線形式で表せると仮定 ln µ(t, s) =αs + βsκt • 生年グループ固有要因の αs, βs と, 時間要因の κt に分解 • データから, 各グループの人口 et,s と死亡者数 dt,s から, mt,s :=dt,s/et,s として, µ(t, s) に当てはめる • dt,s がゼロ件のグループがある場合は注意 12

Lee-Carter モデルの例 • 単純な式だが, いろいろな曲線を表現できる. • 死亡率を合計特殊出生率に置き換えモデリングも可 (右下) • 他の事象にも当てはめられそう? 0 25 50 75 0 25 50 75 1990 1995 2000 2005 2010 0 25 50 75 age age year age αx βx κt lnμx (t) 1990 1995 2000 2005 2010 year Lee-Carter Model (Mortality) <18 22 27 32 37 40< <18 22 27 32 37 40< 1990 1995 2000 2005 2010 <18 22 27 32 37 40< age age year age αx βx κt lnμx (t) 1990 1995 2000 2005 2010 year Lee-Carter Model (Fertility) Figure 3: 日本の人口統計に LC モデルを当てはめた例 13

Lee-Carter モデルの計算アルゴリズム 1/2 • 単回帰モデルと似ている • 最小二乗法推定量の定義は以下の通り. min ˆ αs,ˆ βs,ˆ κt ∑ t,s (ln mt,s − ˆ αs − ˆ βsˆ κt)2 • 2 次関数の極値を求めるのは簡単 • αs の推定量は, ˆ αs := 1 T T ∑ t=1 ln mt,s 14

Lee-Carter モデルの計算アルゴリズム 2/2 • βs, κt を決めるのは最小二乗法では不可能 (双線形モデル) min ∑ t,s [ (ln mt,s − ˆ αs) − ˆ βsˆ κt ] 2 • 成分が at,s := ln mt,s − ˆ αs の T × S 行列 A を特異値分解す ると, A =UΣV⊤ • 最小化の条件は, 特異ベクトルで ˆ βsˆ κt =u⊤ 1,s v1,t • 解が一意にならないため, 以下の制約を満たすよう正規化 ∑ s βs = 1, ∑ t κt = 0 15

Lee-Carter モデルによる予測 • 求めた時系列 {ˆ κ1, ˆ κ2, · · · , ˆ κT} に何らかのモデルを当ては め, ˆ κT+1 , ˆ κT+2 , · · · を外挿し将来の予測を計算 ˆ µT+h,s = exp(ˆ αs + ˆ βsˆ κT+h) • ARIMA モデルを当てはめることが多い • s 方向を延⻑するのは難い? 16

Lee-Carter モデルの計算方法 • R では demography パッケージの lca() 関数 • 入力データを年齢‧年の行列として渡す • tidy like じゃないので少しめんどくさい 1 require(demography) 2 # events はコホートごとのイベント発生件数 3 # population は人口 4 d <- demogdata(events, population, ages, 5 years, type, label, name, lambda) 6 fit <- lca(d) 7 plot(d) # a, b, k の要素をプロット 8 plot(forecast(result_mort, h=50)) # 50期間先将来予測 17

補足: lca() の使い方 • 同一コホートが斜めに並んでいる • 森下 (2019) の図 2 と似ている age|year 1989 1990 · · · 2019 0 x1989,0 x1990,0 · · · x2019,0 1 x1989,1 x1990,1 · · · x2019,1 2 x1989,2 x1990,2 · · · x2019,2 . . . . . . . . . ... . . . 100+ x1989,100+ x1990,100+ · · · x2019,100+ Table 1: LC モデルの入力のイメージ 18

Lee-Carter モデルの拡張 • LC モデルは単回帰モデルと似ていた. • 実際正規分布のモデルともみなせる • カウントデータならポアソン回帰のようにできないか? • できる (Brouhns et al., 2002) 19

ポアソン対数双線形モデルの定義 • et,s × µt,s がイベント発生件数の期待値なのでポアソン分 布に当てはめる dt,s ∼ Po(et,sµt,s), µt,s = exp(αs + βsκt) • 対数尤度は ln L(α, β, κ) = ∑ t,s [dt,s(αs + βsκt) − et,s exp(αs + βsκt)] + C 20

補足: 一般化双線形モデル • 多くのバリエーション (一般化双線形モデル) • 線形回帰モデルに対する一般化線形モデル (GLM) 的存在 • Renshaw and Haberman (2003): 二重双線形モデル (αs + β(1) s κ(1) t + β(2) s κ(2) t ) • Renshaw and Haberman (2003): 過分散ポアソン分布 • 井川 (2014): 自己回帰項の追加, 負の二項分布の適用 • 小暮‧⻑谷川 (2005) も参考になる 21

一般化双線形モデルの計算方法 • 元の LC モデル同様に, パラメータが一意にならない. • Brouhns et al. (2002) では Goodman (1979) のアルゴリズ ムを利用 • ニュートン法の各イテレーションで正規化してるだけ • 後から正規化する方法もできなくはない • Wiśniowski et al. (2015) はモデルが複雑なので MCMC 法 で計算 22

補足: ニュートン法とは • Newton-Raphson 法とも • 非線形連立方程式 f(x) = 0 を解くアルゴリズム • 最尤法は勾配ベクトル =0 (∇ ln L(θ) = 0) を解く • 以下のような反復計算を繰り返すと解に近似できる x(i+1) =x(i) − J−1 f (x(i))f(x(i)) • (⻑所) ブラウワーの不動点定理により, 必ず解に辿り着く • ∵ ヤコビアン = 対数尤度のフィッシャー情報行列 • しかも速い (二次収束) • (短所) 初期値に一番近い解の 1 つしか求めてくれない • 今回は問題にならない欠点 23

離脱予測への応用 • 人口 et,s を t − s 期間前に初訪問したユーザー数, イベント 発生 dt,s を離脱件数に置き換えればそのまま使えそう • 訪問間隔ではなく訪問回数を扱うモデル • αs, βs は初訪問から s 期間経過したユーザー特有の継続要 因 • κt は絶対時間でみた継続要因の推移 • 施策の影響はここで拾えそう • 残差はモデルの当てはまりの確認にしか使わない • 今回はBrouhns et al. (2002) のポアソン対数双線形モデル • モデルの形状は議論の余地あり • 冪分布と指数分布のどっちか • あるいは負の二項分布か 24

実験

使用するデータ • Abema TV のデータは使えない • 別の似たような状況のデータを用意 • 森下 (2019) の方法とどっちが良いかの決定的な証拠には ならない • rfm パッケージのデータが似たような状況なので代用 25

プログラム i 小暮‧⻑谷川 (2005) による R での実装例 (コメントを一部 省略) 1 dl<-1 #dl は対数尤度の増加量。 初期値として 1 を採用する 2 h<-0#h は繰り返し回数を入れる変数(最初は 0) 3 #死亡者数 D を計算する関数の作成 4 d.keisan<-function(e,a,b,k){ e*exp(a+b%*%k) } 5 while(dl>10^-10){ 6 ll0<-sum(d*(a+b%*%k)-e*exp(a+b%*%k)) 7 d0<-d.keisan(e,a,b,k) 8 a=a-apply(d-d0,1,sum)/(-apply(d0,1,sum)) 9 d1<-d.keisan(e,a,b,k) 10 k=k-apply((d-d1)*b,2,sum)/(-apply(d0*b^2,2,sum)) 26

プログラム ii 11 k<-k-mean(k) #kの合計は 0 12 d2<-d.keisan(e,a,b,k) 13 b=b-apply((d-d2)*k,1,sum)/(-apply(d2*k^2,1,sum)) 14 b<-b/sum(b) #βの合計は 1 15 ll1<-sum(d*(a+b%*%k)-e*exp(a+b%*%k)) 16 dl<-ll1-ll0 #対数尤度の増分 17 h<-h+1 #繰り返しの回数のカウント 18 } 27

28

29

まとめ

今回のまとめ/反省 • 森下 (2019) と代案の優劣の決定的証拠はない • アイディアの堀り下げが必要 & データがない • Lee-Carter モデルは単純だが自由度のあるモデル • demography パッケージで計算できる • 一般化双線形モデルも実装が比較的簡単 • コホート単位の分析なら他の課題にも応用できるか? • 阿部 (2011) のような切断分布を使った潜在変数モデルの 応用を考えたほうが良い? • ようやく Tokyo.R っぽい発表ができた? • LC モデルは 5 年くらい前に調べたが, 自分の残したメモを 見返すまでほとんど忘れていた • 時間に余裕を持って作業しましょう 30

参考文献

参考文献 i Brouhns, Natacha, Michael Denuit, and Jeroen K. Vermunt (2002) “A Poisson Log-Bilinear Regression Approach to the Construction of Projected Lifetables,” Insurance: Mathematics and Economics, Vol. 31, No. 3, pp. 373–393, DOI: 10.1016/S0167-6687(02)00185-3. Goodman, Leo A. (1979) “Simple Models for the Analysis of Association in Cross-Classifications Having Ordered Categories,” Journal of the American Statistical Association, Vol. 74, No. 367, pp. 537–552, DOI: 10.2307/2286971. 31

参考文献 ii Lee, Ronald D. (1993) “Modeling and Forecasting the Time Series of US Fertility: Age Distribution, Range, and Ultimate Level.,” International journal of forecasting, Vol. 9, No. 2, pp. 187–202, DOI: 10.1016/0169-2070(93)90004-7. Lee, Ronald D. and Lawrence R. Carter (1992) “Modeling and Forecasting U.S. Mortality,” Journal of the American Statistical Association, Vol. 87, No. 419, pp. 659–671, September, DOI: 10.1080/01621459.1992.10475265. Renshaw, A. E. and Steven Haberman (2003) “Lee-Carter Mortality Forecasting with Age-Specific Enhancement,” Insurance: Mathematics and Economics, Vol. 33, No. 2, pp. 255–272, DOI: 10.1016/S0167-6687(03)00138-0. 32

参考文献 iii Renshaw, Arthur and Steven Haberman (2003) “Lee-Carter Mortality Forecasting: A Parallel Generalized Linear Modelling Approach for England and Wales Mortality Projections,” Journal of the Royal Statistical Society: Series C (Applied Statistics), Vol. 52, No. 1, pp. 119–137, January, DOI: 10.1111/1467-9876.00393. Wiśniowski, Arkadiusz, Peter W. F. Smith, Jakub Bijak, James Raymer, and Jonathan J. Forster (2015) “Bayesian Population Forecasting: Extending the Lee-Carter Method,” Demography, Vol. 52, No. 3, pp. 1035–1059, DOI: 10.1007/s13524-015-0389-y. 33

参考文献 iv 阿部誠 (2011) 「RFM 指標と顧客生涯価値: 階層ベイズモデル を使った非契約型顧客関係管理における消費者行動の分析」 ， 『日本統計学会誌』 ，第 41 巻，第 1 号，51–81 頁, retrieved from here． 井川孝之 (2014) 「残差構造解析による Lee-Carter モデルの拡 張と年金負債評価」 ，博士論文，総合研究大学院大学, retrieved from here． 小暮厚之‧⻑谷川智弘 (2005) 「将来生命表の統計モデリング: Lee-Carter 法とその拡張 – ヒューマンセキュリティへの基 盤研究」 ，総合政策学ワーキングペーパー 71． 34

参考文献 v 森下壮一郎 (2019) 「メディアサービスにおけるユーザの継続 の冪分布に基づくモデル化」 ， 『2019 年度人工知能学会全国 大会（第 33 回） 』 ，2 頁，DOI： 10.11517/pjsai.JSAI2019.0_1J4J301． 35