Update fmda_rnn_spatial.ipynb
[notebooks.git] / fmda / fmda_rnn_spatial.ipynb
blob3622d65b151ebf5eb473317ab44b616636612f50
2  "cells": [
3   {
4    "cell_type": "markdown",
5    "id": "83b774b3-ef55-480a-b999-506676e49145",
6    "metadata": {},
7    "source": [
8     "# v2.1 run RNN with Spatial Training\n",
9     "\n",
10     "This notebook is intended to set up a test where the RNN is run serial by location and compared to the spatial training scheme. Additionally, the ODE model with the augmented KF will be run as a comparison, but note that the RNN models will be predicting entirely without knowledge of the heldout locations, while the augmented KF will be run directly on the test locations.\n"
11    ]
12   },
13   {
14    "cell_type": "markdown",
15    "id": "bbd84d61-a9cd-47b4-b538-4986fb10b98d",
16    "metadata": {},
17    "source": [
18     "## Environment Setup"
19    ]
20   },
21   {
22    "cell_type": "code",
23    "execution_count": null,
24    "id": "83cc1dc4-3dcb-4325-9263-58101a3dc378",
25    "metadata": {},
26    "outputs": [],
27    "source": [
28     "import numpy as np\n",
29     "from utils import print_dict_summary, print_first, str2time, logging_setup\n",
30     "import pickle\n",
31     "import logging\n",
32     "import os.path as osp\n",
33     "from moisture_rnn_pkl import pkl2train\n",
34     "from moisture_rnn import RNNParams, RNNData, RNN \n",
35     "from utils import hash2, read_yml, read_pkl, retrieve_url, Dict\n",
36     "from moisture_rnn import RNN\n",
37     "import reproducibility\n",
38     "from data_funcs import rmse, to_json, combine_nested, process_train_dict\n",
39     "from moisture_models import run_augmented_kf\n",
40     "import copy\n",
41     "import pandas as pd\n",
42     "import matplotlib.pyplot as plt\n",
43     "import yaml\n",
44     "import time"
45    ]
46   },
47   {
48    "cell_type": "code",
49    "execution_count": null,
50    "id": "17db9b90-a931-4674-a447-5b8ffbcdc86a",
51    "metadata": {},
52    "outputs": [],
53    "source": [
54     "logging_setup()"
55    ]
56   },
57   {
58    "cell_type": "code",
59    "execution_count": null,
60    "id": "35319c1c-7849-4b8c-8262-f5aa6656e0c7",
61    "metadata": {},
62    "outputs": [],
63    "source": [
64     "retrieve_url(\n",
65     "    url = \"https://demo.openwfm.org/web/data/fmda/dicts/fmda_nw_202401-05_f05.pkl\", \n",
66     "    dest_path = \"data/fmda_nw_202401-05_f05.pkl\")"
67    ]
68   },
69   {
70    "cell_type": "code",
71    "execution_count": null,
72    "id": "eabdbd9c-07d9-4bae-9851-cca79f321895",
73    "metadata": {},
74    "outputs": [],
75    "source": [
76     "file_paths = ['data/fmda_nw_202401-05_f05.pkl']"
77    ]
78   },
79   {
80    "cell_type": "code",
81    "execution_count": null,
82    "id": "dcca6185-e799-4dd1-8acb-87ad33c411d7",
83    "metadata": {},
84    "outputs": [],
85    "source": [
86     "# read/write control\n",
87     "train_file='data/train.pkl'\n",
88     "train_create=True   # if false, read\n",
89     "train_write=True\n",
90     "train_read=True"
91    ]
92   },
93   {
94    "cell_type": "code",
95    "execution_count": null,
96    "id": "604388de-11ab-45c3-9f0d-80bdff0cca60",
97    "metadata": {},
98    "outputs": [],
99    "source": [
100     "# Params used for data filtering\n",
101     "params_data = read_yml(\"params_data.yaml\") \n",
102     "params_data"
103    ]
104   },
105   {
106    "cell_type": "code",
107    "execution_count": null,
108    "id": "211a1c2f-ba8d-40b8-b29c-daa38af97a26",
109    "metadata": {},
110    "outputs": [],
111    "source": [
112     "# Params used for setting up RNN\n",
113     "params = read_yml(\"params.yaml\", subkey='rnn') \n",
114     "params"
115    ]
116   },
117   {
118    "cell_type": "code",
119    "execution_count": null,
120    "id": "bc0a775b-b587-42ef-8576-e36dc0be3a75",
121    "metadata": {
122     "scrolled": true
123    },
124    "outputs": [],
125    "source": [
126     "if train_create:\n",
127     "    logging.info('creating the training cases from files %s',file_paths)\n",
128     "    # osp.join works on windows too, joins paths using \\ or /\n",
129     "    train = process_train_dict(file_paths, params_data = params_data, verbose=True)\n",
130     "if train_write:\n",
131     "    with open(train_file, 'wb') as file:\n",
132     "        logging.info('Writing the rain cases into file %s',train_file)\n",
133     "        pickle.dump(train, file)\n",
134     "if train_read:\n",
135     "    logging.info('Reading the train cases from file %s',train_file)\n",
136     "    train = read_pkl(train_file)"
137    ]
138   },
139   {
140    "cell_type": "code",
141    "execution_count": null,
142    "id": "23cd60c0-9865-4314-9a96-948c3d400c08",
143    "metadata": {},
144    "outputs": [],
145    "source": [
146     "from itertools import islice\n",
147     "train = {k: train[k] for k in islice(train, 200)}"
148    ]
149   },
150   {
151    "cell_type": "markdown",
152    "id": "efc10cdc-f18b-4781-84da-b8e2eef39981",
153    "metadata": {},
154    "source": [
155     "## Setup Validation Runs"
156    ]
157   },
158   {
159    "cell_type": "code",
160    "execution_count": null,
161    "id": "66f40c9f-c1c2-4b12-bf14-2ada8c26113d",
162    "metadata": {},
163    "outputs": [],
164    "source": [
165     "params = RNNParams(params)\n",
166     "params.update({'epochs': 200, \n",
167     "               'learning_rate': 0.001,\n",
168     "               'activation': ['tanh', 'tanh'], # Activation for RNN Layers, Dense layers respectively.\n",
169     "               'recurrent_layers': 2, 'recurrent_units': 30, \n",
170     "               'dense_layers': 2, 'dense_units': 30,\n",
171     "               'early_stopping_patience': 30, # how many epochs of no validation accuracy gain to wait before stopping\n",
172     "               'batch_schedule_type': 'exp', # Hidden state batch reset schedule\n",
173     "               'bmin': 20, # Lower bound of hidden state batch reset, \n",
174     "               'bmax': params_data['hours'], # Upper bound of hidden state batch reset, using max hours\n",
175     "               'features_list': ['Ed', 'Ew', 'rain', 'elev', 'lon', 'lat', 'solar', 'wind']\n",
176     "              })"
177    ]
178   },
179   {
180    "cell_type": "code",
181    "execution_count": null,
182    "id": "36823193-b93c-421e-b699-8c1ae5719309",
183    "metadata": {},
184    "outputs": [],
185    "source": [
186     "reproducibility.set_seed(123)"
187    ]
188   },
189   {
190    "cell_type": "markdown",
191    "id": "a24d76fc-6c25-43e7-99df-3cd5dbf84fc3",
192    "metadata": {},
193    "source": [
194     "## Spatial Data Training"
195    ]
196   },
197   {
198    "cell_type": "code",
199    "execution_count": null,
200    "id": "3b5371a9-c1e8-4df5-b360-210746f7cd52",
201    "metadata": {},
202    "outputs": [],
203    "source": [
204     "# Start timer for code \n",
205     "start_time = time.time()"
206    ]
207   },
208   {
209    "cell_type": "code",
210    "execution_count": null,
211    "id": "faf93470-b55f-4770-9fa9-3288a2f13fcc",
212    "metadata": {},
213    "outputs": [],
214    "source": [
215     "# Combine Nested Dictionary into Spatial Data\n",
216     "train_sp = Dict(combine_nested(train))"
217    ]
218   },
219   {
220    "cell_type": "code",
221    "execution_count": null,
222    "id": "c0c7f5fb-4c33-45f8-9a2e-38c9ab1cd4e3",
223    "metadata": {},
224    "outputs": [],
225    "source": [
226     "rnn_dat_sp = RNNData(\n",
227     "    train_sp, # input dictionary\n",
228     "    scaler=\"standard\",  # data scaling type\n",
229     "    features_list = params['features_list'] # features for predicting outcome\n",
230     ")\n",
231     "\n",
232     "\n",
233     "rnn_dat_sp.train_test_split(   \n",
234     "    time_fracs = [.8, .1, .1], # Percent of total time steps used for train/val/test\n",
235     "    space_fracs = [.8, .1, .1] # Percent of total timeseries used for train/val/test\n",
236     ")\n",
237     "rnn_dat_sp.scale_data()\n",
238     "\n",
239     "rnn_dat_sp.batch_reshape(\n",
240     "    timesteps = params['timesteps'], # Timesteps aka sequence length for RNN input data. \n",
241     "    batch_size = params['batch_size'] # Number of samples of length timesteps for a single round of grad. descent\n",
242     ")"
243    ]
244   },
245   {
246    "cell_type": "code",
247    "execution_count": null,
248    "id": "7431bc95-d384-40fd-a622-bbc0ee68e5cd",
249    "metadata": {},
250    "outputs": [],
251    "source": [
252     "# Update Params specific to spatial training\n",
253     "params.update({\n",
254     "    'loc_batch_reset': rnn_dat_sp.n_seqs # Used to reset hidden state when location changes for a given batch\n",
255     "})"
256    ]
257   },
258   {
259    "cell_type": "code",
260    "execution_count": null,
261    "id": "4bc11474-fed8-47f2-b9cf-dfdda0d3d3b2",
262    "metadata": {},
263    "outputs": [],
264    "source": [
265     "rnn_sp = RNN(params)\n",
266     "m, errs = rnn_sp.run_model(rnn_dat_sp)"
267    ]
268   },
269   {
270    "cell_type": "code",
271    "execution_count": null,
272    "id": "704ad662-d81a-488d-be3d-e90bf775a5b8",
273    "metadata": {},
274    "outputs": [],
275    "source": [
276     "errs.mean()"
277    ]
278   },
279   {
280    "cell_type": "code",
281    "execution_count": null,
282    "id": "d53571e3-b6cf-49aa-9848-e3c77053283d",
283    "metadata": {},
284    "outputs": [],
285    "source": [
286     "# End Timer\n",
287     "end_time = time.time()\n",
288     "\n",
289     "# Calculate Code Runtime\n",
290     "elapsed_time_sp = end_time - start_time\n",
291     "print(f\"Spatial Training Elapsed time: {elapsed_time_sp:.4f} seconds\")"
292    ]
293   },
294   {
295    "cell_type": "markdown",
296    "id": "7d8292a2-418c-48ed-aff7-ccbe98b046d3",
297    "metadata": {},
298    "source": [
299     "## Run ODE + KF and Compare"
300    ]
301   },
302   {
303    "cell_type": "code",
304    "execution_count": null,
305    "id": "cca12d8c-c0e1-4df4-b2ca-20440485f2f3",
306    "metadata": {},
307    "outputs": [],
308    "source": [
309     "# Get timeseries IDs from previous RNNData object\n",
310     "test_cases = rnn_dat_sp.loc['test_locs']\n",
311     "print(len(test_cases))"
312    ]
313   },
314   {
315    "cell_type": "code",
316    "execution_count": null,
317    "id": "997f2534-7e77-45b3-93bf-d988837dfc0b",
318    "metadata": {},
319    "outputs": [],
320    "source": [
321     "test_ind = rnn_dat_sp.test_ind # Time index for test period start\n",
322     "print(test_ind)"
323    ]
324   },
325   {
326    "cell_type": "code",
327    "execution_count": null,
328    "id": "1e4ffc68-c775-41c6-ac42-f49c76824b43",
329    "metadata": {
330     "scrolled": true
331    },
332    "outputs": [],
333    "source": [
334     "outputs_kf = {}\n",
335     "for case in test_cases:\n",
336     "    print(\"~\"*50)\n",
337     "    print(case)\n",
338     "    # Run Augmented KF\n",
339     "    print('Running Augmented KF')\n",
340     "    train[case]['h2'] = test_ind\n",
341     "    train[case]['scale_fm'] = 1\n",
342     "    m, Ec = run_augmented_kf(train[case])\n",
343     "    y = train[case]['y']        \n",
344     "    train[case]['m_kf'] = m\n",
345     "    print(f\"KF RMSE: {rmse(m[test_ind:],y[test_ind:])}\")\n",
346     "    outputs_kf[case] = {'case':case, 'errs': rmse(m[test_ind:],y[test_ind:])}"
347    ]
348   },
349   {
350    "cell_type": "code",
351    "execution_count": null,
352    "id": "57b19ec5-23f6-44ec-9f71-16d4d69aec68",
353    "metadata": {},
354    "outputs": [],
355    "source": [
356     "df_kf = pd.DataFrame.from_dict(outputs_kf).transpose()\n",
357     "df_kf.head()"
358    ]
359   },
360   {
361    "cell_type": "code",
362    "execution_count": null,
363    "id": "25a9d2fe-83f7-4ef3-a04b-14c970b6e2ba",
364    "metadata": {},
365    "outputs": [],
366    "source": [
367     "df_kf.errs.mean()"
368    ]
369   },
370   {
371    "cell_type": "markdown",
372    "id": "f616bbf8-d89e-4c5b-9e47-59f02246b6f2",
373    "metadata": {},
374    "source": [
375     "## Serial Training"
376    ]
377   },
378   {
379    "cell_type": "code",
380    "execution_count": null,
381    "id": "6fa20e9f-604a-4938-ab68-b71fbb7326df",
382    "metadata": {},
383    "outputs": [],
384    "source": [
385     "# Start timer for code \n",
386     "start_time = time.time()"
387    ]
388   },
389   {
390    "cell_type": "code",
391    "execution_count": null,
392    "id": "f033e78c-a506-4508-a23c-8e6574014872",
393    "metadata": {},
394    "outputs": [],
395    "source": [
396     "# Update Params specific to Serial training\n",
397     "params.update({\n",
398     "    'loc_batch_reset': None, # Used to reset hidden state when location changes for a given batch\n",
399     "    'epochs': 1 # less epochs since fit will be run multiple times over locations\n",
400     "})"
401    ]
402   },
403   {
404    "cell_type": "code",
405    "execution_count": null,
406    "id": "ff1788ec-081b-403f-bcfa-b625f0e3dbe1",
407    "metadata": {},
408    "outputs": [],
409    "source": [
410     "train_cases = rnn_dat_sp.loc['train_locs']\n",
411     "test_cases = rnn_dat_sp.loc['test_locs']"
412    ]
413   },
414   {
415    "cell_type": "code",
416    "execution_count": null,
417    "id": "8a2af45e-e81b-421f-b940-e8779177dd5d",
418    "metadata": {},
419    "outputs": [],
420    "source": [
421     "# Initialize Model with first train case\n",
422     "rnn_dat = RNNData(train[train_cases[0]], params['scaler'], params['features_list'])\n",
423     "rnn_dat.train_test_split(\n",
424     "    time_fracs = [.8, .1, .1]\n",
425     ")\n",
426     "rnn_dat.scale_data()\n",
427     "rnn_dat.batch_reshape(timesteps = params['timesteps'], batch_size = params['batch_size'])"
428    ]
429   },
430   {
431    "cell_type": "code",
432    "execution_count": null,
433    "id": "ac6fecc2-f614-4506-b5f9-05a6eca3b62e",
434    "metadata": {},
435    "outputs": [],
436    "source": [
437     "reproducibility.set_seed()\n",
438     "rnn = RNN(params)"
439    ]
440   },
441   {
442    "cell_type": "code",
443    "execution_count": null,
444    "id": "79b5af30-7d52-410c-9595-e89e9756fd38",
445    "metadata": {
446     "scrolled": true
447    },
448    "outputs": [],
449    "source": [
450     "# Train\n",
451     "for case in train_cases:\n",
452     "    print(\"~\"*50)\n",
453     "    print(f\"Training with Case {case}\")\n",
454     "    rnn_dat_temp = RNNData(train[case], params['scaler'], params['features_list'])\n",
455     "    rnn_dat_temp.train_test_split(\n",
456     "        time_fracs = [.8, .1, .1]\n",
457     "    )\n",
458     "    rnn_dat_temp.scale_data()\n",
459     "    rnn_dat_temp.batch_reshape(timesteps = params['timesteps'], batch_size = params['batch_size'])\n",
460     "    rnn.fit(rnn_dat_temp['X_train'], rnn_dat_temp['y_train'],\n",
461     "           validation_data=(rnn_dat_temp['X_val'], rnn_dat_temp['y_val']))    "
462    ]
463   },
464   {
465    "cell_type": "code",
466    "execution_count": null,
467    "id": "03d716b4-0ff5-4b80-a241-440543ba9b46",
468    "metadata": {
469     "scrolled": true
470    },
471    "outputs": [],
472    "source": [
473     "# Predict\n",
474     "outputs_rnn_serial = {}\n",
475     "test_ind = rnn_dat.test_ind\n",
476     "for i, case in enumerate(test_cases):\n",
477     "    print(\"~\"*50)\n",
478     "    rnn_dat_temp = RNNData(train[case], params['scaler'], params['features_list'])\n",
479     "    rnn_dat_temp.train_test_split(\n",
480     "        time_fracs = [.8, .1, .1]\n",
481     "    )\n",
482     "    rnn_dat_temp.scale_data()\n",
483     "    rnn_dat_temp.batch_reshape(timesteps = params['timesteps'], batch_size = params['batch_size'])    \n",
484     "    X_temp = rnn_dat_temp.scale_all_X()\n",
485     "    m = rnn.predict(X_temp)\n",
486     "    outputs_rnn_serial[case] = {'case':case, 'errs': rmse(m[test_ind:], rnn_dat.y_test)}"
487    ]
488   },
489   {
490    "cell_type": "code",
491    "execution_count": null,
492    "id": "e5a80bae-fe1a-4ec9-b9ac-31d540eaba40",
493    "metadata": {},
494    "outputs": [],
495    "source": [
496     "df_rnn_serial = pd.DataFrame.from_dict(outputs_rnn_serial).transpose()\n",
497     "df_rnn_serial.head()"
498    ]
499   },
500   {
501    "cell_type": "code",
502    "execution_count": null,
503    "id": "0c5b866e-c2bf-4bc1-8f6f-3ba8a9448d07",
504    "metadata": {},
505    "outputs": [],
506    "source": [
507     "df_rnn_serial.errs.mean()"
508    ]
509   },
510   {
511    "cell_type": "code",
512    "execution_count": null,
513    "id": "f5a364cb-01bf-49ad-a704-5aa3c9564967",
514    "metadata": {},
515    "outputs": [],
516    "source": [
517     "# End Timer\n",
518     "end_time = time.time()\n",
519     "\n",
520     "# Calculate Code Runtime\n",
521     "elapsed_time_ser = end_time - start_time\n",
522     "print(f\"Serial Training Elapsed time: {elapsed_time_ser:.4f} seconds\")"
523    ]
524   },
525   {
526    "cell_type": "markdown",
527    "id": "86795281-f8ea-4141-81ea-c53fae830e80",
528    "metadata": {},
529    "source": [
530     "## Compare"
531    ]
532   },
533   {
534    "cell_type": "code",
535    "execution_count": null,
536    "id": "508a6392-49bc-4471-ad8e-814f60119283",
537    "metadata": {},
538    "outputs": [],
539    "source": [
540     "print(f\"Total Test Cases: {len(test_cases)}\")\n",
541     "print(f\"Total Test Hours: {rnn_dat_temp.y_test.shape[0]}\")"
542    ]
543   },
544   {
545    "cell_type": "code",
546    "execution_count": null,
547    "id": "73e8ca05-d17b-4e72-8def-fa77664e7bb0",
548    "metadata": {},
549    "outputs": [],
550    "source": [
551     "print(f\"Spatial Training RMSE: {errs.mean()}\")\n",
552     "print(f\"Serial Training RMSE: {df_rnn_serial.errs.mean()}\")\n",
553     "print(f\"Augmented KF RMSE: {df_kf.errs.mean()}\")"
554    ]
555   },
556   {
557    "cell_type": "code",
558    "execution_count": null,
559    "id": "a73d22ee-707b-44a3-80ab-ad6e671731cf",
560    "metadata": {},
561    "outputs": [],
562    "source": []
563   },
564   {
565    "cell_type": "code",
566    "execution_count": null,
567    "id": "272bfb32-e8e2-49dd-8f90-4b5b09c3a2a2",
568    "metadata": {},
569    "outputs": [],
570    "source": [
571     "print(f\"Spatial Training Elapsed time: {elapsed_time_sp:.4f} seconds\")\n",
572     "print(f\"Serial Training Elapsed time: {elapsed_time_ser:.4f} seconds\")"
573    ]
574   },
575   {
576    "cell_type": "code",
577    "execution_count": null,
578    "id": "38ab08fb-ac97-45be-8907-6f9cd124243b",
579    "metadata": {},
580    "outputs": [],
581    "source": []
582   }
583  ],
584  "metadata": {
585   "kernelspec": {
586    "display_name": "Python 3 (ipykernel)",
587    "language": "python",
588    "name": "python3"
589   },
590   "language_info": {
591    "codemirror_mode": {
592     "name": "ipython",
593     "version": 3
594    },
595    "file_extension": ".py",
596    "mimetype": "text/x-python",
597    "name": "python",
598    "nbconvert_exporter": "python",
599    "pygments_lexer": "ipython3",
600    "version": "3.12.5"
601   }
602  },
603  "nbformat": 4,
604  "nbformat_minor": 5