確率的購買モデル（BG/NBD）｜マーケティングサイエンス

🎓 レベル：発展　|　重要度：A（必須）

📎 関連：マーケティング予算最適化　|　前提：RFM分析（リセンシー・フリークエンシー・マネタリー）

要点（BLUF）

RFM分析（リセンシー・フリークエンシー・マネタリー）が顧客の過去を記述するだけなのに対し、確率的購買モデルは取引履歴から将来の購買回数と生存確率 $P(\text{alive})$ を予測します。代表が BG/NBDで、**NBD（ガンマ・ポアソン混合）が「購買頻度の異質性」**を、**ベータ幾何（BG）が「離脱（もう買わなくなる）」**を表し、2つを合わせたものです。
個人の購買が Poisson( $\lambda$ )、その強度 $\lambda$ が顧客ごとに Gamma( $r,\theta$ ) で散らばると、集団の取引数は NBD になり過分散します（分散 > 平均、 $\mathrm{Var}/E = 1+\theta > 1$ ）。顧客5000人（ $r=0.8,\theta=2.5$ ）の合成データで標本平均 2.04・標本分散 7.06・分散/平均 3.47——ポアソンの $1$ をはるかに超えます。モーメント法で $\theta=2.47,\ r=0.83$ とパラメータも回復できます。
過去 $x$ 回買った顧客の将来期待は、ガンマ・ポアソンのベイズ更新（事後 Gamma( $r+x,\dots$ )）で出ます。これは観測 $x$ をそのまま将来へ伸ばす素朴外挿より集団平均へ縮約され（平均回帰）、多買客の予測は少し下がります。顧客生涯価値（LTV）の「継続率一定」を個人差で精緻化し、将来取引 × 粗利で LTV を予測する土台になります。⚠️ 離脱まで含む完全な BG/NBD・Pareto/NBD の推定はライブラリ（lifetimes・PyMC-Marketing 等）で要最新確認。

1. RFMは記述・BG/NBDは予測：何を予測するのか

RFM分析（リセンシー・フリークエンシー・マネタリー）は、顧客を R・F・M の3軸で過去のとおりに並べる手法でした。 $555$ （最近・高頻度・高額）の優良客や、F・M は高いのに R が悪い離反予兆——どれも「これまでどう買ったか」の記述です。「では、この顧客は次の半年で何回買うか」「いままだ生きている（離脱していない）確率はどれくらいか」には、RFM 単独では答えられません。ここを確率モデルで予測するのが、本ノートの確率的購買モデルです。

発想はこうです。顧客の購買を「コインを振るような確率過程」と見なし、その確率の個人差を分布で表す。すると、ある顧客の過去の買い方（何回買ったか＝F、最後に買ってからどれだけ経つか＝R）を証拠として、ベイズ的に「将来の期待購買回数」と「生存確率」を計算できます。代表が BG/NBD（Beta-Geometric / Negative-Binomial-Distribution）で、2つの部品からなります。

NBD（購買頻度の異質性）：生きている間、各顧客は自分の強度 $\lambda$ でポアソン的に買う。 $\lambda$ は顧客ごとに違い、Gamma 分布で散らばる。集団で見ると**負の二項分布（NBD）**になる。
ベータ幾何（離脱）：各取引のあと、顧客は確率 $p$ で「離脱（もう二度と買わない）」する。 $p$ も顧客ごとに違い、Beta 分布で散らばる。

flowchart LR
  LOG["取引履歴（誰が・いつ・何回買ったか）"] --> NBD["NBD：購買頻度の異質性（ガンマ・ポアソン混合）"]
  LOG --> BG["ベータ幾何：離脱（各取引後に確率pで離脱）"]
  NBD --> PRED["将来の期待取引数"]
  BG --> ALIVE["生存確率 P(alive)"]
  PRED --> CLV["将来取引 × 粗利 = LTV予測"]
  ALIVE --> CLV

本ノートは、この核のうち**「NBD（過分散）＋条件付き予測（縮約）」を最後まで実装**します。これだけで「個人差を確率分布で表す」という確率的購買モデルの心臓部が体感できるからです。離脱（ベータ幾何）を加えた完全な BG/NBD と、連続時間版の Pareto/NBD は概念を §4 で説明し、推定はライブラリへ誘導します（要最新確認）。

この章は顧客生涯価値（LTV）とリテンションとチャーンの続きです。02-01 は「継続率一定」、02-03 は「顧客の異質性で継続率は一定でない」と進みました。BG/NBD はその異質性を確率分布で明示的に書き下し、個々の顧客の将来まで予測する到達点です。

2. NBD：ガンマ・ポアソン混合と過分散（数式）

まず生きている間の購買頻度だけをモデル化します（離脱は §4）。1人の顧客 $i$ は、期間あたり平均 $\lambda_i$ 回というポアソン過程で買うとします。

X_i \mid \lambda_i \sim \mathrm{Poisson}(\lambda_i)

ポアソン単体なら分散＝平均ですが、現実の顧客は頻度がまちまちです。そこで強度 $\lambda_i$ 自体をガンマ分布で散らします（ $\lambda>0$ の連続量を表すのに自然で、ポアソンと共役）。

\lambda_i \sim \mathrm{Gamma}(\text{shape}=r,\ \text{scale}=\theta), \qquad E[\lambda]=r\theta,\ \ \mathrm{Var}[\lambda]=r\theta^{2}

この $\lambda$ を積分消去（混合）すると、集団の取引数 $X$ は**負の二項分布（NBD）**になります。

P(X=x) = \frac{\Gamma(r+x)}{\Gamma(r)\,x!}\left(\frac{1}{1+\theta}\right)^{r}\left(\frac{\theta}{1+\theta}\right)^{x}

平均と分散は、全分散の法則（条件付き期待・分散の分解）で出ます。

E[X] = E[\lambda] = r\theta

\mathrm{Var}[X] = \underbrace{E[\mathrm{Var}(X\mid\lambda)]}_{E[\lambda]} + \underbrace{\mathrm{Var}[E(X\mid\lambda)]}_{\mathrm{Var}[\lambda]} = r\theta + r\theta^{2} = r\theta(1+\theta)

要点は分散/平均比です。

\frac{\mathrm{Var}[X]}{E[X]} = \frac{r\theta(1+\theta)}{r\theta} = 1+\theta > 1

ポアソンなら $1$ 。NBD は $1+\theta$ で、 $1$ を超える分だけ「過分散（overdispersion）」——頻度の個人差が、集団のばらつきを平均以上に膨らませます。逆に言えば、観測した過分散の大きさから $\theta$ が読めます。これがモーメント法です。標本平均 $m$ 、標本分散 $v$ から、

\hat\theta = \frac{v}{m} - 1, \qquad \hat r = \frac{m}{\hat\theta}

条件付き予測（ベイズ更新）。ここが確率モデルの真価です。ある顧客を1期観測して $x$ 回買ったとわかると、その人の $\lambda$ について事後分布が手に入ります。ガンマ事前 × ポアソン尤度は共役で、事後もガンマです。

\lambda \mid x \ \sim\ \mathrm{Gamma}\!\left(\text{shape}=r+x,\ \text{scale}=\frac{\theta}{1+\theta}\right)

すると次の1期の期待購買数は、この事後平均です。

E[\lambda \mid x] = (r+x)\,\frac{\theta}{1+\theta} = \underbrace{\frac{\theta}{1+\theta}}_{w}\,x + \underbrace{\frac{1}{1+\theta}}_{1-w}\,\underbrace{r\theta}_{\text{集団平均}}

右辺の形が肝心です。予測は、観測 $x$ （重み $w=\theta/(1+\theta)$ ）と集団平均 $r\theta$ （重み $1-w$ ）の加重平均——つまり観測値を集団平均の方へ縮約（shrinkage）します。データが少ない（ $\theta$ 小）ほど集団平均を信じ、データが効く（ $\theta$ 大）ほど観測を信じる。これが信用理論（credibility）・平均回帰の数理で、素朴に「 $x$ 回買ったから将来も $x$ 回」と外挿するのとは違います。

3. 過分散・パラメータ回復・条件付き予測（コード）

§2の3点を順に確かめます。(a) ガンマ・ポアソンで生成した取引数が過分散（分散/平均 > 1）になること、(b) モーメント法で $r,\theta$ を回復できること、(c) 過去 $x$ 回の条件付き将来期待が集団平均へ縮約し、素朴外挿と食い違うこと。数式との対応： $E[X]=r\theta$ 、 $\mathrm{Var}/E=1+\theta$ 、 $\hat\theta=v/m-1$ 、 $\hat r=m/\hat\theta$ 、 $E[\lambda\mid x]=(r+x)\theta/(1+\theta)$ 。

import numpy as np
import pandas as pd
from scipy import stats

# ガンマ・ポアソン混合（NBD）：顧客5000人
rng = np.random.default_rng(0)
n = 5000
r_true, theta_true = 0.8, 2.5

# 個人の購買強度 λ_i ~ Gamma(shape=r, scale=θ)、観測取引数 X_i ~ Poisson(λ_i)
lam = rng.gamma(shape=r_true, scale=theta_true, size=n)
X = rng.poisson(lam)

# (a) NBD の過分散（分散 > 平均）
m = X.mean()
v = X.var(ddof=1)
print("=== (a) NBD の過分散（分散 > 平均）===")
print(f"標本平均 m   = {m:.4f}   （理論 rθ = {r_true*theta_true:.4f}）")
print(f"標本分散 v   = {v:.4f}   （理論 rθ(1+θ) = {r_true*theta_true*(1+theta_true):.4f}）")
print(f"分散/平均    = {v/m:.4f}   （理論 1+θ = {1+theta_true:.4f}、ポアソンなら 1）")

# (b) モーメント法でパラメータを回復：θ = v/m - 1、r = m/θ
theta_hat = v/m - 1
r_hat = m/theta_hat
print("\n=== (b) モーメント法でパラメータ回復 ===")
print(f"θ（推定）= v/m - 1 = {theta_hat:.4f}   （真値 {theta_true}）")
print(f"r（推定）= m/θ     = {r_hat:.4f}   （真値 {r_true}）")

# (c) 条件付き期待将来取引：ガンマ・ポアソンのベイズ更新
# 過去1期で x 回観測 → 事後 λ|x ~ Gamma(shape=r+x, scale=θ/(1+θ))
# 次の1期の期待購買数 = 事後平均 = (r+x)·θ/(1+θ)。ここでは真の r,θ を使う
w = theta_true/(1+theta_true)            # データ側の重み（縮約の信頼度）
prior_mean = r_true*theta_true           # 事前平均 = 集団平均 = 2.0
xs = np.array([0, 1, 2, 3, 5, 10])
naive = xs.astype(float)                 # 素朴外挿：観測 x をそのまま将来へ
bayes = (r_true + xs)*theta_true/(1+theta_true)   # 事後平均（縮約後）
tbl = pd.DataFrame({
    "過去の取引数x": xs,
    "素朴外挿": naive,
    "ベイズ予測": bayes,
    "縮約幅": naive - bayes,
})
print("\n=== (c) 条件付き期待将来取引（次の1期）：素朴外挿 vs ベイズ縮約 ===")
print(f"データ重み w = θ/(1+θ) = {w:.4f}、 集団平均（縮約先）= {prior_mean:.2f}")
print(tbl.to_string(index=False, formatters={
    "素朴外挿": "{:.2f}".format,
    "ベイズ予測": "{:.4f}".format,
    "縮約幅": "{:+.4f}".format}))

出力：

=== (a) NBD の過分散（分散 > 平均）===
標本平均 m   = 2.0374   （理論 rθ = 2.0000）
標本分散 v   = 7.0606   （理論 rθ(1+θ) = 7.0000）
分散/平均    = 3.4655   （理論 1+θ = 3.5000、ポアソンなら 1）

=== (b) モーメント法でパラメータ回復 ===
θ（推定）= v/m - 1 = 2.4655   （真値 2.5）
r（推定）= m/θ     = 0.8264   （真値 0.8）

=== (c) 条件付き期待将来取引（次の1期）：素朴外挿 vs ベイズ縮約 ===
データ重み w = θ/(1+θ) = 0.7143、 集団平均（縮約先）= 2.00
 過去の取引数x  素朴外挿  ベイズ予測     縮約幅
       0  0.00 0.5714 -0.5714
       1  1.00 1.2857 -0.2857
       2  2.00 2.0000 +0.0000
       3  3.00 2.7143 +0.2857
       5  5.00 4.1429 +0.8571
      10 10.00 7.7143 +2.2857

出力の意味：

(a) 過分散。生成した取引数は、標本平均 $2.04$ （理論 $r\theta=2.0$ ）に対し標本分散は $7.06$ （理論 $r\theta(1+\theta)=7.0$ ）。分散が平均の3倍以上で、分散/平均は $3.47$ ——理論値 $1+\theta=3.5$ にぴたりで、ポアソンの $1$ をはるかに超えます。頻度の個人差（ $\lambda$ のばらつき）が、集団のばらつきを平均以上に膨らませた証拠です。もしここでポアソン（分散＝平均と仮定）で将来を見積もれば、ばらつきを大幅に過小評価します。

(b) パラメータ回復。その過分散の大きさ $v/m=3.47$ から、 $\hat\theta=v/m-1=2.47$ （真値 $2.5$ ）、 $\hat r=m/\hat\theta=0.83$ （真値 $0.8$ ）と、生成に使った $r,\theta$ をほぼ言い当てました。観測された平均と分散だけから、背後の異質性分布のパラメータが復元できる——これがモーメント法です（最尤推定はより効率的ですが、過分散の直感はモーメント法が明快）。

(c) 条件付き予測の縮約。表が確率モデルの真価です。素朴外挿は「過去 $x$ 回 → 将来も $x$ 回」と $45$ 度に伸ばすだけ。対してベイズ予測 $(r+x)\theta/(1+\theta)$ は、観測 $x$ （重み $w=0.71$ ）と集団平均 $2.0$ （重み $0.29$ ）の加重平均で、集団平均 $2.0$ の方へ縮約します。読み方は2方向です。過去0回の客は素朴なら将来0回ですが、ベイズは $0.57$ ——「たまたま買わなかっただけで強度がゼロとは限らない」と少し上に見ます。逆に過去10回の多買客は、素朴なら将来10回ですが、ベイズは $7.71$ ——「その多さは一部まぐれ」と** $2.29$ も下に**引きます（縮約幅＋）。境目は集団平均と一致する $x=2$ で、ここだけ縮約ゼロ。高頻度客ほど将来も多い（単調）が、素朴外挿よりは控えめ——これが平均回帰で、確率的購買モデルが素朴な外挿に勝る核心です。

分布の過分散と、条件付き予測の縮約を1枚にします。

import numpy as np
import pandas as pd
from scipy import stats
import matplotlib.pyplot as plt
import japanize_matplotlib

rng = np.random.default_rng(0)
n = 5000
r_true, theta_true = 0.8, 2.5
lam = rng.gamma(shape=r_true, scale=theta_true, size=n)
X = rng.poisson(lam)

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(11, 4.3))

# 左：取引数分布 NBD（観測）vs ポアソン（同じ平均2.0）
ks = np.arange(0, 13)
emp = np.array([(X == k).mean() for k in ks])
nbd_pmf = stats.nbinom(n=r_true, p=1/(1+theta_true)).pmf(ks)   # ガンマ・ポアソン＝NBD
pois_pmf = stats.poisson(mu=r_true*theta_true).pmf(ks)          # 同じ平均のポアソン
ax1.bar(ks, emp, width=0.5, color="C0", alpha=0.6, label="観測度数（5000人）")
ax1.plot(ks, nbd_pmf, "o-", color="C3", lw=1.8, label="NBD（ガンマ・ポアソン混合）")
ax1.plot(ks, pois_pmf, "s--", color="gray", lw=1.6, label="ポアソン（同じ平均2.0）")
ax1.set_xlabel("取引数 x")
ax1.set_ylabel("割合")
ax1.set_title("取引数分布：NBD は過分散（厚いゼロと裾）")
ax1.legend()
ax1.grid(alpha=0.3)

# 右：過去x → 予測将来（素朴外挿 vs ベイズ縮約）
xs = np.arange(0, 11)
naive = xs.astype(float)
bayes = (r_true + xs)*theta_true/(1+theta_true)
ax2.plot(xs, naive, "o--", color="gray", lw=1.6, label="素朴外挿（=x、45度線）")
ax2.plot(xs, bayes, "o-", color="C2", lw=1.8, label="ベイズ予測（縮約）")
ax2.axhline(r_true*theta_true, ls=":", color="C3",
            label=f"集団平均 {r_true*theta_true:.1f}（縮約先）")
ax2.set_xlabel("過去の取引数 x（過去1期）")
ax2.set_ylabel("次の1期の期待取引数")
ax2.set_title("条件付き予測は集団平均へ縮約（平均回帰）")
ax2.legend()
ax2.grid(alpha=0.3)

fig.tight_layout()
plt.show()

左図は取引数の分布です。観測度数（青）に NBD（赤）がぴたりと重なる一方、同じ平均 $2.0$ のポアソン（灰）は、ゼロの山を大幅に低く見積もり（ $0.14$ 対実際 $0.37$ ）、裾（ $x\ge5$ ）も薄く読みます。「買わない人が多く、たまに大量に買う人もいる」という現実の偏りは、ポアソンでは表せず NBD でこそ捉えられる——過分散の正体が一目で分かります。右図は条件付き予測で、素朴外挿（灰の $45$ 度線）に対しベイズ予測（緑）が集団平均 $2.0$ （赤点線）へ寝ていきます。両者が交わるのは $x=2$ （集団平均）。多買客ほど両線の差が開く＝縮約が強いさまが見て取れます。

4. 離脱を入れる：BG/NBD と Pareto/NBD（要最新確認）

§3 までは「顧客はずっと生きている」前提で、頻度だけを扱いました。しかし実際には、顧客はいつか離脱します（解約・乗り換え・自然消滅）。RFM の R（最終購入からの経過日数）が長い顧客は、「強度が低いだけで生きている」のか「もう離脱して死んでいる」のか——これを切り分けるのが、離脱モデルです。

BG/NBD は、§2 の NBD（頻度）に、ベータ幾何（BG）の離脱を足したものです。直感は単純で、各取引のあと、顧客は確率 $p$ で離脱する（次はもう買わない）。この $p$ も顧客ごとに違い、Beta 分布で散らばるとします。

p_i \sim \mathrm{Beta}(a,b)

すると、ある顧客の 生存確率 $P(\text{alive})$ は、その人の頻度 F と直近性 R から計算できます。よく買っていた（F 大）のに最近ぱったり来ない（R 大）顧客は「離脱した」と判断され $P(\text{alive})$ が下がり、もともと低頻度の顧客が少し間隔を空けただけなら「まだ生きている」と判断される。§3 の縮約予測に、この $P(\text{alive})$ を掛け合わせて将来取引を割り引くのが、完全な BG/NBD の予測です。連続時間で離脱を指数分布、その率を Gamma で散らした版が Pareto/NBD（BG/NBD の元祖）です。

そして、こうして得た個々の顧客の期待将来取引数に粗利 $m$ を掛け、割引いて足せば、顧客生涯価値（LTV）の「継続率一定」を脱した、個人ごとの LTV 予測になります。RFM（記述）→ BG/NBD（将来取引と生存）→ LTV（価値）という流れが、ここでつながります。

⚠️ 要最新確認：完全な BG/NBD・Pareto/NBD の推定（尤度の最尤化、 $P(\text{alive})$ や条件付き期待取引の閉形式）は実装が込み入るため、本ノートでは数理の核（過分散と条件付き縮約）の実証に留めました。実務では専用ライブラリを使います。Python で長く定番だった lifetimes は現在メンテナンスモード（積極開発は終了）で、後継として PyMC-Marketing（PyMC ベースの marketing toolbox、BG/NBD・MBG/NBD・Pareto/NBD・MMM を含む）への移行が進んでいます。どのライブラリ・API が現役かは変化が速いので、使う前に必ず最新情報を確認してください。ベータ・幾何やガンマ・ポアソンの共役と推定論そのものは、統計／ベイズ統計のテキストに土台があります。

⚠️ よくある誤解

RFM（記述）と BG/NBD（予測）を混同する：RFM分析（リセンシー・フリークエンシー・マネタリー）は「これまでどう買ったか」の要約で、 $555$ でも将来を保証しません。BG/NBD は同じ取引履歴から「次の期に何回買うか・まだ生きているか」を確率で予測します。RFM で振り分け、BG/NBD で価値を見積もる——役割が違う道具です。
過分散を無視してポアソンで予測する：顧客の取引数は頻度の個人差で過分散（本ノートで分散/平均 $3.47$ ）。ポアソン（分散＝平均）で済ませると、ばらつきも裾も過小評価し、「滅多にいないはずの超多買客」が現に居ることを説明できません。頻度の異質性は Gamma で明示的に表すべきです。
条件付き予測を平均回帰なしで外挿する：「過去10回だから将来も10回」は素朴外挿です。確率モデルの予測は集団平均へ縮約し（本ノートで $7.71$ ）、多買客の将来は少し下がります。これは「その多さの一部はまぐれ」という平均回帰で、無視すると優良客の将来を過大評価します。逆に休眠客（過去0回）も、縮約でわずかに上に見る点に注意。
離脱（生存）を入れないと「死んだ客」を生きていると誤る：§3 の縮約予測は頻度だけで、R（直近性）も離脱も使っていません。よく買っていたのに長く来ない顧客を、頻度だけ見て「将来も多買」と予測すると、実は離脱した客を生きていると誤ります。 $P(\text{alive})$ で割り引く BG/NBD・Pareto/NBD が要る所以です。
モデルの値は前提しだいで動く： $r,\theta$ （と離脱の $a,b$ ）の推定、観測期間の長さ、購買の定常性——前提が外れれば予測も動きます。点推定の一発値を信じすぎず、観測期間やパラメータを振った感度で見るのが安全です。完全な推定は要最新確認のライブラリで。