{"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":"# SARIMAとは？\n## ARモデル(自己回帰モデル)\n「過去の自分のデータ」を説明変数とする。\nΦの値によっては定常モデルになる\n$$\n    y_t = c + \\phi_1 y_{t-1} + \\phi_2 y_{t-2} +･･･ + \\phi_p y_{t-p} + \\epsilon_t\n$$\n## MAモデル(移動平均モデル)\nたとえば「今日の売上が予想より多ければ明日の売上は少なくなる」というような現象をモデリングできる。\n定常なモデル  \n\n\n$$\n    y_t = c + \\theta_1 \\epsilon_{t-1} + \\theta_2 \\epsilon_{t-2} +･･･ + \\theta_q \\epsilon_{t-q} + \\epsilon_t\n$$\n\n## ARMAモデル(ARモデル+MAモデル)\n定常性を有する時系列データに対して有効\n$$\n    y = c + \\sum_{i=1}^{p}\\phi_iy_{t-i} + \\sum_{j=1}^q\\theta_j\\epsilon_{j-q} + \\epsilon_t\\\\\n$$\n\n## ARIMAモデル\nARMAを定常性がないモデルにまで拡張したモデル。差分をとることで定常モデルとして扱えるようデータを加工する。\n$$\\Delta y = y_t-y_{t-1}$$\n一階差分で定常性が見れないときは二階差分、三階差分...d階差分と拡張する。\n$$\\Delta^2 y = \\Delta y_t-\\Delta y_{t-1} \\quad \\Delta^3 y = \\Delta^2 y_t-\\Delta^2 y_{t-1}$$\n\n## SARIMAモデル (Seasonal ARIMA model)\nARIMAモデルに周期成分を取り入れたモデル。  \n時系列方向の説明にARIMA(p,d,q) モデルを使うだけでなく、周期方向の説明にもARIMA(P,D,Q)モデルを利用する。  \nSARIMAモデルでは合計7個の次数がある。 時系列方向のARIMA(p,d,q)に加え季節差分方向のARIMA(P,D,Q)と周期 s がある。  \nこれをSARIMA(p,d,q)(P,D,Q)\\[s\\]と表記する。  \n次数の組み合わせ爆発の問題を回避するため、周期 s は作図や自己相関関数をもとに決め打ちし、季節差分の P,D,Qは低く抑えることがよく行われる。\n\n参考:  \nhttps://ai-trend.jp/basic-study/time-series-analysis/sarima_model/  \nhttp://by-oneself.com/r_ar_ma/\n","metadata":{"_kg_hide-output":true}},{"cell_type":"markdown","source":"# statsmodels: 線形回帰、時系列問題を解くライブラリ  \n参考1: [statsmodelsによる線形回帰 入門](https://qiita.com/innovation1005/items/b712ce54a7a697a9bf03)  \n参考2: [Statsmodelsによる時系列分析入門](https://qiita.com/innovation1005/items/6c5263d79ccc67263b2c)  ","metadata":{}},{"cell_type":"code","source":"import statsmodels.api as sm\nfrom statsmodels.tsa.arima_model import ARIMA\nfrom statsmodels.tsa.statespace.sarimax import SARIMAX\nfrom statsmodels.tsa.seasonal import seasonal_decompose\n\nfrom scipy import stats","metadata":{"execution":{"iopub.status.busy":"2022-02-24T14:14:51.326799Z","iopub.execute_input":"2022-02-24T14:14:51.327725Z","iopub.status.idle":"2022-02-24T14:14:53.343823Z","shell.execute_reply.started":"2022-02-24T14:14:51.327575Z","shell.execute_reply":"2022-02-24T14:14:53.342617Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#基本ライブラリ\nimport os\nfrom datetime import datetime\nimport time\nfrom itertools import product\nimport warnings\nwarnings.filterwarnings(\"ignore\")\nimport numpy as np\nimport pandas as pd\n\n# グラフ\nimport plotly.express as px\nimport plotly.graph_objs as go\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nimport math","metadata":{"execution":{"iopub.status.busy":"2022-02-24T14:14:53.346108Z","iopub.execute_input":"2022-02-24T14:14:53.346375Z","iopub.status.idle":"2022-02-24T14:14:55.087733Z","shell.execute_reply.started":"2022-02-24T14:14:53.346347Z","shell.execute_reply":"2022-02-24T14:14:55.086635Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import os\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n    for filename in filenames:\n        print(os.path.join(dirname, filename))\n\ninfo = pd.read_csv(\"/kaggle/input/g-research-crypto-forecasting/asset_details.csv\")\nctrain = pd.read_csv(\"/kaggle/input/g-research-crypto-forecasting/train.csv\")","metadata":{"execution":{"iopub.status.busy":"2022-02-24T14:14:55.088987Z","iopub.execute_input":"2022-02-24T14:14:55.089231Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Data Selection\ndef crypto_sub(asset_id ,data= ctrain ):\n    df = data[data[\"Asset_ID\"]==asset_id].reset_index(drop = True)\n    df['timestamp'] = pd.to_datetime(df['timestamp'], unit='s')\n    df = df.set_index('timestamp')\n    return df\n\n#Subgroup Bitcoin\nbtc1=crypto_sub(asset_id=1)\n\n#Subgroup Eth\neth1=crypto_sub(asset_id=6)\n\n#Subgroup Ada\nada1=crypto_sub(asset_id=3)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Function to Subset the Price variables\ndef mini_data(data):\n    df=data[[\"Close\",\"Low\",\"High\",\"Open\"]]\n    return df\n\nbtc1_mini=mini_data(btc1)\neth1_mini=mini_data(eth1)\nada1_mini=mini_data(ada1)\nbtc1.head()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 各通貨の周期性をチェック\nどの通貨においても周期性は確認できなかった=ARIMAモデルが適している","metadata":{}},{"cell_type":"code","source":"# Function to Plot\nplt.rcParams[\"figure.figsize\"]=(15,7)\n\ndef season_day(data):\n    df=data.resample(\"D\").mean()\n    seasonal_decompose(df.Close).plot()\n    return plt.show()\n\ndef season_month(data):\n    df=data.resample(\"M\").mean()\n    seasonal_decompose(df.Close).plot()\n    return plt.show()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"season_day(data=btc1_mini)\nseason_month(data=btc1_mini)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"season_day(data=eth1_mini)\nseason_month(data=eth1_mini)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"season_day(data=ada1_mini)\nseason_month(data=ada1_mini)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Box-Cox変換\n変数を正規分布に近しい形に変換する\n$$\n    y^{(\\lambda)} = (x^\\lambda - 1)/\\lambda \\quad (\\lambda \\neq 0) \\\\\n    y^{(\\lambda)} = ln(x) \\quad (\\lambda = 0)\n$$\n","metadata":{}},{"cell_type":"code","source":"btc_month=btc1_mini.resample(\"M\").mean()\nbtc_month[\"close_box\"], lmbda=stats.boxcox(btc_month.Close)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# p,qの値を0~3で設定し、一番良いmodelを探す\n# dは1で固定\nqs = range(0, 3)\nps = range(0, 3)\nd=1\nparameters_list = list(product(ps, qs))\n\nresults = []\nbest_aic = float(\"inf\")\nwarnings.filterwarnings('ignore')\nfor param in parameters_list:\n    # close_boxを説明変数にしてmodel学習\n    try:\n        model = SARIMAX(\n            btc_month.close_box,\n            order=(param[0], d, param[1]),\n            # seasonal_order = ()   # 周期性は見られなかったので入力しない=ARIMAになる\n        ).fit(disp=-1)\n    except ValueError:\n        print('bad parameter combination:', param)\n        continue\n    aic = model.aic\n    # 最も良いmodelとaicの値をbest_modelとbest_aicに格納する\n    # (値が小さいほど良い)\n    if aic < best_aic:\n        best_model = model\n        best_aic = aic\n        best_param = param\n    results.append([param, model.aic])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# modelを結果が良い順に並べる\nresult_table = pd.DataFrame(results)\nresult_table.columns = ['parameters', 'aic']\nprint(result_table.sort_values(by = 'aic', ascending=True).head())","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# 一番良い結果\nprint(best_model.summary())","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"best_param","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"best_model.plot_diagnostics(figsize=(15, 12))\nplt.show()\n# https://blog.amedama.jp/entry/sm-decompose-series","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# 逆Box-Cox変換\ndef invboxcox(y,lmbda):\n    if lmbda == 0:\n        return(np.exp(y))\n    else:\n        return(np.exp(np.log(lmbda*y+1)/lmbda))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Prediction\nbtc_month_pred = btc_month[['Close']]\ndate_list = [datetime(2021, 6, 20), datetime(2021, 7, 20), datetime(2021, 8, 21), datetime(2021, 9, 21)]\nfuture = pd.DataFrame(index=date_list, columns= btc_month.columns)\nbtc_month_pred = pd.concat([btc_month_pred, future])\nbtc_month_pred['forecast'] = invboxcox(best_model.predict(start=datetime(2018, 1, 31), end=datetime(2021, 5, 31)), lmbda)\nplt.figure(figsize=(18,10))\nbtc_month_pred.Close.plot()\nbtc_month_pred.forecast.plot(color='r', ls='--', label='Predicted Close')\nplt.legend()\nplt.title('Bitcoin monthly forecast')\nplt.ylabel('USD')\nplt.show()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# # 提出\n# env = gresearch_crypto.make_env()\n# iter_test = env.iter_test()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# df_test_entire = []\n\n# for i, (df_test, sample_prediction_df) in enumerate(iter_test):\n#     # Asset_IDでモデルを分割しているため、1データ毎で予測\n#     for j, test_row in df_test.iterrows():\n#         # test_rowのCloseにBox-Cox変換を掛ける\n#         test_row[\"close_box\"], lmbda=stats.boxcox(test_row[\"Close\"])\n        \n#         # close_boxに対してmodelを作成する\n#         model = SARIMAX(\n#             test_row[\"close_box\"],\n#             order=(best_param[0], d, best_param[1]),\n#         ).fit(disp=-1)\n#         model.predict(test_row)\n#         test_row[]\n    \n    \n#         test_row\n#         env.predict([])\n\n# import gresearch_crypto\n# env = gresearch_crypto.make_env()   # initialize the environment\n# iter_test = env.iter_test()    # an iterator which loops over the test set and sample submission\n# for (test_df, sample_prediction_df) in iter_test:\n#     sample_prediction_df['Target'] = 0  # make your predictions here\n#     env.predict(sample_prediction_df)   # register your predictions","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}