导读：本文是“数据拾光者”专栏的第五十一篇文章，这个系列将介绍在广告行业中自然语言处理和推荐系统实践。本篇分享了kaggle比赛《Corporación Favorita Grocery Sales Forecasting》冠军方案,对商品销量预测相关问题感兴趣的小伙伴可以一起沟通交流。
摘要：本篇分享了kaggle比赛《Corporación Favorita Grocery Sales Forecasting》冠军方案。因为业务需要所以调研了商品销量预测比赛，重点学习了冠军方案的特征工程和模型构建，其中关于时间滑动窗口特征的构建非常巧妙，受益匪浅。对商品销量预测相关问题感兴趣的小伙伴可以一起沟通交流。
下面主要按照如下思维导图进行学习分享：
01
比赛介绍及数据理解
最近因为工作原因需要调研下kaggle比赛《Corporación Favorita Grocery Sales Forecasting》top方案的特征和模型工作，可以借鉴并应用到实际业务中。很多时候我们的任务可能与kaggle中某个比赛是类似的，想又快又好的完成目标其中一条有效的方法就是参考大牛分享的方案。因为很多大牛在比赛打完之后会分享自己的源码，这样相比于我们自己去从0到1构建模型效率会提升很多。不仅如此，参考大牛的方案可以让我们了解当前业界对于此类问题优秀的方案，快速得到很好的baseline，然后快速迭代更新更容易出成果。
该比赛kaggle地址如下：
https://www.kaggle.com/c/favorita-grocery-sales-forecasting/overview
整体来看该比赛就是预测商品的销量，官方提供了2013-2017年各商店商品的销量，参赛队伍需要根据已有数据预测未来一段时间商店商品的销量，下面是每年的训练样本量：
图1 每年的训练样本量
每年各月的训练样本分布如下：
图2 每年各月的训练样本
对官方提供的数据进行整理，下面是数据说明和示例：
图3 数据说明和示例
其中id代表唯一key值，实际无用；date代表日期，store_nbr代表商店id，item_nbr代表商品id，预测的粒度也就是某个商店中的某个商品在某一天的销量；onpromotion代表当天商店的该商品是否在促销；紫色的四个字段都是商店的特征，其中city代表市，state代表州，tpye代表商店层级，一共有A-E五个等级，cluster代表相似商店分组，一共有17个分组；红色的三个字段是商品的特征，其中family是商品分类，总共有33个分类，class代表商品小类，有337个小类别，perishable代表商品是否容易变质；oil_price代表油价，day_type代表假日类型。
02
详解冠军方案
冠军方案介绍地址如下：
https://www.kaggle.com/c/favorita-grocery-sales-forecasting/discussion/47582
2.1 样本选择
虽然官方提供了2013-2017的训练样本，但是作者仅使用了2017年的训练样本。详细如下：
训练样本：20170531 - 20170719 or 20170614 - 20170719
验证样本：20170726 - 20170810
测试集：20170816 - 20170820
2.2 特征工程
整体来看特征包括两大块，第一块是基本特征，第二块是时间滑动窗口特征。
(1)基本特征
基本特征主要包括item_nbr、family、class、perishable、store_nbr、city、state、type、cluster。
(2) 时间滑动窗口特征
这里重点研究时间滑动窗口特征。作者使用num_day个滑动窗口，分别统计item-store粒度、item粒度和store-class粒度的时间滑动窗口特征，关于时间滑动窗口特征介绍如下：
图4 时间滑动窗口特征介绍如下
时间滑动窗口具体特征如下：
- 时间窗口内(最近14/60/140天)：促销天数汇总；
- 时间窗口内(后3/7/14天)：促销天数汇总。这里需要介绍下为什么可以使用之后的促销天数数据，因为在测试集中官方已经给出了未来一段时间某商店某商品是否会进行促销，所以我们可以用未来几天促销的数据；
- 时间窗口内(最近3/7/14/30/60/140天)：
- 销量差值的均值，以时间窗口最近3天为例，用第二天的销量减去第一天的销量，再用第三天的销量减去第二天的销量，将两者取均值就可以得到销量差值的均值，这个特征可以理解为想查看每天的销量增长率；
- 销量每天按0.9衰减之后汇总，以时间窗口3天为例，最近一天销量不变，最近第二天的销量乘以衰减系数0.9，最近第三天的销量乘以衰减系数0.81，然后将三天衰减之后的销量相加；
- 均值、中位数、最小值、最大值和标准差；
- 时间窗口内(上一周最近3/7/14/30/60/140天)：和前一天销量差值的均值、销量每天按0.9衰减之后汇总、均值、中位数、最小值、最大值和标准偏差。这个特征和上一个特征是一样的，只不过计算的是上一周各个特征值，作者想查看前一周的销量各个特征；
- 时间窗口内(最近7/14/30/60/140天)：
- 有销量/促销的天数，分别查看时间窗口内有销量和促销的天数，以时间窗口3天为例，如果这三天都有销量，那么为3；
- 距离上次有销量/促销的天数，以时间窗口3天为例，上一次有销量是昨天，那么该值为1。这个特征主要是查看上一次有销量或者促销对未来商品销量的影响，以促销为例，有些商品近期才做过促销，可能未来几天的销量就会受影响；
- 距离最早有销量/促销的天数，以时间窗口3天为例，最早有销量是最近第三天，那么该值为3；
- 时间窗口内(后15天)促销的天数、距离上次促销的天数、距离最早促销的天数，这个特征和上一个特征类似，只不过查看未来15天各个特征情况；
- 时间窗口内(最近15天)当天的销量；
- 最近4周时间窗口为(每周1-每周日)的销量均值，比如最近4周每周1的销量均值；
- 最近20周时间窗口为(每周1-每周日)的销量均值，比如最近20周每周1的销量均值；
- 时间窗口内(前16到后15天)每天是否促销。
特征加工代码如下：
代码语言：javascript复制# 计算不同时间窗口的特征
def get_timespan(df, dt, minus, periods,freq='D'):
df_result = df[pd.date_range(dt - timedelta(days=minus),periods=periods, freq=freq)]
return df_result
def prepare_dataset(df, promo_df, t2017,is_train=True, name_prefix=None):
X= {
# 以t2017为起点，最近14/60/140天促销汇总
"promo_14_2017": get_timespan(promo_df, t2017, 14,14).sum(axis=1).values,
"promo_60_2017": get_timespan(promo_df, t2017, 60,60).sum(axis=1).values,
"promo_140_2017": get_timespan(promo_df, t2017, 140,140).sum(axis=1).values,
# 以t2017为起点，后3/7/14天促销汇总
"promo_3_2017_aft": get_timespan(promo_df, t2017 timedelta(days=16), 15, 3).sum(axis=1).values,
"promo_7_2017_aft": get_timespan(promo_df, t2017 timedelta(days=16), 15, 7).sum(axis=1).values,
"promo_14_2017_aft": get_timespan(promo_df, t2017 timedelta(days=16), 15, 14).sum(axis=1).values,
}
#t2017为起点
for i in [3, 7, 14, 30, 60, 140]:
tmp = get_timespan(df, t2017, i, i)
# 最近i天里和前一天销量差值的均值
X['diff_%s_mean' % i] = tmp.diff(axis=1).mean(axis=1).values
# 最近i天里销量每天按0.9衰减之后汇总 *******************************
X['mean_%s_decay' % i] = (tmp * np.power(0.9,np.arange(i)[::-1])).sum(axis=1).values
# 最近i天里均值、中位数、最小值、最大值和标准偏差
X['mean_%s' % i] = tmp.mean(axis=1).values
X['median_%s' % i] = tmp.median(axis=1).values
X['min_%s' % i] = tmp.min(axis=1).values
X['max_%s' % i] = tmp.max(axis=1).values
X['std_%s' % i] = tmp.std(axis=1).values
#t2017上一周，前i天各指标值，和上面是一样的
for i in [3, 7, 14, 30, 60, 140]:
tmp = get_timespan(df, t2017 timedelta(days=-7), i, i)
X['diff_%s_mean_2' % i] = tmp.diff(axis=1).mean(axis=1).values
X['mean_%s_decay_2' % i] = (tmp * np.power(0.9,np.arange(i)[::-1])).sum(axis=1).values
X['mean_%s_2' % i] = tmp.mean(axis=1).values
X['median_%s_2' % i] =tmp.median(axis=1).values
X['min_%s_2' % i] = tmp.min(axis=1).values
X['max_%s_2' % i] = tmp.max(axis=1).values
X['std_%s_2' % i] = tmp.std(axis=1).values
#t2017为起点，最近i天内有销量/促销的天数、距离上次有销量的天数、距离最早有销量的天数
for i in [7, 14, 30, 60, 140]:
tmp = get_timespan(df, t2017, i, i)
# 最近i天内有销量的天数
X['has_sales_days_in_last_%s' % i] = (tmp > 0).sum(axis=1).values
# 最近i天内距离上次有销量的天数，如果都没有销量则为i
X['last_has_sales_day_in_last_%s' % i] = i - ((tmp > 0) *np.arange(i)).max(axis=1).values
# 最近i天内距离最早有销量的天数
X['first_has_sales_day_in_last_%s' % i] = ((tmp > 0) * np.arange(i,0, -1)).max(axis=1).values
tmp = get_timespan(promo_df, t2017, i, i)
X['has_promo_days_in_last_%s' % i] = (tmp > 0).sum(axis=1).values
X['last_has_promo_day_in_last_%s' % i] = i - ((tmp > 0) *np.arange(i)).max(axis=1).values
X['first_has_promo_day_in_last_%s' % i] = ((tmp > 0) * np.arange(i,0, -1)).max(axis=1).values
#t2017为起点，未来15天内有促销的天数、距离上次有促销的天数、距离最早有促销的天数
tmp = get_timespan(promo_df, t2017 timedelta(days=16), 15, 15)
X['has_promo_days_in_after_15_days'] = (tmp > 0).sum(axis=1).values
X['last_has_promo_day_in_after_15_days'] = i - ((tmp > 0) *np.arange(15)).max(axis=1).values
X['first_has_promo_day_in_after_15_days'] = ((tmp > 0) *np.arange(15, 0, -1)).max(axis=1).values
#t2017为起点，前i天当天的销量
for i in range(1, 16):
X['day_%s_2017' % i] = get_timespan(df, t2017, i, 1).values.ravel()
#t2017为起点，最近4/20周时间窗口为(每周1-每周日)的销量均值，比如最近4周每周周1的均值；
for i in range(7):
X['mean_4_dow{}_2017'.format(i)] = get_timespan(df, t2017, 28-i, 4,freq='7D').mean(axis=1).values
X['mean_20_dow{}_2017'.format(i)] = get_timespan(df, t2017, 140-i, 20,freq='7D').mean(axis=1).values
#t2017为起点，前后16天当天促销
for i in range(-16, 16):
# 需要把t2017 timedelta(days=i) 转化成str格式，否则会报错
X["promo_{}".format(i)] = promo_df[str(t2017 timedelta(days=i))].values.astype(np.uint8)
X= pd.DataFrame(X)
#y是未来16天当天销量
if is_train:
y = df[pd.date_range(t2017, periods=16)].values
return X, y
if name_prefix is not None:
X.columns = ['%s_%s' % (name_prefix, c) forc in X.columns]
return X
print("Preparing dataset...")
#num_days = 8
num_days = 2
t2017 = date(2017, 5, 31)
X_l, y_l = [], []
for i in range(num_days):
delta = timedelta(days=7 * i)
#store_nbr-item_nbr粒度
X_tmp, y_tmp = prepare_dataset(df_2017, promo_2017, t2017 delta)
#item_nbr粒度
X_tmp2 = prepare_dataset(df_2017_item, promo_2017_item, t2017 delta,is_train=False, name_prefix='item')
X_tmp2.index = df_2017_item.index
X_tmp2 =X_tmp2.reindex(df_2017.index.get_level_values(1)).reset_index(drop=True)
#store-class粒度
X_tmp3 = prepare_dataset(df_2017_store_class, df_2017_promo_store_class,t2017 delta, is_train=False, name_prefix='store_class')
X_tmp3.index = df_2017_store_class.index
#构建多重索引必须从这里pd.MultiIndex.from_frame，源代码会报错
X_tmp3 =X_tmp3.reindex(pd.MultiIndex.from_frame(df_2017_store_class_index)).reset_index(drop=True)
X_tmp3
#将不同粒度的训练数据合并
X_tmp = pd.concat([X_tmp, X_tmp2, X_tmp3, items.reset_index(),stores.reset_index()], axis=1)
X_l.append(X_tmp)
y_l.append(y_tmp)2.3 模型构建
作者分别使用lgb和nn构建模型，最后通过加权求和的方式得到最终结果。
(1) 单模型效果
- model_1 : 0.506 / 0.511 , 16 lgb modelstrained for each day source code；
- model_2 : 0.507 / 0.513 , 16 nn modelstrained for each day source code；
- model3 : 0.512 / 0.515，1 lgb model for 16 days with almost same features as model1；
- model_4 : 0.517 / 0.519，1 nn model based on @sjv's code
其中mode1和model3使用的是传统lgb模型，model2和model4使用的是神经网络模型，下面是神经网络模型结构：
图5 神经网络模型结构
作者使用LSTM作为特征抽取器，后面再加全连接层。因为当时比赛时间比较早，Transformer还没被使用，如果现在要应用到实际业务中，将LSTM替换为Transformer可能会提升模型效果。
神经网络模型构建源码如下：
代码语言：javascript复制def build_model():
model = Sequential()
model.add(LSTM(512, input_shape=(X_train.shape[1],X_train.shape[2])))
model.add(BatchNormalization())
model.add(Dropout(.2))
model.add(Dense(256))
model.add(PReLU())
model.add(BatchNormalization())
model.add(Dropout(.1))
model.add(Dense(256))
model.add(PReLU())
model.add(BatchNormalization())
model.add(Dropout(.1))
model.add(Dense(128))
model.add(PReLU())
model.add(BatchNormalization())
model.add(Dropout(.05))
model.add(Dense(64))
model.add(PReLU())
model.add(BatchNormalization())
model.add(Dropout(.05))
model.add(Dense(32))
model.add(PReLU())
model.add(BatchNormalization())
model.add(Dropout(.05))
model.add(Dense(16))
model.add(PReLU())
model.add(BatchNormalization())
model.add(Dropout(.05))
model.add(Dense(1))
return model(2) 模型融合
作者通过加权求和的方式将多模型结果进行融合，这也是kaggle比赛提分的套路了，最终提交的结果是：
finalmodel=0.42*model1 0.28 * model2 0.18 * model3 0.12 * model4
03
其他top方案
整理了该比赛其他top2-top6的方案，感兴趣的小伙伴可以好好学习下：
- 2st方案：https://www.kaggle.com/c/favorita-grocery-sales-forecasting/discussion/47568#latest-278474
- 3st方案：https://www.kaggle.com/c/favorita-grocery-sales-forecasting/discussion/47560#latest-302253
- 4st方案：https://www.kaggle.com/c/favorita-grocery-sales-forecasting/discussion/47529#latest-271077
- 5st 方案：https://github.com/LenzDu/Kaggle-Competition-Favorita
- 6st方案：https://www.kaggle.com/c/favorita-grocery-sales-forecasting/discussion/47575#latest-269568
04
总结及反思
本篇分享了kaggle比赛《Corporación Favorita Grocery Sales Forecasting》冠军方案。因为业务需要所以调研了商品销量预测比赛，重点学习了冠军方案的特征工程和模型构建，其中关于时间滑动窗口特征的构建非常巧妙，受益匪浅。对商品销量预测相关问题感兴趣的小伙伴可以一起沟通交流。
