<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Kaggle Restaurant Revenue Prediction
Predict annual restaurant sales based on objective measurements
Import libraries and data; explore the data
Let's begin by importing the Python libraries and data that we'll need
Step1:  Now for a bit of exploratory data analysis so we can get to know our data
Step2:  Plot the data
Step3:  I'm sure there are more creative and informative ways to plot the data, but for now it's time to move on.
Massage, munge, preprocess, and visualize the data
Step4:  Find the relevant features
Step5:  Try various machine learning algorithms
Let's try on some algorithms and see how they fit and predict
Step6:  How to format the data for the Kaggle contest submission based on the sampleSubmission.csv file
Step7:  sklearn.neighbors.KNeighborsRegressor
Step8:  sklearn.linear_model.LinearRegression
Step9:  sklearn.neural_network.MLPClassifier
Step10:  Restaurant Revenue Prediction Kaggle solution
This next exercise is from a blog post (linked above) by Bikash Agrawal.
Step11:  Concatenate the train and test data together into a single dataframe to pre-process and featurize both consistently
Step12:  Now is the time for us to impute values for the rare restaurant types (DT and MB).
Instead of trying to predict with values that appear only 1 or 0 times in the training set, we will replace them with one of the other commonly appearing categories by fitting a model that predicts which common category they "should" be.
Step13:  Here we can define and train a model to impute restaurant type.
The grid below just has a range of values that the author has found to work well with random forest type models (of which ExtraTrees is one).
Step14:  Now we can binarize the "P" columns with dummy variables
Step15:  To finish up our data preprocessing, we need to scale all input features to between 0 and 1 (this is especially important for KNN or SVM(SVR) models.
However, we don't want to scale the output, so we'll temporarily 'drop' it.
Step16:  Now we can define and train a Ridge Regression model.
The author tested others from the sklearn library, including SVR, RandomForest, K-nearest Neighbors, but found that Ridge consistently gave the strongest leaderboard results.
One takeaway -- when the training data is small, simplest is often best.
Step17:  So, now we're ready for our final submission to Kaggle
Step18:  One last quick comparison.
Note the x-axis scale change
<ASSISTANT_TASK:>
Python Code:
import pandas as pd
import numpy as np
import sklearn
import matplotlib.pyplot as plt
import seaborn as sns
from IPython.display import display
%matplotlib inline
train_data = pd.read_csv("train.csv")
train_data = train_data.drop('Id', axis=1)
test_data = pd.read_csv("test.csv")
test_data = test_data.drop('Id', axis=1)
Explanation: Kaggle Restaurant Revenue Prediction
Predict annual restaurant sales based on objective measurements
Import libraries and data; explore the data
Let's begin by importing the Python libraries and data that we'll need:
End of explanation
display(train_data[:10])
display(test_data[:10])
train_data.describe()
test_data.describe()
train_data.head()
train_data.tail()
train_data.sample(5)
train_data.keys()
test_data.keys()
test_data.keys()
Explanation: Now for a bit of exploratory data analysis so we can get to know our data:
End of explanation
feature_columns = train_data[['P1', 'P2', 'P3', 'P4',
       'P5', 'P6', 'P7', 'P8', 'P9', 'P10', 'P11', 'P12', 'P13', 'P14', 'P15',
       'P16', 'P17', 'P18', 'P19', 'P20', 'P21', 'P22', 'P23', 'P24', 'P25',
       'P26', 'P27', 'P28', 'P29', 'P30', 'P31', 'P32', 'P33', 'P34', 'P35',
       'P36', 'P37']]
feature_columns.plot.box(figsize=(20, 20))
Explanation: Plot the data
End of explanation
# Cribbed from https://www.kaggle.com/ani310/restaurant-revenue-prediction/restaurant-revenue
# Format the data so that dates are easier to work with.
# Create a column that contains data about the number of days the restaurant has been open.
# Remove the column that has the restaurant's opening date.
train_data['Open Date'] = pd.to_datetime(train_data['Open Date'], format='%m/%d/%Y')
test_data['Open Date'] = pd.to_datetime(test_data['Open Date'], format='%m/%d/%Y')
train_data['OpenDays'] = ""
test_data['OpenDays'] = ""
date_last_train = pd.DataFrame({'Date':np.repeat(['01/01/2015'], [len(train_data)])})
date_last_test = pd.DataFrame({'Date':np.repeat(['01/01/2015'], [len(test_data)])})
date_last_train['Date'] = pd.to_datetime(date_last_train['Date'], format='%m/%d/%Y')
date_last_test['Date'] = pd.to_datetime(date_last_test['Date'], format='%m/%d/%Y')
train_data['OpenDays'] = date_last_train['Date'] - train_data['Open Date']
test_data['OpenDays'] = date_last_test['Date'] - test_data['Open Date']
train_data['OpenDays'] = train_data['OpenDays'].astype('timedelta64[D]').astype(int)
test_data['OpenDays'] = test_data['OpenDays'].astype('timedelta64[D]').astype(int)
train_data = train_data.drop('Open Date', axis=1)
test_data = test_data.drop('Open Date', axis=1)
# Compare the revenue generated by the restaurants in Big Cities vs Other:
city_perc = train_data [["City Group", "revenue"]].groupby(['City Group'], as_index=False).mean()
sns.barplot(x='City Group', y='revenue', data=city_perc)
plt.title("Revenue by city size")
# Convert data from 'City Group' and create columns of indicator variables for 'Big Cities' or 'Other':
city_group_dummy = pd.get_dummies(train_data['City Group'])
train_data = train_data.join(city_group_dummy)
city_group_dummy_test = pd.get_dummies(test_data['City Group'])
test_data = test_data.join(city_group_dummy_test)
train_data = train_data.drop('City Group', axis=1)
test_data = test_data.drop('City Group', axis=1)
# Create scatterplot showing how long a restaurant has been open impacts revenue.
# This will also show any outliers.
plt.scatter(train_data['OpenDays'], train_data['revenue'])
plt.xlabel("Days Open")
plt.ylabel("Revenue")
plt.title("Restaurant revenue by location age")
Explanation: I'm sure there are more creative and informative ways to plot the data, but for now it's time to move on.
Massage, munge, preprocess, and visualize the data
End of explanation
from sklearn.feature_selection import SelectFromModel
# from sklearn.linear_model import LassoCV
from sklearn.ensemble import ExtraTreesClassifier
X_train = train_data.iloc[:, 2:]
y = train_data['revenue']
print("X_train.shape: {}".format(X_train.shape))
clf = ExtraTreesClassifier()
clf = clf.fit(X_train, y)
print("clf.feature_.importances_: \n{}".format(clf.feature_importances_))
model = SelectFromModel(clf, prefit=True)
print(model)
X_train_new = model.transform(X_train)
print("X_train_new.shape: {}".format(X_train_new.shape))
X_train_new = pd.DataFrame(X_train_new)
print(X_train_new[:5])
Explanation: Find the relevant features
End of explanation
from sklearn.ensemble import RandomForestRegressor
# Tweak seaborn visualizations and adapt to Jupyter notebooks:
sns.set_context("notebook", font_scale=1.1)
sns.set_style("ticks")
# Make dataframes for train and test:
X_train = pd.DataFrame({'OpenDaysLog':train_data['OpenDays'].apply(np.log),
                       'Big Cities':train_data['Big Cities'], 'Other':train_data['Other'],
                       'P2':train_data['P2'], 'P8':train_data['P8'], 'P22':train_data['P22'],
                       'P24':train_data['P24'], 'P28':train_data['P28'], 'P26':train_data['P26']})
y_train = train_data['revenue'].apply(np.log)
X_test = pd.DataFrame({'OpenDaysLog':test_data['OpenDays'].apply(np.log),
                       'Big Cities':test_data['Big Cities'], 'Other':test_data['Other'],
                       'P2':test_data['P2'], 'P8':test_data['P8'], 'P22':test_data['P22'],
                       'P24':test_data['P24'], 'P28':test_data['P28'], 'P26':test_data['P26']})
# Time to build the models and make some predictions:
from sklearn import linear_model
cls = RandomForestRegressor(n_estimators=150)
cls.fit(X_train, y_train)
pred = cls.predict(X_test)
pred = np.exp(pred)
pred
cls.score(X_train, y_train)
Explanation: Try various machine learning algorithms
Let's try on some algorithms and see how they fit and predict:
sklearn.ensemble.RandomForestRegressor
Also cribbed from Kaggle.
End of explanation
test_data = pd.read_csv("test.csv")
submission = pd.DataFrame({
    "Id": test_data["Id"],
    "Prediction": pred
})
# submission.to_csv('RandomForestSimple.csv', header=True, index=False)
Explanation: How to format the data for the Kaggle contest submission based on the sampleSubmission.csv file:
End of explanation
from sklearn.neighbors import KNeighborsRegressor
# Use dataframes from sklearn.ensemble.RandomForestRegressor example above.
knn_cls = KNeighborsRegressor(n_neighbors=2)
knn_cls.fit(X_train, y_train)
knn_pred = knn_cls.predict(X_test)
knn_pred = np.exp(knn_pred)
knn_cls.score(X_train, y_train)
Explanation: sklearn.neighbors.KNeighborsRegressor
End of explanation
from sklearn.linear_model import LinearRegression
# Use dataframes from sklearn.ensemble.RandomForestRegressor example above.
lr_cls = LinearRegression()
lr_cls.fit(X_train, y_train)
lr_pred = lr_cls.predict(X_test)
lr_pred = np.exp(lr_pred)
lr_cls.score(X_train, y_train)
Explanation: sklearn.linear_model.LinearRegression
End of explanation
from sklearn.neural_network import MLPRegressor
mlp_cls = MLPRegressor(solver='lbfgs')
mlp_cls.fit(X_train, y_train)
mlp_pred = mlp_cls.predict(X_test)
mlp_pred = np.exp(mlp_pred)
mlp_cls.score(X_train, y_train)
Explanation: sklearn.neural_network.MLPClassifier
End of explanation
import datetime
%pylab inline
from sklearn.model_selection import LeaveOneOut
from sklearn.grid_search import GridSearchCV, RandomizedSearchCV
from sklearn.metrics import mean_squared_error, mean_absolute_error
from sklearn.preprocessing import MinMaxScaler
from sklearn.ensemble import ExtraTreesClassifier 
# Regressors considered:
from sklearn.svm import SVR
from sklearn.neighbors import KNeighborsRegressor
# Regressor chosen by the author for final submission:
from sklearn.linear_model import Ridge
# Kaggle added ~311.5 "fake" data points to the test for each real data point.
# Dividing by this number gives more accurate counts of the "real" data in the test set.
FAKE_DATA_RATIO = 311.5
# Set a random seed:
SEED = 0
# Read in the data provided by Kaggle:
train = pd.read_csv('train.csv', index_col=0, parse_dates=[1])
test = pd.read_csv('test.csv', index_col=0, parse_dates=[1])
print("Training data dimensions: \n{}".format(train.shape))
print("Test data dimensions: \n{}".format(test.shape))
Explanation: Restaurant Revenue Prediction Kaggle solution
This next exercise is from a blog post (linked above) by Bikash Agrawal.
End of explanation
df = pd.concat((test, train), ignore_index=True)
df.describe()
# Convert date strings to "days open" numerical value:
df["Open Date"] = df["Open Date"].apply(pd.to_datetime)
last_date = df["Open Date"].max()
# Create a datetime delta object:
df["Open Date"] = last_date - df["Open Date"]
# Convert the delta object to an int:
df["Open Date"] = df["Open Date"].dt.days + 1
# Scale "days since opening" so that the marginal impact decreases over time.
# This and the similar log transform of City Count below are the modifications
# that were not in the official competition submission.
df["Log Days Opened"] = df["Open Date"].apply(np.log)
df = df.drop(["Open Date"], axis=1)
# Resize plots:
pylab.rcParams['figure.figsize'] = (8, 6)
df[["Log Days Opened", "revenue"]].plot(x="Log Days Opened", y="revenue",
                                       kind='scatter', title="Log (Days Opened) vs Revenue")
# There is a certain set of columns that are either all zero or all non-zero.
# We have added a feature to mark this -- the 'zeros' feature will be 17 for
# these rows and 0 or 1 for the rows which are rarely or never zero.
# Here are the features with the notable zero behavior:
zero_cols = ['P14', 'P15', 'P16', 'P17', 'P18', 'P24', 'P25', 'P26', 'P27',
             'P30', 'P31', 'P32', 'P33', 'P34', 'P35', 'P36', 'P37']
# We make a feature that holds this count of zero columns in the above list:
df['zeros'] = (df[zero_cols] == 0).sum(1)
pylab.rcParams['figure.figsize'] = (20, 8)
fig, axs = plt.subplots(1,2)
fig.suptitle("Distribution of new Zeros features:", fontsize=18)
# There is only one row with a zero count between 0 and 17 in the training set:
df['zeros'].ix[pd.notnull(df.revenue)].value_counts().plot(
                title="Training Set", kind='bar', ax=axs[0])
# In the test set, however, there are many rows with an intermediate count of zeros.
# This is probably an artifact of how the fake test data was generated, and might
# indicate that conditional dependence between columns was not preserved.
df['zeros'].ix[pd.isnull(df.revenue)].value_counts().plot(
                title="Test Set", kind='bar', ax=axs[1], color='red')
# Here we convert two categorical variables, "Restaurant Type", and "City
# Group (Size)" to dummy variables:
pylab.rcParams['figure.figsize'] = (6, 4)
# The two categories of City Group both appear very frequently:
train["City Group"].value_counts().plot(
    title="City Group Distribution in the Training Set", kind='bar')
# Two of the four Restaurant Types (DT and MB) are very rare:
train["Type"].value_counts().plot(
    title="Restaurant Type Distribution in the Training Set", kind='bar')
(test["Type"].value_counts() / FAKE_DATA_RATIO).plot(
    title="Approximate Restaurant Type Distribution in True Test Set",
    kind='bar', color='red')
df = df.join(pd.get_dummies(df['City Group'], prefix="CG"))
df = df.join(pd.get_dummies(df['Type'], prefix="T"))
# Since only n-1 columns are needed to binarize n categories, drop one
# of the new columns and drop the original columns.
# In addition, drop the rare restaurant types.
df = df.drop(["City Group", "Type", "CG_Other", "T_MB", "T_DT"], axis=1)
print(df.shape)
df.describe(include='all')
# Replace city names with the count of their frequency in the training +
# estimated frequency in the test set.
city_counts = (test["City"].value_counts() /
              FAKE_DATA_RATIO).add(train["City"].value_counts(), fill_value=0)
df["City"] = df["City"].replace(city_counts)
print("Some example estimated counts of restaurants per city: \n{}".format(
                                            city_counts.head()))
# Take the natural logarithm of city count so that the marginal effect decreases:
df["Log City Count"] = df["City"].apply(np.log)
df = df.drop(["City"], axis=1)
# The last vertical spread of points below are restaurants in Istanbul.
pylab.rcParams['figure.figsize'] = (8, 6)
df[["Log City Count", "revenue"]].plot(x="Log City Count", y="revenue",
                kind='scatter', title="Log City Count vs Revenue")
Explanation: Concatenate the train and test data together into a single dataframe to pre-process and featurize both consistently:
End of explanation
# tofit are the rows in the training set that belong to one of the common restaurant types:
tofit = df.ix[((df.T_FC==1) | (df.T_IL==1)) & (pd.notnull(df.revenue))]
# tofill are rows in either train or test that belong to one of the rare types:
tofill = df.ix[((df.T_FC==0) & (df.T_IL==0))]
print("Type training set shape: \n{}".format(tofit.shape))
print("Data to impute: \n{}".format(tofill.shape))
# Restaurants with type FC are labeled 1, those with type IL are labeled 0.
y = tofit.T_FC
# Drop the label columns and revenue (which is not in the test set):
X = tofit.drop(["T_FC", "T_IL", "revenue"], axis=1)
Explanation: Now is the time for us to impute values for the rare restaurant types (DT and MB).
Instead of trying to predict with values that appear only 1 or 0 times in the training set, we will replace them with one of the other commonly appearing categories by fitting a model that predicts which common category they "should" be.
End of explanation
model_grid = {'max_depth': [None, 8], 'min_samples_split': [4,9,16],
             'min_samples_leaf': [1,4], 'max_features': ['sqrt', 0.5, None]}
type_model = ExtraTreesClassifier(n_estimators=25, random_state=SEED)
grid = RandomizedSearchCV(type_model, model_grid, n_iter=10, cv=5, scoring="roc_auc")
grid.fit(X, y)
print("Best parameters for Type Model: \n{}".format(grid.best_params_))
type_model.set_params(**grid.best_params_)
type_model.fit(X, y)
imputations = type_model.predict(tofill.drop(["T_FC", "T_IL", "revenue"], axis=1))
df.loc[(df.T_FC==0) & (df.T_IL==0), "T_FC"] = imputations
df = df.drop(["T_IL"], axis=1)
df[:7]
print("% labeled FC in the training set: \n{}".format(df.T_FC.mean()))
print("% of imputed values labeled FC: \n{}".format(np.mean(imputations)))
Explanation: Here we can define and train a model to impute restaurant type.
The grid below just has a range of values that the author has found to work well with random forest type models (of which ExtraTrees is one).
End of explanation
print("Pre-binarizing columns: {}".format(len(df.columns)))
for col in df.columns:
    if col[0] == 'P':
        print(col, len(df[col].unique()), "Unique Values")
        df = df.join(pd.get_dummies(df[col], prefix=col))
        df = df.drop([col, df.columns[-1]], axis=1)
print("Post-binarizing columns: {}".format(len(df.columns)))
Explanation: Now we can binarize the "P" columns with dummy variables:
End of explanation
min_max_scaler = MinMaxScaler()
rev = df.revenue
df = df.drop(['revenue'], axis=1)
df = pd.DataFrame(data=min_max_scaler.fit_transform(df), columns=df.columns, index=df.index)
df = df.join(rev)
# Now that preprocessing is finished, let's have a look at the data before modeling with it:
df.describe()
# Recover the original train and test rows based on revenue (which is null for test rows)
train = df.ix[pd.notnull(df.revenue)]
test = df.ix[pd.isnull(df.revenue)].drop(['revenue'], axis=1)
# Scale revenue by sqrt.
# The reason is to decrease the influence of the few very large revenue values.
y = train.revenue.apply(np.sqrt)
X = train.drop(["revenue"], axis=1)
Explanation: To finish up our data preprocessing, we need to scale all input features to between 0 and 1 (this is especially important for KNN or SVM(SVR) models.
However, we don't want to scale the output, so we'll temporarily 'drop' it.
End of explanation
model_grid = [{'normalize': [True, False], 'alpha': np.logspace(0,10)}]
model = Ridge()
# Use a grid search and leave-one-out CV on the train set to find the best regularization parameter to use.
grid = GridSearchCV(model, model_grid, scoring='neg_mean_squared_error') 
grid.fit(X, y)
print("Best parameters set found on development set: \n{}".format(
                                               grid.best_params_))
# Retrain model on the full training set using the best parameters found in the last step:
model.set_params(**grid.best_params_)
model.fit(X, y)
# Predict on the test set using the trained model:
submission = pd.DataFrame(columns=['Prediction'], index=test.index,
                         data=model.predict(test))
# Convert back to revenue from sqrt(revenue):
submission.Prediction = submission.Prediction.apply(np.square)
submission.Prediction[:7]
Explanation: Now we can define and train a Ridge Regression model.
The author tested others from the sklearn library, including SVR, RandomForest, K-nearest Neighbors, but found that Ridge consistently gave the strongest leaderboard results.
One takeaway -- when the training data is small, simplest is often best.
End of explanation
# Add required column name for Kaggle's submission parser:
submission.index.name='Id'
# Write out the submission:
# submission.to_csv("TFI_Ridge.csv")
# Quick sanity check on the submission:
submission.describe().astype(int)
# Revenue from training set for comparison:
train[['revenue']].describe().astype(int)
Explanation: So, now we're ready for our final submission to Kaggle:
End of explanation
train[['revenue']].plot(kind='kde', title="Training Set Revenue Distribution")
submission.columns = ["predicted revenue"]
submission.plot(kind='kde', title="Prediction Revenue Distribution", color='red')
Explanation: One last quick comparison.
Note the x-axis scale change: the predictions are more conservative and tend to be closer to the mean than the real revenues.
This is pretty standard behavior when using RMSE -- there are big penalties for being very wrong, so the model will tend towards more moderate predictions.
End of explanation
<END_TASK>