第9・10回では、個人別の価格感応度を推定し、1回の購買機会での施策を比較した。今回は、その先に何回買ってくれそうかを考える。
非契約型ビジネスでは、購買がないだけでは離脱したと断定できない。購買頻度の個人差と、観測されない離脱を一緒にモデル化する。
BG/NBDでは、購買回数・最終購買時点・観測期間 \((x,t_x,T)\) から、個人の活動状態と将来購買を予測する。
将来の購買を利益と割引率で金額に直したものがCLVである。CLVが高いことだけでは、施策による増分利益が大きいとは言えない。
前回とのつながり・今回のゴール
前回は、個人別の需要を使ってクーポンを配る相手や値引き額を決めた。同じ値引きでも、購買確率の変化と失うマージンが顧客によって違うため、一律に配る場合と個別に配る場合で期待利益が変わった。
ただし、そこで計算したのは1回の購買機会から得られる利益だった。ネットスーパーであれば、今日の注文だけで顧客との関係が終わるとは限らない。毎週買う人もいれば、数か月に1度だけ買う人もいる。今後の注文を含めると、この顧客からどれくらいの利益が見込めるだろうか。
この問いでも、第8・9回の階層モデルが役に立つ。個人の購買率 \(\lambda_i\) に集団の分布を置き、本人の購買履歴でその分布を更新する。第9回では価格感応度 \(\alpha_i=\exp(\theta_i)\) の事後分布から次の選択を予測した。今回は購買率と離脱確率の条件付き分布から、将来の購買回数を予測する。
今回のゴールは、Poisson-Gammaモデルに離脱を加えたBG/NBDを推定し、活動を続けている確率、将来購買数、顧客生涯価値（customer lifetime value、CLV）まで計算することである。ここでのCLVは、現在の事業環境が続くもとでの予測である。第10回で考えた施策配分に使うなら、施策によってその価値がどれくらい変わるかを、さらに調べる必要がある。
非契約型ビジネスの難しさ
具体例から始めよう。あなたはネットスーパーのCRM担当だとする。手元には顧客ごとの購買履歴（いつ、何を、いくら買ったか）がある。ある顧客Aさんは、3ヶ月間何も買っていない。この事実だけから、次のどちらが正しいか判定できるだろうか。
Aさんはこのネットスーパーを使うのをやめた（離脱した）。
Aさんはまだ顧客だが、たまたま3ヶ月間何も注文しなかっただけ（休眠中）。
サブスクリプションサービス（動画配信、音楽配信など）であれば、この問題は起きない。解約ボタンを押した瞬間に「離脱」が確定的に観測される。これを契約型（contractual）ビジネス と呼ぶ。一方、ECサイト、小売、レストランのような非契約型（non-contractual）ビジネス では、顧客は「もうこのサービスを使わない」と誰にも宣言しない。ただ、静かに来なくなるだけである。この「離脱が直接観測されない」という問題こそが、今回の講義の要点である。
「最後の購買から90日経ったら休眠顧客とみなす」というルールを考えよう。同じ90日でも、普段は毎月買っていた顧客と、半年に1回ほど買う顧客では意味が違う。前者では離脱を疑う理由になるが、後者では普段の購買間隔の範囲かもしれない。活動状態を予測するには、最後に買ってからの期間と、それまでの購買頻度を合わせて見る必要がある。P(alive)は、この2つの情報をモデルのもとで組み合わせ、まだ活動している確率を計算する。
ここで実務でよく使われる語彙を導入しておく。RFM分析 と呼ばれるフレームワークで、顧客を次の3つの軸で特徴づける。
Recency（R） ：最後にいつ買ったか。Frequency（F） ：どれくらいの頻度・回数買っているか。Monetary（M） ：1回あたりいくら使うか。
RFMは古くから実務で使われてきた素朴なセグメンテーション手法（例えば「RFMそれぞれを5段階にスコアリングして125通りのセグメントに分ける」といった使い方）だが、今回学ぶモデルを通して、次の性質に気づくことになる。購買回数 \(x\) 、最終購買時点 \(t_x\) 、観測期間 \(T\) だけで、離脱を含めたモデルの尤度が書ける。 全員を同じ長さだけ観測する本例では \(T\) は共通なので、顧客間で変わる十分な要約はRFMのRとFになる。Mは顧客単価の計算には必要だが、購買頻度・離脱の推定そのものには不要になる。これは後ほど具体的に確認する。
モデルの積み上げ：Poisson-GammaからBG/NBDへ
ステップ1：購買頻度のモデル化（Poisson-Gamma = NBD）
まず離脱を無視して、「顧客がどれくらいの頻度で買うか」だけを考える。個人 \(i\) は、時間あたり一定の率 \(\lambda_i\) でランダムに購買するとする。これは購買イベントがPoisson過程に従うということであり、観測期間 \(T\) における購買回数 \(x_i\) は
\[
x_i \mid \lambda_i \sim \text{Poisson}(\lambda_i T)
\]
に従う。ここで問題が起きる。もし全顧客の \(\lambda_i\) が同じ値 \(\lambda\) だと仮定すると、購買回数の分散は平均に等しくなる（Poissonの性質）。しかし実際の購買データでは、分散が平均よりずっと大きい （過分散、overdispersion）ことがほとんどである。理由は単純で、顧客ごとに購買頻度がそもそも違う からだ。ヘビーユーザーとライトユーザーが混在すれば、集計した購買回数の分布は、単一のPoissonよりずっと裾が重くなる。
そこで第9回と同じ発想を使う。個人の \(\lambda_i\) を、集団の分布からのドローとしてモデル化する。
\[
\lambda_i \sim \text{Gamma}(r, \alpha)
\]
（\(r\) は形状パラメータ、\(\alpha\) はrateパラメータ であり、scaleではない。Rのrgamma(shape = r, rate = alpha)に対応する。）
第2〜10回では \(\alpha\) を需要推定における価格係数として使ってきた（\(\alpha > 0\) 、効用は \(-\alpha p\) ）。今回のBG/NBDの \(\alpha\) はGamma分布のrateパラメータであり、scaleでも価格係数でもない。 この回に限り、\(\lambda, r, \alpha, a, b\) はCounting Your Customers系の文献の標準記法をそのまま踏襲する。文脈で混同しないように注意してほしい。
\(\lambda_i\) を積分消去すると（Poisson-Gammaミクスチャ）、\(x_i\) の周辺分布は負の二項分布（Negative Binomial Distribution, NBD）になることが知られている。これが「NBDモデル」の名前の由来である。このモデルでは、個人差があるぶん購買回数の分散は平均より大きくなる。ただし、実データの過分散がすべて固定的な個人差で説明できるとは限らない。季節性なども、別途確かめる必要がある。
第8・9回と同じように、全員を同じ値とする代わりに個人差を分布で表す。ここではその分布に、正の値をとるGamma分布を使う。
購買が少ない顧客の予測はどう決まるか
離脱がないPoisson-Gammaモデルでは、\(x_i\) 回を \(T_i\) 週間で観測した後の分布は
\[
\lambda_i\mid x_i,T_i\sim\mathrm{Gamma}(r+x_i,\mathrm{rate}=\alpha+T_i)
\]
となる。その平均は
\[
E[\lambda_i\mid x_i,T_i]
=\frac{r+x_i}{\alpha+T_i}
=\frac{\alpha}{\alpha+T_i}\frac r\alpha
+\frac{T_i}{\alpha+T_i}\frac{x_i}{T_i}.
\]
集団平均 \(r/\alpha\) と本人の購買率 \(x_i/T_i\) の加重平均である。観測期間が短ければ集団の情報が強く効き、長ければ本人の実績が強く効く。第8回のpartial poolingと同じ構造が見える。BG/NBDでは活動状態も観測されないので、後でこの更新を「観測終了時点で活動している場合」と組み合わせる。
ステップ2：離脱を加える（BG/NBD）
ここまでは、顧客が購買する可能性はずっと残ると仮定していた。次に離脱を加える。BG/NBD（Beta-Geometric/NBD）モデルの離脱の仮定は次の通りである。
各購買の直後に、確率 \(p_i\) でこの顧客は離脱する（二度と戻ってこない）。
これは幾何分布（geometric distribution）の考え方そのものである。「何回目の購買の後に離脱するか」が幾何分布に従う。そして、離脱確率 \(p_i\) もまた個人ごとに異質だと考え、
\[
p_i \sim \text{Beta}(a, b)
\]
とする。\(p_i\) を積分消去すると、周辺分布はBeta-Geometric分布になる。これがBG/NBDの「BG」の部分である。
BG/NBDより前に提案されたモデルとして、Schmittlein, Morrison, and Colombo (1987) のPareto/NBDモデルがある。Pareto/NBDでは、離脱は「各購買の直後」ではなく、連続時間で確率的に発生する と仮定する（顧客の「寿命」が指数分布に従い、その率がGamma分布で異質という設定。これがPareto分布の名前の由来）。この定式化は理論的にはBG/NBDより自然に見えるかもしれないが、尤度の中にガウス超幾何関数が現れ、実装が格段に難しくなる。
Fader, Hardie, and Lee (2005) が提案したBG/NBDの最大の貢献は、「各購買直後に離脱するかどうかを決める」という少し不自然に見える仮定を置くことで、尤度がべき乗・ガンマ関数・ベータ関数で書ける という点にある。実務家からすれば、Excelでも実装できるほど扱いやすくなった。この講義でBG/NBDを主役に据えるのも、まさにこの実装のしやすさゆえである。原論文のCDNOWデータでは、両モデルは近い予測を示している。ただし、どの業種でも同じ精度になるとは限らない。
まとめると、BG/NBDモデルは次の4つの前提から成る条件付き分布の束 である。
顧客が生きている間、購買はPoisson過程に従う。個人の購買率 \(\lambda_i\) は \(\text{Gamma}(r, \alpha)\) に従って異質。
各購買の直後、確率 \(p_i\) で離脱する。個人の離脱確率 \(p_i\) は \(\text{Beta}(a, b)\) に従って異質。
\(\lambda_i\) と \(p_i\) は互いに独立。離脱後は二度と購買しない。
ここでは、初回購買の直後を時点0とし、その時点では全員が活動中だとする。初回購買自体は回数 \(x\) に含めず、その後の反復購買の直後だけに離脱が起きる。この仮定のもとでは、反復購買が0回の顧客はまだ離脱していない。後で \(x=0\) のP(alive)が1になるのは、このモデルの約束によるものであり、現実の顧客が必ず活動中だと確かめられたわけではない。
購買率の低下と離脱をどう区別するか
長い間買わないという1つの履歴だけから、離脱したかどうかを直接確かめることはできない。ここでは、活動中の購買率が一定であること、離脱が購買直後に起きること、個人差の分布、購買率と離脱確率の独立性などの制約を置いて、2つの経路の相対的な重みを計算する。季節によって購買率が下がる顧客を、離脱した顧客と取り違える可能性もある。得られた活動確率は、このような仮定のもとでの予測として読む必要がある。
データ構造：RFMのRとFで尤度が書ける
顧客 \(i\) について、観測期間 \([0, T]\) における購買履歴から、次の3つの数字だけを取り出す。
\(x\) ：観測期間中の反復購買回数（初回購買を除く）。\(t_x\) ：最後の購買が起きた時点（\(0 \le t_x \le T\) ）。\(T\) ：観測期間の長さ（この顧客が何週間観測されているか）。
この \((x, t_x, T)\) という3つの数字だけで、BG/NBDモデルの尤度が書ける 。個々の購買が正確にいつ起きたか（購買間隔の詳細）は不要である。全員で観測期間 \(T\) が共通なら、顧客ごとに必要なのは「何回買ったか」と「最後にいつ買ったか」だけになる。これがRFMのR（recency、\(t_x\) やそこから作る「最終購買からの経過時間」に対応）とF（frequency、\(x\) に対応）が十分な要約になる、という先ほど予告した性質である。
まず、顧客の購買タイムラインを実際に描いて、「見た目だけでは離脱したかどうか分からない」ことを体感しておこう。
<- function (lambda, p_drop, T) {<- 0 <- c ()repeat {<- rexp (1 , rate = lambda)<- t + waitif (t_next > T) break <- t_next<- c (events, t)if (runif (1 ) < p_drop) break # 購買直後に確率p_dropで離脱 <- 20 <- 52 <- 0.8 ; alpha_demo <- 6 ; a_demo <- 0.7 ; b_demo <- 2.5 <- rgamma (n_demo, shape = r_demo, rate = alpha_demo)<- rbeta (n_demo, a_demo, b_demo)<- lapply (1 : n_demo, function (i) {<- simulate_customer_events (lambda_demo[i], p_demo[i], T_demo)# 購買0回の顧客は rep(i, 0) が長さ0になり、0行のdata.frameになる（安全） data.frame (customer = rep (i, length (ev)), event_time = ev)<- do.call (rbind, timeline_list)ggplot (timeline_df, aes (x = event_time, y = customer)) + annotate ("segment" , x = 0 , xend = T_demo, y = 1 : n_demo, yend = 1 : n_demo,color = "gray85" , linewidth = 3 ) + geom_point (color = "#4C72B0" , size = 2.2 ) + labs (x = "週" , y = "顧客ID" , title = "20人の購買タイムライン（観測期間52週）" ) + scale_y_continuous (breaks = 1 : n_demo)
この図を見て、「この人は離脱している」「この人はまだ生きている」を目で判定してみてほしい。かなり難しいはずだ。例えば、観測期間の後半にぽつんと点が1つあって、その後何もない顧客がいる。これは「たまたま最近買って、その後まだ次の購買のタイミングが来ていないだけ」なのか、「その購買を最後に離脱した」のか、見た目だけでは判別できない。この曖昧さを統計的に解消するのが、これから作るモデルの仕事である。
尤度とMLE
BG/NBDの対数尤度
顧客の \((x, t_x, T)\) が観測されたときの尤度を導出する（導出の詳細はFader, Hardie, and Lee (2005) に譲り、ここでは結果と実装に集中する）。推定に必要な、パラメータ \((r, \alpha, a, b)\) に依存する尤度カーネルは
\[
\mathcal K(r,\alpha,a,b \mid x, t_x, T) = \frac{B(a, b+x)}{B(a,b)} \cdot \frac{\Gamma(r+x)\,\alpha^r}{\Gamma(r)} \cdot \frac{1}{(\alpha+T)^{r+x}}
\;+\; \mathbb{1}[x>0] \cdot \frac{B(a+1, b+x-1)}{B(a,b)} \cdot \frac{\Gamma(r+x)\,\alpha^r}{\Gamma(r)} \cdot \frac{1}{(\alpha+t_x)^{r+x}}
\]
である。\(x>0\) のとき、要約統計 \((x,t_x,T)\) の完全な同時密度には、さらにデータだけに依存する共通因子 \(t_x^{x-1}/(x-1)!\) が掛かる。この因子は \((r,\alpha,a,b)\) に依存しないためMLEでは落としてよく、後で計算する \(A_0/(A_0+A_1)\) でも相殺される。ここで \(B(\cdot,\cdot)\) はベータ関数である。カーネルの構造は次のように読める。
第1項 ：「観測期間の終わり \(T\) まで、この顧客はまだ生きている」という経路のカーネルへの寄与。第2項 （\(x>0\) のときだけ存在）：「最後の購買 \(t_x\) の直後にこの顧客は離脱していた」という経路のカーネルへの寄与。
つまり、このカーネルは「まだ生きている」経路と「\(t_x\) で離脱していた」経路の2つの寄与の和 になっている。\(x=0\) （反復購買がない）の顧客には、そもそも「離脱するタイミング」自体が存在しない（購買が1回もないので、購買直後に離脱するという幾何分布の仮定が働く場面がない）ため、第2項は存在せず、第1項だけになる。
上の \(\mathcal K\) は、パラメータ推定とP(alive)の比を計算するには十分である。しかし、\(t_x^{x-1}/(x-1)!\) を省いたまま、異なる \(t_x\) を持つ履歴どうしの「起こりやすさ」を比較してはいけない。そうした密度比較では、省略したデータ依存因子を戻す必要がある。
この尤度の実装で、初学者がハマりやすい罠が2つある。
\(x=0\) の場合分けを忘れる\(x=0\) のときは第2項を計算してはいけない（\(B(a+1, b-1)\) のような、そもそも意味をなさない項を作ってしまう）。if (x > 0)で明示的に分岐すること。対数の中の和でオーバーフローする ：尤度は「2つの項の和」なので、対数尤度を計算するには \(\ln(A_0 + A_1)\) の形になる。\(A_0, A_1\) をそのまま計算してから対数を取ると、パラメータの値次第で \(A_0\) や \(A_1\) が極端に大きい・小さい数値になり、桁あふれ（overflow/underflow）を起こす。第2回でsoftmaxのオーバーフロー対策をやったのと同じ発想で、logsumexpトリック を使う。\(\ln(e^{a} + e^{b}) = m + \ln(e^{a-m} + e^{b-m})\) （\(m = \max(a,b)\) ）という恒等式を使えば、指数の引数が0以下に収まり、数値的に安定する。
Rで実装する前に、この尤度カーネルの実装ロジックを検証しておく（本講義の方針：数値計算のロジックはまずPython/numpyでプロトタイプしてから、Rに書き起こす）。同じ購買回数 \(x\) のもとで \(t_x\) を大きくすると、表示した \(\mathcal K\) は単調に減少する。これは第1項が \(t_x\) に依存せず、第2項が \((\alpha+t_x)^{-(r+x)}\) に比例するというカーネルの代数的性質 である。省略した \(t_x^{x-1}/(x-1)!\) を含む完全な要約統計の密度まで単調減少する、あるいは後の履歴ほど起こりにくい、という意味ではない。一方、P(alive)は2経路の相対比なので共通因子が相殺され、\(t_x\) が大きいほど増加する。カーネルの値と顧客状態の事後確率を区別しよう。
Rシミュ2：DGP→自作対数尤度→optimでの真値回収
N=2,000人の顧客について、52週間のcalibration期間のデータを生成する。真のパラメータは \(r=0.8, \alpha=6, a=0.7, b=2.5\) とする（このパラメータのもとで、週あたりの平均購買率は \(E[\lambda] = r/\alpha \approx 0.133\) 回、平均的な離脱確率は \(E[p] = a/(a+b) \approx 0.22\) 程度になる）。
<- 2000 <- 52 # calibration期間（週） <- 26 # holdout期間（週） <- T_cal + T_holdout<- 0.8 ; alpha_true <- 6 ; a_true <- 0.7 ; b_true <- 2.5 <- rgamma (N, shape = r_true, rate = alpha_true)<- rbeta (N, a_true, b_true)# 各顧客の購買イベントを観測期間全体(calibration+holdout)でシミュレーションする <- function (lambda, p_drop, T) {<- 0 <- c ()repeat {<- rexp (1 , rate = lambda)<- t + waitif (t_next > T) break <- t_next<- c (events, t)if (runif (1 ) < p_drop) break <- vector ("list" , N)for (i in 1 : N) {<- simulate_full_history (lambda_true[i], p_true[i], T_total)# calibration期間(52週)だけを見て (x, t_x) を作る(holdoutの情報は使わない) <- sapply (all_events, function (ev) sum (ev <= T_cal))<- sapply (all_events, function (ev) {<- ev[ev <= T_cal]if (length (ev_cal) == 0 ) 0 else max (ev_cal)cat ("calibration期間の購買回数xの分布: \n " )table (cut (x_cal, breaks = c (- 1 , 0 , 1 , 2 , 3 , 5 , 10 , Inf ),labels = c ("0" , "1" , "2" , "3" , "4-5" , "6-10" , "11+" )))
0 1 2 3 4-5 6-10 11+
324 534 329 188 248 249 128
<- data.frame (x = x_cal, tx = tx_cal)ggplot (scatter_df, aes (x = tx, y = x)) + geom_jitter (alpha = 0.3 , width = 0.5 , height = 0.15 , color = "#4C72B0" ) + labs (x = "t_x（最終購買時点、週）" , y = "x（購買回数）" ,title = "顧客ごとの (t_x, x) 散布図" )
購買回数0の顧客が一定割合いること、購買回数が多い顧客ほど最終購買が観測期間の後半に集中しがちであること（頻繁に買うのでいつ観測を打ち切っても直近に購買がある）が見て取れる。
対数尤度を実装する。
# logsumexp: log(exp(a) + exp(b)) を数値的に安定に計算する <- function (log_a, log_b) {<- pmax (log_a, log_b)+ log (exp (log_a - m) + exp (log_b - m))# BG/NBDの対数尤度（ベクトル化: x, tx, Tは同じ長さのベクトル） <- function (r, alpha, a, b, x, tx, T) {<- lbeta (a, b) # 定数(全顧客共通) # 第1項(常に生きている経路) <- lbeta (a, b + x) - ln_Bab + lgamma (r + x) - lgamma (r) + * log (alpha) - (r + x) * log (alpha + T)# 第2項(t_xで離脱していた経路)。x=0の行では式自体を評価しない <- rep (- Inf , length (x))<- x > 0 <- lbeta (a + 1 , b + x[idx] - 1 ) - ln_Bab + lgamma (r + x[idx]) - lgamma (r) + * log (alpha) - (r + x[idx]) * log (alpha + tx[idx])logsumexp2 (ln_A0, ln_A1)# optimに渡す負の対数尤度。thetaはlog変換した(r,alpha,a,b)を受け取る(正値制約を外すため) <- function (theta_log, x, tx, T) {<- exp (theta_log)<- params[1 ]; alpha <- params[2 ]; a <- params[3 ]; b <- params[4 ]- sum (bgnbd_loglik (r, alpha, a, b, x, tx, T))
\(r, \alpha, a, b\) はいずれも正の値でなければならない。制約付き最適化を避けるため、\(\theta = \log(r, \alpha, a, b)\) というパラメータに変換し、optim()には制約のない\(\theta\) を探索させる。最適化が終わったらexp(theta)で元のスケールに戻す。これは第2回でロジットモデルの推定をしたときには不要だった工夫だが（効用のパラメータには符号の制約がなかった）、正値制約があるパラメータでは頻出のテクニックである。
<- log (c (r = 1 , alpha = 1 , a = 1 , b = 1 ))<- optim (par = theta_init,fn = neg_loglik_bgnbd,x = x_cal, tx = tx_cal, T = rep (T_cal, N),method = "L-BFGS-B" ,hessian = TRUE <- exp (opt_result$ par)names (theta_hat) <- c ("r" , "alpha" , "a" , "b" )
r alpha a b
0.8317126 6.5066501 0.7598365 2.8325799
<- c (r = r_true, alpha = alpha_true, a = a_true, b = b_true)<- data.frame (parameter = names (true_params),true_value = true_params,estimate = round (theta_hat, 4 )kable (recovery_df, caption = "真値とMLE推定値の比較（N=2,000、calibration 52週）" )
真値とMLE推定値の比較（N=2,000、calibration 52週）
r
r
0.8
0.8317
alpha
alpha
6.0
6.5067
a
a
0.7
0.7598
b
b
2.5
2.8326
このseedでは推定値は真値に近い。\(a,b\) の誤差が\(r,\alpha\) よりやや大きいのは、離脱が直接観測されず、最終購買後の沈黙と購買回数の組み合わせから間接的に推定されるためである。ただし、1回のシミュレーションで常にこの順序になるわけではない。推定精度は反復シミュレーションのbias・RMSE・被覆率で評価するのが原則である。
推論：HessianからのSE
# Hessianはlog-scaleのthetaについてのものなので、delta法で元のスケールのSEに変換する <- solve (opt_result$ hessian)<- sqrt (diag (vcov_log))# delta法: g(theta) = exp(theta) のとき se(g(theta)) ~ g(theta) * se(theta) <- theta_hat * se_log<- theta_hat - 1.96 * se_orig<- theta_hat + 1.96 * se_orig<- data.frame (parameter = names (true_params),true_value = true_params,estimate = round (theta_hat, 4 ),se = round (se_orig, 4 ),ci_lower = round (ci_lower, 4 ),ci_upper = round (ci_upper, 4 )kable (inference_df, caption = "推定値・標準誤差・95%信頼区間（真値との比較）" )
推定値・標準誤差・95%信頼区間（真値との比較）
r
r
0.8
0.8317
0.0349
0.7633
0.9001
alpha
alpha
6.0
6.5067
0.3711
5.7792
7.2341
a
a
0.7
0.7598
0.0953
0.5730
0.9467
b
b
2.5
2.8326
0.4710
1.9094
3.7558
推定値のずれを標準誤差と比べてみよう。1回のシミュレーションで真値が区間に入るかどうかだけでは、信頼区間の被覆率は判断できない。Hessianの逆行列から得るのは、正しく特定された独立同分布の尤度を前提とするモデルベース標準誤差 である。実データで顧客間依存やモデルの誤特定を疑う場合には、サンドイッチ分散やクラスタリングなどを別途検討する必要がある。
P(alive)：この顧客はまだ生きているか
推定したパラメータを使えば、「この顧客が現時点でまだ生きている確率」P(alive)を計算できる。これは先ほどの尤度の中に現れた「第1項」と「第2項」の相対的な大きさから直接求まる。
\[
P(\text{alive} \mid x, t_x, T) = \frac{A_0}{A_0 + A_1}
\]
ここで \(A_0, A_1\) は先ほどの尤度式の第1項・第2項（\(B(\cdot,\cdot)\) を含む部分全体）である。\(x=0\) のときは常に\(P(\text{alive})=1\) になる（一度も買っていない顧客について「離脱した」という証拠は存在しないため、モデル上は「まだ最初の購買機会を待っている」という扱いになる）。
Fader, Hardie, and Lee (2005) はこれを整理して、次の閉形式を導いている。
\[
P(\text{alive} \mid x, t_x, T) = \left[1 + \mathbb{1}[x>0]\cdot\frac{a}{b+x-1}\left(\frac{\alpha+T}{\alpha+t_x}\right)^{r+x}\right]^{-1}
\]
<- function (r, alpha, a, b, x, tx, T) {<- (alpha + T) / (alpha + tx)<- ifelse (x > 0 , (a / (b + x - 1 )) * ratio^ (r + x), 0 )1 / (1 + term)# 全顧客にP(alive)を付与する <- palive_bgnbd (theta_hat["r" ], theta_hat["alpha" ], theta_hat["a" ], theta_hat["b" ],cat ("P(alive)の要約統計: \n " )
Min. 1st Qu. Median Mean 3rd Qu. Max.
0.00000196 0.18600454 0.61761399 0.54810663 0.86335947 1.00000000
ggplot (data.frame (palive = palive_all), aes (x = palive)) + geom_histogram (bins = 30 , fill = "#4C72B0" , alpha = 0.8 ) + labs (x = "P(alive)" , y = "顧客数" , title = "生存確率P(alive)の分布" )
生存確率が0付近と1付近に山ができる、両極端に分かれた分布になる。これは活動状態を実際に観測した結果ではなく、モデルに基づく条件付き確率である。極端な確率が出ていても、仮定が妥当かどうかはholdoutの購買予測などで確認する必要がある。
Rシミュ3：(x, t_x)平面でのP(alive)ヒートマップ
x（購買回数）とt_x（最終購買時点）の全ての組み合わせについてP(alive)を計算し、ヒートマップにする。
<- 0 : 20 <- seq (0 , 52 , by = 2 )<- expand.grid (x = x_grid, tx = tx_grid)# t_x > T になる組み合わせ、あるいは x=0なのにtx>0という不可能な組み合わせは除く <- heatmap_df[heatmap_df$ tx <= T_cal & ! (heatmap_df$ x == 0 & heatmap_df$ tx > 0 ), ]$ palive <- palive_bgnbd ("r" ], theta_hat["alpha" ], theta_hat["a" ], theta_hat["b" ],$ x, heatmap_df$ tx, T_calggplot (heatmap_df, aes (x = tx, y = x, fill = palive)) + geom_tile () + scale_fill_gradient2 (low = "#C44E52" , mid = "white" , high = "#4C72B0" ,midpoint = 0.5 , limits = c (0 , 1 )) + labs (x = "t_x（最終購買時点、週）" , y = "x（購買回数）" , fill = "P(alive)" ,title = "P(alive)ヒートマップ：recencyとfrequencyの交互作用（T=52週）" )
このヒートマップから、非自明な交互作用 が読み取れる。単純に「最近買っていれば生きている、買っていなければ離脱している」という単調な関係ではない。同じ最終購買時点（\(t_x\) 、横軸の同じ位置）でも、購買回数（\(x\) 、縦軸）によってP(alive)の判定が大きく変わる。\(x\) を増やすとP(alive)は横ばい〜微増した後に低下することがあり、本例では沈黙が十分長く、購買回数が多い領域で「離脱している」と判定されやすくなる 。これは全範囲での単調関係ではない。
同じ「最後の購買から8週間、observation windowの終わりまで沈黙」という状況でも、購買回数によってP(alive)が劇的に変わる例を見てみよう。
# T=52週、最終購買が44週目(=最後の購買から8週間の沈黙)で固定し、xを変える <- 44 <- data.frame (x = c (1 , 2 , 3 , 5 , 8 , 12 , 20 ))$ palive <- palive_bgnbd ("r" ], theta_hat["alpha" ], theta_hat["a" ], theta_hat["b" ],$ x, tx_fixed, T_calkable (example_df, digits = 4 ,caption = "recency固定(最終購買から8週間沈黙)でも、frequencyでP(alive)が変わる" )
recency固定(最終購買から8週間沈黙)でも、frequencyでP(alive)が変わる
1
0.7401
2
0.7688
3
0.7836
5
0.7923
8
0.7793
12
0.7340
20
0.5732
購買回数が少ないうちはP(alive)が横ばい〜微増するが、ある水準を超えると急激に低下していく。\(x\) が多い履歴は「この顧客の購買率 \(\lambda\) は高い」という情報と、「過去の購買直後には離脱しなかったので \(p\) は低そうだ」という情報を同時に与える。固定した8週間の沈黙に対しては、やがて前者が支配的になる。高い購買率でまだ生存しているなら8週間も購買がない事象は起こりにくいため、直近の購買後に離脱した経路の相対尤度が高くなる のである。逆に、もともと購買頻度が低い顧客にとっては、8週間の沈黙は通常の購買間隔の範囲内かもしれない。したがってP(alive)はrecencyとfrequencyを別々の単調スコアとして足すのではなく、両者をモデル内で組み合わせて判定する。
将来購買の予測とholdout検証
事後分布から将来をシミュレーションする
BG/NBDでは、将来の購買回数の条件付き期待値を公式で計算できる。その公式にはガウス超幾何関数（\({}_2F_1\) ）が現れる。今回は、購買と離脱がどう起きるかを確認しながら、シミュレーションで計算しよう。
「現在は活動中か」「次にいつ買うか」「その購買後に離脱するか」を順に引けば、将来の購買履歴を1本作れる。この計算を繰り返し、購買回数を平均する。第9回と同じく、個人パラメータを1つの値に決めず、その条件付き分布を使って予測することになる。
具体的な手順は次の通りである。ここでのGamma分布は第2引数がrate であり、\(\alpha\) と\(T\) はともに「週」を基準にした量でなければ足せない。
各顧客の \((x, t_x)\) を所与として、\(P(\text{alive})\) を計算する（すでに実装済み）。
calibration終了時点で生存しているという経路のもとでは、Poisson-Gamma共役性から \[\lambda \mid x,T,\text{alive at }T \sim \text{Gamma}(r+x,\ \text{rate}=\alpha+T)\] となる。
同じ生存経路では、観測された\(x\) 回の購買直後にすべて離脱しなかったので、Beta-Bernoulli共役性から \[p \mid x,\text{alive at }T \sim \text{Beta}(a,b+x)\] となる。
P(alive)でcalibration終了時点の活動状態を引く。離脱済みなら将来購買は0。生存なら事後分布から\((\lambda,p)\) を引き、次の待ち時間を\(\text{Exp}(\lambda)\) から引く。
予測期間内に購買が起きたら1回を記録し、その購買直後に確率\(p\) で離脱させる 。生存したときだけ次の待ち時間へ進む。
これを多数回繰り返して平均を取れば、予測期間中の追加離脱まで含む条件付き期待購買数のモンテカルロ近似が得られる。
ここではMLE \(\hat\theta=(\hat r,\hat\alpha,\hat a,\hat b)\) を固定したplug-in予測を行う。将来パスの乱数を増やせばモンテカルロ誤差は小さくなるが、パラメータ推定誤差まで消えるわけではない。予測区間やCLVの不確実性を報告する実務では、顧客単位のbootstrapやパラメータの近似分布・事後分布からも乱数を引き、両方の不確実性を伝播させる。
calibration終了 T
├─ 離脱済み（確率 1 - P(alive)）→ 将来購買 0
└─ 生存（確率 P(alive)）
├─ λ ~ Gamma(r+x, rate=α+T)
├─ p ~ Beta(a, b+x)
└─ 待ち時間 Δ ~ Exp(λ)
├─ Δ が予測期間外 → 終了
└─ Δ が予測期間内 → 購買を1回記録
├─ 確率 p で離脱 → 終了
└─ 確率 1-p で生存 → 次の Δ を引く（繰り返し）重要なのは、「Tで生存」と判定した後も永久に生きるわけではないことだ。BG/NBDでは将来の各購買も新しい離脱機会 になる。
# 生存している1顧客を、予測開始時点0からhorizon_weeks週だけ前向きに走らせる <- function (lambda, p_drop, horizon_weeks) {<- 0 <- numeric (0 )repeat {<- t + rexp (1 , rate = lambda)if (t_next > horizon_weeks) break <- t_next<- c (events, t)if (runif (1 ) < p_drop) break # 将来の各購買直後にも離脱判定 # 1顧客の予測期間合計購買数をS本の将来パスで近似する <- function (r_hat, alpha_hat, a_hat, b_hat,S = 30 ) {<- rgamma (S, shape = r_hat + x, rate = alpha_hat + T_cal)<- rbeta (S, shape1 = a_hat, shape2 = b_hat + x)<- runif (S) < p_alive<- numeric (S)for (s in seq_len (S)) {if (! is_alive[s]) next <- length (simulate_future_path (lambda_post[s], p_post[s], H))mean (n_events)
離脱がなければ、長さ\(H\) の区間の購買数は\(\text{Poisson}(\lambda H)\) なので、rpois(1, lambda * H)でよい。しかしBG/NBDでは、予測期間中の1回目の購買直後に離脱すれば2回目以降は起きない。将来の各購買が次の購買の有無を変えるため、今回は到着時刻を順番に生成して、購買直後の離脱判定を挟む必要がある。
# holdout期間を連続時間（単位=週）で前向きに走らせ、(w-1,w]ごとに集計する <- function (r_hat, alpha_hat, a_hat, b_hat,S = 20 ) {<- rgamma (S, shape = r_hat + x, rate = alpha_hat + T_cal)<- rbeta (S, shape1 = a_hat, shape2 = b_hat + x)<- runif (S) < p_alive<- matrix (0 , nrow = S, ncol = n_weeks)for (s in seq_len (S)) {if (! is_alive[s]) next <- simulate_future_path (lambda_post[s], p_post[s], n_weeks)if (length (event_times) > 0 ) {# ceiling(t)により、時刻tを実績と同じ区間(w-1,w]へ割り当てる <- tabulate (ceiling (event_times), nbins = n_weeks)colMeans (weekly_counts)# 全顧客についてholdoutの週次予測購買数を積み上げる(計算量を抑えるためS=20) <- matrix (0 , nrow = N, ncol = T_holdout)for (i in 1 : N) {<- simulate_holdout_weekly ("r" ], theta_hat["alpha" ], theta_hat["a" ], theta_hat["b" ],S = 20 <- colSums (pred_weekly_matrix)<- rowSums (pred_weekly_matrix) # 個人ごとの予測合計(デシル分析用) # 実測: holdout期間(52-78週)の週次購買数 <- sapply (1 : T_holdout, function (w) {sum (sapply (all_events, function (ev) sum (ev > T_cal + w - 1 & ev <= T_cal + w)))<- sapply (all_events, function (ev) sum (ev > T_cal))<- data.frame (series = c ("実績" , "BG/NBD予測" ),purchases_26_weeks = c (sum (actual_weekly_total), sum (pred_weekly_total))kable (tracking_summary, digits = 1 , caption = "holdout 26週の総購買数（再計算結果）" )
holdout 26週の総購買数（再計算結果）
実績
1940.0
BG/NBD予測
1933.3
(i) 週次累積購買のトラッキングプロット
<- data.frame (week = 1 : T_holdout,actual_cum = cumsum (actual_weekly_total),predicted_cum = cumsum (pred_weekly_total)%>% pivot_longer (cols = c (actual_cum, predicted_cum), names_to = "type" , values_to = "cum_purchases" )ggplot (tracking_df, aes (x = week, y = cum_purchases, color = type, linetype = type)) + geom_line (linewidth = 1 ) + scale_color_manual (values = c ("actual_cum" = "black" , "predicted_cum" = "#C44E52" ),labels = c ("実績" , "予測" )) + scale_linetype_manual (values = c ("actual_cum" = "solid" , "predicted_cum" = "dashed" ),labels = c ("実績" , "予測" )) + labs (x = "holdout期間の週" , y = "累積購買数（全顧客合計）" , color = "" , linetype = "" ,title = "累積購買数のトラッキング：実績 vs BG/NBD予測" )
実績と予測の累積曲線は近い軌跡を描く。BG/NBDのようなモデルで、この種のトラッキングプロット（実務ではCDNOWデータでの検証が有名）を描き、モデルが集計レベルで妥当な予測をしているかを確認するのが標準的な検証手順である。
(ii) 個人レベル：予測デシル別の実績平均
<- data.frame (pred = pred_indiv_total, actual = actual_indiv_total) %>% mutate (decile = ntile (pred, 10 )) %>% group_by (decile) %>% summarise (pred_mean = mean (pred), actual_mean = mean (actual), n = n ())<- decile_df %>% select (decile, pred_mean, actual_mean) %>% pivot_longer (cols = c (pred_mean, actual_mean), names_to = "type" , values_to = "value" )ggplot (decile_long, aes (x = decile, y = value, fill = type)) + geom_col (position = position_dodge (width = 0.7 ), width = 0.6 ) + scale_fill_manual (values = c ("pred_mean" = "#4C72B0" , "actual_mean" = "#DD8452" ),labels = c ("実績平均" , "予測平均" )) + scale_x_continuous (breaks = 1 : 10 ) + labs (x = "予測購買数デシル（1=最も低い予測、10=最も高い予測）" ,y = "holdout期間の購買数（平均）" , fill = "" ,title = "予測デシル別の実績平均：校正確認" )
予測が高いデシルほど実績もおおむね高くなるかを確認する。有限標本やモンテカルロ誤差で局所的な逆転は起こりうるので、完全な単調性を合格条件にはしない。これは順位付けだけでなく、各デシルで予測平均と実績平均の水準が合っているかという校正の確認でもある。
CLV（顧客生涯価値）へ
購買回数の期待値が予測できれば、そこに1回あたりの粗利（マージン）と割引率を掛け合わせて、顧客価値を計算できる。以下は無限期間の厳密な「生涯」ではなく、2年間で打ち切った予測ホライズンCLV である。
\[
\text{CLV}_i = \text{margin} \times \sum_{w=1}^{\text{horizon}} \frac{E[\text{購買数（$w$週目）} \mid \text{パラメータ}, x_i, t_{x,i}]}{(1+d)^w}
\]
\(d\) は割引率（週次に換算したもの）である。将来のキャッシュフローほど現在価値は小さくなる、というファイナンスの基本的な考え方をそのまま輸入している。
<- 1500 # 1回の購買あたりの粗利(円) <- 0.10 # 実効年率 <- (1 + annual_discount)^ (1 / 52 ) - 1 # 実効週率へ変換 <- 104 # 2年先まで <- function (r_hat, alpha_hat, a_hat, b_hat, x, T_cal, p_alive,S = 200 ) {<- rgamma (S, shape = r_hat + x, rate = alpha_hat + T_cal)<- rbeta (S, shape1 = a_hat, shape2 = b_hat + x)<- runif (S) < p_alive<- numeric (S)for (s in seq_len (S)) {if (! is_alive[s]) next <- simulate_future_path (lambda_post[s], p_post[s], horizon_weeks)if (length (event_times) > 0 ) {<- ceiling (event_times)<- sum (margin / (1 + weekly_d)^ event_weeks)mean (path_value)# 上位顧客の例(頻繁に買い、最近も購買がある)と平均的顧客の例を比較する <- compute_clv (theta_hat["r" ], theta_hat["alpha" ], theta_hat["a" ], theta_hat["b" ],x = 15 , T_cal = T_cal,p_alive = palive_bgnbd (theta_hat["r" ], theta_hat["alpha" ],"a" ], theta_hat["b" ], 15 , 50 , T_cal),horizon_weeks = horizon_weeks, weekly_d = weekly_discount, margin = margin)<- compute_clv (theta_hat["r" ], theta_hat["alpha" ], theta_hat["a" ], theta_hat["b" ],x = 3 , T_cal = T_cal,p_alive = palive_bgnbd (theta_hat["r" ], theta_hat["alpha" ],"a" ], theta_hat["b" ], 3 , 30 , T_cal),horizon_weeks = horizon_weeks, weekly_d = weekly_discount, margin = margin)cat ("上位顧客例(x=15,t_x=50)のCLV:" , round (clv_top), "円 \n " )
上位顧客例(x=15,t_x=50)のCLV: 24956 円
cat ("平均的顧客例(x=3,t_x=30)のCLV:" , round (clv_avg), "円 \n " )
平均的顧客例(x=3,t_x=30)のCLV: 3583 円
cat ("CLV比(上位/平均):" , round (clv_top / clv_avg, 2 ), "倍 \n " )
履歴が違うと、同じマージンを仮定してもCLVは変わる。この例の金額は、今後104週間の利益の現在価値であり、無限の将来を足した値ではない。また、高いCLVだけを理由にサポートやクーポンの予算を増やすことはできない。
ここまでで計算してきたCLVは、あくまで「この顧客が将来どれだけ買うと予測されるか 」である。これは重要だが、次の推論はまったく別の主張であり、CLVの数字からは何の保証も得られない 。
「CLVが高い顧客に広告・クーポン・特別オファーを出せば、売上が増える。」
CLVが高い顧客は、放っておいても買う顧客かもしれない。その顧客に広告を打っても、広告のおかげで追加的に生まれる売上（incrementality）はゼロに近いかもしれない——広告費が完全に無駄になる。逆に、CLVが中程度で「放っておいたら離脱しそうな」顧客の方が、施策によって行動が大きく変わる（つまり因果効果が大きい）可能性もある。
CLVは「この顧客が今後どれだけの価値を生むかという予測」であり、「この顧客に施策を打ったらどれだけ価値が変わるかという因果効果」ではない。 顧客への資源配分を最適化するには、CLVの予測だけでなく、施策の因果効果（第10回で扱った増分利益やCATE）を別途推定する必要がある。
Venkatesan and Kumar (2004) はCLVに基づく資源配分のフレームワークを提案した代表的な研究である。Ascarza (2018) は2つのフィールド実験と機械学習を組み合わせ、離脱リスクが最も高い顧客が、引き止め施策の最良の対象とは限らない ことを示した。その研究設定では、離脱リスクではなく施策への反応の異質性に基づくターゲティングの方が有効だった。これは「高リスク顧客への施策は常に無効」という一般命題ではなく、予測リスクをそのまま因果的な施策配分ルールにしてはいけない、という警告である。
実証研究コーナー
この論文は、非契約型ビジネスの購買と離脱を、推定しやすいモデルで表すことを目指した。実証ではオンライン音楽小売店CDNOWの2,357人 を使い、78週間のデータを前半39週間の推定期間と後半39週間の予測期間に分けている。初回購買の時点が顧客ごとに違うため、推定期間内の観測長 \(T_i\) も異なる。購買直後に離脱するBG/NBDを提案し、このデータではPareto/NBDに近い予測が得られることを示した。これは、すべての業種で両モデルの精度が同じになるという結果ではない。今回も、推定に使わなかった期間で予測を確かめる点を引き継いでいる。著者公開論文
AbeはPareto/NBDを階層ベイズに拡張し、購買率と離脱率の相関や顧客属性をモデルに取り入れた。実証にはCDNOW、日本の百貨店、音楽CDチェーンの3種類 のデータを使っている。百貨店の例では400人を52週間 観測し、前後26週間を推定と予測検証に分けている。MCMCによって個人パラメータの分布を得るため、将来購買や活動確率もその分布に沿って計算できる。個人の履歴が短い場合には、不確実性が大きいことも論文で示されている。第9回の個人別需要と同じく、個人パラメータを1つの確定した数値として扱わないという接点がある。著者公開稿、§4
この論文は、保険会社のデータで複数のCLV予測モデルを比較した。最初の1年でモデルを推定し、続く4年間 の利益を予測する。割引率は年率10% とし、予測時点で既に観測されている利益を評価対象から外している。現状が続くとする単純なモデル、Tobitモデル、商品別の選択モデルや期間モデルを比べると、複雑化が大きな精度改善につながるとは限らなかった。これは契約が観測できる保険の分析であり、BG/NBDそのものの性能比較ではない。今回のモデルでも、理論の細かさだけで判断せず、簡単な予測基準とholdoutで比較することが大切になる。論文、§4.1
ビジネスの現場で
多くの企業では「最終購買から90日経過で休眠」のような一律のルールで顧客をセグメント分けしている。しかし今回学んだP(alive)を使えば、顧客ごとの購買頻度の違いを反映した、統計的に一貫性のある活動状態の推定 ができる。月1回買う客の90日沈黙と、半年に1回しか買わない客の90日沈黙では、含意がまったく違う。一律ルールをP(alive)に置き換えることは、実装コストはさほど高くない一方で、将来の購買を予測する精度が改善する可能性がある。ただし、これはあくまで「予測」の改善であり、離脱防止施策自体が効果を持つかどうかは別途検証が必要である点は、CLVのwarningと同様に留意すべきである。
ある広告によって追加的に獲得できる顧客の、将来利益の現在価値が平均5,000円だとする。広告費をその追加顧客数で割った金額が5,000円を上回るなら、この利益だけでは広告費を回収できない。ここで必要なのは、既存顧客全体の平均CLVではなく、その広告で増える顧客の価値である。広告を見た人が全員、広告のおかげで加入したわけでもない。第10回のランダム化実験などで追加顧客数を調べ、対応する顧客集団のCLVと合わせて費用を比べる。
今回のBG/NBDは非契約型ビジネス（EC、小売）を想定したモデルだった。もしあなたの分析対象がサブスクリプション（契約型）ビジネスであれば、離脱は直接観測できるので、より標準的な生存分析（survival analysis）の手法がそのまま使える。さらに、顧客が「今月解約するかどうか」を将来を見据えて（forward-lookingに）判断しているとモデル化したい場合は、第13回で学ぶ動的離散選択モデルの枠組みが必要になる。今回学んだ「観測されない状態（離脱）を統計的に推定する」という発想は、契約型・非契約型を問わず、顧客行動のモデリング全般に通底する考え方である。
非契約型ビジネスでは離脱が直接観測されない。\((x,t_x,T)\) という最小限の情報から活動状態・将来購買を推定でき、\(T\) が共通なら顧客間で変わる要約はRFMのR・Fだけになる。
BG/NBDは、購買頻度をGamma異質性のPoisson過程（NBD）、離脱を各購買直後のBeta異質性の幾何分布としてモデル化する。これは第9回の階層モデル（個人パラメータに分布を与えるshrinkageの発想）の応用である。
BG/NBDの対数尤度はガンマ関数・ベータ関数だけで書け、optim()で自作MLEにより真値を回収できる。x=0の場合分けとlogsumexpによる数値安定化が実装の勘所である。
P(alive)は閉形式で計算でき、購買頻度と最終購買時点の非自明な交互作用を生む。本例ではP(alive)は購買回数とともに横ばい〜微増した後、沈黙が十分長い高頻度領域で低下する。
将来購買数を予測するには、まずP(alive)で観測終了時点の活動状態を引く。活動中なら、\(\lambda\sim\mathrm{Gamma}(r+x,\text{rate}=\alpha+T)\) 、\(p\sim\mathrm{Beta}(a,b+x)\) を引き、将来の各購買直後にも離脱判定を行う。離脱済みなら将来購買は0とする。予測はholdout期間の実績と比べて確認する。
CLVはマージン・割引率・予測購買数から計算できるが、CLVは予測であり因果効果ではない 。「CLVが高い顧客に施策を打つ」ことと「施策が効く顧客に施策を打つ」ことはまったく別の問題であり、後者には第10回のように施策を実施した場合としない場合の比較が必要になる。
宿題
ブラウザ実験室（playground11.html ）で、今回のシミュレーションの数値を自分でいじって遊んでみよう。インストール不要、ブラウザだけでRが動く（初回ロードのみ30秒〜1分かかる）。
ChatGPTなどの生成AIに「BG/NBDモデルとPareto/NBDモデルの違いは何か、なぜBG/NBDの方が実務で使われやすいのか」と聞いてみて、今日の講義の説明と比べてみよう。
今日のpalive_bgnbd()関数のヒートマップで、\(a, b\) （離脱確率のBeta分布パラメータ）を変えると、ヒートマップの形がどう変わるか試してみよう。\(a/(a+b)\) （平均離脱確率）を大きくすると、どのあたりの色が変わるだろうか。
coding課題は assignment11.qmd を参照。ネットスーパーの顧客ベースを自分でシミュレーションし、真値回収・P(alive)・holdout検証・CLV試算までの一連の流れを自分の手で実装してもらう。
この先の講義とのつながり
第12回のゲスト講義を挟み、第13回から消費者の動学的な選択を扱う。今回も将来の購買を予測したが、顧客自身が将来を考えて選んだわけではない。購買率や離脱確率に従ってイベントが発生するモデルだった。
例えば、今日は買わずに来週の値下げを待つ、定期購入を始めて将来の買い物を変える、といった行動では、今日の選択と将来の状態を一緒に考える必要がある。第13回では、分析者が将来を予測することに加えて、消費者も将来を見て意思決定するモデルへ進む。
参考文献
Fader, Peter S. and Bruce G. S. Hardie (2009) “Probability Models for Customer-Base Analysis,” Journal of Interactive Marketing , 23(1), 61-69.
Schmittlein, David C., Donald G. Morrison, and Richard Colombo (1987) “Counting Your Customers: Who Are They and What Will They Do Next?” Management Science , 33(1), 1-24.
Fader, Peter S., Bruce G. S. Hardie, and Ka Lok Lee (2005) “Counting Your Customers the Easy Way: An Alternative to the Pareto/NBD Model,” Marketing Science , 24(2), 275-284.
Abe, Makoto (2009) “Counting Your Customers One by One: A Hierarchical Bayes Extension to the Pareto/NBD Model,” Marketing Science , 28(3), 541-553.
Donkers, Bas, Peter C. Verhoef, and Martijn G. de Jong (2007) “Modeling CLV: A Test of Competing Models in the Insurance Industry,” Quantitative Marketing and Economics , 5(2), 163-190.
Venkatesan, Rajkumar and V. Kumar (2004) “A Customer Lifetime Value Framework for Customer Selection and Resource Allocation Strategy,” Journal of Marketing , 68(4), 106-125.
Ascarza, Eva (2018) “Retention Futility: Targeting High-Risk Customers Might Be Ineffective,” Journal of Marketing Research , 55(1), 80-98.
