Задача: дано порядка 70000 данных на train - время и различные срезы данных по продажам для каждого товара item_id. Нужно научиться предсказывать продажи на 3 недели вперед.
Подключим все нужное и считаем данные.
In [1]:
import pandas as pd
from sklearn import model_selection, metrics
import numpy as np
import matplotlib.pyplot as plt
import seaborn
import xgboost
import os
%pylab inline
train = pd.read_csv("train.tsv")
test = pd.read_csv("test.tsv")
sample_submission = pd.read_csv("sample_submission.tsv")
sample_submission_a = pd.read_csv("boost_submission_a.csv")
sample_submission_b = pd.read_csv("boost_submission_b.csv")
In [2]:
def score_func(y_true, y_pred):
y_true = np.array(y_true)
y_pred = np.array(y_pred)
return (np.abs(y_true - y_pred) / (np.abs(y_true) + np.abs(y_pred))).mean() * 200.
return s
scorer = metrics.make_scorer(score_func=score_func, greater_is_better=False)
Будем обучаться на всех данных.
In [3]:
frac = 1.
train = train.sample(frac=frac, random_state=42)
X = train.drop(['Num','y'], axis=1)
y = train['y']
pd.set_option('max_columns', 64)
In [4]:
X.head(20)
Out[4]:
Можно заметить, что фичи f31-f60 полная копия f1-f30. Проверим этом и избавимся от них. Сделаем 3 модели xgboost. Первая - со всеми данными на всех фичах. Вторая - без повторяющихся f31-f60. Последняя - еще меньше фичей.
Теперь будем последовательно подбирать параметры xgboost. Перебирать все параметры сразу - слишком долго считается. Поэтому пойдем от основных параметров, как число деревьев или глубина, к менее значимым (возможно). Сначала будем искать параметр с большим шагом, потом уменьшать шаг.
Так как данные зависят от времени, то будем использовать специальное разбиение учитывающее временную зависимость (нельзя брать куски данных из "будущего").
Считаются модели долго, поэтому перезапускать не стоит. Перед созданием модели подогнанные параметры сохранены. Как показала практика - первая самая простая модель в комбинации с работой над данными дала наилучшие результаты, поэтому код подбора ее параметров приведен, остальных - нет.
In [5]:
params_height = {
'learning_rate': [0.1],
'n_estimators': [70, 100, 120, 150, 175, 200, 225]
}
zero = model_selection.GridSearchCV(xgboost.XGBRegressor(silent=False), params_height,
scoring=scorer, cv=model_selection.TimeSeriesSplit(),
fit_params={'eval_metric' : 'mae'})
In [6]:
%%time
zero.fit(X, y)
Out[6]:
In [7]:
# to produce 'beep' sound when finished
length = 0.4
frequency = 1000
os.system('play --no-show-progress --null --channels 1 synth %s sine %f' % (length, frequency))
Out[7]:
In [8]:
zero.best_params_
Out[8]:
In [9]:
params_height = {
'learning_rate': [0.1],
'n_estimators': [115, 120, 125, 130, 135]
}
zero = model_selection.GridSearchCV(xgboost.XGBRegressor(silent=False), params_height,
scoring=scorer, cv=model_selection.TimeSeriesSplit(),
fit_params={'eval_metric' : 'mae'})
In [10]:
%%time
zero.fit(X, y)
Out[10]:
In [11]:
length = 0.4
frequency = 1000
os.system('play --no-show-progress --null --channels 1 synth %s sine %f' % (length, frequency))
Out[11]:
In [12]:
zero.best_params_
Out[12]:
In [13]:
params2 = {
'n_estimators' : [115],
'learning_rate': [0.1],
'max_depth': [2, 6, 10, 14, 18],
'min_child_weight' : [1, 3, 5]
}
zero = model_selection.GridSearchCV(xgboost.XGBRegressor(silent=False), params2,
scoring=scorer, cv=model_selection.TimeSeriesSplit(),
fit_params={'eval_metric' : 'mae'})
In [14]:
%%time
zero.fit(X, y)
Out[14]:
In [15]:
length = 0.4
frequency = 1000
os.system('play --no-show-progress --null --channels 1 synth %s sine %f' % (length, frequency))
Out[15]:
In [16]:
zero.best_params_
Out[16]:
In [17]:
params3 = {
'n_estimators' : [115],
'learning_rate': [0.1],
'max_depth': [18],
'min_child_weight' : [1],
'gamma': [i / 10.0 for i in range(0, 5)]
}
zero = model_selection.GridSearchCV(xgboost.XGBRegressor(silent=False), params3,
scoring=scorer, cv=model_selection.TimeSeriesSplit(),
fit_params={'eval_metric' : 'mae'})
In [18]:
%%time
zero.fit(X, y)
Out[18]:
In [19]:
length = 0.4
frequency = 1000
os.system('play --no-show-progress --null --channels 1 synth %s sine %f' % (length, frequency))
Out[19]:
In [20]:
zero.best_params_
Out[20]:
In [21]:
params4 = {
'gamma': [0.3],
'n_estimators' : [115],
'learning_rate': [0.1],
'max_depth': [18],
'subsample': [i / 10.0 for i in range(6, 10)],
'colsample_bytree': [i / 10.0 for i in range(6, 10)]
}
zero = model_selection.GridSearchCV(xgboost.XGBRegressor(silent=False), params4,
scoring=scorer, cv=model_selection.TimeSeriesSplit(),
fit_params={'eval_metric' : 'mae'})
In [22]:
%%time
zero.fit(X, y)
Out[22]:
In [23]:
length = 0.4
frequency = 1000
os.system('play --no-show-progress --null --channels 1 synth %s sine %f' % (length, frequency))
Out[23]:
In [24]:
zero.best_params_
Out[24]:
Так как считается все долго, введем сохраним параметры выше.
In [12]:
zero_best_params = {
'colsample_bytree': 0.9,
'gamma': 0.3,
'learning_rate': 0.1,
'max_depth': 18,
'n_estimators': 115,
'subsample': 0.9
}
In [13]:
best_zero = xgboost.XGBRegressor(**zero_best_params)
In [14]:
%%time
model = best_zero
model.fit(X, y, eval_metric='mae')
test_drop = test.drop(['Num'], axis=1)
preds = model.predict(test_drop)
print len(preds)
print len(sample_submission)
In [15]:
sample_submission['y'] = preds
In [16]:
sample_submission['y'] = sample_submission['y'].map(lambda x: x if x > 0 else 0.0)
In [17]:
sample_submission.to_csv("boost_submission_zero.csv", sep=',', index=False)
In [18]:
model_selection.cross_val_score(model, X, y, scoring=scorer)
Out[18]:
Данное решение обеспечивает SMAPE на public leaderboard примерно 23, что похоже на полученную ранее оценку.
Проверим, что f31-f60 полная копия f1-f30.
In [19]:
sum(sum(X[X.columns[i + 4]] != X[X.columns[i + 34]] for i in xrange(30)))
Out[19]:
Удалим эти копии и подгоним модель.
In [20]:
X = X.drop(X.columns[34:], axis=1)
Сохраним параметры
In [21]:
a_best_params = {
'colsample_bytree': 0.9,
'gamma': 0.3,
'learning_rate': 0.01,
'max_depth': 18,
'n_estimators': 115*10,
'subsample': 0.9
}
In [22]:
best = xgboost.XGBRegressor(**a_best_params)
In [23]:
%%time
model = best
model.fit(X, y, eval_metric='mae')
test_drop = test.drop(['Num'], axis=1)
test_drop = test_drop.drop(test_drop.columns[34:], axis=1)
preds = model.predict(test_drop)
print len(preds)
print len(sample_submission_a)
In [24]:
sample_submission_a['y'] = preds
In [25]:
sample_submission_a['y'] = sample_submission_a['y'].map(lambda x: x if x > 0 else 0.0)
In [26]:
sample_submission_a.to_csv("boost_submission_a.csv", sep=',', index=False)
In [27]:
model_selection.cross_val_score(model, X, y, scoring=scorer)
Out[27]:
Теперь, учитывая скоррелированность фичей f, попробуем убрать часть из них и посмотреть на результаты.
Посчитаем коэффициенты корреляции между оставшимися признаками. Те, которые наиболее коррелируют стоит убрать (один из двух).
In [28]:
X.corr('kendall')
Out[28]:
Можно заметить, что фичи f коррелируют между собой очень сильно, причем чем дальше (номера) - тем меньше. Поэтому попробуем выбрать несколько из них.
In [29]:
columns = X.columns[[0, 1, 2, 3, 4, 14, 24, -1]]
columns
Out[29]:
In [30]:
X = X[columns]
Сохраним параметры
In [31]:
b_best_params = {
'colsample_bytree': 0.9,
'gamma': 0.2,
'learning_rate': 0.01,
'max_depth': 21,
'n_estimators': 85*10,
'subsample': 0.9
}
In [32]:
best_b = xgboost.XGBRegressor(**b_best_params)
In [33]:
%%time
model = best_b
model.fit(X, y, eval_metric='mae')
preds = model.predict(test_drop[test_drop.columns[[0, 1, 2, 3, 4, 14, 24, -1]]])
print len(preds)
print len(sample_submission_b)
In [34]:
sample_submission_b['y'] = preds
In [35]:
sample_submission_b['y'] = sample_submission_b['y'].map(lambda x: x if x > 0 else 0.0)
In [36]:
sample_submission_b.to_csv("boost_submission_b.csv", sep=',', index=False)
In [37]:
model_selection.cross_val_score(model, X, y, scoring=scorer)
Out[37]:
Хоть уменьшение кол-ва скоррелированных признаков должно было положительно сказаться на результатах, заметного улучшения нет.
In [38]:
train = pd.read_csv("train.tsv")
test = pd.read_csv("test.tsv")
train = train.drop(train.columns[36:], axis=1)
sorted_train = train.sort_values(['item_id', 'year', 'week', 'shift'])
check_id = sorted_train['item_id'][0]
sorted_train[sorted_train['item_id']==check_id].head(20)
Out[38]:
Можно заметить, что значения фичей f для одного item_id очень похожи. Сделаем срез по shift=1.
In [39]:
sorted_train[sorted_train['shift']==1][sorted_train['item_id']==check_id].head(20)
Out[39]:
Видно, что для одного item_id при shift=1 значения фичей по побочным диагоналям совпадают. Запомним это.
Построим график зависимости target переменной от, например, последней f30, помня, что все остальные по диагонали совпдают с ней.
In [40]:
first_id_data = sorted_train[sorted_train['item_id']==check_id]
first_id_data = first_id_data[first_id_data['shift']==1]
time = np.arange(len(first_id_data))
In [41]:
plt.figure(figsize=(16, 10))
plt.plot(time, first_id_data['y'], label='y')
plt.plot(time, first_id_data['f30'], label='f30')
plt.xlabel('time')
plt.grid(True)
plt.legend(loc='best')
plt.show()
Видно, что значения совпадают с точностью до множителя и сдвига на 1 шаг по времени.
Поэтому ответом будут значения f30, умноженные на посчитанный коэффициент между y и f30. Но для последних недель для каждого item_id соответствующего значения f30 нет (из-за сдвига по времени). Поэтому для них будем использовать предсказания полученные ранее.
In [48]:
boost_sample = pd.read_csv('boost_submission_zero.csv')
boost_sample = dict(boost_sample.values)
In [49]:
uids = test['item_id'].unique()
train = train[['Num', 'item_id', 'shift', 'week', 'year', 'f30', 'y']]
test = test[['Num', 'item_id', 'shift', 'week', 'year', 'f30']]
def get_table(item_id):
id_table_train = train[train['item_id'] == item_id]
id_table_test = test[test['item_id'] == item_id]
id_table_train = id_table_train.sort_values(['year', 'week', 'shift'])
id_table_train = id_table_train[id_table_train['shift'] == 1]
id_table_test_full = id_table_test.sort_values(['year', 'week', 'shift'])
id_table_test = id_table_test_full[id_table_test_full['shift'] == 1]
id_table_test_shift = id_table_test_full[id_table_test_full['shift'] == 2]
id_table_test['Extra_Num_2'] = id_table_test_shift['Num'].values
id_table_test_shift = id_table_test_full[id_table_test_full['shift'] == 3]
id_table_test['Extra_Num_3'] = id_table_test_shift['Num'].values
full_table = id_table_train.append(id_table_test)
full_table = full_table.drop(['year', 'week', 'shift', 'item_id'], axis=1)
full_table['f30'] = (list(full_table['f30']) + [None])[1:]
return full_table
In [50]:
with open('submission.csv', 'w') as f:
f.write('Num,y\n')
for item_id in uids:
table = get_table(item_id=item_id)
w = []
for row in table[table['Extra_Num_2'].isnull()].values:
try:
w.append(row[-2] / row[-1])
except TypeError:
print table
w0 = np.median(w)
table = table.drop(['y'], axis=1)[table['Extra_Num_2'].notnull()]
for row in table.values:
if not np.isnan(row[-1]):
y = row[-1] / w0
else:
y = np.median([boost_sample[row[0]], boost_sample[row[1]], boost_sample[row[2]]])
f.write('{},{}\n'.format(int(row[0]), y))
f.write('{},{}\n'.format(int(row[1]), y))
f.write('{},{}\n'.format(int(row[2]), y))
Попробовав 3 варианта xgboost, оказывается, что первый работает получше. SMAPE на public liderboard дали 8.06, 13, 8.37. Таким образом лучшие результаты дали простая модель с хорошими параметрами + работа с данными.
In [2]:
train['y'].mean()
Out[2]:
In [5]:
y = pd.read_csv('sample_submission.tsv')
y.head()
Out[5]:
Таким образом sample_submission.tsv - среднее по y на train.
In [ ]: