几行代码估计你的客户的价值并制定相应的策略
β-几何负二项分布(BG-NBD)模型是一个有影响力的概率模型，用于描述客户行为和预测客户终身价值(CLV)。在该系列的[上一篇文章]如何预测客户终身价值CLV中，我们已经探讨了这个模型的感知、假设和数学推导。建议回顾一下!
正如你可能从那篇文章中推断的那样，手动编码所有这些BG-NBD方程并不容易。幸运的是，有 lifetimes 的 Python 库可以为我们完成所有的流程重任。lifetimes 抽象化了BG-NBD的内部复杂性，使我们能够专注于从模型中获得洞察力和商业价值。
这个很棒的库是由Cam-Davidson Pillon创建的，他也是Probabilistic Programming & Bayesian Methods for Hackers一书的作者，这是所有贝叶斯爱好者的必读书。
本文探讨了 lifetimes 的来龙去脉。以下是我们的议程：
-
首先，我们要看一下我们的数据集在建模之前需要有的必要的表格格式。 -
第二，我们将看到如何使用 lifetimes 来训练和评估 BG-NBD模型。 -
第三，我们将了解 lifetimes 有价值的各种见解和分析。 -
最后，我们将讨论几个由 lifetimes 实现的实际商业应用。
数据准备
数据转换
通常情况下，企业在 "交易" 表中记录他们的交易。这样的表有代表单个交易的行和代表交易不同方面的列，如涉及的客户、交易的时间戳及其价值。
一个交易表的例子是由 lifetimes 提供的。这个表记录了CDNow的交易，这是一家网络公司，销售CD和音乐产品。它包含了从1997年1月1日到1998年6月30日（总共78周）来自2,357位客户的6,919笔交易。我们将在本文中使用它。以下是如何加载它。
import pandas as pd
from lifetimes.datasets import load_dataset
transactions = load_dataset(
filename='CDNOW_sample.txt',
header=None,
delim_whitespace=True,
names=['customer_id', 'customer_index', 'date', 'quantity', 'amount'],
converters={'date': lambda x: pd.to_datetime(x, format="%Y%m%d")}
)
transactions.head(3)
在拟合BG-NBD模型之前，我们需要首先将这个表重组为典型的 "RFM"格式。RFM格式的表有代表单个客户的行和以下列。
-
frequency 表示客户重复购买的数量。 -
recency 代表客户最近一次购买时的 "年龄"。 -
T 代表观察期结束时顾客的 "年龄"。 -
monetary_value , 是一个可选的列，代表一个给定客户购买的平均价值，不包括第一次购买。
lifetimes 提供了一个实用函数，方便从 "事务 "格式转换为RFM格式。
from lifetimes.utils import summary_data_from_transaction_data
rfm = summary_data_from_transaction_data(transactions=transactions,
customer_id_col='customer_id',
datetime_col='date',
monetary_value_col = 'amount',
observation_period_end=pd.to_datetime('1998-06-30'),
freq='W')
rfm.head(3)
我们对RFM数据框的第一行的解释如下。
-
在观察期内，顾客 4总共进行了4次购买（因此有3次重复购买）。 -
他第一次和最后一次购买之间的时间跨度是49周。 -
从他第一次购买到观察期结束的时间跨度为78周。 -
他的交易的平均值，不包括他的第一次交易，是23.73美元。
请注意，当 frequency、recency 和 monetary_value 列的数值为零时（如客户18），表示该客户在观察期内只进行了一次购买。
校准-观察分割
当把我们的表格改成RFM格式时，我们可以选择进行校准-观察分割，这与我们所熟悉的训练-测试分割程序的做法相似。校准-观察分割将我们的交易分为两部分，取决于它们是属于 校准期 还是 观察期 。发生在校准期的交易被用来训练模型，而发生在观察期的观察数据，被用来验证模型。
在这个例子中，我们将校准期和观察期分别设置为52周和26周。下面是我们如何使用 lifetimes 包进行分割。
from lifetimes.utils import calibration_and_holdout_data
rfm_cal_holdout = calibration_and_holdout_data(transactions=transactions,
customer_id_col='customer_id',
datetime_col='date',
monetary_value_col = 'amount',
freq='W',
calibration_period_end='1998-01-01',
observation_period_end='1998-06-30' )
我们看到，现在每个客户都与用于模型训练的校准（calibration）特征和用于模型验证的保持期holdout（观察期）特征相关。
训练、预测和评估
拟合
在 lifetimes 中，模型拟合遵循熟悉的 scikit-learn 步骤，即实例化一个模型对象（可选择包括超参数设置）并将其拟合到校准（训练）数据。
from lifetimes import BetaGeoFitter
# instantiation of BG-NBD model
bgf = BetaGeoFitter(penalizer_coef=0.0)
# fitting of BG-NBD model
bgf.fit(frequency=rfm_cal_holdout['frequency_cal'],
recency=rfm_cal_holdout['recency_cal'],
T=rfm_cal_holdout['T_cal'])
就是这样--两行字，我们就有了一个可以工作的BG-NBD模型!
你可能还记得第一篇文章**，BG-NBD模型由一个Gamma分布和一个Beta分布组成。Gamma分布的参数 r 和 α 以及Beta分布的参数 a 和 b 可以通过拟合模型对象的 .summary 属性获得。
bgf.summary
拟合模型的参数
在下一节中，我们将从这些参数中得出一些商业洞察力。
通过比较模型生成的人工数据和真实数据来评估模型拟合度
俗话说："垃圾进，垃圾出"，这句话适用于BG-NBD建模：依赖一个拟合不好的模型会导致灾难性的商业决策。**因此，在使用我们的模型进行预测和解释之前，我们应该首先评估模型的性能。
lifetimes 提供了几种评估我们模型的方法。第一个方法是比较我们的真实校准数据和从BG-NBD模型生成的分布中采样的人工数据之间的频率。函数 plot_period_transactions 只用一行就可以完成这个任务。
from lifetimes.plotting import plot_period_transactions
_ = plot_period_transactions(bgf)
该图显示，训练数据中的频率分布与我们的拟合模型人工生成的频率分布基本一致。**这表明，我们假设产生数据的物理过程被我们的模型很好地捕捉到了。
作出预测
拟合BG-NBD模型的主要目标之一是进行预测。在本节中，我们将看到如何生成预测，并随后用于进一步评估模型的性能。
当给定一个具有特定 frequency 、recency 和 T 的客户时，一个拟合的模型对象可以产生两种类型的预测。
-
他在未来 k 个时期的购买数量 -
他在观察期结束时处于活跃状态的概率。
让我们看一个例子来说明如何做到这一点。首先，我们将选择一个样本客户，检查他在校准期和观察期的 frequency 、recency 和 T。
# 首先，我们选择一个样本客户。
sample_customer = rfm_cal_holdout.iloc[20]
# 让我们检查一下他在校准期和观察期的频率、频度和T。
sample_customer
我们看到，在观察期内，客户进行了2+1=3次交易（"+ 1 "来自于将不包括第一次交易的 重复 交易的数量转换为 总 交易）。
现在让我们把这个 "真实 "的数字与BG-NBD产生的预测进行比较。下面的代码得出了客户在未来26周内（观察期的长度）的预期交易数量。
# 这个函数计算了在给定时间长度内的条件性预期交易数
# we set it to 26 weeks (the length of the observation period)
n_transactions_pred = bgf.predict(t=26,
frequency=sample_customer['frequency_cal'],
recency=sample_customer['recency_cal'],
T=sample_customer['T_cal'])
n_transactions_pred # = 0.7647440846242359
我们看到，预测的交易数（0.76次）低于实际交易数（3次）。
我们同样可以预测某个客户在校准期结束/观察期开始时仍然活跃/活着的概率。
alive_prob = bgf.conditional_probability_alive(frequency=sample_customer['frequency_cal'],
recency=sample_customer['recency_cal'],
T=sample_customer['T_cal'])
alive_prob # = 0.57089896
由于该客户在观察期内确实进行了一些购买，我们可以绝对肯定地知道，在校准期结束时他_是活跃的。因此，预测的概率为0.57是相对来说偏低。
真实和预测的交易数量的比较
在了解了如何对一个人的交易数量进行预测后，我们可以将该程序扩展到整个客户群。然后可以将预测结果与真实交易数据进行比较，以衡量我们模型的准确性。
# 观察期内的真实交易数，等于frequency_holdout + 1
rfm_cal_holdout["n_transactions_holdout_real"] = rfm_cal_holdout["frequency_holdout"] + 1
# the predicted number of transactions in the next 26 weeks (length of the observation period)
rfm_cal_holdout["n_transactions_holdout_pred"] = bgf.predict(t=26,
frequency=rfm_cal_holdout['frequency_cal'],
recency=rfm_cal_holdout['recency_cal'],
T=rfm_cal_holdout['T_cal'])
# comparison of the real and predicted transactions
rfm_cal_holdout[["n_transactions_holdout_real", "n_transactions_holdout_pred"]].head()
如果需要更严格的评估，可以对这两栏进行典型的回归指标，如RMSE。
from sklearn.metrics import mean_squared_error
RMSE = mean_squared_error(y_true = rfm_cal_holdout["n_transactions_holdout_real"],
y_pred = rfm_cal_holdout["n_transactions_holdout_pred"],
squared = False)
RMSE # = 1.3536793286521
洞察力和解释 将MLE的Gamma和Beta分布可视化 一旦我们对我们的模型的性能感到满意，我们就可以继续从它那里获得洞察力了。
一个很好的起点是对估计的Gamma和Beta参数进行提取和可视化。正如前文所阐述的，Gamma和Beta分布具有重要的商业意义。Gamma分布告诉我们客户群的交易率分布，而Beta分布反映了停用概率的分布。
mport numpy as np
import scipy.stats as stats
from matplotlib import pyplot as plt
fig, (ax_gamma, ax_beta) = plt.subplots(ncols=1, nrows =2, figsize = (20, 16))
x_gamma = np.linspace(0, 10, 1000)
y_gamma = stats.gamma.pdf(x_gamma, a=bgf.params_["alpha"], scale=bgf.params_["r"])
ax_gamma.plot(x_gamma, y_gamma, "-")
ax_gamma.set_title(f'Gamma distribution (alpha = {bgf.params_["alpha"]:.2f}, scale (r) = {bgf.params_["r"]:.2f})')
ax_gamma.set_xlabel(r'$\lambda$')
ax_gamma.set_ylabel(r'$P(\lambda)$')
x_beta = np.linspace(0, 1, 100)
y_beta = stats.beta.pdf(x_beta, a=bgf.params_["a"], b=bgf.params_["b"])
ax_beta.plot(x_beta, y_beta, "-")
ax_beta.set_title(f'Beta distribution (a = {bgf.params_["a"]:.2f}, b = {bgf.params_["b"]:.2f})')
ax_beta.set_xlabel('p')
ax_beta.set_ylabel('P(p)')
plt.show()
下面的代码显示了如何从模型中提取Gamma和Beta参数并绘制分布图。
从图中，我们可以看到Gamma分布是相对健康的，大部分𝜆在2附近。这意味着我们的客户预计将以每周2次交易的速度购物。
贝塔分布似乎也很健康，大部分 p 发现在0附近。这意味着我们的客户不太可能很快停用。
Frequency/Recency/Future Purchases 矩阵
除了所述的模型拟合评估外，拟合的模型对象还可用于一些解释性分析。
例如，我们可以看一下它的频率/时间/未来购买矩阵。这个矩阵显示了不同的频率-频率组合是如何产生不同的未来购买的预期数量的。
from lifetimes.plotting import plot_frequency_recency_matrix
_ = plot_frequency_recency_matrix(bgf)
Frequency/Recency/Future Purchases Matrix
我们可以看到，我们最好的客户在右下角--这些人最近经常购买，所以我们对他们回来的期望很高。
相反，我们最有希望的顾客在右上角--他们经常购买，然后不再来了，我们已经几个月没有看到他们了。很可能他们已经找到了另一家商店（即他们已经停用了）。.
Frequency/Recency/Probability Alive 矩阵
同样地，我们也可以产生 Frequency/Recency/Probability Alive 的概率矩阵。这个矩阵与前一个矩阵的轴线相同，但现在每个单元格的阴影告诉我们不同 Frequency/Recency 组合的客户的生存概率。
from lifetimes.plotting import plot_probability_alive_matrix
_ = plot_probability_alive_matrix(bgf)
Frequency/Recency/Probability Alive Matrix
不出所料，我们可以看到类似的模式，其中最好和最差的客户分别在矩阵的右上方和右下方。
商业应用
使用BG-NBD进行CLV估算
我们已经看到，我们如何使用BG-NBD模型来预测生存的概率 p ，以及在接下来的 k 期的购买数量。这两个预测可以反过来用来计算客户在未来 k 个时期的价值的粗略估计。该估计值可按以下方式计算。
这一估计依赖于以下两个简单的假设。
-
活着的概率 p 在下一个 k 时期保持不变 -
在接下来的 k 个时期，购买的平均价值等于观察期的平均购买价值。
请注意，这些天真的假设具有薄弱的理论基础--在现实中，p 通常在每次购买后都会发生变化，未来的购买值可能与过去的值有很大的出入。因此，得出的估计是一个粗略的估计。如果你有兴趣，存在更复杂的概率模型，如Gamma-Gamma模型，可以以更严格的方式计算未来的购买价值。
这就是我们在Pandas中进行上述计算的方法。
# the predicted number of transactions in the next 10 weeks
rfm_cal_holdout["n_transactions_10_pred"] = bgf.predict(t=10,
frequency=rfm_cal_holdout['frequency_cal'],
recency=rfm_cal_holdout['recency_cal'],
T=rfm_cal_holdout['T_cal'])
# the probability of being alive
rfm_cal_holdout["alive_prob"] = bgf.conditional_probability_alive(frequency=rfm_cal_holdout['frequency_cal'],
recency=rfm_cal_holdout['recency_cal'],
T=rfm_cal_holdout['T_cal'])
# multiplication of alive probability x number of purchases x average past purchase
rfm_cal_holdout["value_10_pred"] = rfm_cal_holdout["alive_prob"]* \
rfm_cal_holdout["n_transactions_10_pred"]*\
rfm_cal_holdout["monetary_value_cal"]
rfm_cal_holdout[["value_10_pred", "alive_prob", "n_transactions_10_pred", "monetary_value_cal"]].head()
由此产生的 value_10_pred 列封装了我们的客户在未来10周的估计价值。然后我们可以用这个估计值来衡量我们业务的健康状况。例如，我们可以获得 value_10_pred 的汇总统计。
rfm_cal_holdout["value_10_pred"].describe()
我们看到，平均而言，我们希望客户在未来10周内花费约5.18美元。
我们还可以绘制柱状图。
我们看到，我们的很大一部分客户的 value_10_pred 接近于0。结果是5.18美元的平均数是由一些异常值（少数客户的 value_10_pred 非常高）驱动的。在下一节，我们将看到如何进一步利用这一结果。
主动干预
这三个新栏目，alive/prob，n_transactions_10_pred 和 value10_pred ，将我们的客户以前看不见的方面清晰化和量化。除了提供信息之外，这些新特征还可以用来推动具体的行动。例如，我们可以利用这些新特征有选择地、主动地接触不同的客户群，以提高他们的预期CLV。
我们之前已经得出结论，我们的整体客户价值主要归功于少数高价值客户。由于这些客户是我们的收入驱动力，我们会想特别关注他们。一个想法是偶尔向他们发送祝贺信息/优惠券，以鼓励他们持续的忠诚度。
另一个想法是与高流失概率的客户（即那些活着的概率低的客户）联系，劝阻他们放弃我们的业务。
需求预测
正如我们之前看到的，value10_pred 估计了我们的客户在未来10个时期的购买量。另一种利用方式是将此值视为来自客户的预测需求。这种预测可以反过来支持我们的供应链决策（例如，库存补充或制造计划）。
注意事项：BG-NBD不是作为一个时间序列模型来设想的；它不具备时间序列的特征，如季节性和趋势的计算。因此，当涉及到预测时，谨慎的做法是不要只依赖BG-NBD的预测，而是将其与ARIMA等时间序列模型相结合。
总结
我们已经看到 lifetimes 是如何让我们方便地了解我们客户群的购买行为，并从中获得强大的、可操作的洞察力。
虽然这很了不起，但请记住，来自 lifetimes 的估计是最大似然估计，是纯粹从我们提供的数据推断出来的。
在本系列的第三篇文章中，我将向你展示实现BG-NBD的另一种方式，这次是使用由PyMC3支持的贝叶斯层次结构建模。我们将看到，这种实现方式提供了更多的灵活性，因为它允许我们在建模步骤中引入领域知识。
我希望能在那里见到你!