{"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":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n# It is defined by the kaggle/python Docker image: https://github.com/kaggle/docker-python\n# For example, here's several helpful packages to load\n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\n\n# Input data files are available in the read-only \"../input/\" directory\n# For example, running this (by clicking run or pressing Shift+Enter) will list all files under the input directory\n\nimport os\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n for filename in filenames:\n print(os.path.join(dirname, filename))\n\n# You can write up to 20GB to the current directory (/kaggle/working/) that gets preserved as output when you create a version using \"Save & Run All\" \n# You can also write temporary files to /kaggle/temp/, but they won't be saved outside of the current session","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2022-08-09T08:24:31.052135Z","iopub.execute_input":"2022-08-09T08:24:31.052472Z","iopub.status.idle":"2022-08-09T08:24:31.076526Z","shell.execute_reply.started":"2022-08-09T08:24:31.052438Z","shell.execute_reply":"2022-08-09T08:24:31.075824Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\n\nimport statsmodels.api as sm\nfrom statsmodels.tsa.deterministic import DeterministicProcess,CalendarFourier\nfrom statsmodels.graphics.tsaplots import plot_pacf,plot_acf\nfrom statsmodels.tsa.stattools import adfuller\n\nfrom matplotlib import pyplot as plt, style\nstyle.use('seaborn-darkgrid')\nimport seaborn as sns\nsns.set_style('darkgrid')\nimport plotly.express as px\nfrom tqdm import tqdm_notebook\n\nimport gc\ngc.enable()\nfrom warnings import filterwarnings, simplefilter\nfilterwarnings('ignore')\nsimplefilter('ignore')","metadata":{"execution":{"iopub.status.busy":"2022-08-09T08:24:32.931535Z","iopub.execute_input":"2022-08-09T08:24:32.932683Z","iopub.status.idle":"2022-08-09T08:24:36.114781Z","shell.execute_reply.started":"2022-08-09T08:24:32.932635Z","shell.execute_reply":"2022-08-09T08:24:36.113802Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Favorita stores forecasting\n\n\n\n\n* Place Ecuador\n* 54 Stores\n\n\n* 33 Family products\n* 1,782 Distinct pairs for modeling!","metadata":{}},{"cell_type":"code","source":"df_train=pd.read_csv(\"/kaggle/input/store-sales-time-series-forecasting/train.csv\",index_col=0,parse_dates=[\"date\"],infer_datetime_format=True)\ndf_train.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:35:21.586034Z","iopub.execute_input":"2022-08-08T22:35:21.586584Z","iopub.status.idle":"2022-08-08T22:35:28.436693Z","shell.execute_reply.started":"2022-08-08T22:35:21.586542Z","shell.execute_reply":"2022-08-08T22:35:28.435090Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_train[\"date\"]=df_train.date.dt.to_period(\"D\")\ndf_train=df_train.set_index(\"date\").to_timestamp().reset_index()\ndf_train=df_train.set_index([\"date\",\"store_nbr\",\"family\"]).sort_index()\ndf_train.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:35:41.216691Z","iopub.execute_input":"2022-08-08T22:35:41.217170Z","iopub.status.idle":"2022-08-08T22:35:43.708207Z","shell.execute_reply.started":"2022-08-08T22:35:41.217136Z","shell.execute_reply":"2022-08-08T22:35:43.706642Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Example MultiIndex selection","metadata":{}},{"cell_type":"code","source":"pd.set_option('display.max_columns', 500)\nfrom itertools import product\ndate=df_train.reset_index()[\"date\"].unique()\nstore=df_train.reset_index()[\"store_nbr\"].unique()\nparameters = product(date, store)\nparameters_list = list(parameters)\nprint(parameters_list[100])\ndf_train.loc[parameters_list[100],[\"sales\"]].head()\n","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:36:20.015099Z","iopub.execute_input":"2022-08-08T22:36:20.015683Z","iopub.status.idle":"2022-08-08T22:36:20.439533Z","shell.execute_reply.started":"2022-08-08T22:36:20.015643Z","shell.execute_reply":"2022-08-08T22:36:20.438144Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Total family product sales per store","metadata":{}},{"cell_type":"code","source":"df_date_str=df_train.groupby([\"date\",\"store_nbr\"]).sales.sum().reset_index()\npx.line(df_date_str.sort_values([\"store_nbr\", \"date\"]), x='date', y='sales', color='store_nbr',title = \"Total Sales\" )","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:36:47.533775Z","iopub.execute_input":"2022-08-08T22:36:47.534240Z","iopub.status.idle":"2022-08-08T22:36:51.661475Z","shell.execute_reply.started":"2022-08-08T22:36:47.534207Z","shell.execute_reply":"2022-08-08T22:36:51.659774Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Summary**\n\nIn every store there was a similar pattern till **june 2015**, from that moment on sales were more stable. Maybe we should use that stable time frame for model development.\n\n**Closed**:\n\n**12** didn't work from 30.03.2015-28.05.2015. **14** didn't work from 4.08.2014-10.09.2014. **18** didn't work from 15.08.2016-2.12.2016. **24** didn't work from 14.04.2014-23.07.2014. **25** didn't work from 22.04.2016-26.10.2016.\n\n**Opened**:\n\n**20** opened 12.02.2015. **21** opened 23.7.2015. **22** opened 8.10. 2015. **29** opened 19.03.2015. **36** opened 8.05.2013. **42** opened 20.8.2015. **52** opened 19.4.2017 god. **53** opened 28.3.2014\n\nWe also noticed that there is a trend in total sales in almost every store it could be people are bying more family products, maybe more items are introduced in each store as time went on or stores had more and more customers over time. We will inspect this latter.","metadata":{}},{"cell_type":"code","source":"df_test=pd.read_csv(\"/kaggle/input/store-sales-time-series-forecasting/test.csv\",parse_dates=[\"date\"],infer_datetime_format=True)\ndf_test[\"date\"]=df_test.date.dt.to_period(\"D\")\ndf_test=df_test.set_index(\"date\").to_timestamp().reset_index()\ndf_test=df_test.set_index([\"date\",\"store_nbr\",\"family\"]).sort_index()\ndf_test.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:37:57.836981Z","iopub.execute_input":"2022-08-08T22:37:57.837496Z","iopub.status.idle":"2022-08-08T22:37:57.920295Z","shell.execute_reply.started":"2022-08-08T22:37:57.837457Z","shell.execute_reply":"2022-08-08T22:37:57.918913Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_oil = pd.read_csv(\"/kaggle/input/store-sales-time-series-forecasting/oil.csv\",parse_dates=[\"date\"],infer_datetime_format=True,\n index_col = 'date').to_period('D')\ndf_oil.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-09T08:24:58.088532Z","iopub.execute_input":"2022-08-09T08:24:58.088872Z","iopub.status.idle":"2022-08-09T08:24:58.129700Z","shell.execute_reply.started":"2022-08-09T08:24:58.088838Z","shell.execute_reply":"2022-08-09T08:24:58.128685Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Changing to a more suitable column name \"oil_price\"","metadata":{}},{"cell_type":"code","source":"df_oil=df_oil.rename(columns={\"dcoilwtico\":\"oil_price\"})\ndf_oil.isna().sum()","metadata":{"execution":{"iopub.status.busy":"2022-08-09T08:25:02.856685Z","iopub.execute_input":"2022-08-09T08:25:02.856963Z","iopub.status.idle":"2022-08-09T08:25:02.867351Z","shell.execute_reply.started":"2022-08-09T08:25:02.856934Z","shell.execute_reply":"2022-08-09T08:25:02.866602Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"43 missing values we had to umpute using Intepolation and first value using \"bfill\" copy of the next value.\n\nOil price time series before imputation\n","metadata":{}},{"cell_type":"code","source":"df_oil.plot(figsize=(15,8))","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:39:35.649376Z","iopub.execute_input":"2022-08-08T22:39:35.649885Z","iopub.status.idle":"2022-08-08T22:39:36.061688Z","shell.execute_reply.started":"2022-08-08T22:39:35.649845Z","shell.execute_reply":"2022-08-08T22:39:36.060431Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_oil[\"oil_price\"]=df_oil.oil_price.interpolate(method=\"linear\")\ndf_oil[\"oil_price\"]= df_oil.oil_price.fillna(method ='bfill')\ndf_oil.isna().sum()","metadata":{"execution":{"iopub.status.busy":"2022-08-09T08:25:07.723927Z","iopub.execute_input":"2022-08-09T08:25:07.724298Z","iopub.status.idle":"2022-08-09T08:25:07.744310Z","shell.execute_reply.started":"2022-08-09T08:25:07.724255Z","shell.execute_reply":"2022-08-09T08:25:07.743615Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_oil.plot(figsize=(15,8))","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:40:02.785005Z","iopub.execute_input":"2022-08-08T22:40:02.785881Z","iopub.status.idle":"2022-08-08T22:40:03.214018Z","shell.execute_reply.started":"2022-08-08T22:40:02.785830Z","shell.execute_reply":"2022-08-08T22:40:03.212551Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**ADF** stationary test for Oil TS","metadata":{}},{"cell_type":"code","source":"def adfuller_test(time_series):\n result=adfuller(time_series)\n labels = ['ADF Test Statistic','p-value','#Lags Used','Number of Observations Used']\n for value,label in zip(result,labels):\n print(label+' : '+str(value) )\n if result[1] <= 0.05:\n print(\"Time series has no unit root and is stationary\")\n else:\n print(\"Time series has a unit root, indicating it is non-stationary \")","metadata":{"execution":{"iopub.status.busy":"2022-08-09T08:25:14.317039Z","iopub.execute_input":"2022-08-09T08:25:14.317653Z","iopub.status.idle":"2022-08-09T08:25:14.324226Z","shell.execute_reply.started":"2022-08-09T08:25:14.317612Z","shell.execute_reply":"2022-08-09T08:25:14.323029Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"adfuller_test(df_oil.oil_price)","metadata":{"execution":{"iopub.status.busy":"2022-08-09T08:25:15.380020Z","iopub.execute_input":"2022-08-09T08:25:15.380726Z","iopub.status.idle":"2022-08-09T08:25:15.477327Z","shell.execute_reply.started":"2022-08-09T08:25:15.380682Z","shell.execute_reply":"2022-08-09T08:25:15.476122Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"As we can see Oil TS is not statinary maybe we can make TS stationary by doing first lag differencing ","metadata":{}},{"cell_type":"code","source":"df_oil_diff=df_oil.diff()\nadfuller_test(df_oil_diff.oil_price[1:])","metadata":{"execution":{"iopub.status.busy":"2022-08-09T08:25:20.491651Z","iopub.execute_input":"2022-08-09T08:25:20.491962Z","iopub.status.idle":"2022-08-09T08:25:20.583241Z","shell.execute_reply.started":"2022-08-09T08:25:20.491933Z","shell.execute_reply":"2022-08-09T08:25:20.581999Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Great now that Oil TS is statinary we can plot partial_autocorelation and see if maybe there is AR(q) proces-(auto regresion by order of q).\n\nBut first lets plot Differenced Oil TS and see how it looks","metadata":{}},{"cell_type":"code","source":"df_oil_diff.plot(figsize=(15,8))","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:41:15.207394Z","iopub.execute_input":"2022-08-08T22:41:15.207787Z","iopub.status.idle":"2022-08-08T22:41:15.583436Z","shell.execute_reply.started":"2022-08-08T22:41:15.207757Z","shell.execute_reply":"2022-08-08T22:41:15.582119Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"_=plot_pacf(df_oil_diff.oil_price[1:], lags = 10)","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:41:26.993551Z","iopub.execute_input":"2022-08-08T22:41:26.993960Z","iopub.status.idle":"2022-08-08T22:41:27.263584Z","shell.execute_reply.started":"2022-08-08T22:41:26.993928Z","shell.execute_reply":"2022-08-08T22:41:27.262692Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"As we can see from partial_autocorelation plot there is no AR(q) proces. Future oil prices are independet from it's past values.\n\nWe call this type of Time Series **Random Walk**, it can't be modeled becouse future values are random from it's past values! If you want to model **Random Walk** TS type use **Naive Forecasting**.","metadata":{}},{"cell_type":"markdown","source":"# Indicators \nThere are many ways in which we can model TS. Like Box-Jenkins ARIMA models, Exponential smooting, using libraries like Darts and FbProphet etc... and with ML aproach by treating TS as a regression problem. In the future notebooks I will be be posting ML forecasting aproach with hybrid model, one regression model and with parallel multiprocessing trainig with FbProphet if your using cloud computing this will train very fast even with 1782 models!. This notebook will be used for creating indicator and EDA analysis.","metadata":{}},{"cell_type":"markdown","source":"## Oil indicator\nIm not sure if current oil_price is data lakege, but i will treat it as such. Since we also have oil_price in the test set guess we don't have to run a prediction, as we saw the only prediction we can use for future oil_prices is Naive forecasting.\n\nWe will be using MA Oil filter with weekly window size","metadata":{}},{"cell_type":"code","source":"indicators = pd.DataFrame(index = pd.date_range('2013-01-01', '2017-08-31')).to_period('D')\nindicators.index.rename(\"date\",inplace=True)\nindicators.index=indicators.index.to_timestamp()\ndf_oil['avg_oil'] = df_oil['oil_price'].rolling(7).mean()\ndf_oil.index=df_oil.index.to_timestamp()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:43:33.379320Z","iopub.execute_input":"2022-08-08T22:43:33.379852Z","iopub.status.idle":"2022-08-08T22:43:33.399707Z","shell.execute_reply.started":"2022-08-08T22:43:33.379814Z","shell.execute_reply":"2022-08-08T22:43:33.398403Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Relative to lag values, what is the coorelation between oil_price lags or MA oil_price lags and total product sales?\n","metadata":{}},{"cell_type":"code","source":"def oil_MA_lags(lag=0):\n lags=[]\n for i in range(lag):\n a1 = pd.merge(df_train.groupby([\"date\", \"family\"]).sales.sum().reset_index(), df_oil.avg_oil.shift(i+1).reset_index(),on=[\"date\"], how = \"inner\")\n c1 = a1.groupby(\"family\").corr(\"spearman\").reset_index()\n c1 = c1[c1.level_1==\"avg_oil\"][[\"family\", \"sales\"]].sort_values(\"sales\",ascending=False)\n lags.append(np.mean(abs(c1.sales)))\n return lags\n\ndef oil_price_lags(lag=0):\n lags_price=[]\n for i in range(lag):\n a1 = pd.merge(df_train.groupby([\"date\", \"family\"]).sales.sum().reset_index(), df_oil.oil_price.shift(i+1).reset_index(),on=[\"date\"], how = \"inner\")\n c1 = a1.groupby(\"family\").corr(\"spearman\").reset_index()\n c1 = c1[c1.level_1==\"oil_price\"][[\"family\", \"sales\"]].sort_values(\"sales\",ascending=False)\n lags_price.append(np.mean(abs(c1.sales)))\n return lags_price\n\noil_price_lag=oil_price_lags(100)\noil_MA_lag=oil_MA_lags(100)","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:43:56.144014Z","iopub.execute_input":"2022-08-08T22:43:56.145168Z","iopub.status.idle":"2022-08-08T22:44:56.998386Z","shell.execute_reply.started":"2022-08-08T22:43:56.145123Z","shell.execute_reply":"2022-08-08T22:44:56.997305Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig,ax=plt.subplots(1,1,figsize=(12,7))\nplt.plot(range(1,101),oil_price_lag)\nax.set_title(\"Oil_price and Sales corelation lags\")\nplt.xlabel(\"Lags\")\nplt.ylabel(\"Oil_price and Sales Corelation\")","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:44:57.000235Z","iopub.execute_input":"2022-08-08T22:44:57.000640Z","iopub.status.idle":"2022-08-08T22:44:57.261253Z","shell.execute_reply.started":"2022-08-08T22:44:57.000609Z","shell.execute_reply":"2022-08-08T22:44:57.260161Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig,ax=plt.subplots(1,1,figsize=(12,7))\nplt.plot(range(1,101),oil_MA_lag)\nax.set_title(\"MA_Oil and Sales corelation lags\")\nplt.xlabel(\"Lags\")\nplt.ylabel(\"MA_Oil and Sales Corelation\")","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:44:57.262602Z","iopub.execute_input":"2022-08-08T22:44:57.262928Z","iopub.status.idle":"2022-08-08T22:44:57.530019Z","shell.execute_reply.started":"2022-08-08T22:44:57.262898Z","shell.execute_reply":"2022-08-08T22:44:57.528824Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We clearly see bumps from oil_price and sales correlation plot. So it's better to use MA oil_price as a indicator for prediction. How many lags depends from your inference. Let's pick 3 MA lags but remove current because of data leakage. Maybe we don't know what will be tomorow oil_price if we model it. Lucky future oil_price is given so we don't have to do naive forecasting for test set?","metadata":{}},{"cell_type":"code","source":"n_lags = 3\nfor l in range(1, n_lags + 1): #Note that weekend data from oil is missing \n indicators[f'oil_lags{l}'] = df_oil.avg_oil.shift(l)\n indicators[f'oil_lags{l}']=indicators[f'oil_lags{l}'].interpolate(method=\"linear\")\n indicators[f'oil_lags{l}']= indicators[f'oil_lags{l}'].fillna(method ='bfill')\nindicators.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:44:57.532099Z","iopub.execute_input":"2022-08-08T22:44:57.532467Z","iopub.status.idle":"2022-08-08T22:44:57.555548Z","shell.execute_reply.started":"2022-08-08T22:44:57.532436Z","shell.execute_reply":"2022-08-08T22:44:57.554503Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Transaction indicator\nTransactions can be understood as a number of customers in one store at one point in time making a purchase. Transaction is also **data leakage**!","metadata":{}},{"cell_type":"code","source":"df_transaction=pd.read_csv(\"/kaggle/input/store-sales-time-series-forecasting/transactions.csv\",parse_dates=[\"date\"],index_col=[\"date\"],infer_datetime_format=True).to_period(\"D\")\ndf_transaction.index=df_transaction.index.to_timestamp()\ndf_trans_sales=pd.merge(df_train.groupby([\"date\", \"store_nbr\"]).sales.sum().reset_index(),df_transaction.reset_index(),on=[\"date\",\"store_nbr\"])\ndf_trans_sales.tail(10)","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:45:20.610706Z","iopub.execute_input":"2022-08-08T22:45:20.611137Z","iopub.status.idle":"2022-08-08T22:45:20.963870Z","shell.execute_reply.started":"2022-08-08T22:45:20.611103Z","shell.execute_reply":"2022-08-08T22:45:20.962673Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"cor_trans_sales=df_trans_sales.groupby(\"store_nbr\").corr(\"spearman\").reset_index()\ncor_trans_sales=cor_trans_sales[cor_trans_sales.level_1==\"transactions\"][[\"store_nbr\",\"sales\"]].sort_values(\"sales\",ascending=False)\ncor_trans_sales.store_nbr=cor_trans_sales.store_nbr.astype(\"category\")\npx.bar(cor_trans_sales,x=\"sales\",y=\"store_nbr\",color=\"store_nbr\",title='Transaction and Sales Corelation:{:,.4f}'.format(np.mean(abs(cor_trans_sales.sales))))","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:45:35.148388Z","iopub.execute_input":"2022-08-08T22:45:35.148883Z","iopub.status.idle":"2022-08-08T22:45:35.775675Z","shell.execute_reply.started":"2022-08-08T22:45:35.148843Z","shell.execute_reply":"2022-08-08T22:45:35.774467Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"px.scatter(df_trans_sales.sort_values([\"store_nbr\", \"date\"]), x='sales', y='transactions', color='store_nbr',trendline = \"ols\", trendline_color_override = \"red\",title = \"Sales and transactions\" )","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:45:45.426298Z","iopub.execute_input":"2022-08-08T22:45:45.426808Z","iopub.status.idle":"2022-08-08T22:45:45.725204Z","shell.execute_reply.started":"2022-08-08T22:45:45.426759Z","shell.execute_reply":"2022-08-08T22:45:45.723731Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_trans_sales.store_nbr=df_trans_sales.store_nbr.astype(\"category\")\npx.scatter(df_trans_sales.sort_values([\"store_nbr\", \"date\"]), x='sales', y='transactions', color='store_nbr',hover_data=cor_trans_sales,trendline = \"ols\", trendline_color_override = \"red\",title = \"Sales and transactions per Store\" )","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:46:02.244414Z","iopub.execute_input":"2022-08-08T22:46:02.245009Z","iopub.status.idle":"2022-08-08T22:46:03.041834Z","shell.execute_reply.started":"2022-08-08T22:46:02.244954Z","shell.execute_reply":"2022-08-08T22:46:03.040275Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"You can inspect how much each store sales is corelated with transactions or simply put corelation \"**How much total items in store is sold if there are X nbr of customers (transactions) that day**\"","metadata":{}},{"cell_type":"code","source":"px.line(df_trans_sales.sort_values([\"store_nbr\", \"date\"]), x='date', y='transactions', color='store_nbr',title = \"Transactions\" )","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:46:30.538899Z","iopub.execute_input":"2022-08-08T22:46:30.539427Z","iopub.status.idle":"2022-08-08T22:46:32.690793Z","shell.execute_reply.started":"2022-08-08T22:46:30.539389Z","shell.execute_reply":"2022-08-08T22:46:32.689472Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Observing every transaction TS per store, we notice that it doesn't have upward trend and instability till june 2015 as total sales have. We can infer that instability is a **product supply problem**. While total sales trend is maybe influenced by **inflation** or rise in **income**. We are not predicting total revenue but a number of selled family items in stores. So in order for rise in sales **inflation** must be lower.\n\nObserving stores 20,21,22,29,36,52,53 which are opened later, we infer their opening effect, more customers in a first few days. After effect transactions are more stable. We can also observe some missing data, the first day after New Year except store nbr 25, however stores are closed that day so there are no transactions. There are also missing data when store is not yet opened or it is closed for a small amount of time.\nWe have to predict transactions if we want to use them as a indicator.For transaction forecasting **DirRec** strategy will be used with 16 days horison.\n![](\"https://i.imgur.com/B7KAvAO.png\")
","metadata":{}},{"cell_type":"code","source":"def make_lags(ts, lags, lead_time=1):\n return pd.concat(\n {\n f'trans_lag_{i}': ts.shift(i)\n for i in range(lead_time, lags + lead_time)\n },\n axis=1)\ndef make_multistep_target(ts, steps):\n return pd.concat(\n {f'trans_step_{i + 1}': ts.shift(-i)\n for i in range(steps)},\n axis=1)\n\nX_trans=df_transaction.reset_index().set_index([\"date\",\"store_nbr\"]).unstack(\"store_nbr\").loc[\"2015-6-1\":,\"transactions\"] # for forecasting we use data since 2015-6-1\nY_trans=df_transaction.reset_index().set_index([\"date\",\"store_nbr\"]).unstack(\"store_nbr\").loc[\"2015-6-1\":,\"transactions\"]\nX_trans=make_lags(X_trans,7) #using weekly lags \nY_trans=make_multistep_target(Y_trans,17)#17 lags because first one is transaction for current date\nX_trans=X_trans.loc[\"2015-6-8\":].fillna(0)\nY_trans=Y_trans.loc[:\"2017-7-30\"].fillna(0)\nY_trans,X_trans=Y_trans.align(X_trans,join=\"inner\",axis=0)\nX_trans.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:47:14.889641Z","iopub.execute_input":"2022-08-08T22:47:14.890089Z","iopub.status.idle":"2022-08-08T22:47:15.404878Z","shell.execute_reply.started":"2022-08-08T22:47:14.890039Z","shell.execute_reply":"2022-08-08T22:47:15.403867Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"Y_trans.tail()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:47:25.471983Z","iopub.execute_input":"2022-08-08T22:47:25.472625Z","iopub.status.idle":"2022-08-08T22:47:26.096855Z","shell.execute_reply.started":"2022-08-08T22:47:25.472571Z","shell.execute_reply":"2022-08-08T22:47:26.095712Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_trans=X_trans.stack(\"store_nbr\").sort_index()\nY_trans=Y_trans.stack(\"store_nbr\").sort_index()\nX_trans.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:47:41.484992Z","iopub.execute_input":"2022-08-08T22:47:41.485402Z","iopub.status.idle":"2022-08-08T22:47:41.549924Z","shell.execute_reply.started":"2022-08-08T22:47:41.485368Z","shell.execute_reply":"2022-08-08T22:47:41.548774Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from sklearn.multioutput import RegressorChain\nfrom lightgbm import LGBMRegressor\nfrom sklearn.model_selection import GridSearchCV\nfrom sklearn.metrics import mean_squared_log_error as msle\nfrom sklearn.metrics import make_scorer\n\ndef rmsle_scorer(y,y_pred):\n return np.sqrt(msle(y,y_pred))\n\nRMSLE= make_scorer(rmsle_scorer,greater_is_better=False)\n\nmodel = RegressorChain(base_estimator=LGBMRegressor())\nmodel.get_params()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:47:57.383742Z","iopub.execute_input":"2022-08-08T22:47:57.384383Z","iopub.status.idle":"2022-08-08T22:47:57.714528Z","shell.execute_reply.started":"2022-08-08T22:47:57.384346Z","shell.execute_reply":"2022-08-08T22:47:57.713528Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import warnings\nwarnings.filterwarnings(\"ignore\")\n\nparams={\"base_estimator__n_estimators\":[70,100,120],'base_estimator__reg_lambda':[2,1,0], #training will take a little while on kaggle \n 'base_estimator__min_child_samples':[5,10,15], #sorry for that\n \"base_estimator__max_depth\":[3,5]}\n\ntrans_estimators=[]\nstore_nb=X_trans.reset_index(\"store_nbr\").store_nbr.unique()\na=X_trans.reset_index(\"store_nbr\").copy()\nb=Y_trans.reset_index(\"store_nbr\").copy()\n\nfor i in store_nb:\n clf=GridSearchCV(model,params,cv=10,n_jobs=-1,verbose=0,scoring=\"neg_mean_squared_error\")# for some reason RMSLE won't work here\n clf.fit(a.loc[a.store_nbr==i,a.columns.difference([\"store_nbr\"])],b.loc[b.store_nbr==i,b.columns.difference([\"store_nbr\"])])# don't mind the errors\n trans_estimators.append(clf)","metadata":{"execution":{"iopub.status.busy":"2022-08-08T22:48:51.569055Z","iopub.execute_input":"2022-08-08T22:48:51.569920Z","iopub.status.idle":"2022-08-09T02:03:57.955726Z","shell.execute_reply.started":"2022-08-08T22:48:51.569865Z","shell.execute_reply":"2022-08-09T02:03:57.954275Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_trans=df_transaction.reset_index().set_index([\"date\",\"store_nbr\"]).unstack(\"store_nbr\").loc[\"2015-6-1\":,\"transactions\"]\nX_trans=make_lags(X_trans,7)\n\nX_trans=X_trans.loc[\"2015-6-8\":].fillna(0)\nX_trans=X_trans.stack(\"store_nbr\")\nY_trans_predictions= pd.DataFrame(index=X_trans.reset_index()[\"date\"].unique())\nY_trans_predictions[\"store_nbr\"]=None\nindex=Y_trans_predictions.index\n\na=X_trans.reset_index(\"store_nbr\").copy()\nfor i,model in enumerate(trans_estimators):\n Y_pred=pd.DataFrame(model.predict(a.loc[a.store_nbr==i+1,a.columns.difference([\"store_nbr\"])]),index=index,columns=Y_trans.columns)\n Y_pred[\"store_nbr\"]=i+1\n Y_trans_predictions=pd.concat([Y_trans_predictions,Y_pred])\nY_trans_predictions=Y_trans_predictions.reset_index().rename(columns={\"index\":\"date\"}).set_index([\"date\",\"store_nbr\"]).sort_index()\n\nY_trans_predictions.dropna(inplace=True)\nY_trans_predictions.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-09T02:03:57.958600Z","iopub.execute_input":"2022-08-09T02:03:57.959514Z","iopub.status.idle":"2022-08-09T02:03:59.697781Z","shell.execute_reply.started":"2022-08-09T02:03:57.959462Z","shell.execute_reply":"2022-08-09T02:03:59.696704Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import datetime\ntrans_indic=Y_trans_predictions.loc[:,\"trans_step_1\"] #taking trans_step_1 as a Current day Transaction prediction \ntrans_indic.name=\"transactions\"\ntrans_indic=trans_indic.reset_index()\ntrans_indic.set_index([\"date\",\"store_nbr\"],inplace=True)\nfor step in range(2,18):\n a=Y_trans_predictions.loc[\"2017-08-15\",\"trans_step_{}\".format(step)] # Our test set has 16 day horison so we take data from the last day of training set\n a.name=\"transactions\" # up to trans_step_17 which is our 16 day horison\n a=a.reset_index()\n a.date=a[\"date\"]+datetime.timedelta(days=step-1) # adding one day at the time\n a.set_index([\"date\",\"store_nbr\"],inplace=True)\n trans_indic=pd.concat([trans_indic,a]).sort_index()\n\ntrans_indic.tail()","metadata":{"execution":{"iopub.status.busy":"2022-08-09T02:03:59.699612Z","iopub.execute_input":"2022-08-09T02:03:59.700385Z","iopub.status.idle":"2022-08-09T02:04:00.024045Z","shell.execute_reply.started":"2022-08-09T02:03:59.700343Z","shell.execute_reply":"2022-08-09T02:04:00.022819Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Since we took training data from **2015-6-8** these are the stores since then that are opened or closed for a short time.\n\n**Closed**\n\n**18** from 15.8.2016-2.12.2016. **25** from 22.4.2016- 26.10.2016.\n\n\n---\n\n**Opened**:\n\n**21** from 23.7.2015. **22** from 8.10. 2015. **42** from 20.8.2015. **52** from 19.4.2017.\n\nHow well did we forecast stores that were closed at any point in time.","metadata":{}},{"cell_type":"code","source":"import plotly.graph_objs as go\n\nstores=[18,21,22,25,42,52]\na=trans_indic.reset_index()\nb=df_transaction.reset_index()\nplt.figure(figsize=(30,25))\nfor i,store in enumerate(stores):\n fig = go.Figure([\n go.Scatter(\n name='Transactions',\n x=b.loc[b.store_nbr==store,\"date\"],\n y=b.loc[b.store_nbr==store,\"transactions\"],\n mode='lines',\n line=dict(color='rgb(31, 119, 180)')\n \n ),\n go.Scatter(\n \n name='Predicted Transactions',\n x=a.loc[a.store_nbr==store,\"date\"],\n y=a.loc[a.store_nbr==store,\"transactions\"],\n mode='lines',\n line=dict(color='rgb(179, 129, 21)')\n \n )])\n fig.update_layout(\n yaxis_title='Transactions',\n title=f'Transactions for store {store}',\n hovermode=\"x\"\n )\n fig.show()","metadata":{"execution":{"iopub.status.busy":"2022-08-09T02:04:00.026830Z","iopub.execute_input":"2022-08-09T02:04:00.027228Z","iopub.status.idle":"2022-08-09T02:04:00.364571Z","shell.execute_reply.started":"2022-08-09T02:04:00.027193Z","shell.execute_reply":"2022-08-09T02:04:00.363351Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"indicators=pd.merge(trans_indic.reset_index(),indicators.reset_index(),on=[\"date\"]).set_index(\"date\")\nindicators.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-09T02:04:00.366611Z","iopub.execute_input":"2022-08-09T02:04:00.367080Z","iopub.status.idle":"2022-08-09T02:04:00.398083Z","shell.execute_reply.started":"2022-08-09T02:04:00.367041Z","shell.execute_reply":"2022-08-09T02:04:00.397138Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Removing rows when stores were closed\n","metadata":{}},{"cell_type":"code","source":"indicators=indicators.reset_index().set_index([\"date\",\"store_nbr\"])\nb=indicators.reset_index()\na=df_train.reset_index().groupby([\"date\",\"store_nbr\"]).sales.sum()\na=a.loc[a.values==0].reset_index().drop(columns=[\"sales\"])\nb=b.merge(a,how = 'outer', indicator = True).set_index([\"date\",\"store_nbr\"])\nclosed_stores_index=b.loc[b._merge==\"both\"].index\n\nindicators.drop(closed_stores_index,inplace=True)\n\ndf_train_index=a.set_index([\"date\",\"store_nbr\"]).index\ndf_train.drop(df_train_index,inplace=True)\nfor store in stores:\n print(f'{store}:{np.min(indicators.loc[indicators.index.get_level_values(1)==store,\"transactions\"])}')","metadata":{"execution":{"iopub.status.busy":"2022-08-09T02:04:00.399866Z","iopub.execute_input":"2022-08-09T02:04:00.400247Z","iopub.status.idle":"2022-08-09T02:06:30.294687Z","shell.execute_reply.started":"2022-08-09T02:04:00.400214Z","shell.execute_reply":"2022-08-09T02:06:30.293510Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Above are transactions on the first day when stores opened. So to remove first day opening effect we remove from data any transactions below 79. Since we can't realy predict transactions on the opening day because we have no prior data.","metadata":{}},{"cell_type":"code","source":"indicators=indicators.loc[indicators.transactions>79,:]","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Detecting and removing anomalies","metadata":{}},{"cell_type":"code","source":"df_trans_sales[\"sale/trans\"]=df_trans_sales[\"sales\"]/df_trans_sales[\"transactions\"]\npx.line(df_trans_sales,x=\"date\",y=\"sale/trans\",color=\"store_nbr\",title = \"Avg items sold per Customer\")","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"It looks like there was an error in data or some kind of transaction anomaly (less likey) in data captured 2016-1-2 and 2016-1-4. Best is to delete data from that period so it can't influence model.\n\n### 2016-1-2 data anomaly","metadata":{}},{"cell_type":"code","source":"df_trans_sales.loc[df_trans_sales.date==\"2016-1-2\",:].head()","metadata":{"execution":{"iopub.status.busy":"2022-08-09T02:06:30.296081Z","iopub.execute_input":"2022-08-09T02:06:30.296478Z","iopub.status.idle":"2022-08-09T02:06:30.314266Z","shell.execute_reply.started":"2022-08-09T02:06:30.296443Z","shell.execute_reply":"2022-08-09T02:06:30.313004Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### 2016-1-4 data anomaly","metadata":{}},{"cell_type":"code","source":"df_trans_sales.loc[df_trans_sales.date==\"2016-1-4\",:].sort_values([\"store_nbr\"]).head()","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"It does seems extremely unlikely that customers on average bought 2770.951000 and 1461.933902 items respectively. We also see a lot of missing values from stores these days.It could be some kind of an error in collecting data.","metadata":{}},{"cell_type":"code","source":"df_trans_sales.set_index([\"date\",\"store_nbr\"],inplace=True)\nindicators.reset_index(inplace=True)\nindicators.set_index([\"date\",\"store_nbr\"],inplace=True)\ndf_transaction.reset_index(inplace=True)\ndf_transaction.set_index([\"date\",\"store_nbr\"],inplace=True)\n\ndf_train.drop([\"2016-1-2\",\"2016-1-4\"],level=0,inplace=True)\ndf_trans_sales.drop(index=[\"2016-1-2\",\"2016-1-4\"],level=0,inplace=True)\ndf_transaction.drop(index=[\"2016-1-2\",\"2016-1-4\"],level=0,inplace=True)\nindicators.drop(index=[\"2016-1-2\",\"2016-1-4\"],level=0,inplace=True)","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_trans_sales.reset_index(inplace=True)\npx.line(df_trans_sales.loc[df_trans_sales.date>\"2015-6-1\",:],x=\"date\",y=\"sale/trans\",color=\"store_nbr\",title = \"Sales per Customer\")","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"You can inspect those impulses and figure out why did customers buy more items that day. Where there a promotions? Maybe a holiday or something else?\n\nYou can also check if instability exist til june 2015 as we saw in sales, and indeed it does which means for sure that was because of product supply issue. Also check if in some stores there was upward trend, and indeed some stores have trend, that on average customers are **buying more products!** Which at the end is also result of Ecuador economic boom as we will see down below.\nSo trend in sales is not only because Favorita opened more stores it was also influence by customers bying more products. Maybe Favorita could open even more stores?! But that would require some geospatial analysis, involving store places.\n### Inflation and income in Ecuador\n\nInflation is lower over the years. That could have an affect in bying more items which could reflect in customers on average bying more products, but we didn't see that indicator in Sales per customer graph.\nEcuador Economic Growth 1960-2022\n\n","metadata":{}},{"cell_type":"markdown","source":"From the data above we clearly see that Ecudor economy is booming. So why was there trend in sales? As more stores are opened more sales wiil it be and also people are bying more products in some stores. Positive economic trend lift many Ecuadorians, maybe thats why Favorita has decided to open more stores over time as people in EC are getting wealthier.","metadata":{}},{"cell_type":"markdown","source":"### Holiday, Events indicators","metadata":{}},{"cell_type":"code","source":"a = df_transaction.reset_index().copy()\na[\"year\"] = a.date.dt.year\na[\"dayofweek\"] = a.date.dt.dayofweek+1\na = a.groupby([\"year\", \"dayofweek\"]).transactions.mean().reset_index()\npx.line(a, x=\"dayofweek\", y=\"transactions\" , color = \"year\", title = \"Avg Transactions (Day of the Week)\")","metadata":{"execution":{"iopub.status.busy":"2022-08-09T02:06:30.315769Z","iopub.execute_input":"2022-08-09T02:06:30.316227Z","iopub.status.idle":"2022-08-09T02:06:30.437642Z","shell.execute_reply.started":"2022-08-09T02:06:30.316193Z","shell.execute_reply":"2022-08-09T02:06:30.436449Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"a = df_transaction.reset_index().copy()\na[\"year\"] = a.date.dt.year\na[\"dayofweek\"] = a.date.dt.dayofweek+1\na = a.groupby([\"year\", \"dayofweek\"]).transactions.sum().reset_index()\npx.line(a, x=\"dayofweek\", y=\"transactions\" , color = \"year\", title = \"Total Transactions (Day of the Week)\")","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We can infer weekly seasonality, stores have the most customers on saturday. And lowest on thursday. \n\nInteresting conclusion from first graph (Avg Transactions) is that in 2014 on average Favorita had most customers and lowest in 2016 which is simmilar to 2017 prediction year. However we know not every store is opened in 2013,2014 and 2015 so we take these results with a grain of salt. As in later years when more stores are opened more customers will be spread out.\n\nLooking at Total transaction graph, Favorita has indeed more customers in total over the years (2017 is ignored data is not complete) as more stores are opened.","metadata":{}},{"cell_type":"code","source":"a = df_train.reset_index().copy()\na[\"year\"] = a.date.dt.year\na[\"dayofweek\"] = a.date.dt.dayofweek+1\na = a.groupby([\"year\", \"dayofweek\"]).sales.mean().reset_index()\npx.line(a, x=\"dayofweek\", y=\"sales\" , color = \"year\", title = \"Avg Sales (day of the week)\")","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_holiday=pd.read_csv(\"/kaggle/input/store-sales-time-series-forecasting/holidays_events.csv\",parse_dates=[\"date\"],index_col=[\"date\"],infer_datetime_format=True).to_period(\"D\")\nprint(\"Unique locale names:\",df_holiday.locale_name.unique())\nprint(\"Unique locale:\",df_holiday.locale.unique())\nprint(\"Unique types:\",df_holiday.type.unique())","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_holiday=df_holiday[df_holiday.locale==\"National\"]\ndf_holiday=df_holiday.groupby(df_holiday.index).first()\ndf_holiday.index=df_holiday.index.to_timestamp()\ndf_holiday[df_holiday.type==\"Event\"]","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"indicators=pd.merge(indicators.reset_index(),df_holiday.reset_index(), on=\"date\",how=\"left\")\nindicators[\"Day_of_week\"]=indicators.date.dt.dayofweek+1\nindicators[\"wrkDay\"]=1\nindicators=indicators.set_index([\"date\",\"store_nbr\"])\nindicators.loc[indicators.Day_of_week>5,\"wrkDay\"]=0\nindicators.loc[(indicators.type==\"Holiday\") & (indicators.transferred==False),\"wrkDay\"]=0\nindicators.loc[indicators.type==\"Transfer\",\"wrkDay\"]=0\nindicators.loc[indicators.type==\"Bridge\",\"wrkDay\"]=0\nindicators.loc[indicators.type==\"Work Day\",\"wrkDay\"]=1\n\nindicators[\"disaster\"]=0\nindicators[\"disaster_step_lag\"]=0\nindicators[\"Black Friday\"]=0\nindicators[\"Cyber Monday\"]=0\nindicators[\"Mothers Day\"]=0\n\ndf_holiday.loc[df_holiday.type==\"Additional\",:]","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"indicators.loc[(indicators.type==\"Event\") & (indicators.description.str.contains(pat=\"Terremoto\")),\"disaster\"]=1\nindicators.loc[(indicators.type==\"Event\") & (indicators.description.str.contains(pat=\"Black Friday\")),\"Black Friday\"]=1\nindicators.loc[(indicators.type==\"Event\") & (indicators.description.str.contains(pat=\"Cyber Monday\")),\"Cyber Monday\"]=1\nindicators.loc[(indicators.type==\"Event\") & (indicators.description.str.contains(pat=\"Dia de la Madre\")),\"Mothers Day\"]=1\n\nfor i,date in enumerate(pd.date_range(\"2016-04-17\",\"2016-05-16\")):\n indicators.loc[(indicators.index.get_level_values(0)==date) & (indicators.type==\"Event\") \n & (indicators.description.str.contains(pat=\"\\+\")),\"disaster_step_lag\"]=i+1\n\nindicators[\"Christmass_lag\"]=0\nfor y in indicators.index.get_level_values(0).year.unique():\n lag=0\n for d in range(19,25):\n lag=lag+1\n indicators.loc[(indicators.index.get_level_values(0).year==y) & \n (indicators.index.get_level_values(0).month==12) &\n (indicators.index.get_level_values(0).day==d),\"Christmass_lag\"]=lag\n\nindicators = pd.get_dummies(indicators, columns = ['Day_of_week'])\nindicators.drop(columns=[\"locale\",\"locale_name\",\"description\",\"transferred\",\"type\"],inplace=True)\nindicators['pay_day']=0\nindicators.loc[(indicators.reset_index(\"store_nbr\").index.is_month_end) \n | (indicators.reset_index(\"store_nbr\").index.day == 15), 'pay_day'] = 1\n\nindicators.head()","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Family products analysis","metadata":{}},{"cell_type":"code","source":"a=df_train.loc[\"2015-06-08\":].reset_index().copy()\na=a.groupby([\"date\",\"family\"]).sales.sum().reset_index()\npx.line(a.sort_values([\"date\",\"family\"]),x=\"date\",y=\"sales\",color=\"family\",title=\"Total sales per family Products\")","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### School season indicator\n\nSchool and Office supplies have seasonality in april, may, august and september.","metadata":{}},{"cell_type":"code","source":"indicators['school&Office']=0\nindicators.loc[indicators.index.get_level_values(0).month.isin([4, 5, 8, 9]),'school&Office']=1","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Checking family products stationarity","metadata":{}},{"cell_type":"code","source":"for family in a.family.unique():\n print(f'\\nStacionarity for {family}:\\n')\n adfuller_test(a.loc[a.family==family,\"sales\"])\n","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"***NonStationary Time Series***\n\n\n---\n**Books,Dairy,Eggs,HOME APPLIANCES,LAWN AND GARDEN,PET SUPPLIES,PREPARED FOOD**S","metadata":{}},{"cell_type":"code","source":"df_train=df_train.loc[\"2015-06-08\":] #selecting date from which we will train our model","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**BOOKS** were seling since 10.10.2016 could be a trial product, after some time sales are almost zero.\n\n**EGGS** after 20.06.2016 has regular bigger imulses than before maybe from that date data should be used for modeling.\n\n**FROZEN FOODS** has yearly impulses in december (christmas new year)\n\n**GROCERY I** has impulse shortly after earthquake\n\n**HOME CARE** has impulse shortly after earthquake\n\n**LAWN AND GARDEN** is introduced in more stores after 3.12.2016\n\n**MAGASINES** after 10.10.2015 jump in sales\n\n**PERSONAL CARE** has impulses shortly after earthquake\n\n\n**PREPARED FOODS** has downward trend in 15.04.2017 could be lack of selling items. Maybe after that date data should be used for training.\n\n**SCHOOL AND OFFICE** is already mention.\n\nRest of Family products are more or less stable with a few impulses.","metadata":{}},{"cell_type":"code","source":"a=df_train.reset_index()\nfor fam in a.family.unique():\n fig=px.line(a.loc[a.family==fam,:],x=\"date\",y=\"sales\",color=\"store_nbr\",title=f'{fam}')\n fig.show()","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Sezonality for individual family products- Periodogram","metadata":{}},{"cell_type":"code","source":"def plot_periodogram(ts, detrend='linear', ax=None,family=None):\n from scipy.signal import periodogram\n fs = pd.Timedelta(\"1Y\") / pd.Timedelta(\"1D\")\n freqencies, spectrum = periodogram(\n ts,\n fs=fs,\n detrend=detrend,\n window=\"boxcar\",\n scaling='spectrum',\n )\n if ax is None:\n fig, ax = plt.subplots(figsize=(12,5))\n ax.step(freqencies, spectrum, color=\"purple\")\n ax.set_xscale(\"log\")\n ax.set_xticks([1, 2, 4, 6, 12, 26, 52, 104])\n ax.set_xticklabels(\n [\n \"Annual (1)\",\n \"Semiannual (2)\",\n \"Quarterly (4)\",\n \"Bimonthly (6)\",\n \"Monthly (12)\",\n \"Biweekly (26)\",\n \"Weekly (52)\",\n \"Semiweekly (104)\",\n ],\n rotation=30,\n )\n ax.ticklabel_format(axis=\"y\", style=\"sci\", scilimits=(0, 0))\n ax.set_ylabel(\"Variance\")\n ax.set_title(family)\n return ax\n\n%matplotlib inline\n\nt=df_train.reset_index().groupby([\"date\",\"family\"]).sales.sum().reset_index(\"family\")\nfor i,fam in enumerate(t.family.unique()):\n plot_periodogram(t.loc[t.family==fam,\"sales\"],family=fam)","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"y = indicators.unstack(['store_nbr'])\nfourier = CalendarFourier(freq = 'W', order = 3)\nfourier1 = CalendarFourier(freq = 'M', order = 10)\nfourier2=CalendarFourier(freq = 'Y', order = 15)\ndp = DeterministicProcess(index = y.index,\n order = 1,\n seasonal = False,\n constant = True,\n additional_terms = [fourier,fourier1,fourier2],\n drop = True)\nx = dp.in_sample()\nindicators = pd.merge(x.reset_index().rename(columns={\"index\":\"date\"}),indicators.reset_index(),on=\"date\")\nindicators.set_index([\"date\",\"store_nbr\"],inplace=True)\nindicators.head()","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Note: We didn't add stores nad places indicators becouse they are static features and we have 1782 separate models.","metadata":{}},{"cell_type":"code","source":"X_train,Y_train=indicators.align(df_train,join=\"inner\",axis=0)\nX_train=X_train.join(Y_train.onpromotion)\nY_train.drop(columns=\"onpromotion\",inplace=True)\n\nX_test,Y_test=indicators.align(df_test,join=\"inner\",axis=0)\nX_test=X_test.join(Y_test.onpromotion)","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_train.to_csv(\"./X_train.csv\")\nY_train.to_csv(\"./Y_train.csv\")\n\nX_test.to_csv(\"./X_test.csv\")\nY_test.to_csv(\"./Y_test.csv\")","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Thank you for reading this Notebook. :)\n\nI will be posting in upcoming Notebooks how to forecast with FbProphet using parallel computing.\nForecasting with hybrid models and with LGBMRegressor. Will also be adding more indicators to help with forecasting problems. ","metadata":{}}]}