{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"**This notebook is an exercise in the [Time Series](https://www.kaggle.com/learn/time-series) course.  You can reference the tutorial at [this link](https://www.kaggle.com/ryanholbrook/time-series-as-features).**\n\n---\n","metadata":{}},{"cell_type":"markdown","source":"# Introduction #\n\nRun this cell to set everything up!","metadata":{}},{"cell_type":"code","source":"# Setup feedback system\nfrom learntools.core import binder\nbinder.bind(globals())\nfrom learntools.time_series.ex4 import *\n\n# Setup notebook\nfrom pathlib import Path\nfrom learntools.time_series.style import *  # plot style settings\nfrom learntools.time_series.utils import plot_lags, make_lags, make_leads\n\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport numpy as np\nimport seaborn as sns\nfrom sklearn.linear_model import LinearRegression\nfrom sklearn.metrics import mean_squared_log_error\nfrom statsmodels.graphics.tsaplots import plot_pacf\nfrom statsmodels.tsa.deterministic import CalendarFourier, DeterministicProcess\n\n\ncomp_dir = Path('../input/store-sales-time-series-forecasting')\n\nstore_sales = pd.read_csv(\n    comp_dir / 'train.csv',\n    usecols=['store_nbr', 'family', 'date', 'sales', 'onpromotion'],\n    dtype={\n        'store_nbr': 'category',\n        'family': 'category',\n        'sales': 'float32',\n        'onpromotion': 'uint32',\n    },\n    parse_dates=['date'],\n    infer_datetime_format=True,\n)\nstore_sales['date'] = store_sales.date.dt.to_period('D')\nstore_sales = store_sales.set_index(['store_nbr', 'family', 'date']).sort_index()\n\nfamily_sales = (\n    store_sales\n    .groupby(['family', 'date'])\n    .mean() \n    .unstack('family')\n    .loc['2017', ['sales', 'onpromotion']]\n)","metadata":{"execution":{"iopub.status.busy":"2022-07-28T16:18:42.412071Z","iopub.execute_input":"2022-07-28T16:18:42.412526Z","iopub.status.idle":"2022-07-28T16:19:00.106378Z","shell.execute_reply.started":"2022-07-28T16:18:42.412434Z","shell.execute_reply":"2022-07-28T16:19:00.105220Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"----------------------------------------------------------------------------\n\nNot every product family has sales showing cyclic behavior, and neither does the series of average sales. Sales of school and office supplies, however, show patterns of growth and decay not well characterized by trend or seasons. In this question and the next, you'll model cycles in sales of school and office supplies using lag features.\n\nTrend and seasonality will both create serial dependence that shows up in correlograms and lag plots. To isolate any purely *cyclic* behavior, we'll start by deseasonalizing the series. Use the code in the next cell to deseasonalize *Supply Sales*. We'll store the result in a variable `y_deseason`.","metadata":{}},{"cell_type":"code","source":"supply_sales = family_sales.loc(axis=1)[:, 'SCHOOL AND OFFICE SUPPLIES']\ny = supply_sales.loc[:, 'sales'].squeeze()\n\nfourier = CalendarFourier(freq='M', order=4)\ndp = DeterministicProcess(\n    constant=True,\n    index=y.index,\n    order=1,\n    seasonal=True,\n    drop=True,\n    additional_terms=[fourier],\n)\nX_time = dp.in_sample()\nX_time['NewYearsDay'] = (X_time.index.dayofyear == 1)\n\nmodel = LinearRegression(fit_intercept=False)\nmodel.fit(X_time, y)\ny_deseason = y - model.predict(X_time)\ny_deseason.name = 'sales_deseasoned'\n\nax = y_deseason.plot()\nax.set_title(\"Sales of School and Office Supplies (deseasonalized)\");","metadata":{"execution":{"iopub.status.busy":"2022-07-28T16:19:00.108954Z","iopub.execute_input":"2022-07-28T16:19:00.109403Z","iopub.status.idle":"2022-07-28T16:19:00.532286Z","shell.execute_reply.started":"2022-07-28T16:19:00.109359Z","shell.execute_reply":"2022-07-28T16:19:00.531145Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Does this deseasonalized series show cyclic patterns? To confirm our intuition, we can try to isolate cyclic behavior using a moving-average plot just like we did with trend. The idea is to choose a window long enough to smooth over short-term seasonality, but short enough to still preserve the cycles.\n\n# 1) Plotting cycles\n\nCreate a seven-day moving average from `y`, the series of supply sales. Use a centered window, but don't set the `min_periods` argument.","metadata":{}},{"cell_type":"code","source":"# YOUR CODE HERE\ny_ma = y.rolling(7, center=True).mean()\n\n\n# Plot\nax = y_ma.plot()\nax.set_title(\"Seven-Day Moving Average\");\n\n# Check your answer\nq_1.check()","metadata":{"execution":{"iopub.status.busy":"2022-07-28T16:22:59.873897Z","iopub.execute_input":"2022-07-28T16:22:59.874312Z","iopub.status.idle":"2022-07-28T16:23:00.201713Z","shell.execute_reply.started":"2022-07-28T16:22:59.874280Z","shell.execute_reply":"2022-07-28T16:23:00.200670Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Lines below will give you a hint or solution code\n#q_1.hint()\n# q_1.solution()","metadata":{"execution":{"iopub.status.busy":"2022-07-28T16:23:08.945192Z","iopub.execute_input":"2022-07-28T16:23:08.945623Z","iopub.status.idle":"2022-07-28T16:23:08.950643Z","shell.execute_reply.started":"2022-07-28T16:23:08.945590Z","shell.execute_reply":"2022-07-28T16:23:08.949288Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Do you see how the moving average plot resembles the plot of the deseasonalized series? In both, we can see cyclic behavior indicated.\n\n-------------------------------------------------------------------------------\n\nLet's examine our deseasonalized series for serial dependence. Take a look at the partial autocorrelation correlogram and lag plot.","metadata":{}},{"cell_type":"code","source":"plot_pacf(y_deseason, lags=8);\nplot_lags(y_deseason, lags=8, nrows=2);","metadata":{"execution":{"iopub.status.busy":"2022-07-28T16:23:10.651883Z","iopub.execute_input":"2022-07-28T16:23:10.652287Z","iopub.status.idle":"2022-07-28T16:23:12.039161Z","shell.execute_reply.started":"2022-07-28T16:23:10.652255Z","shell.execute_reply":"2022-07-28T16:23:12.038015Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 2) Examine serial dependence in *Store Sales*\n\nAre any of the lags significant according to the correlogram? Does the lag plot suggest any relationships that weren't apparent from the correlogram?\n\nAfter you've thought about your answer, run the next cell.","metadata":{}},{"cell_type":"code","source":"# View the solution (Run this cell to receive credit!)\nq_2.check()","metadata":{"execution":{"iopub.status.busy":"2022-07-28T16:23:38.253826Z","iopub.execute_input":"2022-07-28T16:23:38.254237Z","iopub.status.idle":"2022-07-28T16:23:38.263690Z","shell.execute_reply.started":"2022-07-28T16:23:38.254207Z","shell.execute_reply":"2022-07-28T16:23:38.262686Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"-------------------------------------------------------------------------------\n\nRecall from the tutorial that a *leading indicator* is a series whose values at one time can be used to predict the target at a future time -- a leading indicator provides \"advance notice\" of changes in the target.\n\nThe competition dataset includes a time series that could potentially be useful as a leading indicator -- the `onpromotion` series, which contains the number of items on a special promotion that day. Since the company itself decides when to do a promotion, there's no worry about \"lookahead leakage\"; we could use Tuesday's `onpromotion` value to forecast sales on Monday, for instance.\n\nUse the next cell to examine leading and lagging values for `onpromotion` plotted against sales of school and office supplies.","metadata":{}},{"cell_type":"code","source":"onpromotion = supply_sales.loc[:, 'onpromotion'].squeeze().rename('onpromotion')\n\n# Drop days without promotions\nplot_lags(x=onpromotion.loc[onpromotion > 1], y=y_deseason.loc[onpromotion > 1], lags=3, leads=3, nrows=1);","metadata":{"execution":{"iopub.status.busy":"2022-07-28T16:24:37.820162Z","iopub.execute_input":"2022-07-28T16:24:37.820826Z","iopub.status.idle":"2022-07-28T16:24:38.876410Z","shell.execute_reply.started":"2022-07-28T16:24:37.820788Z","shell.execute_reply":"2022-07-28T16:24:38.875118Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 3) Examine time series features\n\nDoes it appear that either leading or lagging values of `onpromotion` could be useful as a feature?","metadata":{}},{"cell_type":"code","source":"q_3.check()","metadata":{"execution":{"iopub.status.busy":"2022-07-28T16:25:33.381294Z","iopub.execute_input":"2022-07-28T16:25:33.382298Z","iopub.status.idle":"2022-07-28T16:25:33.391475Z","shell.execute_reply.started":"2022-07-28T16:25:33.382260Z","shell.execute_reply":"2022-07-28T16:25:33.390362Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"-------------------------------------------------------------------------------\n\n# 4) Create time series features\n\nCreate the features indicated in the solution to Question 3. If no features from that series would be useful, use an empty dataframe `pd.DataFrame()` as your answer.","metadata":{}},{"cell_type":"code","source":"# YOUR CODE HERE: Make features from `y_deseason`\nX_lags = make_lags(y_deseason, lags=1)\n\n# YOUR CODE HERE: Make features from `onpromotion`\n# You may want to use `pd.concat`\nX_promo = pd.concat([\n    make_lags(onpromotion, lags=1),\n    onpromotion,\n    make_leads(onpromotion, leads=1),\n], axis=1)\n\nX = pd.concat([X_lags, X_promo], axis=1)\ny, X = y.align(X, join='inner')\n\n# Check your answer\nq_4.check()","metadata":{"execution":{"iopub.status.busy":"2022-07-28T16:29:08.995598Z","iopub.execute_input":"2022-07-28T16:29:08.995995Z","iopub.status.idle":"2022-07-28T16:29:09.014757Z","shell.execute_reply.started":"2022-07-28T16:29:08.995964Z","shell.execute_reply":"2022-07-28T16:29:09.013628Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Lines below will give you a hint or solution code\n#q_4.hint()\n# q_4.solution()","metadata":{"execution":{"iopub.status.busy":"2022-07-28T16:28:53.523829Z","iopub.execute_input":"2022-07-28T16:28:53.524248Z","iopub.status.idle":"2022-07-28T16:28:53.533577Z","shell.execute_reply.started":"2022-07-28T16:28:53.524216Z","shell.execute_reply":"2022-07-28T16:28:53.532425Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Use the code in the next cell if you'd like to see predictions from the resulting model.","metadata":{}},{"cell_type":"code","source":"# from sklearn.model_selection import train_test_split\n\n# X_train, X_valid, y_train, y_valid = train_test_split(X, y, test_size=30, shuffle=False)\n\n# model = LinearRegression(fit_intercept=False).fit(X_train, y_train)\n# y_fit = pd.Series(model.predict(X_train), index=X_train.index).clip(0.0)\n# y_pred = pd.Series(model.predict(X_valid), index=X_valid.index).clip(0.0)\n\n# rmsle_train = mean_squared_log_error(y_train, y_fit) ** 0.5\n# rmsle_valid = mean_squared_log_error(y_valid, y_pred) ** 0.5\n# print(f'Training RMSLE: {rmsle_train:.5f}')\n# print(f'Validation RMSLE: {rmsle_valid:.5f}')\n\n# ax = y.plot(**plot_params, alpha=0.5, title=\"Average Sales\", ylabel=\"items sold\")\n# ax = y_fit.plot(ax=ax, label=\"Fitted\", color='C0')\n# ax = y_pred.plot(ax=ax, label=\"Forecast\", color='C3')\n# ax.legend();","metadata":{"execution":{"iopub.status.busy":"2022-07-28T16:29:25.504418Z","iopub.execute_input":"2022-07-28T16:29:25.505893Z","iopub.status.idle":"2022-07-28T16:29:25.513108Z","shell.execute_reply.started":"2022-07-28T16:29:25.505841Z","shell.execute_reply":"2022-07-28T16:29:25.511512Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"-------------------------------------------------------------------------------\n\nWinners of Kaggle forecasting competitions have often included moving averages and other rolling statistics in their feature sets. Such features seem to be especially useful when used with GBDT algorithms like XGBoost.\n\nIn Lesson 2 you learned how to compute moving averages to estimate trends. Computing rolling statistics to be used as features is similar except we need to take care to avoid lookahead leakage. First, the result should be set at the right end of the window instead of the center -- that is, we should use `center=False` (the default) in the `rolling` method. Second, the target should be lagged a step.","metadata":{}},{"cell_type":"markdown","source":"# 5) Create statistical features\n\nEdit the code in the next cell to create the following features:\n- 14-day rolling median (`median`) of lagged target\n- 7-day rolling standard deviation (`std`) of lagged target\n- 7-day sum (`sum`) of items \"on promotion\", with centered window","metadata":{}},{"cell_type":"code","source":"y_lag = supply_sales.loc[:, 'sales'].shift(1)\nonpromo = supply_sales.loc[:, 'onpromotion']\n\n# 28-day mean of lagged target\nmean_7 = y_lag.rolling(7).mean()\n# YOUR CODE HERE: 14-day median of lagged target\nmedian_14 = y_lag.rolling(14).median()\n# YOUR CODE HERE: 7-day rolling standard deviation of lagged target\nstd_7 = y_lag.rolling(7).std()\n# YOUR CODE HERE: 7-day sum of promotions with centered window\npromo_7 = onpromo.rolling(7, center=True).sum()\n\n\n# Check your answer\nq_5.check()","metadata":{"execution":{"iopub.status.busy":"2022-07-28T16:30:53.883068Z","iopub.execute_input":"2022-07-28T16:30:53.883454Z","iopub.status.idle":"2022-07-28T16:30:53.908869Z","shell.execute_reply.started":"2022-07-28T16:30:53.883423Z","shell.execute_reply":"2022-07-28T16:30:53.908085Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Lines below will give you a hint or solution code\n#q_5.hint()\n# q_5.solution()","metadata":{"execution":{"iopub.status.busy":"2022-07-28T16:30:56.812887Z","iopub.execute_input":"2022-07-28T16:30:56.813307Z","iopub.status.idle":"2022-07-28T16:30:56.817932Z","shell.execute_reply.started":"2022-07-28T16:30:56.813274Z","shell.execute_reply":"2022-07-28T16:30:56.816780Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Check out the Pandas [`Window` documentation](https://pandas.pydata.org/pandas-docs/stable/reference/window.html) for more statistics you can compute. Also try \"exponential weighted\" windows by using `ewm` in place of `rolling`; exponential decay is often a more realistic representation of how effects propagate over time.","metadata":{}},{"cell_type":"markdown","source":"# Keep Going #\n\n[**Create hybrid forecasters**](https://www.kaggle.com/ryanholbrook/hybrid-models) and combine the strengths of two machine learning algorithms.","metadata":{}},{"cell_type":"markdown","source":"---\n\n\n\n\n*Have questions or comments? Visit the [course discussion forum](https://www.kaggle.com/learn/time-series/discussion) to chat with other learners.*","metadata":{}}]}