基于BG/NBD概率模型的用户CLV预测

# pip install lifetimes

import pandas as pd
import numpy as np
import warnings
import matplotlib.pyplot as plt
import seaborn as sns
import lifetimes
import toad 
from lifetimes.utils import *
from lifetimes import BetaGeoFitter
from lifetimes import GammaGammaFitter
from lifetimes.plotting import *

# 初始化设置
%matplotlib inline
pd.set_option('display.max_columns', None) # 显示所有列
sns.set(style="ticks")

# 读取数据
raw_data = pd.read_excel('Online Retail.xlsx')
raw_data.head()

# 查看数据信息
toad.detector.detect(raw_data)

# 数据处理
df=raw_data.copy()
# 剔除缺失的CustomerID
df.dropna(subset=['CustomerID'],inplace=True)
# 剔除Quantity<=0 or UnitPrice<=0的数据
df.query('Quantity > 0 & UnitPrice >0', inplace=True)
# 生成结果数据
raw_result = df.copy()
# 再次查看数据信息
toad.detector.detect(raw_result)

# 函数方式-通过lifetimes的summary_data_from_transaction_data
df_model=raw_result.copy()
# 生成价格数据
df_model['Sales']=df_model['Quantity']*df_model['UnitPrice']
# 生成观察日期
last_date=df_model.InvoiceDate.max()
df_temp = summary_data_from_transaction_data(df_model, 'CustomerID', 'InvoiceDate',
                         monetary_value_col='Sales',
                         observation_period_end=last_date)
df_temp.sort_values('monetary_value', ascending=True).head()

# 构造模型数据
df_model_finall=df_temp.copy()

# 构造模型输入数据RF:T
df_model=raw_result.copy()
# 生成价格数据
df_model['Sales']=df_model['Quantity']*df_model['UnitPrice']
# 生成观察日期
last_date=df_model.InvoiceDate.max()
# 手动计算
df_model.set_index('CustomerID',inplace=True)
df_temp2=df_model.groupby('CustomerID').agg(
 # 计算RF：T
 max_date=('InvoiceDate','max'),
 min_date=('InvoiceDate','min'),
 frequency=('InvoiceDate',lambda x: pd.Series(x.dt.to_period('D')).nunique()), # 按日计算频次
 recency=('InvoiceDate',lambda x: (max(x.dt.to_period('D')) - 
                                   min(x.dt.to_period('D')))/np.timedelta64(1, 'D')), # 按日计算最近时间差
 T=('InvoiceDate',lambda x: (pd.Timestamp(last_date).to_period('D')-
                             min(pd.Series(x.dt.to_period('D'))))/np.timedelta64(1,'D')), # 按日计算观察时间差
 Invoice_count=('InvoiceNo','nunique'), # 计算订单数
 Stock_count=('StockCode','count'), # 计算产品数
 Sales_sum=('Sales','sum'), # 计算销售额总额
)
df_temp2['Sales_mean']=df_temp2['Sales_sum']/df_temp2['frequency']
df_temp2.sort_values('Sales_mean', ascending=True).head()

# summary_data_from_transaction_data可以通过参数设置日期差的方式，是否包含首次购买
df_model=df.copy()
# 生成价格数据
df_model['Sales']=df_model['Quantity']*df_model['UnitPrice']
# 生成观察日期
last_date=df_model.InvoiceDate.max()
df_temp3 = summary_data_from_transaction_data(df_model, 'CustomerID', 'InvoiceDate',
                      monetary_value_col='Sales',
                      observation_period_end=last_date,
                      freq='W',
                      include_first_transaction='True'
                      )
df_temp3.sort_values('frequency', ascending=True).head()

# 模型拟合
bgf = BetaGeoFitter(penalizer_coef=0)
bgf.fit(df_model_finall['frequency'], df_model_finall['recency'], df_model_finall['T'])

<lifetimes.BetaGeoFitter: fitted with 4338 subjects, a: 0.00, alpha: 68.91, b: 6.75, r: 0.83>

# 用户预期交易热力图
fig = plt.figure(figsize=(12,8))
plot_frequency_recency_matrix(bgf, T=1, cmap='cool')
plt.show()

# 用户留存概率热力图
fig = plt.figure(figsize=(12,8))
plot_probability_alive_matrix(bgf, cmap='cool')
plt.show()

# 预测用户下个时期（t）的预期购买量
t = 30
df_model_finall['predicted_purchases'] = bgf.conditional_expected_number_of_purchases_up_to_time(t, 
                                  df_model_finall['frequency'], df_model_finall['recency'], df_model_finall['T'])
df_model_finall.sort_values(by='predicted_purchases', ascending=False).head()

# 我们仅估算至少有一次重复购买的客户
df_gg_model=df_model_finall[df_model_finall['frequency']>0]
df_gg_model.shape

(2790, 5)

# 前提假设：购买频次和购买金额无相关性
df_gg_model[['monetary_value', 'frequency']].corr()

# 模型拟合
ggf = GammaGammaFitter(penalizer_coef = 0)
ggf.fit(df_gg_model['frequency'],
        df_gg_model['monetary_value'])

<lifetimes.GammaGammaFitter: fitted with 2790 subjects, p: 2.10, q: 3.45, v: 485.89>

# 预测每笔交易的预期收益
ggf.conditional_expected_average_profit(
        df_gg_model['frequency'],
        df_gg_model['monetary_value'])

CustomerID
12347.0    569.978836
12348.0    333.784235
12352.0    376.175359
12356.0    324.039419
12358.0    539.907126
              ...    
18272.0    474.368524
18273.0    201.838133
18282.0    260.340479
18283.0    174.532812
18287.0    492.169257
Length: 2790, dtype: float64

# 预测用户未来的CLV
bgf.fit(df_gg_model['frequency'], df_gg_model['recency'], df_gg_model['T'])
ggf.customer_lifetime_value(
    bgf, # bgf预测用户未来的预期购买量
    df_gg_model['frequency'],
    df_gg_model['recency'],
    df_gg_model['T'],
    df_gg_model['monetary_value'],
    time=12, # 用户的预期寿命，以月为单位
    freq='D', # T的单位,默认为'D'
    discount_rate=0.01 # monthly discount rate ~ 12.7% annually
)

CustomerID
12347.0    3194.349694
12348.0    1175.905597
12352.0    2472.058198
12356.0     979.177238
12358.0    1940.851493
              ...     
18272.0    3111.075131
18273.0     723.867524
18282.0    1026.433547
18283.0    1966.557687
18287.0    2070.391699
Name: clv, Length: 2790, dtype: float64

# 预期购买次数校验
plot_period_transactions(bgf)
plt.show()

# 分数据集校准

# 设置校准结束时间，观察结束时间，对数据集划分
summary_cal_holdout = calibration_and_holdout_data(raw_result, 'CustomerID', 'InvoiceDate',
                                        calibration_period_end='2011-05-01',
                                        observation_period_end='2011-12-9' )  
# 拟合查看校准结果
bgf = BetaGeoFitter(penalizer_coef=0.1)
bgf.fit(summary_cal_holdout['frequency_cal'], summary_cal_holdout['recency_cal'], 
        summary_cal_holdout['T_cal'])
plot_calibration_purchases_vs_holdout_purchases(bgf, summary_cal_holdout)
plt.show()

/Users/heinrich/opt/anaconda3/lib/python3.8/site-packages/pandas/core/series.py:726: RuntimeWarning: invalid value encountered in sqrt
  result = getattr(ufunc, method)(*inputs, **kwargs)
/Users/heinrich/opt/anaconda3/lib/python3.8/site-packages/pandas/core/series.py:726: RuntimeWarning: invalid value encountered in log
  result = getattr(ufunc, method)(*inputs, **kwargs)