{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.12","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"gpu","dataSources":[{"sourceId":30894,"databundleVersionId":3147400,"sourceType":"competition"}],"dockerImageVersionId":30527,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"### CryptoCurrency (BTC-USD) Forecasting","metadata":{}},{"cell_type":"code","source":"import tensorflow as tf\nimport os\nimport warnings \nimport time\nimport pandas as pd\nimport random\nimport gc\nimport numpy as np\nimport time\nfrom sklearn.preprocessing import StandardScaler, Normalizer\nimport matplotlib.pyplot as plt\nfrom statsmodels.tools.eval_measures import meanabs, rmse\nfrom statsmodels.tools.eval_measures import meanabs\nfrom sklearn.preprocessing import MinMaxScaler\nimport gc\nimport math as mt\nfrom IPython.display import display\nimport tensorflow_probability as tfp\n\nfrom tensorflow_probability import bijectors as tfb\nfrom tensorflow_probability import distributions as tfd\n\nwarnings.filterwarnings('ignore')","metadata":{"execution":{"iopub.status.busy":"2024-10-15T09:31:42.136045Z","iopub.execute_input":"2024-10-15T09:31:42.136760Z","iopub.status.idle":"2024-10-15T09:31:52.120491Z","shell.execute_reply.started":"2024-10-15T09:31:42.136729Z","shell.execute_reply":"2024-10-15T09:31:52.119688Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Data Extraction and Preprocessing | With 1 month data (1 minute unit-stick)  ","metadata":{}},{"cell_type":"code","source":"df_asset_details = pd.read_csv('../input/g-research-crypto-forecasting/asset_details.csv')\nprint('df_asset_details: ',df_asset_details.shape)\ndf_train = pd.read_csv('../input/g-research-crypto-forecasting/train.csv').iloc[:350000] #~3 months worth of data-points\nprint('train_shape: ', df_train.shape)","metadata":{"execution":{"iopub.status.busy":"2024-10-15T09:31:56.323212Z","iopub.execute_input":"2024-10-15T09:31:56.323882Z","iopub.status.idle":"2024-10-15T09:33:07.608260Z","shell.execute_reply.started":"2024-10-15T09:31:56.323850Z","shell.execute_reply":"2024-10-15T09:33:07.607330Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"display(df_asset_details)\n\ntarget_coins = ['Bitcoin', 'Ethereum', 'Binance Coin']\n\ndf_asset_dict = {k : v for k, v in zip(df_asset_details['Asset_Name'].values, df_asset_details['Asset_ID'].values)}\ndf_btc = df_train[df_train['Asset_ID'] == df_asset_dict[target_coins[0]]]\ndf_eth = df_train[df_train['Asset_ID'] == df_asset_dict[target_coins[1]]]\ndf_bnb = df_train[df_train['Asset_ID'] == df_asset_dict[target_coins[-1]]]\n\nprint('BTC-Len: ', len(df_btc), ', ETH-Len: ', len(df_eth), ', BNB-Len: ', len(df_bnb))\n\nprint()\ndisplay(df_btc.head(5))","metadata":{"execution":{"iopub.status.busy":"2024-10-15T09:33:13.039702Z","iopub.execute_input":"2024-10-15T09:33:13.040064Z","iopub.status.idle":"2024-10-15T09:33:13.094183Z","shell.execute_reply.started":"2024-10-15T09:33:13.040033Z","shell.execute_reply":"2024-10-15T09:33:13.093194Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig, axs = plt.subplots(1, 1, figsize = (15, 4))\n\naxs.plot(df_btc.Close, ',r', label = 'Close')\naxs.plot(df_btc.High, ',b', alpha = 0.2, label = 'High')\naxs.plot(df_btc.Low, ',y', alpha = 0.2, label = 'Low')\n\naxs.set_title(target_coins[0] + ' (H/C/L) | 1min-stick ')\n\nplt.legend()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-10-15T09:33:15.669623Z","iopub.execute_input":"2024-10-15T09:33:15.669966Z","iopub.status.idle":"2024-10-15T09:33:16.088480Z","shell.execute_reply.started":"2024-10-15T09:33:15.669939Z","shell.execute_reply":"2024-10-15T09:33:16.087561Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The prior slicing of the total sample size allows us to compare predictions with OOD close prices.","metadata":{}},{"cell_type":"code","source":"total_sample_size = len(df_btc)\nprint('total_sample_size: ', total_sample_size)\n\npred_sample_size = 5000\ntrain_sample_size = total_sample_size - pred_sample_size\n\nprint('train_sample_size: ', train_sample_size, ', pred_sample_size: ', pred_sample_size)","metadata":{"execution":{"iopub.status.busy":"2024-10-15T09:33:19.698869Z","iopub.execute_input":"2024-10-15T09:33:19.700099Z","iopub.status.idle":"2024-10-15T09:33:19.705836Z","shell.execute_reply.started":"2024-10-15T09:33:19.700053Z","shell.execute_reply":"2024-10-15T09:33:19.704872Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_btc_train, df_btc_test = df_btc.iloc[:train_sample_size], df_btc.iloc[train_sample_size+1:]\nprint(df_btc_train.shape, df_btc_test.shape)\n\nfig, axs = plt.subplots(1, 1, figsize = (8, 3))\n\naxs.plot(df_btc_train.index, df_btc_train.Close, ',b', label = 'train')\naxs.plot(df_btc_test.index, df_btc_test.Close, ',r', label = 'OOD test-set')\naxs.set_title('BTC-Close-Price')\n\nplt.legend()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-10-15T09:33:21.924824Z","iopub.execute_input":"2024-10-15T09:33:21.925600Z","iopub.status.idle":"2024-10-15T09:33:22.222485Z","shell.execute_reply.started":"2024-10-15T09:33:21.925562Z","shell.execute_reply":"2024-10-15T09:33:22.221571Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#### Custom-scaling","metadata":{}},{"cell_type":"code","source":"def min_max_scale(data, new):\n    \n    new_min, new_max = new\n    \n    min_values = np.min(data, axis=0)\n    max_values = np.max(data, axis=0)\n\n    scaled_data = (data - min_values) * (new_max - new_min) / (max_values - min_values) + new_min\n\n    return scaled_data, min_values, max_values\n\ndef descale(scaled_data, new_min, new_max):\n    \n    min_values = np.min(scaled_data, axis=0)\n    max_values = np.max(scaled_data, axis=0)\n\n    original_data = (scaled_data - min_values) * (new_max - new_min) / (max_values - min_values) + new_min\n    \n    return original_data","metadata":{"execution":{"iopub.status.busy":"2024-10-13T20:28:49.206181Z","iopub.execute_input":"2024-10-13T20:28:49.206511Z","iopub.status.idle":"2024-10-13T20:28:49.213343Z","shell.execute_reply.started":"2024-10-13T20:28:49.206486Z","shell.execute_reply":"2024-10-13T20:28:49.212278Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"x_obs_close = df_btc_train.Close.to_numpy()\nx_obs_close_test = df_btc_test.Close.to_numpy()\n\nx_obs_all = df_btc_train[['Count', 'Open', 'High', 'Low', 'Volume', 'Close']].to_numpy()\nx_obs_all_test = df_btc_test[['Count', 'Open', 'High', 'Low', 'Volume', 'Close']].to_numpy()\n\nprint('Close-Price: ', x_obs_close.shape, x_obs_close_test.shape)\nprint('Multi-Feature-Space: ', x_obs_all.shape, x_obs_all_test.shape)\n\ntrain_split, test_split = 0.70, 0.30\ntrain_max = round(train_split * x_obs_close.shape[0])\nval_max = (1 - train_split) * x_obs_close.shape[0]\ntest_max = round(test_split * val_max)\nval_max = round(val_max - test_max)","metadata":{"execution":{"iopub.status.busy":"2024-10-15T09:33:34.691973Z","iopub.execute_input":"2024-10-15T09:33:34.692685Z","iopub.status.idle":"2024-10-15T09:33:34.703094Z","shell.execute_reply.started":"2024-10-15T09:33:34.692652Z","shell.execute_reply":"2024-10-15T09:33:34.702130Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#### PreProcessing","metadata":{}},{"cell_type":"code","source":"scaler_close = StandardScaler()\nscaler_all = StandardScaler()\n\nnorm_close = Normalizer()\nnorm_all = Normalizer()","metadata":{"execution":{"iopub.status.busy":"2024-10-15T09:33:36.718105Z","iopub.execute_input":"2024-10-15T09:33:36.718918Z","iopub.status.idle":"2024-10-15T09:33:36.723441Z","shell.execute_reply.started":"2024-10-15T09:33:36.718870Z","shell.execute_reply":"2024-10-15T09:33:36.722334Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"### Need to perform min-max scaling on the entire dataspace with the train_set distribution\n### information only\n### Assuming Lower Volatitly Levels for BTC\n\nnew = 0, 1\n\nprint(f'min: {min(x_obs_close)}, max: {max(x_obs_close)}')\n\n# x_obs_close_scaled, min_og, max_og = min_max_scale(x_obs_close, new)\n\nx_obs_close_scaled = scaler_close.fit_transform(x_obs_close.reshape(-1, 1))\nx_obs_close_scaled = x_obs_close_scaled.reshape(-1, 1)\n\n# x_obs_close_scaled = norm_close.fit_transform(x_obs_close.reshape(-1, 1))\n# x_obs_close_scaled = x_obs_close_scaled.reshape(-1, 1)\n\n# x_obs_close_test_scaled, min_og, max_og = min_max_scale(x_obs_close_test, new)\n# x_obs_close_test_scaled = x_obs_close_test_scaled.reshape(-1, 1)\n\nprint('close-train-shape: ', x_obs_close_scaled.shape)\nprint(f'scaled-min: {min(x_obs_close_scaled)}, scaled-max: {max(x_obs_close_scaled)}')\n\n# x_obs_descaled = descale(x_obs_close_scaled, min_og, max_og)\n\nx_obs_descaled = scaler_close.inverse_transform(x_obs_close_scaled)\nprint('Descaled: ', np.min(x_obs_descaled), ', ', np.max(x_obs_descaled))\ndel x_obs_descaled\ngc.collect()\n\n'----------------------------------------------------------------------'\n\n''' For Multi-Feature dataspace might need to use event-based normalization - So need to add a inside layer instead of preprocessing '''\n\n# x_obs_all_scaled, min_og, max_og = min_max_scale(x_obs_all, new)\n\n# x_obs_all_scaled = scaler_all.fit_transform(x_obs_all)\n# x_obs_all_scaled = x_obs_all_scaled.reshape(-1, x_obs_all.shape[-1])\n\n# x_obs_all_scaled = norm_all.fit_transform(x_obs_all)\n# x_obs_all_scaled = x_obs_all_scaled.reshape(-1, x_obs_all.shape[-1])\n\n# x_obs_all_test_scaled, min_og, max_og = min_max_scale(x_obs_all_test, new)\n# x_obs_all_test_scaled = x_obs_all_test_scaled.reshape(-1, x_obs_all.shape[-1])\n\nx_obs_all_scaled = x_obs_all\n\nprint()\nprint('multi-train-shape: ', x_obs_all_scaled.shape)\nprint(f'scaled-min: {min(np.ravel(x_obs_all_scaled[:,-1]))}, scaled-max: {max(np.ravel(x_obs_all_scaled[:, -1]))}')\n\n# x_obs_descaled = descale(x_obs_all_scaled, min_og, max_og)\n\n# x_obs_all_descaled = scaler_all.inverse_transform(x_obs_all_scaled)\n# print('Descaled: ', np.min(x_obs_all_descaled), ', ', np.max(x_obs_all_descaled))\n# del x_obs_all_descaled\n# gc.collect()\n\nfig, axs = plt.subplots(1, 2, figsize = (10, 3))\n\naxs[0].plot(x_obs_close_scaled, ',b')\naxs[0].set_title('Min-Max Normalized (ONLY)')\n\naxs[1].plot(x_obs_all_scaled[:, -1], ',r')\naxs[1].set_title('Min-Max Normalized (ALL)')\n\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-10-15T09:33:42.540051Z","iopub.execute_input":"2024-10-15T09:33:42.540405Z","iopub.status.idle":"2024-10-15T09:33:43.260944Z","shell.execute_reply.started":"2024-10-15T09:33:42.540375Z","shell.execute_reply":"2024-10-15T09:33:43.260015Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"seed = 28\nnp.random.seed(seed)\n\n\ndef Data_func(data, past, delay, minn, maxx = None, shuffle=False, batch_size=128, rate=6):\n    if maxx is None:\n        maxx = len(data) - delay -1\n    i = minn + past\n    batch_no = 1\n\n    while True:\n\n        if shuffle:\n            rows = np.random.randint(minn + past, maxx, size = batch_size)\n        else:\n            rows = np.arange(i, min(i + batch_size, maxx))\n        \n        i+=len(rows)\n        if i + batch_size > maxx: ## when you increase the past too much this statement gets executed at the first batch creation // add more data to overcome\n            break\n        \n        samples = np.zeros((len(rows), past // rate, data.shape[-1]))\n        targets = np.zeros((len(rows),))\n        cbase_preds = np.zeros((len(rows),))\n    \n        for j, row in enumerate(rows):\n            \n            indices = range(rows[j] - past, rows[j], rate) ## take every point\n            samples[j] = data[indices]\n            targets[j] = data[rows[j] + delay - 1][-1] ### acessing Close price HERE!\n       \n        if batch_no == 1:\n            new_sample = samples\n            new_targets = targets\n            print('Ini_samples: ', samples.shape, 'Ini_targets: ', targets.shape)\n        else:\n            new_sample = np.concatenate([new_sample, samples], axis=0)\n            new_targets = np.concatenate([new_targets, targets], axis=0)\n        batch_no += 1\n        \n    return new_sample, new_targets\n\nbatch_size = 32\npast = 360 #  6 hours in past {Localised}\ndelay = 0 # After 15 minute into the future\nsample_rate = 5\n\n## sampling data every (sample_rate) minute | Information Flow into the <val> and <test> set | Transforming from unit-minute chart to sample_rate minute chart\ntrain_data_close, train_target_close = Data_func(x_obs_close_scaled, past, delay, minn = 0, maxx = train_max - 1, batch_size = batch_size, shuffle=True, rate = sample_rate)\nval_data_close, val_target_close = Data_func(x_obs_close_scaled, past, delay, minn = train_max + 1, maxx = train_max + val_max - 1, batch_size = batch_size, rate = sample_rate)\ntest_data_close, test_target_close = Data_func(x_obs_close_scaled, past, delay, minn = train_max + val_max + 1, maxx = train_max + val_max + test_max, batch_size = batch_size, rate = sample_rate)\n\n## sampling data every (sample_rate) minute | Information Flow into the <val> and <test> set | Transforming from unit-minute chart to sample_rate minute chart\ntrain_data_all, train_target_all = Data_func(x_obs_all_scaled, past, delay, minn = 0, maxx = train_max - 1, batch_size = batch_size, shuffle=True, rate = sample_rate)\nval_data_all, val_target_all = Data_func(x_obs_all_scaled, past, delay, minn = train_max + 1, maxx = train_max + val_max - 1, batch_size = batch_size, rate = sample_rate)\ntest_data_all, test_target_all = Data_func(x_obs_all_scaled, past, delay, minn = train_max + val_max + 1, maxx = train_max + val_max + test_max, batch_size = batch_size, rate = sample_rate)\n\nprint()\nprint(train_data_close.shape, train_target_close.shape, val_data_close.shape, val_target_close.shape, test_data_close.shape, test_target_close.shape)\nprint('Min-Max-Range')\nprint(min(np.ravel(train_data_close)), max(np.ravel(train_target_close)))\nprint(min(np.ravel(val_data_close)), max(np.ravel(val_target_close)))\nprint(min(np.ravel(test_data_close)), max(np.ravel(test_target_close)))\n\nprint()\nprint(train_data_all.shape, train_target_all.shape, val_data_all.shape, val_target_all.shape, test_data_all.shape, test_target_all.shape)\nprint('Min-Max-Range')\nprint(min(np.ravel(train_data_all)), max(np.ravel(train_target_all)))\nprint(min(np.ravel(val_data_all)), max(np.ravel(val_target_all)))\nprint(min(np.ravel(test_data_all)), max(np.ravel(test_target_all)))","metadata":{"execution":{"iopub.status.busy":"2024-10-15T09:33:52.579458Z","iopub.execute_input":"2024-10-15T09:33:52.579814Z","iopub.status.idle":"2024-10-15T09:34:12.844023Z","shell.execute_reply.started":"2024-10-15T09:33:52.579783Z","shell.execute_reply":"2024-10-15T09:34:12.843078Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#### 2-Step PreProcessing - Normalization","metadata":{}},{"cell_type":"code","source":"def inv_norm(test_target):\n    test_target = np.ravel(test_target)\n    test_target = test_target * std[-1]\n    test_target = test_target + mean[-1]\n    return test_target","metadata":{"execution":{"iopub.status.busy":"2023-07-27T19:23:11.153969Z","iopub.execute_input":"2023-07-27T19:23:11.154825Z","iopub.status.idle":"2023-07-27T19:23:11.161847Z","shell.execute_reply.started":"2023-07-27T19:23:11.154786Z","shell.execute_reply":"2023-07-27T19:23:11.161068Z"},"jupyter":{"source_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sample_seed = 5\n\nfig, axs = plt.subplots(1, 3, figsize = (15, 6))\n\naxs[1].plot(np.arange(0, len(train_data[sample_seed])), np.squeeze(train_data[sample_seed]), 'k')\naxs[1].set_title(' - Scaling (1st PreProcess) - ')\n\naxs[0].plot(np.arange(0, len(train_data[sample_seed])), inv_norm(np.squeeze(train_data[sample_seed])), 'k')\naxs[0].set_title(' - Normalized (2nd PreProcess) - ')\n\naxs[2].plot(np.arange(0, len(train_data[sample_seed])), descale(inv_norm(np.squeeze(train_data[sample_seed])), min_og, max_og), 'k')\naxs[2].set_title(' - De-Scaled - ')\n\nplt.tight_layout(pad = 2)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-10-12T14:19:02.849204Z","iopub.execute_input":"2024-10-12T14:19:02.849670Z","iopub.status.idle":"2024-10-12T14:19:03.592835Z","shell.execute_reply.started":"2024-10-12T14:19:02.849636Z","shell.execute_reply":"2024-10-12T14:19:03.590982Z"},"jupyter":{"source_hidden":true,"outputs_hidden":true},"collapsed":true,"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Statistical-TS Workflow <br>","metadata":{}},{"cell_type":"code","source":"print('Train-Len: ', len(df_btc_train.Close))\nprint('Test-Len: ', len(df_btc_test.Close))","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:39:30.198934Z","iopub.execute_input":"2024-10-15T11:39:30.199297Z","iopub.status.idle":"2024-10-15T11:39:30.205016Z","shell.execute_reply.started":"2024-10-15T11:39:30.199269Z","shell.execute_reply":"2024-10-15T11:39:30.203821Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#### Local-Linear Trend Model (statsmodels) --v1","metadata":{}},{"cell_type":"code","source":"x_sample_train = df_btc_train.Close.iloc[-1500:] \nx_sample_test = df_btc_test.Close.iloc[:50]\n\nprint(x_sample_train.shape, x_sample_test.shape)","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:39:31.795833Z","iopub.execute_input":"2024-10-15T11:39:31.796781Z","iopub.status.idle":"2024-10-15T11:39:31.802063Z","shell.execute_reply.started":"2024-10-15T11:39:31.796744Z","shell.execute_reply":"2024-10-15T11:39:31.801105Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import pandas as pd\nimport statsmodels.api as sm\nimport matplotlib.pyplot as plt\n\nlocal_linear_model = sm.tsa.UnobservedComponents(x_sample_train, level='local linear trend')\nlocal_linear_result = local_linear_model.fit()\n\nprint(local_linear_result.summary())\n\nforecast_local_linear = local_linear_result.get_forecast(steps = len(x_sample_test))\nforecast_local_linear_values = forecast_local_linear.predicted_mean\n\nplt.figure(figsize=(8, 4))\n\nxlin_1 = np.arange(0, len(x_sample_train))\nxlin_2 = np.arange(len(x_sample_train), len(x_sample_train) + len(x_sample_test))\n\nplt.plot(xlin_1, x_sample_train, label='Observed Close Prices', color = 'black')\nplt.plot(xlin_2, forecast_local_linear_values, label='Forecasted Prices (Local Linear Trend)', color='red')\nplt.plot(xlin_2, x_sample_test, label = 'True Forecast', color = 'blue')\nplt.title('Local Linear Trend Model - BTC Close Price Prediction')\nplt.legend()\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:39:33.516318Z","iopub.execute_input":"2024-10-15T11:39:33.517051Z","iopub.status.idle":"2024-10-15T11:39:35.853165Z","shell.execute_reply.started":"2024-10-15T11:39:33.517016Z","shell.execute_reply":"2024-10-15T11:39:35.852195Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":" #### Local-Linear Trend Model --v2","metadata":{}},{"cell_type":"code","source":"x_obs_scaled = x_obs_scaled[:200].reshape(-1, 1)\nprint(x_obs_scaled.shape)","metadata":{"execution":{"iopub.status.busy":"2023-07-27T21:02:04.232433Z","iopub.execute_input":"2023-07-27T21:02:04.233251Z","iopub.status.idle":"2023-07-27T21:02:04.240129Z","shell.execute_reply.started":"2023-07-27T21:02:04.233205Z","shell.execute_reply":"2023-07-27T21:02:04.239267Z"},"jupyter":{"source_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"num_timesteps = len(x_obs_scaled)\n\ndef build_model(approximate_unconstrained_rates):\n    trend = tfp.sts.LocalLinearTrend(\n      observed_time_series=approximate_unconstrained_rates)\n    return tfp.sts.Sum([trend],\n                     observed_time_series=approximate_unconstrained_rates)\n\npositive_bijector = tfb.Softplus()  \n\napproximate_unconstrained_rates = positive_bijector.inverse(tf.convert_to_tensor(x_obs_scaled))\nsts_model = build_model(approximate_unconstrained_rates)\n\ndef sts_with_poisson_likelihood_model():\n    # Encode the parameters of the STS model as random variables.\n    param_vals = []\n    for param in sts_model.parameters:\n        param_val = yield param.prior\n        param_vals.append(param_val)\n\n    ''' WE yield only those R.V whose Joint-distrib we want to find'''    \n        \n    ## instantiate this model as a distribution over the given time-steps\n    unconstrained_rate = yield sts_model.make_state_space_model(num_timesteps, param_vals)\n    rate = positive_bijector.forward(unconstrained_rate[..., 0])\n    obs = yield tfd.Poisson(rate, name='observed_counts')\n\nmodel = tfd.JointDistributionCoroutineAutoBatched(sts_with_poisson_likelihood_model)","metadata":{"execution":{"iopub.status.busy":"2023-07-27T21:02:31.535421Z","iopub.execute_input":"2023-07-27T21:02:31.535713Z","iopub.status.idle":"2023-07-27T21:02:31.569551Z","shell.execute_reply.started":"2023-07-27T21:02:31.535681Z","shell.execute_reply":"2023-07-27T21:02:31.568818Z"},"jupyter":{"source_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"pinned_model = model.experimental_pin(observed_counts = x_obs_scaled)\nconstraining_bijector = pinned_model.experimental_default_event_space_bijector()","metadata":{"execution":{"iopub.status.busy":"2023-07-27T21:02:34.005159Z","iopub.execute_input":"2023-07-27T21:02:34.005448Z","iopub.status.idle":"2023-07-27T21:02:44.232328Z","shell.execute_reply.started":"2023-07-27T21:02:34.005416Z","shell.execute_reply":"2023-07-27T21:02:44.231512Z"},"jupyter":{"source_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Allow external control of sampling to reduce test runtimes.\nnum_results = 500 \nnum_results = int(num_results)\n\nnum_burnin_steps = 100 \nnum_burnin_steps = int(num_burnin_steps)\n\nsampler = tfp.mcmc.TransformedTransitionKernel(\n    tfp.mcmc.NoUTurnSampler(\n        target_log_prob_fn=pinned_model.unnormalized_log_prob,\n        step_size=0.1),\n    bijector=constraining_bijector)\n\nadaptive_sampler = tfp.mcmc.DualAveragingStepSizeAdaptation(\n    inner_kernel=sampler,\n    num_adaptation_steps=int(0.8 * num_burnin_steps),\n    target_accept_prob=0.75)\n\ninitial_state = constraining_bijector.forward(\n    type(pinned_model.event_shape)(\n        *(tf.random.normal(part_shape)\n          for part_shape in constraining_bijector.inverse_event_shape(\n              pinned_model.event_shape))))\n\n# Speed up sampling by tracing with `tf.function`.\n# @tf.function(autograph=False, jit_compile=True)\ndef do_sampling():\n    return tfp.mcmc.sample_chain(\n      kernel=adaptive_sampler,\n      current_state=initial_state,\n      num_results=num_results,\n      num_burnin_steps=num_burnin_steps,\n      trace_fn=None)\n\nt0 = time.time()\nsamples = do_sampling()\nt1 = time.time()\nprint(\"Inference ran in {:.2f}s.\".format(t1-t0))\n","metadata":{"execution":{"iopub.status.busy":"2023-07-27T21:05:58.846328Z","iopub.execute_input":"2023-07-27T21:05:58.846632Z","iopub.status.idle":"2023-07-27T21:05:58.852046Z","shell.execute_reply.started":"2023-07-27T21:05:58.846597Z","shell.execute_reply":"2023-07-27T21:05:58.851023Z"},"jupyter":{"source_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"f = plt.figure(figsize=(12, 4))\nfor i, param in enumerate(sts_model.parameters):\n    ax = f.add_subplot(1, len(sts_model.parameters), i + 1)\n    ax.plot(samples[i])\n    ax.set_title(\"{} samples\".format(param.name))","metadata":{"execution":{"iopub.status.busy":"2023-07-27T21:03:03.925497Z","iopub.status.idle":"2023-07-27T21:03:03.926462Z","shell.execute_reply.started":"2023-07-27T21:03:03.926208Z","shell.execute_reply":"2023-07-27T21:03:03.926234Z"},"jupyter":{"source_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"param_samples = samples[:-1]\nunconstrained_rate_samples = samples[-1][..., 0]\nrate_samples = positive_bijector.forward(unconstrained_rate_samples)\n\nplt.figure(figsize=(10, 4))\nmean_lower, mean_upper = np.percentile(rate_samples, [10, 90], axis=0)\npred_lower, pred_upper = np.percentile(np.random.poisson(rate_samples), [10, 90], axis=0)\n\n_ = plt.plot(x_obs_scaled, color=\"blue\", ls='-', label='observed', alpha=0.7)\n_ = plt.plot(np.mean(rate_samples, axis=0), label='rate', color=\"green\", ls='dashed', lw=2, alpha=0.7)\n_ = plt.fill_between(np.arange(0, 150), mean_lower, mean_upper, color='green', alpha=0.2)\n_ = plt.fill_between(np.arange(0, 150), pred_lower, pred_upper, color='grey', label='obs', alpha=0.2)\nplt.xlabel(\"Day\")\nplt.ylabel(\"Daily Sample Size\")\nplt.title(\"Posterior Mean\")\nplt.legend()","metadata":{"execution":{"iopub.status.busy":"2023-07-27T20:55:59.268679Z","iopub.status.idle":"2023-07-27T20:55:59.269298Z","shell.execute_reply.started":"2023-07-27T20:55:59.269020Z","shell.execute_reply":"2023-07-27T20:55:59.269047Z"},"jupyter":{"source_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Forecast","metadata":{}},{"cell_type":"code","source":"def sample_forecasted_counts(sts_model, posterior_latent_rates,\n                             posterior_params, num_steps_forecast,\n                             num_sampled_forecasts):\n\n    # Forecast the future latent unconstrained rates, given the inferred latent\n    # unconstrained rates and parameters.\n    unconstrained_rates_forecast_dist = tfp.sts.forecast(sts_model,\n    observed_time_series=unconstrained_rate_samples,\n    parameter_samples=posterior_params,\n    num_steps_forecast=num_steps_forecast)\n\n    # Transform the forecast to positive-valued Poisson rates.\n    rates_forecast_dist = tfd.TransformedDistribution(\n      unconstrained_rates_forecast_dist,\n      positive_bijector)\n\n    # Sample from the forecast model following the chain rule:\n    # P(counts) = P(counts | latent_rates)P(latent_rates)\n    sampled_latent_rates = rates_forecast_dist.sample(num_sampled_forecasts)\n    sampled_forecast_counts = tfd.Poisson(rate=sampled_latent_rates).sample()\n\n    return sampled_forecast_counts, sampled_latent_rates\n\nforecast_samples, rate_samples = sample_forecasted_counts(\n   sts_model,\n   posterior_latent_rates=unconstrained_rate_samples,\n   posterior_params=param_samples,\n   # Days to forecast:\n   num_steps_forecast=30,\n   num_sampled_forecasts=100)\n\nforecast_samples = np.squeeze(forecast_samples)","metadata":{"execution":{"iopub.status.busy":"2023-07-27T20:55:59.270536Z","iopub.status.idle":"2023-07-27T20:55:59.271171Z","shell.execute_reply.started":"2023-07-27T20:55:59.270894Z","shell.execute_reply":"2023-07-27T20:55:59.270933Z"},"jupyter":{"source_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def plot_forecast_helper(data, forecast_samples, CI=90):\n    \n    plt.figure(figsize=(10, 4))\n    forecast_median = np.median(forecast_samples, axis=0)\n\n    num_steps = len(data)\n    num_steps_forecast = forecast_median.shape[-1]\n\n    plt.plot(np.arange(num_steps), data, lw=2, color='blue', linestyle='-', label='Observed Data', alpha=0.7)\n\n    forecast_steps = np.arange(num_steps, num_steps+num_steps_forecast)\n\n    CI_interval = [(100 - CI)/2, 100 - (100 - CI)/2]\n    lower, upper = np.percentile(forecast_samples, CI_interval, axis=0)\n\n    plt.plot(forecast_steps, forecast_median, lw=2, ls='-', color='orange', label=str(CI) + '% Forecast Interval', alpha=0.7)\n    plt.fill_between(forecast_steps,\n                       lower,\n                       upper, color='orange', alpha=0.2)\n\n    plt.xlim([0, num_steps+num_steps_forecast])\n    ymin, ymax = min(np.min(forecast_samples), np.min(data)), max(np.max(forecast_samples), np.max(data))\n    yrange = ymax-ymin\n    plt.title(\"{}\".format('Observed time series with ' + str(num_steps_forecast) + ' Day Forecast'))\n    plt.xlabel('Day')\n    plt.ylabel('Daily Sample Size')\n    plt.ylim([0, max(x_obs_scaled) + 50])\n    ","metadata":{"execution":{"iopub.status.busy":"2023-07-27T20:55:59.272354Z","iopub.status.idle":"2023-07-27T20:55:59.272943Z","shell.execute_reply.started":"2023-07-27T20:55:59.272688Z","shell.execute_reply":"2023-07-27T20:55:59.272720Z"},"jupyter":{"source_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_forecast_helper(x_obs_scaled, forecast_samples, CI=20) # 20% Confidence-Interval for forecast\nplt.plot((len(x_obs_scaled)) + np.arange(0, len(x_obs_test_scaled)), x_obs_test_scaled, 'r--', label = 'actual-forecast')\nplt.legend()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-07-27T20:55:59.274102Z","iopub.status.idle":"2023-07-27T20:55:59.274717Z","shell.execute_reply.started":"2023-07-27T20:55:59.274471Z","shell.execute_reply":"2023-07-27T20:55:59.274498Z"},"jupyter":{"source_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#### Local Linear Trend Model --v3","metadata":{}},{"cell_type":"code","source":"import tensorflow_probability as tfp\nimport tensorflow as tf\n\n# Define a Local Linear Trend model\nlocal_linear_trend = tfp.sts.LocalLinearTrend(\n    observed_time_series=time_series_data)  # Pass your time-series data here\n\n# Build the model\nmodel = tfp.sts.Sum([local_linear_trend], observed_time_series=time_series_data)\n\n# Fit the model using variational inference (to find the posterior distribution of parameters)\nvariational_posteriors = tfp.sts.build_factored_surrogate_posterior(model=model)\n\n# Train the model using variational inference\nnum_variational_steps = 200  # Number of iterations for optimization\noptimizer = tf.optimizers.Adam(learning_rate=0.1)\n\n@tf.function\ndef train():\n    elbo_loss_curve = tfp.vi.fit_surrogate_posterior(\n        target_log_prob_fn=model.joint_log_prob(observed_time_series=time_series_data),\n        surrogate_posterior=variational_posteriors,\n        optimizer=optimizer,\n        num_steps=num_variational_steps)\n    return elbo_loss_curve\n\nelbo_loss = train()\n\n# Forecast future values\nnum_forecast_steps = 30  # Predict 30 future steps\nforecast_dist = tfp.sts.forecast(\n    model=model,\n    observed_time_series=time_series_data,\n    parameter_samples=variational_posteriors.sample(50),\n    num_steps_forecast=num_forecast_steps)\n\nforecast_mean = forecast_dist.mean().numpy()  # Mean of forecast\nforecast_stddev = forecast_dist.stddev().numpy()  # Standard deviation for uncertainty\n","metadata":{"jupyter":{"source_hidden":true}},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#### ARIMA model","metadata":{}},{"cell_type":"code","source":"import pandas as pd\nfrom statsmodels.tsa.arima.model import ARIMA\nimport matplotlib.pyplot as plt\n\narima_model = ARIMA(x_sample_train, order = (2, 1, 2))\narima_result = arima_model.fit()\n\nprint(arima_result.summary())\n\nforecast_arima = arima_result.get_forecast(steps = len(x_sample_test))\nforecast_arima_values = forecast_arima.predicted_mean\n\nxlin_1 = np.arange(0, len(x_sample_train))\nxlin_2 = np.arange(len(x_sample_train), len(x_sample_train) + len(x_sample_test))\n\nplt.figure(figsize=(8, 4))\nplt.plot(xlin_1, x_sample_train, label='Observed Close Prices', color = 'black')\nplt.plot(xlin_2, forecast_arima_values, label='Forecasted Prices (ARIMA)', color='red')\nplt.plot(xlin_2, x_sample_test, label = 'True Forecast', color = 'blue')\nplt.title('ARIMA Model - BTC Close Price Prediction')\nplt.legend()\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:39:41.192635Z","iopub.execute_input":"2024-10-15T11:39:41.193414Z","iopub.status.idle":"2024-10-15T11:39:41.882910Z","shell.execute_reply.started":"2024-10-15T11:39:41.193377Z","shell.execute_reply":"2024-10-15T11:39:41.881930Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"----","metadata":{}},{"cell_type":"markdown","source":"### Machine Learning Models (ANN and RNN)","metadata":{}},{"cell_type":"markdown","source":"#### Artificial/Dense Neural Nets (A/DNN) ","metadata":{}},{"cell_type":"code","source":"print(test_data_close.shape, test_target_close.shape)\nprint(test_data_all.shape, test_target_all.shape)\n\nprint(' -- One-Sample-Test -- ')\n\nrn = random.randint(0, len(test_data_all))\nprint('Random: ', rn)\n\nfig, axs = plt.subplots(1, 2, figsize = (10, 5))\n\nidx_range = np.arange(len(test_data_close[rn]))\naxs[0].plot(idx_range, test_data_close[rn], 'k')\naxs[0].plot(idx_range[-1] + delay + 1, test_target_close[rn], '.r')\naxs[0].set_title('Close (only)')\n\nprint(test_target_all[rn])\nidx_range = np.arange(len(test_data_all[rn][:, -1]))\naxs[1].plot(idx_range, test_data_all[rn][:, -1], 'k')\naxs[1].plot(idx_range[-1] + delay + 1, test_target_all[rn], '.r')\naxs[1].set_title('Close (All)')\n\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:39:45.232458Z","iopub.execute_input":"2024-10-15T11:39:45.232941Z","iopub.status.idle":"2024-10-15T11:39:45.731392Z","shell.execute_reply.started":"2024-10-15T11:39:45.232908Z","shell.execute_reply":"2024-10-15T11:39:45.730331Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from tensorflow.keras.layers import Normalization\n\nkeras_normalizer = Normalization()\nkeras_normalizer.adapt(train_data_all)\n\ndnn_model = tf.keras.models.Sequential([\n    tf.keras.layers.Input(shape=(72, train_data_all.shape[-1]), name = 'Input'),\n    keras_normalizer,\n    tf.keras.layers.Dense(64, activation = 'linear'),\n    tf.keras.layers.Dense(128, activation='relu'),\n    tf.keras.layers.Dense(256, activation = 'relu'),\n    tf.keras.layers.Dense(128, activation = 'relu'),\n    tf.keras.layers.Dense(64, activation = 'relu'),\n    tf.keras.layers.Dense(32, activation = 'relu'),\n    tf.keras.layers.Flatten(),\n    tf.keras.layers.Dense(1, activation = 'linear')\n    \n])\n# es = tf.keras.callbacks.EarlyStopping(monitor='val_loss', patience = 10)\ndnn_model.compile(optimizer='Adam', loss='mse', metrics=['mae'])\n\ndnn_model.summary()","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:39:48.553129Z","iopub.execute_input":"2024-10-15T11:39:48.553517Z","iopub.status.idle":"2024-10-15T11:39:50.067855Z","shell.execute_reply.started":"2024-10-15T11:39:48.553482Z","shell.execute_reply":"2024-10-15T11:39:50.066950Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ft = time.time()\ndnn_hist = dnn_model.fit(train_data_all, train_target_all, epochs=50, validation_data = (val_data_all, val_target_all))#, callbacks=[es])\nlt = time.time()\n\nprint('Training_Done ...  ')\nprint('Execution Time: ', (lt - ft)/60)","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:39:53.801504Z","iopub.execute_input":"2024-10-15T11:39:53.801906Z","iopub.status.idle":"2024-10-15T11:42:41.395939Z","shell.execute_reply.started":"2024-10-15T11:39:53.801871Z","shell.execute_reply":"2024-10-15T11:42:41.395019Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"dnn_history = dnn_hist.history\nepochs_arr = np.arange(len(dnn_history['loss']))\n\nplt.plot(epochs_arr[5:], dnn_history['loss'][5:], label = 'train_loss', color = 'red', marker = '.')\nplt.plot(epochs_arr[5:], dnn_history['val_loss'][5:], label = 'val_loss', color = 'blue', marker = '.')\nplt.title('Loss-Curve')\n\nplt.legend()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:43:04.514874Z","iopub.execute_input":"2024-10-15T11:43:04.515276Z","iopub.status.idle":"2024-10-15T11:43:04.728514Z","shell.execute_reply.started":"2024-10-15T11:43:04.515238Z","shell.execute_reply":"2024-10-15T11:43:04.727524Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test_mae_err = dnn_model.evaluate(test_data_all, test_target_all, steps = 512)[-1]\nprint(f'The Absolute Error in Prediction of Close Price of BTC is {test_mae_err} USD, which is more than std-dev of Close Price in {sample_rate} minute time-frame chart')","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:43:07.552878Z","iopub.execute_input":"2024-10-15T11:43:07.553258Z","iopub.status.idle":"2024-10-15T11:43:10.147618Z","shell.execute_reply.started":"2024-10-15T11:43:07.553230Z","shell.execute_reply":"2024-10-15T11:43:10.146586Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"rn = random.randint(0, len(test_data_all))\nprint('random-seed: ', rn)\n# rn = 2306\n\ntest_sample, test_target = test_data_all[rn].reshape(1, -1, test_data_all.shape[-1]), test_target_all[rn].reshape(-1, 1) \nprint(test_sample.shape, test_target.shape)\n\ntest_pred = dnn_model.predict(test_sample)\n\n# test_sample_descaled = np.squeeze(descale(test_sample, min_og, max_og))[:, -2]\n# test_pred_descaled = descale(test_pred, min_og, max_og)[0][-2]\n# test_target_descaled = descale(test_target, min_og, max_og)[0][-2]\n\n# print(' --- descaled-Estimation --- ')\n# print(test_pred_descaled.shape, test_target_descaled.shape)\n\nfig = plt.figure(figsize=(6, 3))\n\nx_lin = np.arange(0, len(np.squeeze(test_sample)[:, -2]))\n\nplt.plot(x_lin, np.squeeze(test_sample)[:, -1], color = 'black', label = 'sample')\nplt.plot(x_lin[-1] + delay + 1, test_pred[0][0], '.b', label = 'pred')\nplt.plot(x_lin[-1] + delay + 1, test_target[0][0], '.r',  label = 'test')\nplt.title('ANN - Point. Estimation Prediction')\n\nplt.legend()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:43:12.460786Z","iopub.execute_input":"2024-10-15T11:43:12.461691Z","iopub.status.idle":"2024-10-15T11:43:13.044923Z","shell.execute_reply.started":"2024-10-15T11:43:12.461641Z","shell.execute_reply":"2024-10-15T11:43:13.043952Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#### Recurrent Neural Nets (RNNs) ","metadata":{}},{"cell_type":"code","source":"keras_normalizer_rnn = Normalization()\nkeras_normalizer_rnn.adapt(train_data_all)\n\ncnn_rnn_model = tf.keras.models.Sequential([\n    tf.keras.layers.Input(shape=(72, train_data_all.shape[-1])),\n    keras_normalizer_rnn,\n    tf.keras.layers.Conv1D(64, 3),\n    tf.keras.layers.Dense(128, activation = 'relu'),\n    tf.keras.layers.GRU(256, return_sequences=True),\n    tf.keras.layers.GRU(256, return_sequences=True),\n    tf.keras.layers.GRU(128),\n    tf.keras.layers.Dense(64, activation = 'relu'),\n    tf.keras.layers.Dense(1, activation = 'linear')\n]) \n\ncnn_rnn_model.compile(optimizer='Adam', loss='mse', metrics=['mae'])\ncnn_rnn_model.summary()","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:43:17.826684Z","iopub.execute_input":"2024-10-15T11:43:17.827408Z","iopub.status.idle":"2024-10-15T11:43:19.829090Z","shell.execute_reply.started":"2024-10-15T11:43:17.827373Z","shell.execute_reply":"2024-10-15T11:43:19.827973Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ft = time.time()\ncnn_rnn_hist = cnn_rnn_model.fit(train_data_all, train_target_all, epochs = 20, validation_data = (val_data_all, val_target_all))\nlt = time.time()\n\nprint('Training_Done...')\nprint('Execution Time: ', (lt - ft)/60)","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:43:23.273170Z","iopub.execute_input":"2024-10-15T11:43:23.273953Z","iopub.status.idle":"2024-10-15T11:52:08.365140Z","shell.execute_reply.started":"2024-10-15T11:43:23.273873Z","shell.execute_reply":"2024-10-15T11:52:08.363552Z"},"collapsed":true,"jupyter":{"outputs_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"cnn_rnn_hist_history = cnn_rnn_hist.history \nepochs_arr = np.arange(len(cnn_rnn_hist_history['loss']))\n\nplt.plot(epochs_arr[5:], cnn_rnn_hist_history['loss'][5:], label = 'train_loss', color = 'red', marker = '.')\nplt.plot(epochs_arr[5:], cnn_rnn_hist_history['val_loss'][5:], label = 'val_loss', color = 'blue', marker = '.')\nplt.title('Loss-Curve')\n\nplt.legend()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:52:11.618833Z","iopub.execute_input":"2024-10-15T11:52:11.619229Z","iopub.status.idle":"2024-10-15T11:52:11.668655Z","shell.execute_reply.started":"2024-10-15T11:52:11.619197Z","shell.execute_reply":"2024-10-15T11:52:11.667366Z"},"collapsed":true,"jupyter":{"outputs_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test_mae_err = cnn_rnn_model.evaluate(test_data_all, test_target_all, steps = 512)[-1]\nprint(f'The Absolute Error in Prediction of Close Price of BTC is {test_mae_err} USD, which is more than std-dev of Close Price in {sample_rate} minute time-frame chart')","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:52:35.222606Z","iopub.execute_input":"2024-10-15T11:52:35.222993Z","iopub.status.idle":"2024-10-15T11:52:39.336950Z","shell.execute_reply.started":"2024-10-15T11:52:35.222961Z","shell.execute_reply":"2024-10-15T11:52:39.335943Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"rn = random.randint(0, len(test_data_all))\nprint('random-seed: ', rn)\n# rn = 2306\n\ntest_sample, test_target = test_data_all[rn].reshape(1, -1, test_data_all.shape[-1]), test_target_all[rn].reshape(-1, 1) \nprint(test_sample.shape, test_target.shape)\n\ntest_pred = dnn_model.predict(test_sample)\n\nfig = plt.figure(figsize=(6, 3))\n\nx_lin = np.arange(0, len(np.squeeze(test_sample)[:, -2]))\n\nplt.plot(x_lin, np.squeeze(test_sample)[:, -1], color = 'black', label = 'sample')\nplt.plot(x_lin[-1] + delay + 1, test_pred[0][0], '.b', label = 'pred')\nplt.plot(x_lin[-1] + delay + 1, test_target[0][0], '.r',  label = 'test')\nplt.title('RNN - Point. Estimation Prediction')\n\nplt.legend()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:52:45.707782Z","iopub.execute_input":"2024-10-15T11:52:45.708169Z","iopub.status.idle":"2024-10-15T11:52:46.015395Z","shell.execute_reply.started":"2024-10-15T11:52:45.708136Z","shell.execute_reply":"2024-10-15T11:52:46.014415Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"------------------------","metadata":{}},{"cell_type":"markdown","source":"### Bayesian Neural Networks (BNN) Workflow","metadata":{}},{"cell_type":"code","source":"hidden_units = [32, 16]\nlearning_rate = 0.001  \ndataset_size = len(x_obs_close_scaled)\nbatch_size = 32","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:15:01.658797Z","iopub.execute_input":"2024-10-15T11:15:01.659786Z","iopub.status.idle":"2024-10-15T11:15:01.664176Z","shell.execute_reply.started":"2024-10-15T11:15:01.659750Z","shell.execute_reply":"2024-10-15T11:15:01.663163Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We are assuming that the posterior will follow a multivariate normal distribution, where the distribution of the learnable parameters is parameterized by the mean and covariance.\n\nTo represent the covariance matrix, we use a lower triangular matrix, which we then learn, instead of a full <code>n × n</code> covariance matrix. This is because of <b> Cholesky decomposition </b>, which states that any symmetric and positive definite covariance matrix can be written as the product of a lower triangular matrix and its transpose.\n\nBy learning the lower triangular matrix, we reduce computational complexity compared to learning the full covariance matrix.","metadata":{}},{"cell_type":"code","source":"def prior(kernel_size, bias_size, dtype=None):\n    n = kernel_size + bias_size\n    prior_model = tf.keras.Sequential(\n        [\n            tfp.layers.DistributionLambda(lambda t: tfp.distributions.MultivariateNormalDiag(loc=tf.zeros(n), scale_diag=tf.ones(n)))\n        ]\n    )\n    return prior_model\n\n\ndef posterior(kernel_size, bias_size, dtype=None):\n    n = kernel_size + bias_size\n    posterior_model = tf.keras.Sequential(\n        [\n            ## Student-t might be a better choice\n            tfp.layers.VariableLayer(tfp.layers.MultivariateNormalTriL.params_size(n), dtype=dtype), ## Creates trainable parameeters\n            tfp.layers.MultivariateNormalTriL(n), ## Uses the trainable parameters to -> to create distrib. of weights\n        ]\n    )\n    return posterior_model","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:15:03.470959Z","iopub.execute_input":"2024-10-15T11:15:03.471339Z","iopub.status.idle":"2024-10-15T11:15:03.478483Z","shell.execute_reply.started":"2024-10-15T11:15:03.471308Z","shell.execute_reply":"2024-10-15T11:15:03.477591Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def NLL(targets, est_distrib):\n    return - est_distrib.log_prob(targets)","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:15:05.006621Z","iopub.execute_input":"2024-10-15T11:15:05.007300Z","iopub.status.idle":"2024-10-15T11:15:05.011323Z","shell.execute_reply.started":"2024-10-15T11:15:05.007266Z","shell.execute_reply":"2024-10-15T11:15:05.010405Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"det_units = [128, 64, 64, 32]\n\n#### Only Supplying the BTC Closing Price - NOT multi-feature Space\ndef create_bnn_model(train_size):\n\n    inputs = tf.keras.layers.Input(shape = train_data_close.shape[1:])\n    x = tf.keras.layers.Dense(det_units[0], activation = 'tanh')(inputs)\n    x = tf.keras.layers.BatchNormalization()(x)\n    x = tf.keras.layers.GRU(det_units[1], activation = 'tanh', return_sequences = True)(x)\n    x = tf.keras.layers.LayerNormalization()(x)    \n    x = tf.keras.layers.GRU(det_units[2], activation = 'tanh', return_sequences = False)(x)\n    x = tf.keras.layers.LayerNormalization()(x)\n    x = tf.keras.layers.Dense(det_units[-1], activation = 'sigmoid')(x)\n    x = tf.keras.layers.BatchNormalization()(x)\n\n    # Create hidden layers with weight uncertainty using the DenseVariational layer.\n    for units in hidden_units:\n        x = tfp.layers.DenseVariational(\n            units=units,\n            make_prior_fn=prior,\n            make_posterior_fn=posterior,\n            kl_weight = 1 / train_size,\n            activation=\"sigmoid\",\n        )(x)\n\n    # The output is deterministic: a single point estimate | Distribution.\n    ''' Distrib. Output | Independent Normal (Gaussian) distributions for each element in the data (Probabilistic) '''\n    x = tf.keras.layers.Dense(units=2)(x) # |meu, Sigma| --> N_distrib \n    outputs = tfp.layers.IndependentNormal(1)(x)\n\n    ''' Deterministic Output (Point Estimation) '''\n#     outputs = tf.keras.layers.Dense(units=1)(x)\n    \n    model = tf.keras.models.Model(inputs=inputs, outputs=outputs)\n    return model\n\n\ntrain_dataset = tf.data.Dataset.from_tensor_slices((train_data_close, train_target_close)).batch(batch_size)\nval_dataset = tf.data.Dataset.from_tensor_slices((val_data_close, val_target_close)).batch(batch_size)\ntf.print(train_dataset, val_dataset)\n\nloss = tf.keras.losses.MSE\n\nbnn_model = create_bnn_model(dataset_size)\nbnn_model.compile(optimizer = tf.keras.optimizers.RMSprop(learning_rate=learning_rate), loss = loss, metrics=[tf.keras.metrics.RootMeanSquaredError()])\n\nbnn_model.summary()","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:27:18.897489Z","iopub.execute_input":"2024-10-15T11:27:18.898244Z","iopub.status.idle":"2024-10-15T11:27:20.223076Z","shell.execute_reply.started":"2024-10-15T11:27:18.898203Z","shell.execute_reply":"2024-10-15T11:27:20.222041Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"num_epochs = 30\n\nes = tf.keras.callbacks.EarlyStopping(monitor = 'val_loss', patience = 5)\n\nprint(\"Start training the model...\")\n\nft = time.time()\n\nhist_values = bnn_model.fit(train_dataset, epochs=num_epochs, validation_data = val_dataset, callbacks = [es])\n\nlt = time.time()\n\nprint(f\"Model training finished in {(lt - ft) / 60} mins\")","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:27:22.740033Z","iopub.execute_input":"2024-10-15T11:27:22.740421Z","iopub.status.idle":"2024-10-15T11:34:28.051385Z","shell.execute_reply.started":"2024-10-15T11:27:22.740390Z","shell.execute_reply":"2024-10-15T11:34:28.050568Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Point-Estimate","metadata":{}},{"cell_type":"code","source":"data_dict = {}\n\ndata_dict['train'] = [train_data_close, train_target_close]\ndata_dict['val'] = [val_data_close, val_target_close]\ndata_dict['test'] = [test_data_close, test_target_close]","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:35:11.490880Z","iopub.execute_input":"2024-10-15T11:35:11.491759Z","iopub.status.idle":"2024-10-15T11:35:11.497623Z","shell.execute_reply.started":"2024-10-15T11:35:11.491712Z","shell.execute_reply":"2024-10-15T11:35:11.496527Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"_rate = 300\n\ndef Pred_Plot(pred_data, pred_target):\n\n    _ran = min(len(pred_data), 300)\n    \n    iterations = 10\n    \n    x_lin = np.arange(0, _ran) \n\n    predicted = []\n    \n    for _ in range(iterations):\n        predicted.append(bnn_model(pred_data[:_ran]).numpy())\n    \n    predicted = np.concatenate(predicted, axis=1)\n    \n    ## intermediate stochastic components\n    _st = 3\n    predicted_three = predicted[:, _st : _st + (_rate//100)] \n\n    prediction_mean = np.mean(predicted, axis=1).tolist()\n    prediction_min = np.min(predicted, axis=1).tolist()\n    prediction_max = np.max(predicted, axis=1).tolist()\n    prediction_range = (np.max(predicted, axis=1) - np.min(predicted, axis=1)).tolist()\n    \n    inv_test_pred = prediction_mean\n    inv_test_target = pred_target[:_ran]\n    lower = prediction_min\n    upper = prediction_max\n    predicted_three = predicted_three\n\n    return x_lin, inv_test_pred, inv_test_target, lower, upper, predicted_three\n\n\nfig, axs = plt.subplots(1, 1, figsize=(8, 4))\n\nnames = list(data_dict.keys())\n\nft = time.time()\n\n_color = ['blue', 'green', 'yellow']\n\nfor i, nm in enumerate(names[1:2]):\n    \n    k_st = 0\n    x_lin, inv_p, inv_t, lower, upper, p_three = Pred_Plot(data_dict[nm][0], data_dict[nm][-1])\n    \n    if i == 0:\n        print(x_lin.shape, np.array(inv_p).shape, np.array(inv_t).shape, np.array(lower).shape, np.array(upper).shape, p_three.shape)\n    \n    inv_t = scaler_close.inverse_transform(np.array(inv_t).reshape(-1, 1))\n    inv_p = scaler_close.inverse_transform(np.array(inv_p).reshape(-1, 1))\n    lower = scaler_close.inverse_transform(np.array(lower).reshape(-1, 1)).reshape(len(lower), )\n    upper = scaler_close.inverse_transform(np.array(upper).reshape(-1, 1)).reshape(len(upper), )\n    p_three = scaler_close.inverse_transform(p_three)\n    \n    axs.plot(x_lin, inv_t, 'k', label = 'target')\n    axs.plot(x_lin, inv_p, 'r--',  label = 'Prediction (BNN)')\n    \n    for step, c in enumerate(_color):\n        axs.plot(x_lin, p_three[:, step], color = c, alpha = 0.2, label = 'Stoch. Components')\n\n    axs.fill_between(x_lin, lower, upper, color='orange', alpha=0.4, label = 'min-max range')\n    axs.set_title('Prediction')\n    axs.legend()\n    \nlt = time.time()\n \nprint('time-taken: ', (lt - ft)/60, 'mins')\n\nplt.tight_layout(pad = 2)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:26:21.212488Z","iopub.execute_input":"2024-10-15T11:26:21.212858Z","iopub.status.idle":"2024-10-15T11:26:23.295201Z","shell.execute_reply.started":"2024-10-15T11:26:21.212826Z","shell.execute_reply":"2024-10-15T11:26:23.294232Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test_dataset = tf.data.Dataset.from_tensor_slices((test_data_close, test_target_close))\nmse_error, rmse_error = bnn_model.evaluate(test_data_close, test_target_close)\n\nprint(\"MSE-Error (normalized): \", mse_error, \", RMSE-Error (normalized):\",rmse_error)","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:26:47.800959Z","iopub.execute_input":"2024-10-15T11:26:47.801307Z","iopub.status.idle":"2024-10-15T11:26:50.375066Z","shell.execute_reply.started":"2024-10-15T11:26:47.801280Z","shell.execute_reply":"2024-10-15T11:26:50.373890Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Probabilistic","metadata":{}},{"cell_type":"code","source":"pred_data = test_data_close\npred_target = test_target_close\n\nprint(pred_data.shape, pred_target.shape)\n\nci = 0.675\nmax_predict = 500\n\n_ran = min(len(pred_data), max_predict)\n\nx_lin = np.arange(0, _ran) \n\n## Just samples a value at each t from the surrogate posterior distrb.\nprediction_distribution = bnn_model(pred_data[:_ran]) ### Gives stochastic value at each time-step (i,e n random variables for n t-steps) for each run\nprediction_mean = prediction_distribution.mean().numpy().tolist() ### We get a distribution object, giving seperate distributions for each-time step  \nprediction_stdv = prediction_distribution.stddev().numpy()\n\ninv_test_pred = scaler_close.inverse_transform(prediction_mean)\ninv_test_target = scaler_close.inverse_transform(pred_target[:_ran].reshape(-1, 1))\n\n# The 75% CI \nupper = (prediction_mean + (ci * prediction_stdv)).tolist()\nlower = (prediction_mean - (ci * prediction_stdv)).tolist()\n\nupper = np.array(upper).reshape(-1, 1)\nlower = np.array(lower).reshape(-1, 1) \n\nupper = scaler_close.inverse_transform(upper).reshape(len(upper), )\nlower = scaler_close.inverse_transform(lower).reshape(len(lower), )\n\nfig = plt.figure(figsize=(8, 3))\nplt.plot(x_lin, inv_test_target, 'k', label = 'test')\nplt.plot(x_lin, inv_test_pred, 'r--',  label = 'Prediction (BNN)')\n# plt.fill_between(x_lin, lower, upper, color='orange', alpha=0.5, label = '75% CI')\nplt.title('Test-Performance')\n\nplt.legend()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:38:56.140184Z","iopub.execute_input":"2024-10-15T11:38:56.140749Z","iopub.status.idle":"2024-10-15T11:38:56.630840Z","shell.execute_reply.started":"2024-10-15T11:38:56.140710Z","shell.execute_reply":"2024-10-15T11:38:56.629685Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test_dataset = tf.data.Dataset.from_tensor_slices((test_data_close, test_target_close))\nmse_error, rmse_error = bnn_model.evaluate(test_data_close, test_target_close)\n\nprint(\"MSE-Error: \", mse_error, \", RMSE-Error:\",rmse_error)","metadata":{"execution":{"iopub.status.busy":"2024-10-15T11:37:24.974267Z","iopub.execute_input":"2024-10-15T11:37:24.975477Z","iopub.status.idle":"2024-10-15T11:37:26.304687Z","shell.execute_reply.started":"2024-10-15T11:37:24.975429Z","shell.execute_reply":"2024-10-15T11:37:26.303697Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"---","metadata":{}},{"cell_type":"markdown","source":"SESSION STOPPER","metadata":{}},{"cell_type":"code","source":"import time\nfor i in range(50):\n    print('Time: ', i)\n    time.sleep(30*60)","metadata":{"execution":{"iopub.status.busy":"2023-07-23T19:43:41.326701Z","iopub.execute_input":"2023-07-23T19:43:41.327043Z"},"trusted":true},"execution_count":null,"outputs":[]}]}