{"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":"- Configuration","metadata":{}},{"cell_type":"code","source":"!pip install darts","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import warnings\nwarnings.filterwarnings('ignore')\nwarnings.simplefilter(action='ignore', category=FutureWarning)\n\nimport numpy as np\nimport pandas as pd\n\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nimport plotly.graph_objects as go\nfrom plotly.subplots import make_subplots\n\nfrom datetime import datetime\nfrom functools import reduce\n\nfrom darts import TimeSeries\nfrom darts.models import (\n    NaiveSeasonal,\n    NaiveDrift,\n    Prophet,\n    ExponentialSmoothing,\n    ARIMA,\n    AutoARIMA,\n    RegressionEnsembleModel,\n    RegressionModel,\n    FFT\n)\n\nfrom darts.metrics import mape, mase\nfrom darts.utils.statistics import check_seasonality, plot_acf, plot_residuals_analysis, plot_hist\n\nimport warnings\nwarnings.filterwarnings(\"ignore\")\nimport logging\nlogging.disable(logging.CRITICAL)","metadata":{"execution":{"iopub.execute_input":"2021-11-27T22:23:55.878659Z","iopub.status.busy":"2021-11-27T22:23:55.878391Z","iopub.status.idle":"2021-11-27T22:23:58.104276Z","shell.execute_reply":"2021-11-27T22:23:58.103381Z","shell.execute_reply.started":"2021-11-27T22:23:55.878629Z"}},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"- Utils","metadata":{}},{"cell_type":"code","source":"def group_trading_timeseries(df, time_step_size):\n    #Copied from https://www.kaggle.com/alexisalvarez/eda-resampling-the-power-of-technical-analysis\n    df = df.sort_index()\n    df = df.resample('{}'.format(time_step_size)).agg({\n            'Count': lambda df: df.sum(),\n            'Open': lambda s: s.iloc[0] ,\n            'High': lambda df: df.max(), \n            'Low': lambda df: df.min(),\n            'Close': lambda df: df.iloc[-1], \n            'Volume': lambda df: df.sum(),\n            'VWAP' : lambda df: df.mean(),\n            'Target' : lambda df: df.mean()\n            }) \n    return df","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Preprocessing of the data","metadata":{}},{"cell_type":"code","source":"df = pd.read_csv(\"../input/g-research-crypto-forecasting/train.csv\")\nasset = pd.read_csv(\"../input/g-research-crypto-forecasting/asset_details.csv\")\ndf[\"timestamp\"] = df[\"timestamp\"].apply(lambda x: datetime.fromtimestamp(x))\ndf[\"date\"] = df[\"timestamp\"].astype('datetime64[s]')\ndf.head()","metadata":{"execution":{"iopub.execute_input":"2021-11-27T22:23:58.119424Z","iopub.status.busy":"2021-11-27T22:23:58.119154Z","iopub.status.idle":"2021-11-27T22:25:41.915045Z","shell.execute_reply":"2021-11-27T22:25:41.913587Z","shell.execute_reply.started":"2021-11-27T22:23:58.119375Z"}},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df = pd.merge(df, asset, on=\"Asset_ID\")\ndf.set_index(\"date\", inplace=True)\ndf = df[[\"Asset_Name\", \"Weight\", \"Count\",\n         \"Open\", \"High\", \"Low\",\n         \"Close\", \"Volume\", \"VWAP\",\n         \"Target\"]]\ndf.head()","metadata":{"execution":{"iopub.execute_input":"2021-11-27T22:25:41.917904Z","iopub.status.busy":"2021-11-27T22:25:41.917257Z","iopub.status.idle":"2021-11-27T22:25:53.255520Z","shell.execute_reply":"2021-11-27T22:25:53.254795Z","shell.execute_reply.started":"2021-11-27T22:25:41.917853Z"}},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df.isnull().mean() * 100","metadata":{"execution":{"iopub.execute_input":"2021-11-27T22:25:53.257664Z","iopub.status.busy":"2021-11-27T22:25:53.257181Z","iopub.status.idle":"2021-11-27T22:25:56.298199Z","shell.execute_reply":"2021-11-27T22:25:56.297518Z","shell.execute_reply.started":"2021-11-27T22:25:53.257631Z"}},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df.dropna(inplace=True)","metadata":{"execution":{"iopub.execute_input":"2021-11-27T22:25:56.299987Z","iopub.status.busy":"2021-11-27T22:25:56.299495Z","iopub.status.idle":"2021-11-27T22:26:01.291291Z","shell.execute_reply":"2021-11-27T22:26:01.290177Z","shell.execute_reply.started":"2021-11-27T22:25:56.299953Z"}},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## First trials with Bitcoin","metadata":{}},{"cell_type":"code","source":"df_btc = df.loc[df[\"Asset_Name\"] == \"Bitcoin\"]\ndf_btc = group_trading_timeseries(df_btc, time_step_size=\"1d\")\ndf_btc.head()","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig = go.Figure(data=[go.Candlestick(x=df_btc.index,\n                                     open=df_btc['Open'],\n                                     high=df_btc['High'],\n                                     low=df_btc['Low'],\n                                     close=df_btc['Close'])])\nfig.update_layout(xaxis_rangeslider_visible=False)\nfig.update_layout(\n    title='Bitcoin from Janv. 2018 to Sep. 2021',\n    yaxis_title='USD',\n    shapes = [dict(\n        x0='2020-03-11', x1='2020-03-11', y0=0, y1=1, xref='x', yref='paper',\n        line_width=2)],\n    annotations=[dict(\n        x='2020-03-11', y=0.05, xref='x', yref='paper',\n        showarrow=False, xanchor='left', text='Covid Sanitary Crisis Begin')]\n)\nfig.show()","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for i in df_btc.columns:\n    print('\\033[1m'+i+'\\033[0m')\n    print(\"Mean: \", round(df_btc[i].mean(), 3))\n    print(\"Median: \", round(df_btc[i].median(), 3))\n    print(\"Standard Deviation: \", round(df_btc[i].std(), 3))\n    print(\"Variance: \", round(df_btc[i].var(), 3))\n    print(\"-------------------\")","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for i in df_btc.columns:\n    sns.displot(x=i,\n                data=df_btc,\n                kde=True)\n    plt.title(\"Distribution of the variable: \" +  i)\n    plt.show()","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Modelisation in Opening data","metadata":{}},{"cell_type":"code","source":"df_btc_open = df_btc[[\"Open\"]]","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_btc_open = TimeSeries.from_dataframe(df_btc_open)","metadata":{"execution":{"iopub.execute_input":"2021-11-27T22:26:01.293091Z","iopub.status.busy":"2021-11-27T22:26:01.292840Z"}},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_btc_open.plot();","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Creating train/val set","metadata":{}},{"cell_type":"code","source":"train_btc, val_btc = df_btc_open.split_before(pd.Timestamp('20200921'))\ntrain_btc.plot(label='training')\nval_btc.plot(label='validation')\nplt.legend();","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Naive Models","metadata":{}},{"cell_type":"code","source":"naive_model = NaiveSeasonal(K=1)\nnaive_model.fit(train_btc)\nnaive_forecast = naive_model.predict(365)\n\ndf_btc_open.plot(label='actual')\nnaive_forecast.plot(label='naive forecast (K=1)')\nplt.legend();","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"seasonal_model = NaiveSeasonal(K=12)\nseasonal_model.fit(train_btc)\nseasonal_forecast = seasonal_model.predict(365)\n\ndf_btc_open.plot(label='actual')\nseasonal_forecast.plot(label='naive forecast (K=12)')\nplt.legend();","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"drift_model = NaiveDrift()\ndrift_model.fit(train_btc)\ndrift_forecast = drift_model.predict(365)\n\ncombined_forecast = drift_forecast + seasonal_forecast - train_btc.last_value()\n\ndf_btc_open.plot()\ncombined_forecast.plot(label='combined')\ndrift_forecast.plot(label='drift')\nplt.legend();","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(\"Mean absolute percentage error for the combined naive drift + seasonal: {:.2f}%.\".format(\n      mape(df_btc_open, combined_forecast)))","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Probabilistic Methods","metadata":{}},{"cell_type":"code","source":"model_es = ExponentialSmoothing()\nmodel_es.fit(train_btc)\nprobabilistic_forecast = model_es.predict(len(val_btc), num_samples=500)\n\ndf_btc_open.plot(label='actual')\nprobabilistic_forecast.plot(label='probabilistic forecast')\nplt.legend()\nplt.show()","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"probabilistic_forecast.plot(low_quantile=0.01, high_quantile=0.99, label='1-99th percentiles')\nprobabilistic_forecast.plot(low_quantile=0.2, high_quantile=0.8, label='20-80th percentiles')","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Regression approaches","metadata":{}},{"cell_type":"code","source":"models = [NaiveSeasonal(6), NaiveSeasonal(12), NaiveDrift()]\n\nmodel_predictions = [m.historical_forecasts(df_btc_open,\n                                            start=pd.Timestamp('20200921'),\n                                            forecast_horizon=12,\n                                            stride=12,\n                                            last_points_only=False,\n                                            verbose=True)\n                     for m in models]\n\nmodel_predictions = [reduce((lambda a, b: a.append(b)), model_pred) for model_pred in model_predictions]","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"model_predictions_stacked = model_predictions[0]\nfor model_prediction in model_predictions[1:]:\n    model_predictions_stacked = model_predictions_stacked.stack(model_prediction)","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\" We build the regression model, and tell it to use the current predictions\n\"\"\"\nregr_model = RegressionModel(lags=None, lags_future_covariates=[0])\n\n\"\"\" Our target series is what we want to predict (the actual data)\n    It has to have the same time index as the features series:\n\"\"\"\nseries_target = df_btc_open.slice_intersect(model_predictions[0])\n\n\"\"\" Here we backtest our regression model\n\"\"\"\nensemble_pred = regr_model.historical_forecasts(\n    series=series_target, future_covariates=model_predictions_stacked,\n    start=pd.Timestamp('20200922'), forecast_horizon=3, verbose=True\n)","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig, ax = plt.subplots(2,2,figsize=(12,6))\nax = ax.ravel()\n\nfor i, m in enumerate(models):\n    df_btc_open.plot(label='actual', ax=ax[i])\n    model_predictions[i].plot(label=str(m), ax=ax[i])\n\n    # intersect last part, to compare all the methods over the duration of the ensemble forecast\n    model_pred = model_predictions[i].slice_intersect(ensemble_pred)\n\n    mape_model = mape(df_btc_open, model_pred)\n    ax[i].set_title('\\nMAPE: {:.2f}%'.format(mape_model))\n    ax[i].legend()\n\ndf_btc_open.plot(label='actual', ax=ax[3])\nensemble_pred.plot(label='Ensemble', ax=ax[3])\nax[3].set_title('\\nMAPE, ensemble: {:.2f}%'.format(mape(df_btc_open, ensemble_pred)))\nax[3].legend()\n\nprint('\\nRegression coefficients for the individual models:')\nfor i, m in enumerate(models):\n    print('Learned coefficient for {}: {:.2f}'.format(m, regr_model.model.coef_[i]))\nplt.tight_layout();","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ensemble_model = RegressionEnsembleModel(\n    forecasting_models=[NaiveSeasonal(6), NaiveSeasonal(12), NaiveDrift()],\n    regression_train_n_points=12)\n\nensemble_model.fit(train_btc)\nensemble_pred = ensemble_model.predict(365)\n\ndf_btc_open.plot(label='actual')\nensemble_pred.plot(label='Ensemble forecast')\nplt.title('MAPE = {:.2f}%'.format(mape(ensemble_pred, df_btc_open)))\nplt.legend();","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ensemble_pred_hist = ensemble_model.historical_forecasts(df_btc_open,\n                                                    start=pd.Timestamp('20200922'),\n                                                    forecast_horizon=3,\n                                                    verbose=True)\ndf_btc_open.plot(label='actual')\nensemble_pred_hist.plot(label='Ensemble forecast')\nplt.title('Historical forecast: MAPE = {:.2f}%'.format(mape(ensemble_pred_hist, df_btc_open)))\nplt.legend();","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}