<jupyter_start><jupyter_text>Predict future sales We are asking you to predict total sales for every product and store in the next month. By solving this competition you will be able to apply and enhance your data science skills.You are provided with daily historical sales data. The task is to forecast the total amount of products sold in every shop for the test set. Note that the list of shops and products slightly changes every month. Creating a robust model that can handle such situations is part of the challenge. Data fields description:* ID - an Id that represents a (Shop, Item) tuple within the test set* shop_id - unique identifier of a shop* item_id - unique identifier of a product* item_category_id - unique identifier of item category* date_block_num - a consecutive month number, used for convenience. January 2013 is 0, February 2013 is 1,..., October 2015 is 33* date - date in format dd/mm/yyyy* item_cnt_day - number of products sold. You are predicting a monthly amount of this measure* item_price - current price of an item* item_name - name of item* shop_name - name of shop* item_category_name - name of item category Dependencies<jupyter_code>import datetime
import warnings
import numpy as np
import pandas as pd
import seaborn as sns
import matplotlib.pyplot as plt
from xgboost import XGBRegressor
from xgboost import plot_importance
from sklearn.metrics import mean_squared_error
from sklearn.linear_model import LinearRegression
from sklearn.neighbors import KNeighborsRegressor
from sklearn.ensemble import RandomForestRegressor
from sklearn.preprocessing import StandardScaler, MinMaxScaler

%matplotlib inline
sns.set(style="darkgrid")
pd.set_option('display.float_format', lambda x: '%.2f' % x)
warnings.filterwarnings("ignore")<jupyter_output><empty_output><jupyter_text>Loading data<jupyter_code>sales = pd.read_csv('../competitive-data-science-predict-future-sales/sales_train.csv')
test = pd.read_csv('../competitive-data-science-predict-future-sales/test.csv')
submission = pd.read_csv('../competitive-data-science-predict-future-sales/sample_submission.csv')
items = pd.read_csv('../competitive-data-science-predict-future-sales/items.csv')
item_categories = pd.read_csv('../competitive-data-science-predict-future-sales/item_categories.csv')
shops = pd.read_csv('../competitive-data-science-predict-future-sales/shops.csv')<jupyter_output><empty_output><jupyter_text>Join data sets<jupyter_code>train = sales.join(items, on='item_id', rsuffix='_').join(shops, on='shop_id', rsuffix='_').join(item_categories, on='item_category_id', rsuffix='_').drop(['item_id_', 'shop_id_', 'item_category_id_'], axis=1)<jupyter_output><empty_output><jupyter_text>Let's take a look at the raw data<jupyter_code>print('Train rows: ', train.shape[0])
print('Train columns: ', train.shape[1])
train.head().T
train.describe()<jupyter_output><empty_output><jupyter_text>Time period of the dataset<jupyter_code>print('Min date from train set: %s' +train['date'].min())
print('Max date from train set: %s' +train['date'].max())<jupyter_output>Min date from train set: %s01.01.2013
Max date from train set: %s31.12.2014<jupyter_text>Data leakagesAbout data leakages I'll only be using only the "shop_id" and "item_id" that appear on the test set.<jupyter_code>test_shop_ids = test['shop_id'].unique()
test_item_ids = test['item_id'].unique()
# Only shops that exist in test set.
lk_train = train[train['shop_id'].isin(test_shop_ids)]
# Only items that exist in test set.
lk_train = lk_train[lk_train['item_id'].isin(test_item_ids)]
print('Data set size before leaking:', train.shape[0])
print('Data set size after leaking:', lk_train.shape[0])<jupyter_output>Data set size before leaking: 2935849
Data set size after leaking: 1224439<jupyter_text>Data cleaning    Only records with "item_price" > 0.<jupyter_code>train = train.query('item_price > 0')<jupyter_output><empty_output><jupyter_text>Data preprocessing* I'm dropping the text features since I won't be doing anything with them.* We are asked to predict total sales for every product and store in the next month, and our data is given by day, so let's remove unwanted columns and aggregate the data by month.<jupyter_code># Select only useful features.
train_monthly = lk_train[['date', 'date_block_num', 'shop_id', 'item_category_id', 'item_id', 'item_price', 'item_cnt_day']]
# Group by month in this case "date_block_num" and aggregate features.
train_monthly = train_monthly.sort_values('date').groupby(['date_block_num', 'shop_id', 'item_category_id', 'item_id'], as_index=False)
train_monthly = train_monthly.agg({'item_price':['sum', 'mean'], 'item_cnt_day':['sum', 'mean','count']})
# Rename features.
train_monthly.columns = ['date_block_num', 'shop_id', 'item_category_id', 'item_id', 'item_price', 'mean_item_price', 'item_cnt', 'mean_item_cnt', 'transactions']<jupyter_output><empty_output><jupyter_text>To mimic the real behavior of the data we have to create the missing records from the loaded dataset, so for each month we need to create the missing records for each shop and item, since we don't have data for them I'll replace them with 0.<jupyter_code># Build a data set with all the possible combinations of ['date_block_num','shop_id','item_id'] so we won't have missing records.
shop_ids = train_monthly['shop_id'].unique()
item_ids = train_monthly['item_id'].unique()
empty_df = []
for i in range(34):
    for shop in shop_ids:
        for item in item_ids:
            empty_df.append([i, shop, item])
    
empty_df = pd.DataFrame(empty_df, columns=['date_block_num','shop_id','item_id'])
# Merge the train set with the complete set (missing records will be filled with 0).
train_monthly = pd.merge(empty_df, train_monthly, on=['date_block_num','shop_id','item_id'], how='left')
train_monthly.fillna(0, inplace=True)<jupyter_output><empty_output><jupyter_text>New dataset<jupyter_code>train_monthly.head().T
train_monthly.describe().T
# Extract time based features.
train_monthly['year'] = train_monthly['date_block_num'].apply(lambda x: ((x//12) + 2013))
train_monthly['month'] = train_monthly['date_block_num'].apply(lambda x: (x % 12 + 1))
train_monthly['no-of-days-in-month'] = pd.to_datetime(train_monthly['month'],
                                                      format='%m').dt.days_in_month<jupyter_output><empty_output><jupyter_text>EDA<jupyter_code># Grouping data for EDA.
gp_month_mean = train_monthly.groupby(['month'], as_index=False)['item_cnt'].mean()
gp_month_sum = train_monthly.groupby(['month'], as_index=False)['item_cnt'].sum()
gp_category_mean = train_monthly.groupby(['item_category_id'], as_index=False)['item_cnt'].mean()
gp_category_sum = train_monthly.groupby(['item_category_id'], as_index=False)['item_cnt'].sum()
gp_shop_mean = train_monthly.groupby(['shop_id'], as_index=False)['item_cnt'].mean()
gp_shop_sum = train_monthly.groupby(['shop_id'], as_index=False)['item_cnt'].sum()<jupyter_output><empty_output><jupyter_text>How sales behaves along the year?<jupyter_code>f, axes = plt.subplots(2, 1, figsize=(22, 10), sharex=True)
sns.lineplot(x="month", y="item_cnt", data=gp_month_mean, ax=axes[0]).set_title("Monthly mean")
sns.lineplot(x="month", y="item_cnt", data=gp_month_sum, ax=axes[1]).set_title("Monthly sum")
plt.show()<jupyter_output><empty_output><jupyter_text>As we can see we have a trending increase of item sales count (mean) towards the ending of the year.  What category sells more?<jupyter_code>f, axes = plt.subplots(2, 1, figsize=(22, 10), sharex=True)
sns.barplot(x="item_category_id", y="item_cnt", data=gp_category_mean, ax=axes[0], palette="rocket").set_title("Monthly mean")
sns.barplot(x="item_category_id", y="item_cnt", data=gp_category_sum, ax=axes[1], palette="rocket").set_title("Monthly sum")
plt.show()<jupyter_output><empty_output><jupyter_text>Also only few of the categories seems to hold most of the sell count.  What shop sells more?<jupyter_code>f, axes = plt.subplots(2, 1, figsize=(22, 10), sharex=True)
sns.barplot(x="shop_id", y="item_cnt", data=gp_shop_mean, ax=axes[0], palette="rocket").set_title("Monthly mean")
sns.barplot(x="shop_id", y="item_cnt", data=gp_shop_sum, ax=axes[1], palette="rocket").set_title("Monthly sum")
plt.show()<jupyter_output><empty_output><jupyter_text>Most of the shops have a similar sell rate, but 3 of them have a much higher rate, this may be a indicative of the shop size.  Checking for outliers<jupyter_code>sns.jointplot(x="item_cnt", y="item_price", data=train_monthly, height=8)
plt.show()
sns.jointplot(x="item_cnt", y="transactions", data=train_monthly, height=8)
plt.show()<jupyter_output><empty_output><jupyter_text>Feature "item_cnt" distribution.<jupyter_code>plt.subplots(figsize=(22, 8))
sns.boxplot(train_monthly['item_cnt'])
plt.show()<jupyter_output><empty_output><jupyter_text>Removing outliersI'll treat "item_cnt" > 20 and = 400000 as outliers, so I'll remove them.<jupyter_code>train_monthly = train_monthly.query('item_cnt >= 0 and item_cnt <= 20 and item_price < 400000')<jupyter_output><empty_output><jupyter_text>Creating the labelOur label will be the "item_cnt" of the next month, as we are dealing with a forecast problem.<jupyter_code>train_monthly['item_cnt_month'] = train_monthly.sort_values('date_block_num').groupby(['shop_id', 'item_id'])['item_cnt'].shift(-1)<jupyter_output><empty_output><jupyter_text>Feature engineering Unitary item prices.<jupyter_code>train_monthly['item_price_unit'] = train_monthly['item_price'] // train_monthly['item_cnt']
train_monthly['item_price_unit'].fillna(0, inplace=True)<jupyter_output><empty_output><jupyter_text>Group based features.<jupyter_code>gp_item_price = train_monthly.sort_values('date_block_num').groupby(['item_id'], as_index=False).agg({'item_price':[np.min, np.max]})
gp_item_price.columns = ['item_id', 'hist_min_item_price', 'hist_max_item_price']

train_monthly = pd.merge(train_monthly, gp_item_price, on='item_id', how='left')<jupyter_output><empty_output><jupyter_text>How much each item's price changed from its (lowest/highest) historical price.<jupyter_code>train_monthly['price_increase'] = train_monthly['item_price'] - train_monthly['hist_min_item_price']
train_monthly['price_decrease'] = train_monthly['hist_max_item_price'] - train_monthly['item_price']<jupyter_output><empty_output><jupyter_text>Rolling window based features (window = 3 months).<jupyter_code># Min value
f_min = lambda x: x.rolling(window=3, min_periods=1).min()
# Max value
f_max = lambda x: x.rolling(window=3, min_periods=1).max()
# Mean value
f_mean = lambda x: x.rolling(window=3, min_periods=1).mean()
# Standard deviation
f_std = lambda x: x.rolling(window=3, min_periods=1).std()

function_list = [f_min, f_max, f_mean, f_std]
function_name = ['min', 'max', 'mean', 'std']

for i in range(len(function_list)):
    train_monthly[('item_cnt_%s' % function_name[i])] = train_monthly.sort_values('date_block_num').groupby(['shop_id', 'item_category_id', 'item_id'])['item_cnt'].apply(function_list[i])

# Fill the empty std features with 0
train_monthly['item_cnt_std'].fillna(0, inplace=True)<jupyter_output><empty_output><jupyter_text>Lag based features.<jupyter_code>lag_list = [1, 2, 3]

for lag in lag_list:
    ft_name = ('item_cnt_shifted%s' % lag)
    train_monthly[ft_name] = train_monthly.sort_values('date_block_num').groupby(['shop_id', 'item_category_id', 'item_id'])['item_cnt'].shift(lag)
    # Fill the empty shifted features with 0
    train_monthly[ft_name].fillna(0, inplace=True)<jupyter_output><empty_output><jupyter_text>Item sales count trend.<jupyter_code>train_monthly['item_trend'] = train_monthly['item_cnt']

for lag in lag_list:
    ft_name = ('item_cnt_shifted%s' % lag)
    train_monthly['item_trend'] -= train_monthly[ft_name]

train_monthly['item_trend'] /= len(lag_list) + 1
# Add seasonality features
train_monthly['year'] = train_monthly['date_block_num'].apply(lambda x: ((x//12) + 2013))
train_monthly['month'] = train_monthly['date_block_num'].apply(lambda x: (x % 12 + 1))
train_monthly['no-of-days-in-month'] = pd.to_datetime(train_monthly['month'], 
                                                          format='%m').dt.days_in_month<jupyter_output><empty_output><jupyter_text>Dataset after feature engineering<jupyter_code>train_monthly.to_csv("train_monthly.csv",index=False)
train_monthly.head().T
train_monthly.describe().T<jupyter_output><empty_output><jupyter_text>Train/validation split* As we know the test set in on the future, so we should try to simulate the same distribution on our train/validation split.* Our train set will be the first 3~28 blocks, validation will be last 5 blocks (29~32) and test will be block 33.* I'm leaving the first 3 months out because we use a 3 month window to generate features, so these first 3 month won't have really windowed useful features.<jupyter_code>train_set = train_monthly.query('date_block_num >= 3 and date_block_num < 27').copy()
validation_set = train_monthly.query('date_block_num >= 27 and date_block_num < 32').copy()
test_set = train_monthly.query('date_block_num == 33').copy()

train_set.dropna(subset=['item_cnt_month'], inplace=True)
validation_set.dropna(subset=['item_cnt_month'], inplace=True)
test_set.dropna(subset=['item_cnt_month'], inplace=True)

train_set.dropna(inplace=True)
validation_set.dropna(inplace=True)
test_set.dropna(inplace=True)

print('Train set records:', train_set.shape[0])
print('Validation set records:', validation_set.shape[0])
print('Test set records:', test_set.shape[0])

print('Train set records: %s (%.f%% of complete data)' % (train_set.shape[0], ((train_set.shape[0]/train_monthly.shape[0])*100)))
print('Validation set records: %s (%.f%% of complete data)' % (validation_set.shape[0], ((validation_set.shape[0]/train_monthly.shape[0])*100)))<jupyter_output>Train set records: 4748953
Validation set records: 989522
Test set records: 0
Train set records: 4748953 (71% of complete data)
Validation set records: 989522 (15% of complete data)<jupyter_text>Mean encoding.* done after the train/validation split.<jupyter_code># Shop mean encoding.
gp_shop_mean = train_set.groupby(['shop_id']).agg({'item_cnt_month': ['mean']})
gp_shop_mean.columns = ['shop_mean']
gp_shop_mean.reset_index(inplace=True)
# Item mean encoding.
gp_item_mean = train_set.groupby(['item_id']).agg({'item_cnt_month': ['mean']})
gp_item_mean.columns = ['item_mean']
gp_item_mean.reset_index(inplace=True)
# Shop with item mean encoding.
gp_shop_item_mean = train_set.groupby(['shop_id', 'item_id']).agg({'item_cnt_month': ['mean']})
gp_shop_item_mean.columns = ['shop_item_mean']
gp_shop_item_mean.reset_index(inplace=True)
# Year mean encoding.
gp_year_mean = train_set.groupby(['year']).agg({'item_cnt_month': ['mean']})
gp_year_mean.columns = ['year_mean']
gp_year_mean.reset_index(inplace=True)
# Month mean encoding.
gp_month_mean = train_set.groupby(['month']).agg({'item_cnt_month': ['mean']})
gp_month_mean.columns = ['month_mean']
gp_month_mean.reset_index(inplace=True)

# Add meand encoding features to train set.
train_set = pd.merge(train_set, gp_shop_mean, on=['shop_id'], how='left')
train_set = pd.merge(train_set, gp_item_mean, on=['item_id'], how='left')
train_set = pd.merge(train_set, gp_shop_item_mean, on=['shop_id', 'item_id'], how='left')
train_set = pd.merge(train_set, gp_year_mean, on=['year'], how='left')
train_set = pd.merge(train_set, gp_month_mean, on=['month'], how='left')
# Add meand encoding features to validation set.
validation_set = pd.merge(validation_set, gp_shop_mean, on=['shop_id'], how='left')
validation_set = pd.merge(validation_set, gp_item_mean, on=['item_id'], how='left')
validation_set = pd.merge(validation_set, gp_shop_item_mean, on=['shop_id', 'item_id'], how='left')
validation_set = pd.merge(validation_set, gp_year_mean, on=['year'], how='left')
validation_set = pd.merge(validation_set, gp_month_mean, on=['month'], how='left')

test_set = pd.merge(test_set, gp_shop_mean, on=['shop_id'], how='left')
test_set = pd.merge(test_set, gp_item_mean, on=['item_id'], how='left')
test_set = pd.merge(test_set, gp_shop_item_mean, on=['shop_id', 'item_id'], how='left')
test_set = pd.merge(test_set, gp_year_mean, on=['year'], how='left')
test_set = pd.merge(test_set, gp_month_mean, on=['month'], how='left')
# Create train and validation sets and labels. 
X_train = train_set.drop(['item_cnt_month', 'date_block_num'], axis=1)
Y_train = train_set['item_cnt_month'].astype(int)
X_validation = validation_set.drop(['item_cnt_month', 'date_block_num'], axis=1)
Y_validation = validation_set['item_cnt_month'].astype(int)
X_test = test_set.drop(['item_cnt_month', 'date_block_num'], axis=1)
Y_test = test_set['item_cnt_month'].astype(int)
# Integer features (used by catboost model).
int_features = ['shop_id', 'item_id', 'year', 'month','no-of-days-in-month']

X_train[int_features] = X_train[int_features].astype('int32')
X_validation[int_features] = X_validation[int_features].astype('int32')
X_test[int_features] = X_test[int_features].astype('int32')<jupyter_output><empty_output><jupyter_text>Build test setWe want to predict for "date_block_num" 34 so our test set will be block 33 and our predictions should reflect block 34 values. In other words we use block 33 because we want to forecast values for block 34.<jupyter_code>latest_records = pd.concat([train_set, validation_set]).drop_duplicates(subset=['shop_id', 'item_id'], keep='last')
X_test = pd.merge(test, latest_records, on=['shop_id', 'item_id'], how='left', suffixes=['', '_'])
X_test['year'] = 2015
X_test['month'] = 9
X_test.drop('item_cnt_month', axis=1, inplace=True)
X_test[int_features] = X_test[int_features].astype('int32')
X_test = X_test[X_train.columns]<jupyter_output><empty_output><jupyter_text>Replacing missing values.<jupyter_code>sets = [X_train, X_validation, X_test]

# This was taking too long.
# Replace missing values with the median of each item.
# for dataset in sets:
#     for item_id in dataset['item_id'].unique():
#         for column in dataset.columns:
#             item_median = dataset[(dataset['item_id'] == item_id)][column].median()
#             dataset.loc[(dataset[column].isnull()) & (dataset['item_id'] == item_id), column] = item_median

# Replace missing values with the median of each shop.            
for dataset in sets:
    for shop_id in dataset['shop_id'].unique():
        for column in dataset.columns:
            shop_median = dataset[(dataset['shop_id'] == shop_id)][column].median()
            dataset.loc[(dataset[column].isnull()) & (dataset['shop_id'] == shop_id), column] = shop_median
            
# Fill remaining missing values on test set with mean.
X_test.fillna(X_test.mean(), inplace=True)
# I'm dropping "item_category_id", we don't have it on test set and would be a little hard to create categories for items that exist only on test set.
X_train.drop(['item_category_id'], axis=1, inplace=True)
X_validation.drop(['item_category_id'], axis=1, inplace=True)
X_test.drop(['item_category_id'], axis=1, inplace=True)<jupyter_output><empty_output><jupyter_text>Test set<jupyter_code>X_train.to_csv("X_train.csv",index=False)
Y_train.to_csv("Y_train.csv",index=False)
X_validation.to_csv("X_validation.csv",index=False)
Y_validation.to_csv("Y_validation.csv",index=False)
X_test.to_csv("X_test.csv",index=False)
Y_test.to_csv("Y_test.csv",index=False)
X_test.describe().T<jupyter_output><empty_output><jupyter_text>Modeling the data Tree based models Catboost<jupyter_code>cat_features = [0, 1, 7, 8]

catboost_model = CatBoostRegressor(
    iterations=500,
    max_ctr_complexity=4,
    random_seed=0,
    od_type='Iter',
    od_wait=25,
    verbose=50,
    depth=4
)

catboost_model.fit(
    X_train, Y_train,
    cat_features=cat_features,
    eval_set=(X_validation, Y_validation)
)
print('Model params:', catboost_model.get_params())<jupyter_output><empty_output><jupyter_text>Catboost feature importance<jupyter_code>feature_score = pd.DataFrame(list(zip(X_train.dtypes.index, catboost_model.get_feature_importance(Pool(X_train, label=Y_train, cat_features=cat_features)))), columns=['Feature','Score'])
feature_score = feature_score.sort_values(by='Score', ascending=False, inplace=False, kind='quicksort', na_position='last')

plt.rcParams["figure.figsize"] = (19, 6)
ax = feature_score.plot('Feature', 'Score', kind='bar', color='c')
ax.set_title("Catboost Feature Importance Ranking", fontsize = 14)
ax.set_xlabel('')
rects = ax.patches
labels = feature_score['Score'].round(2)

for rect, label in zip(rects, labels):
    height = rect.get_height()
    ax.text(rect.get_x() + rect.get_width()/2, height + 0.35, label, ha='center', va='bottom')

plt.show()
catboost_train_pred = catboost_model.predict(X_train)
catboost_val_pred = catboost_model.predict(X_validation)
catboost_test_pred = catboost_model.predict(X_test)
print('Train rmse:', np.sqrt(mean_squared_error(Y_train, catboost_train_pred)))
print('Validation rmse:', np.sqrt(mean_squared_error(Y_validation, catboost_val_pred)))<jupyter_output><empty_output><jupyter_text>Let's see how catboost performed with this "prediction x label" plot.The closer the points are to the middle dashed line the better are the predictions.<jupyter_code>X_train= pd.read_csv("X_train.csv")
Y_train= pd.read_csv("Y_train.csv")
X_validation= pd.read_csv("X_validation.csv")
Y_validation= pd.read_csv("Y_validation.csv")
X_test= pd.read_csv("X_test.csv")
Y_test= pd.read_csv("Y_test.csv")
def model_performance_sc_plot(predictions, labels, title):
    # Get min and max values of the predictions and labels.
    min_val = max(max(predictions), max(labels))
    max_val = min(min(predictions), min(labels))
    # Create dataframe with predicitons and labels.
    performance_df = pd.DataFrame({"Label":labels})
    performance_df["Prediction"] = predictions
    # Plot data
    sns.jointplot(y="Label", x="Prediction", data=performance_df, kind="reg", height=7)
    plt.plot([min_val, max_val], [min_val, max_val], 'm--')
    plt.title(title, fontsize=9)
    plt.show()
    
# model_performance_sc_plot(catboost_train_pred, Y_train, 'Train')
model_performance_sc_plot(catboost_val_pred, Y_validation, 'Validation')<jupyter_output><empty_output><jupyter_text>XGBoost<jupyter_code># Use only part of features on XGBoost.
xgb_features = ['item_cnt','item_cnt_mean', 'item_cnt_std', 'item_cnt_shifted1', 
                'item_cnt_shifted2', 'item_cnt_shifted3', 'shop_mean', 
                'shop_item_mean', 'item_trend', 'mean_item_cnt']
xgb_train = X_train[xgb_features]
xgb_val = X_validation[xgb_features]
xgb_test = X_test[xgb_features]
xgb_model = XGBRegressor(max_depth=8, 
                         n_estimators=500, 
                         min_child_weight=1000,  
                         colsample_bytree=0.7, 
                         subsample=0.7, 
                         eta=0.3, 
                         seed=0,
                        tree_method = 'exact')
xgb_model.fit(X_train, 
              Y_train, 
              eval_metric="rmse", 
              eval_set=[(X_train, Y_train), (X_validation, Y_validation)], 
              verbose=20, 
              early_stopping_rounds=10)
xgb_model.best_estimator_<jupyter_output><empty_output><jupyter_text>XGBoost feature importance<jupyter_code>plt.rcParams["figure.figsize"] = (15, 6)
plot_importance(xgb_model)
plt.show()
xgb_train_pred = xgb_model.predict(xgb_train)
xgb_val_pred = xgb_model.predict(xgb_val)
xgb_test_pred = xgb_model.predict(xgb_test)
print('Train rmse:', np.sqrt(mean_squared_error(Y_train, xgb_train_pred)))
print('Validation rmse:', np.sqrt(mean_squared_error(Y_validation, xgb_val_pred)))<jupyter_output><empty_output><jupyter_text>Let's see how XGBoosting performed with this "prediction x label" plot.<jupyter_code># model_performance_sc_plot(xgb_train_pred, Y_train, 'Train')
model_performance_sc_plot(xgb_val_pred, Y_validation, 'Validation')<jupyter_output><empty_output><jupyter_text>Random forest<jupyter_code>rf_model = RandomForestRegressor(n_estimators=50, max_depth=7, random_state=0, n_jobs=-1)
rf_model.fit(xgb_train, Y_train)
rf_train_pred = rf_model.predict(rf_train)
rf_val_pred = rf_model.predict(rf_val)
rf_test_pred = rf_model.predict(rf_test)
print('Train rmse:', np.sqrt(mean_squared_error(Y_train, rf_train_pred)))
print('Validation rmse:', np.sqrt(mean_squared_error(Y_validation, rf_val_pred)))<jupyter_output><empty_output><jupyter_text>Let's see how random forest performed with this "prediction x label" plot.<jupyter_code># model_performance_sc_plot(rf_train_pred, Y_train, 'Train')
model_performance_sc_plot(rf_val_pred, Y_validation, 'Validation')<jupyter_output><empty_output><jupyter_text>Linear models Linear Regression<jupyter_code># Use only part of features on linear Regression.
lr_features = ['item_cnt', 'item_cnt_shifted1', 'item_trend', 'mean_item_cnt', 'shop_mean']
lr_train = X_train[lr_features]
lr_val = X_validation[lr_features]
lr_test = X_test[lr_features]<jupyter_output><empty_output><jupyter_text>Normalizing features<jupyter_code>lr_scaler = MinMaxScaler()
lr_scaler.fit(lr_train)
lr_train = lr_scaler.transform(lr_train)
lr_val = lr_scaler.transform(lr_val)
lr_test = lr_scaler.transform(lr_test)
lr_model = LinearRegression(n_jobs=-1)
lr_model.fit(lr_train, Y_train)
lr_train_pred = lr_model.predict(lr_train)
lr_val_pred = lr_model.predict(lr_val)
lr_test_pred = lr_model.predict(lr_test)
print('Train rmse:', np.sqrt(mean_squared_error(Y_train, lr_train_pred)))
print('Validation rmse:', np.sqrt(mean_squared_error(Y_validation, lr_val_pred)))<jupyter_output><empty_output><jupyter_text>Let's see how linear regression performed with this "prediction x label" plot.<jupyter_code># model_performance_sc_plot(lr_train_pred, Y_train, 'Train')
model_performance_sc_plot(lr_val_pred, Y_validation, 'Validation')<jupyter_output><empty_output><jupyter_text>Clustering models KNN Regressor<jupyter_code># Use only part of features on KNN.
knn_features = ['item_cnt', 'item_cnt_mean', 'item_cnt_std', 'item_cnt_shifted1',
                'item_cnt_shifted2', 'shop_mean', 'shop_item_mean', 
                'item_trend', 'mean_item_cnt']

# Subsample train set (using the whole data was taking too long).
X_train_sampled = X_train[:100000]
Y_train_sampled = Y_train[:100000]

knn_train = X_train_sampled[knn_features]
knn_val = X_validation[knn_features]
knn_test = X_test[knn_features]<jupyter_output><empty_output><jupyter_text>Normalizing features<jupyter_code>knn_scaler = MinMaxScaler()
knn_scaler.fit(knn_train)
knn_train = knn_scaler.transform(knn_train)
knn_val = knn_scaler.transform(knn_val)
knn_test = knn_scaler.transform(knn_test)
knn_model = KNeighborsRegressor(n_neighbors=9, leaf_size=13, n_jobs=-1)
knn_model.fit(knn_train, Y_train_sampled)
knn_train_pred = knn_model.predict(knn_train)
knn_val_pred = knn_model.predict(knn_val)
knn_test_pred = knn_model.predict(knn_test)
print('Train rmse:', np.sqrt(mean_squared_error(Y_train_sampled, knn_train_pred)))
print('Validation rmse:', np.sqrt(mean_squared_error(Y_validation, knn_val_pred)))<jupyter_output><empty_output><jupyter_text>Let's see how knn performed with this "prediction x label" plot.<jupyter_code># model_performance_sc_plot(knn_train_pred, Y_train_sampled, 'Train')
model_performance_sc_plot(knn_val_pred, Y_validation, 'Validation')<jupyter_output><empty_output><jupyter_text>Create new datasets with the predictions from first level models.* Here I'll be using a simple ensembling technique, I'll use the 1st level models predictions as the input for the 2nd level model, this way the 2nd level model will basically use the 1st level models predictions as features and learn where to give more weight.* To use this technique I also need to use the 1st level models and make predictions on the test set, so I can use them on the 2nd level model.* I could also pass the complete validation set with extra features (the 1st level models prediction) to the 2nd level model and let it do a little more work on finding the solution.<jupyter_code># Dataset that will be the train set of the ensemble model.
first_level = pd.DataFrame(catboost_val_pred, columns=['catboost'])
first_level['xgbm'] = xgb_val_pred
first_level['random_forest'] = rf_val_pred
first_level['linear_regression'] = lr_val_pred
first_level['knn'] = knn_val_pred
first_level['label'] = Y_validation.values
first_level.head(20)
# Dataset that will be the test set of the ensemble model.
first_level_test = pd.DataFrame(catboost_test_pred, columns=['catboost'])
first_level_test['xgbm'] = xgb_test_pred
first_level_test['random_forest'] = rf_test_pred
first_level_test['linear_regression'] = lr_test_pred
first_level_test['knn'] = knn_test_pred
first_level_test.head()<jupyter_output><empty_output><jupyter_text>Ensembling* To combine the 1st level model predictions, I'll use a simple linear regression.* As I'm only feeding the model with predictions I don't need a complex model. Ensemble architecture:* 1st level:    * Catboost    * XGBM    * Random forest    * Linear Regression    * KNN* 2nd level;    * Linear Regression     Here is an  image to help the understanding    2nd level model as a linear regression* This is the model that will combine the other ones to hopefully make an overall better prediction.* If the inputs to this mode were more complex, could be a good idea to split the data into train and validation again, this way you can check if the metal model is overfitting.<jupyter_code>meta_model = LinearRegression(n_jobs=-1)<jupyter_output><empty_output><jupyter_text>Trained on validation set using the 1st level models predictions as features.<jupyter_code># Drop label from dataset.
first_level.drop('label', axis=1, inplace=True)
meta_model.fit(first_level, Y_validation)<jupyter_output><empty_output><jupyter_text>Make predictions on test set using the 1st level models predictions as features.<jupyter_code>ensemble_pred = meta_model.predict(first_level)
final_predictions = meta_model.predict(first_level_test)<jupyter_output><empty_output><jupyter_text>Ensemble model metrics on validation set.<jupyter_code>print('Train rmse:', np.sqrt(mean_squared_error(ensemble_pred, Y_validation)))<jupyter_output><empty_output><jupyter_text>Let's see how the meta model performed with this "prediction x label" plot.<jupyter_code>model_performance_sc_plot(ensemble_pred, Y_validation, 'Validation')<jupyter_output><empty_output><jupyter_text>Output dataframe.<jupyter_code>prediction_df = pd.DataFrame(test['ID'], columns=['ID'])
prediction_df['item_cnt_month'] = final_predictions.clip(0., 20.)
prediction_df.to_csv('submission.csv', index=False)
prediction_df.head(10)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.12754020289004142, 'fraction_numerical': 0.011090452425220993, 'mean_word_length': 7.010974653775803, 'pattern_counts': {'":': 1, '<': 206, '<?xml version=': 0, '>': 208, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '51', 'dir': 'jupyter-structured-clean-dedup', 'id': '228664', 'license': 'MIT', 'n_tokens_mistral': 10941, 'n_tokens_neox': 10483, 'n_words': 2277, 'path': 'Code/.ipynb_checkpoints/run_model-checkpoint.ipynb', 'repo_name': 'vartikatewari/KaggleTimeseriesPredictSales', 'size': '30657'}