{"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/trend).**\n\n---\n","metadata":{}},{"cell_type":"markdown","source":"\n# Introduction #","metadata":{}},{"cell_type":"markdown","source":"Run 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.ex2 import *\n\n# Setup notebook\nfrom pathlib import Path\nfrom learntools.time_series.style import *  # plot style settings\n\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport numpy as np\nimport seaborn as sns\nfrom sklearn.linear_model import LinearRegression\n\ndata_dir = Path('../input/ts-course-data/')\ncomp_dir = Path('../input/store-sales-time-series-forecasting')\n\nretail_sales = pd.read_csv(\n    data_dir / \"us-retail-sales.csv\",\n    parse_dates=['Month'],\n    index_col='Month',\n).to_period('D')\nfood_sales = retail_sales.loc[:, 'FoodAndBeverage']\nauto_sales = retail_sales.loc[:, 'Automobiles']\n\ndtype = {\n    'store_nbr': 'category',\n    'family': 'category',\n    'sales': 'float32',\n    'onpromotion': 'uint64',\n}\nstore_sales = pd.read_csv(\n    comp_dir / 'train.csv',\n    dtype=dtype,\n    parse_dates=['date'],\n    infer_datetime_format=True,\n)\nstore_sales = store_sales.set_index('date').to_period('D')\nstore_sales = store_sales.set_index(['store_nbr', 'family'], append=True)\naverage_sales = store_sales.groupby('date').mean()['sales']","metadata":{"lines_to_next_cell":0,"execution":{"iopub.status.busy":"2022-07-28T01:54:32.218932Z","iopub.execute_input":"2022-07-28T01:54:32.219388Z","iopub.status.idle":"2022-07-28T01:54:41.232033Z","shell.execute_reply.started":"2022-07-28T01:54:32.219294Z","shell.execute_reply":"2022-07-28T01:54:41.230962Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"-------------------------------------------------------------------------------","metadata":{}},{"cell_type":"markdown","source":"# 1) Determine trend with a moving average plot\n\nThe *US Retail Sales* dataset contains monthly sales data for a number of retail industries in the United States. Run the next cell to see a plot of the *Food and Beverage* series.","metadata":{}},{"cell_type":"code","source":"ax = food_sales.plot(**plot_params)\nax.set(title=\"US Food and Beverage Sales\", ylabel=\"Millions of Dollars\");","metadata":{"execution":{"iopub.status.busy":"2022-07-28T01:54:41.233977Z","iopub.execute_input":"2022-07-28T01:54:41.234416Z","iopub.status.idle":"2022-07-28T01:54:41.589022Z","shell.execute_reply.started":"2022-07-28T01:54:41.234384Z","shell.execute_reply":"2022-07-28T01:54:41.587964Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Now make a moving average plot to estimate the trend for this series.","metadata":{}},{"cell_type":"code","source":"# YOUR CODE HERE: Add methods to `food_sales` to compute a moving\n# average with appropriate parameters for trend estimation.\ntrend = food_sales.rolling(\n    window=12,\n    center=True,\n    min_periods=6,\n).mean()\n\n\n# Check your answer\nq_1.check()\n\n# Make a plot\nax = food_sales.plot(**plot_params, alpha=0.5)\nax = trend.plot(ax=ax, linewidth=3)","metadata":{"lines_to_next_cell":0,"execution":{"iopub.status.busy":"2022-07-28T01:57:46.608360Z","iopub.execute_input":"2022-07-28T01:57:46.608782Z","iopub.status.idle":"2022-07-28T01:57:46.929762Z","shell.execute_reply.started":"2022-07-28T01:57:46.608748Z","shell.execute_reply":"2022-07-28T01:57:46.928529Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Uncomment to get a hint or solution\n#q_1.hint()\n# q_1.solution()","metadata":{"lines_to_next_cell":0,"execution":{"iopub.status.busy":"2022-07-28T01:57:58.726788Z","iopub.execute_input":"2022-07-28T01:57:58.727205Z","iopub.status.idle":"2022-07-28T01:57:58.731766Z","shell.execute_reply.started":"2022-07-28T01:57:58.727163Z","shell.execute_reply":"2022-07-28T01:57:58.730777Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"-------------------------------------------------------------------------------\n\n# 2) Identify trend\n\nWhat order polynomial trend might be appropriate for the *Food and Beverage Sales* series? Can you think of a non-polynomial curve that might work even better?\n\nOnce you've thought about it, run this cell for some discussion.","metadata":{}},{"cell_type":"code","source":"# View the solution (Run this cell to receive credit!)\nq_2.check()","metadata":{"execution":{"iopub.status.busy":"2022-07-28T01:58:01.167914Z","iopub.execute_input":"2022-07-28T01:58:01.168331Z","iopub.status.idle":"2022-07-28T01:58:01.177589Z","shell.execute_reply.started":"2022-07-28T01:58:01.168297Z","shell.execute_reply":"2022-07-28T01:58:01.176435Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"-------------------------------------------------------------------------------\n\nWe'll continue using the time series of average sales in this lesson. Run this cell to see a moving average plot of `average_sales` estimating the trend.","metadata":{}},{"cell_type":"code","source":"trend = average_sales.rolling(\n    window=365,\n    center=True,\n    min_periods=183,\n).mean()\n\nax = average_sales.plot(**plot_params, alpha=0.5)\nax = trend.plot(ax=ax, linewidth=3)","metadata":{"execution":{"iopub.status.busy":"2022-07-28T01:54:41.925769Z","iopub.execute_input":"2022-07-28T01:54:41.926556Z","iopub.status.idle":"2022-07-28T01:54:42.218913Z","shell.execute_reply.started":"2022-07-28T01:54:41.926510Z","shell.execute_reply":"2022-07-28T01:54:42.217807Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 3) Create a Trend Feature\n\nUse `DeterministicProcess` to create a feature set for a cubic trend model. Also create features for a 90-day forecast.","metadata":{}},{"cell_type":"code","source":"from statsmodels.tsa.deterministic import DeterministicProcess\n\ny = average_sales.copy()  # the target\n\n# YOUR CODE HERE: Instantiate `DeterministicProcess` with arguments\n# appropriate for a cubic trend model\ndp = DeterministicProcess(index=y.index, order=3)\n\n# YOUR CODE HERE: Create the feature set for the dates given in y.index\nX = dp.in_sample()\n\n# YOUR CODE HERE: Create features for a 90-day forecast.\nX_fore = dp.out_of_sample(steps=90)\n\n\n# Check your answer\nq_3.check()","metadata":{"execution":{"iopub.status.busy":"2022-07-28T01:59:25.148806Z","iopub.execute_input":"2022-07-28T01:59:25.149230Z","iopub.status.idle":"2022-07-28T01:59:25.166753Z","shell.execute_reply.started":"2022-07-28T01:59:25.149196Z","shell.execute_reply":"2022-07-28T01:59:25.165586Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Lines below will give you a hint or solution code\n#q_3.hint()\n# q_3.solution()","metadata":{"execution":{"iopub.status.busy":"2022-07-28T01:59:30.296465Z","iopub.execute_input":"2022-07-28T01:59:30.296876Z","iopub.status.idle":"2022-07-28T01:59:30.302362Z","shell.execute_reply.started":"2022-07-28T01:59:30.296844Z","shell.execute_reply":"2022-07-28T01:59:30.301072Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"You can see the a plot of the result by running the next cell.","metadata":{}},{"cell_type":"code","source":"model = LinearRegression()\nmodel.fit(X, y)\n\ny_pred = pd.Series(model.predict(X), index=X.index)\ny_fore = pd.Series(model.predict(X_fore), index=X_fore.index)\n\nax = y.plot(**plot_params, alpha=0.5, title=\"Average Sales\", ylabel=\"items sold\")\nax = y_pred.plot(ax=ax, linewidth=3, label=\"Trend\", color='C0')\nax = y_fore.plot(ax=ax, linewidth=3, label=\"Trend Forecast\", color='C3')\nax.legend();","metadata":{"execution":{"iopub.status.busy":"2022-07-28T01:59:33.174292Z","iopub.execute_input":"2022-07-28T01:59:33.174785Z","iopub.status.idle":"2022-07-28T01:59:33.661034Z","shell.execute_reply.started":"2022-07-28T01:59:33.174744Z","shell.execute_reply":"2022-07-28T01:59:33.659905Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"--------------------------------------------------------------------------------\n\nOne way to fit more complicated trends is to increase the order of the polynomial you use. To get a better fit to the somewhat complicated trend in *Store Sales*, we could try using an order 11 polynomial.","metadata":{}},{"cell_type":"code","source":"from statsmodels.tsa.deterministic import DeterministicProcess\n\ndp = DeterministicProcess(index=y.index, order=11)\nX = dp.in_sample()\n\nmodel = LinearRegression()\nmodel.fit(X, y)\n\ny_pred = pd.Series(model.predict(X), index=X.index)\n\nax = y.plot(**plot_params, alpha=0.5, title=\"Average Sales\", ylabel=\"items sold\")\nax = y_pred.plot(ax=ax, linewidth=3, label=\"Trend\", color='C0')\nax.legend();","metadata":{"execution":{"iopub.status.busy":"2022-07-28T01:59:42.181788Z","iopub.execute_input":"2022-07-28T01:59:42.182193Z","iopub.status.idle":"2022-07-28T01:59:42.596588Z","shell.execute_reply.started":"2022-07-28T01:59:42.182160Z","shell.execute_reply":"2022-07-28T01:59:42.595482Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 4) Understand risks of forecasting with high-order polynomials\n\nHigh-order polynomials are generally not well-suited to forecasting, however. Can you guess why?","metadata":{}},{"cell_type":"code","source":"# View the solution (Run this cell to receive credit!)\nq_4.check()","metadata":{"execution":{"iopub.status.busy":"2022-07-28T01:59:55.747075Z","iopub.execute_input":"2022-07-28T01:59:55.747505Z","iopub.status.idle":"2022-07-28T01:59:55.759668Z","shell.execute_reply.started":"2022-07-28T01:59:55.747473Z","shell.execute_reply":"2022-07-28T01:59:55.758198Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Uncomment the next line for a hint\nq_4.hint()","metadata":{"execution":{"iopub.status.busy":"2022-07-28T01:59:53.348827Z","iopub.execute_input":"2022-07-28T01:59:53.349253Z","iopub.status.idle":"2022-07-28T01:59:53.358468Z","shell.execute_reply.started":"2022-07-28T01:59:53.349211Z","shell.execute_reply":"2022-07-28T01:59:53.357599Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Run this cell to see the same 90-day forecast using an order 11 polynomial. Does it confirm your intuition?","metadata":{}},{"cell_type":"code","source":"X_fore = dp.out_of_sample(steps=90)\ny_fore = pd.Series(model.predict(X_fore), index=X_fore.index)\n\nax = y.plot(**plot_params, alpha=0.5, title=\"Average Sales\", ylabel=\"items sold\")\nax = y_pred.plot(ax=ax, linewidth=3, label=\"Trend\", color='C0')\nax = y_fore.plot(ax=ax, linewidth=3, label=\"Trend Forecast\", color='C3')\nax.legend();","metadata":{"execution":{"iopub.status.busy":"2022-07-28T02:00:02.163539Z","iopub.execute_input":"2022-07-28T02:00:02.163935Z","iopub.status.idle":"2022-07-28T02:00:02.693206Z","shell.execute_reply.started":"2022-07-28T02:00:02.163904Z","shell.execute_reply":"2022-07-28T02:00:02.692266Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"--------------------------------------------------------------------------------\n\n# (Optional) Fit trend with splines\n\n*Splines* are a nice alternative to polynomials when you want to fit a trend. The *Multivariate Adaptive Regression Splines* (MARS) algorithm in the `pyearth` library is powerful and easy to use. There are a lot of hyperparameters you may want to investigate.","metadata":{}},{"cell_type":"code","source":"from pyearth import Earth\n\n# Target and features are the same as before\ny = average_sales.copy()\ndp = DeterministicProcess(index=y.index, order=1)\nX = dp.in_sample()\n\n# Fit a MARS model with `Earth`\nmodel = Earth()\nmodel.fit(X, y)\n\ny_pred = pd.Series(model.predict(X), index=X.index)\n\nax = y.plot(**plot_params, title=\"Average Sales\", ylabel=\"items sold\")\nax = y_pred.plot(ax=ax, linewidth=3, label=\"Trend\")","metadata":{"execution":{"iopub.status.busy":"2022-07-28T02:00:11.890499Z","iopub.execute_input":"2022-07-28T02:00:11.890941Z","iopub.status.idle":"2022-07-28T02:00:18.075028Z","shell.execute_reply.started":"2022-07-28T02:00:11.890906Z","shell.execute_reply":"2022-07-28T02:00:18.073915Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Forecasting complicated trends like this will typically be difficult (if not impossible). With historical data, however, you can use splines to isolate other patterns in a time series by *detrending*.","metadata":{}},{"cell_type":"code","source":"y_detrended = y - y_pred   # remove the trend from store_sales\n\ny_detrended.plot(**plot_params, title=\"Detrended Average Sales\");","metadata":{"execution":{"iopub.status.busy":"2022-07-28T02:00:40.015672Z","iopub.execute_input":"2022-07-28T02:00:40.016108Z","iopub.status.idle":"2022-07-28T02:00:40.229702Z","shell.execute_reply.started":"2022-07-28T02:00:40.016073Z","shell.execute_reply":"2022-07-28T02:00:40.228470Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Keep Going #\n\n[**Model seasonality**](https://www.kaggle.com/ryanholbrook/seasonality), another common type of time dependence, with indicators and Fourier features.","metadata":{}},{"cell_type":"markdown","source":"https://machinelearningmastery.com/xgboost-for-time-series-forecasting/","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":{}}]}