6_Bonus_Projects.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "62321d43-d6be-4c35-904d-8baad835d3ff",
   "metadata": {},
   "source": [
    "# Part VI: Bonus projects\n",
    "\n",
    "Authors: Jordi Bolibar & Facundo Sapienza"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "49bf331d-4fda-4b9a-9aed-c239aabc6250",
   "metadata": {
    "user_expressions": []
   },
   "source": [
    "In notebook 5 we re-computed a distributed dataset which allows us to work with 2D gridded data instead of glacier-wide data. In this notebook, we will introduce two different bonus projects which can help you exploit these datasets further. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ff63f385-9125-455b-ad21-f2b43c8917c4",
   "metadata": {
    "user_expressions": []
   },
   "source": [
    "For these two projects, we will use the following libraries:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "478d0888-0204-4ef1-8878-255845ea2880",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "2023-04-20 08:46:47: oggm.cfg: Reading default parameters from the OGGM `params.cfg` configuration file.\n",
      "2023-04-20 08:46:47: oggm.cfg: Multiprocessing switched OFF according to the parameter file.\n",
      "2023-04-20 08:46:47: oggm.cfg: Multiprocessing: using all available processors (N=4)\n",
      "2023-04-20 08:46:49: oggm.cfg: PARAMS['border'] changed from `80` to `10`.\n",
      "2023-04-20 08:46:49: oggm.cfg: Multiprocessing switched ON after user settings.\n"
     ]
    }
   ],
   "source": [
    "## Auxiliary libraries\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import geopandas as gpd\n",
    "import xarray as xr\n",
    "import salem\n",
    "from pathlib import Path\n",
    "import os\n",
    "\n",
    "# OGGM\n",
    "import oggm.cfg as cfg\n",
    "from oggm import utils, workflow, tasks, graphics\n",
    "\n",
    "cfg.initialize(logging_level='WARNING') # initialize OGGM\n",
    "cfg.PARAMS['border'] = 10\n",
    "cfg.PARAMS['use_multiprocessing'] = True "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9aa0fdf7-43df-4def-9139-a8a16e706924",
   "metadata": {
    "user_expressions": []
   },
   "source": [
    "## 6.1. Bonus project No1: Distributed glacier mass balance\n",
    "\n",
    "<img src=\"Figures/eye_logo.png\" width=\"75\"/>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b062142f-d23f-4b31-ba07-f9b258a25765",
   "metadata": {
    "user_expressions": []
   },
   "source": [
    "In this project, we propose to go beyond the mass balance modelling we did in the previous notebooks, and to attempt to learn the spatial distribution of the glacier mass balance directly from the gridded geodetic mass balance data from Hugonnet et al. (2021). Instead of attempting to model a single data point per glacier, we will simulate the spatial distribution of the mass balance (i.e. its altitudinal gradient). "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d6371c38-b00e-492d-a23c-5a623ebb5191",
   "metadata": {
    "user_expressions": []
   },
   "source": [
    "## 6.1.1 Building a distributed dataset"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3a6f1260-352b-4633-bef0-c0ede9bf4bd5",
   "metadata": {
    "user_expressions": []
   },
   "source": [
    "In the previous notebooks we worked with tabular data. As you can see in the plot above, now we will be working with distributed (2D) data. Therefore, we need to reshape and adjust our dataset for this new problem."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "9d33007a-2f00-4ea6-babd-ac99fc11360c",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Index(['Unnamed: 0.1', 'Unnamed: 0', 'RGI_ID', 'PDD_2D', 'rain_2D', 'snow_2D',\n",
       "       'topo', 'aspect', 'slope', 'dis_from_border', 'glacier_mask',\n",
       "       'millan_ice_thickness', 'hugonnet_dhdt'],\n",
       "      dtype='object')"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df_subset = pd.read_csv('Data/df_distributed_subset.csv')\n",
    "df_subset.columns"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1089fd33-d4bf-47c1-82eb-c4fbf328d441",
   "metadata": {},
   "source": [
    "We clean the DataFrame a little bit"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "0ccf7c4a-c762-4f9e-a338-b5ca72ab37c1",
   "metadata": {},
   "outputs": [],
   "source": [
    "df_subset = df_subset.drop(['Unnamed: 0.1', 'Unnamed: 0'], axis=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "43a76e48-ce78-436f-9c3c-287146e9e7bd",
   "metadata": {},
   "source": [
    "This is a small subset of the whole dataset, a lightweight version to start working with."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "54871a99-391e-431d-bd7a-cd458a0142ea",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>RGI_ID</th>\n",
       "      <th>PDD_2D</th>\n",
       "      <th>rain_2D</th>\n",
       "      <th>snow_2D</th>\n",
       "      <th>topo</th>\n",
       "      <th>aspect</th>\n",
       "      <th>slope</th>\n",
       "      <th>dis_from_border</th>\n",
       "      <th>glacier_mask</th>\n",
       "      <th>millan_ice_thickness</th>\n",
       "      <th>hugonnet_dhdt</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>RGI60-08.00005</td>\n",
       "      <td>609.424791</td>\n",
       "      <td>24237.964135</td>\n",
       "      <td>16140.702605</td>\n",
       "      <td>1024.7790</td>\n",
       "      <td>5.331003</td>\n",
       "      <td>0.102807</td>\n",
       "      <td>1997.02900</td>\n",
       "      <td>0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.047467</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>RGI60-08.00005</td>\n",
       "      <td>615.582145</td>\n",
       "      <td>24237.964135</td>\n",
       "      <td>16140.702605</td>\n",
       "      <td>1034.9397</td>\n",
       "      <td>5.075106</td>\n",
       "      <td>0.112442</td>\n",
       "      <td>1955.76070</td>\n",
       "      <td>0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.091581</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>RGI60-08.00005</td>\n",
       "      <td>626.451758</td>\n",
       "      <td>24237.964135</td>\n",
       "      <td>16049.763045</td>\n",
       "      <td>1052.8350</td>\n",
       "      <td>4.781625</td>\n",
       "      <td>0.125466</td>\n",
       "      <td>1916.53980</td>\n",
       "      <td>0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.083521</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>RGI60-08.00005</td>\n",
       "      <td>632.948060</td>\n",
       "      <td>24237.964135</td>\n",
       "      <td>16049.763045</td>\n",
       "      <td>1063.4498</td>\n",
       "      <td>4.490769</td>\n",
       "      <td>0.125620</td>\n",
       "      <td>1879.49460</td>\n",
       "      <td>0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.034988</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>RGI60-08.00005</td>\n",
       "      <td>637.226233</td>\n",
       "      <td>24237.964135</td>\n",
       "      <td>16049.763045</td>\n",
       "      <td>1070.4403</td>\n",
       "      <td>4.222532</td>\n",
       "      <td>0.151456</td>\n",
       "      <td>1844.75610</td>\n",
       "      <td>0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.034173</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>53565</th>\n",
       "      <td>RGI60-08.00014</td>\n",
       "      <td>1053.637775</td>\n",
       "      <td>23169.686477</td>\n",
       "      <td>11568.195539</td>\n",
       "      <td>1352.4573</td>\n",
       "      <td>1.996007</td>\n",
       "      <td>0.119973</td>\n",
       "      <td>357.61713</td>\n",
       "      <td>0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>-0.010636</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>53566</th>\n",
       "      <td>RGI60-08.00014</td>\n",
       "      <td>1051.288168</td>\n",
       "      <td>23169.686477</td>\n",
       "      <td>11721.653822</td>\n",
       "      <td>1349.4329</td>\n",
       "      <td>2.012404</td>\n",
       "      <td>0.114451</td>\n",
       "      <td>373.89438</td>\n",
       "      <td>0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>-0.016246</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>53567</th>\n",
       "      <td>RGI60-08.00014</td>\n",
       "      <td>1049.039072</td>\n",
       "      <td>23169.686477</td>\n",
       "      <td>11721.653822</td>\n",
       "      <td>1346.5271</td>\n",
       "      <td>2.080615</td>\n",
       "      <td>0.099085</td>\n",
       "      <td>390.62260</td>\n",
       "      <td>0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>-0.022051</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>53568</th>\n",
       "      <td>RGI60-08.00014</td>\n",
       "      <td>1046.852415</td>\n",
       "      <td>23169.686477</td>\n",
       "      <td>11721.653822</td>\n",
       "      <td>1343.7019</td>\n",
       "      <td>2.261690</td>\n",
       "      <td>0.075171</td>\n",
       "      <td>407.74625</td>\n",
       "      <td>0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>-0.024491</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>53569</th>\n",
       "      <td>RGI60-08.00014</td>\n",
       "      <td>1044.694706</td>\n",
       "      <td>23169.686477</td>\n",
       "      <td>11721.653822</td>\n",
       "      <td>1340.9142</td>\n",
       "      <td>2.432525</td>\n",
       "      <td>0.062766</td>\n",
       "      <td>425.21760</td>\n",
       "      <td>0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>-0.025489</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>53570 rows × 11 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "               RGI_ID       PDD_2D       rain_2D       snow_2D       topo  \\\n",
       "0      RGI60-08.00005   609.424791  24237.964135  16140.702605  1024.7790   \n",
       "1      RGI60-08.00005   615.582145  24237.964135  16140.702605  1034.9397   \n",
       "2      RGI60-08.00005   626.451758  24237.964135  16049.763045  1052.8350   \n",
       "3      RGI60-08.00005   632.948060  24237.964135  16049.763045  1063.4498   \n",
       "4      RGI60-08.00005   637.226233  24237.964135  16049.763045  1070.4403   \n",
       "...               ...          ...           ...           ...        ...   \n",
       "53565  RGI60-08.00014  1053.637775  23169.686477  11568.195539  1352.4573   \n",
       "53566  RGI60-08.00014  1051.288168  23169.686477  11721.653822  1349.4329   \n",
       "53567  RGI60-08.00014  1049.039072  23169.686477  11721.653822  1346.5271   \n",
       "53568  RGI60-08.00014  1046.852415  23169.686477  11721.653822  1343.7019   \n",
       "53569  RGI60-08.00014  1044.694706  23169.686477  11721.653822  1340.9142   \n",
       "\n",
       "         aspect     slope  dis_from_border  glacier_mask  \\\n",
       "0      5.331003  0.102807       1997.02900             0   \n",
       "1      5.075106  0.112442       1955.76070             0   \n",
       "2      4.781625  0.125466       1916.53980             0   \n",
       "3      4.490769  0.125620       1879.49460             0   \n",
       "4      4.222532  0.151456       1844.75610             0   \n",
       "...         ...       ...              ...           ...   \n",
       "53565  1.996007  0.119973        357.61713             0   \n",
       "53566  2.012404  0.114451        373.89438             0   \n",
       "53567  2.080615  0.099085        390.62260             0   \n",
       "53568  2.261690  0.075171        407.74625             0   \n",
       "53569  2.432525  0.062766        425.21760             0   \n",
       "\n",
       "       millan_ice_thickness  hugonnet_dhdt  \n",
       "0                       0.0       0.047467  \n",
       "1                       0.0       0.091581  \n",
       "2                       0.0       0.083521  \n",
       "3                       0.0       0.034988  \n",
       "4                       0.0       0.034173  \n",
       "...                     ...            ...  \n",
       "53565                   0.0      -0.010636  \n",
       "53566                   0.0      -0.016246  \n",
       "53567                   0.0      -0.022051  \n",
       "53568                   0.0      -0.024491  \n",
       "53569                   0.0      -0.025489  \n",
       "\n",
       "[53570 rows x 11 columns]"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df_subset"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c6eadd17-ee33-4f93-ac06-695405bd8bbb",
   "metadata": {
    "user_expressions": []
   },
   "source": [
    "### TO BE CONTINUED\n",
    "\n",
    "See how to build training matrices with all geodetic MB data, and how to re-structure the training dataset from the previous notebooks to give it the same shape as the flattened matrices of this dataset. For many data (e.g. climate) only available at low resolution, we can just repeat the same value for precipitations and do a rough downscaling for temperature.  "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1d67ffd7-53e0-4e98-bcff-e7c2c81753ac",
   "metadata": {
    "user_expressions": []
   },
   "source": [
    "## 6.2. Bonus project No2: Inferring glacier ice thickness\n",
    "\n",
    "<img src=\"Figures/eye_logo.png\" width=\"75\"/>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9df7a886-355d-4507-81a3-dcc57fc3a2dc",
   "metadata": {
    "user_expressions": []
   },
   "source": [
    "In this project, we propose to simulate something completely different to what we have been doing in the previous slides. Instead of working on glacier mass balance, we will try to infer glacier ice thickness from a data-driven perspective. OGGM provides access to many different glacier datasets, including Glathida, the global dataset of glacier ice thickness observations. These observations cover a reduced number of glaciers (ca 2700), and only partially in terms of glacier surface. Nonetheless, all the data points available can be used to train a machine learning model to infer glacier ice thickness from different input features of interest.\n",
    "\n",
    "These different input features could be:\n",
    "- Glacier surface elevation\n",
    "- Surface slope\n",
    "- Ice surface velocities"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e231f127-0c9a-4c10-aefb-b1bcd9db476f",
   "metadata": {},
   "outputs": [],
   "source": [
    "parent_path = os.path.dirname(Path().resolve())\n",
    "workspace_path = os.path.join(parent_path, 'OGGM_data_Finse_glathida')\n",
    "#workspace_path = '/home/jovyan/shared/glacier-ml-2022/Mass_Balance_ML_Modelling/Data'\n",
    "\n",
    "if not os.path.exists(workspace_path):\n",
    "    os.mkdir(workspace_path)\n",
    "\n",
    "cfg.PATHS['working_dir'] = workspace_path"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "08bfddbe-aa8d-498b-bf5e-56c1ec68b313",
   "metadata": {},
   "outputs": [],
   "source": [
    "gtd_file = utils.file_downloader('https://cluster.klima.uni-bremen.de/~oggm/glathida/glathida-v3.1.0/data/TTT_per_rgi_id.h5')\n",
    "\n",
    "glathida = pd.HDFStore(gtd_file)\n",
    "rgi_ids = glathida.keys()\n",
    "rgi_ids = [id[1:] for id in rgi_ids]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0e30bacc-e594-44aa-939c-3e4002ac29ca",
   "metadata": {},
   "outputs": [],
   "source": [
    "# We use the directories with the shop data in it: \"W5E5_w_data\"\n",
    "base_url = 'https://cluster.klima.uni-bremen.de/~oggm/gdirs/oggm_v1.6/L3-L5_files/2023.1/elev_bands/W5E5_w_data/'\n",
    "gdirs = workflow.init_glacier_directories(rgi_ids, from_prepro_level=3, prepro_base_url=base_url, prepro_border=10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f018694b-021c-44db-9654-7266758f66e9",
   "metadata": {},
   "outputs": [],
   "source": [
    "df = pd.read_hdf(gtd_file, key=rgi_id)\n",
    "df.plot(y='POINT_LAT', x='POINT_LON', c='THICKNESS', kind='scatter', cmap='viridis');"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}