Skip to content
Snippets Groups Projects
CANDOGLIA_t_unit.ipynb 7.66 MiB
Newer Older
Marco Mazzolini's avatar
Marco Mazzolini committed
{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "d2efdfdb",
   "metadata": {},
   "source": [
    "# CANDOGLIA RIVER BASIN\n",
    "\n",
    "13/10/2021\n",
    "\n",
    "In this notebook results of different feature selections are compared for the Candoglia basin (of which we have around 18 years of data)\n",
    "\n",
    "Input data is clipped from ERA5 metereological reanalysis quantile mapped and downscaled.\n",
    "\n",
    "Monthly averages (for the previous year) of pecipitation, temperature and potential evapotranspiration are selected as input.\n",
    "\n",
    "\n",
    "The settings are the following:\n",
    "\n",
    "    A) 180 features are selected with PCA, the same numeriosity as setting C) ;\n",
    "\n",
    "    B) 36 features are selectedwith PCA, the same numeriosity as setting D) ;\n",
    "        \n",
    "    C) metereological inputs spatial statistics are used as input: mean, the 5th, 25th, 75th and 95th quantiles are selected.\n",
    "        \n",
    "    D) metereological inputs are spatially averaged.\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "19d31cfb",
   "metadata": {},
   "source": [
    "import sys\n",
    "sys.path.append('/time_unit')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "efcc49ff",
   "metadata": {},
   "outputs": [],
   "source": [
    "from sf_runoff import daily_climatology, spatial_avg_daily_input, spatial_stats_daily_input, compute_anomalies\n",
    "from nested_CV import SVR_nested_CV_gridsearch, SVR_PCA_nested_CV_gridsearch\n",
    "from test import evaluate_prediction, plot_prediction, plot_anomalies\n",
    "from test import nested_CV_PCA_SVR_predict, nested_CV_SVR_predict\n",
    "from classic_CV_predict import classic_CV_PCA_SVR_predict, classic_CV_SVR_predict\n",
    "\n",
    "\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "from scipy.stats import gaussian_kde\n",
    "\n",
    "from sklearn.svm import SVR\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "from sklearn.pipeline import make_pipeline\n",
    "from sklearn.compose import TransformedTargetRegressor\n",
    "from sklearn.model_selection import GridSearchCV,TimeSeriesSplit\n",
    "from sklearn.metrics import mean_squared_error\n",
    "from sklearn.decomposition import PCA\n",
    "\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "import os\n",
    "\n",
    "import pdb\n",
    "import seaborn as sns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "cd514cbe",
   "metadata": {},
   "outputs": [],
   "source": [
    "path=r'C:\\Users\\mmazzolini\\OneDrive - Scientific Network South Tyrol\\Documents\\conda\\daily_input\\\\'\n",
    "\n",
    "daily_input = pd.read_csv(path+'CANDOGLIA_TOCE_2000_2019.csv', index_col=0, parse_dates=True)\n",
    "\n",
    "daily_input_TPE = spatial_avg_daily_input(daily_input)\n",
    "\n",
    "daily_input_stat = spatial_stats_daily_input(daily_input)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6799ce57",
   "metadata": {},
   "source": [
    "import sys, importlib\n",
    "importlib.reload(sys.modules['test'])\n"
   ]
  },
  {
   "cell_type": "code",
Marco Mazzolini's avatar
Marco Mazzolini committed
   "execution_count": 4,
Marco Mazzolini's avatar
Marco Mazzolini committed
   "id": "a024e5fe",
   "metadata": {},
   "outputs": [],
   "source": [
Marco Mazzolini's avatar
Marco Mazzolini committed
    "t_unit=30\n",
Marco Mazzolini's avatar
Marco Mazzolini committed
    "\n",
    "\n",
    "#define the possible parameters value (where Gridsearch is applied)\n",
    "\n",
    "C_range=np.logspace(-2, 2, 10)\n",
    "epsilon_range=np.logspace(-6, -2, 10)\n",
    "#n_range = [17, 50, 200]\n",
    "components_range = [5*3*12]\n",
    "#do not enlarge t_range for now\n",
Marco Mazzolini's avatar
Marco Mazzolini committed
    "t_range=[12]\n",
Marco Mazzolini's avatar
Marco Mazzolini committed
    "n_splits=5\n",
    "test_size=365"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e7d5c48a",
   "metadata": {},
   "source": [
    "# A) PCA+SVR"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "18861993",
   "metadata": {},
   "source": [
    "### TRAIN A PCA+SVR MODEL "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "aacb3a01",
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 't_unit' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m-----------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mNameError\u001b[0m                       Traceback (most recent call last)",
      "\u001b[1;32mC:\\Users\\MMAZZO~1\\AppData\\Local\\Temp/ipykernel_11492/1336338538.py\u001b[0m in \u001b[0;36m<module>\u001b[1;34m\u001b[0m\n\u001b[0;32m      1\u001b[0m \u001b[0mt_unit\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m      2\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m----> 3\u001b[1;33m \u001b[0mC\u001b[0m\u001b[1;33m,\u001b[0m\u001b[0meps\u001b[0m\u001b[1;33m,\u001b[0m\u001b[0mn\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mSVR_PCA_nested_CV_gridsearch\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mdaily_input\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mC_range\u001b[0m\u001b[1;33m,\u001b[0m  \u001b[0mepsilon_range\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mcomponents_range\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mt_range\u001b[0m\u001b[1;33m,\u001b[0m\u001b[0mn_splits\u001b[0m\u001b[1;33m,\u001b[0m\u001b[0mtest_size\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m      4\u001b[0m \u001b[0mprint\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;34mf'C={C}'\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m      5\u001b[0m \u001b[0mprint\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;34mf'eps={eps}'\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32m~\\OneDrive - Scientific Network South Tyrol\\Documents\\conda\\Runoff_prediction\\nested_CV.py\u001b[0m in \u001b[0;36mSVR_PCA_nested_CV_gridsearch\u001b[1;34m(daily_input, C_range, epsilon_range, components_range, t_range, n_splits, test_size)\u001b[0m\n\u001b[0;32m    104\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    105\u001b[0m     \u001b[1;32mfor\u001b[0m \u001b[0mt_length\u001b[0m \u001b[1;32min\u001b[0m \u001b[0mt_range\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 106\u001b[1;33m         \u001b[0mit_matrix\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mcreate_it_matrix\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mdaily_input\u001b[0m\u001b[1;33m,\u001b[0m\u001b[0mt_length\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    107\u001b[0m         \u001b[0mtscv\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mTimeSeriesSplit\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mgap\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mt_unit\u001b[0m \u001b[1;33m,\u001b[0m\u001b[0mn_splits\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mn_splits\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mtest_size\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mtest_size\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    108\u001b[0m         \u001b[0msets\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mtscv\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0msplit\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mit_matrix\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mindex\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32m~\\OneDrive - Scientific Network South Tyrol\\Documents\\conda\\Runoff_prediction\\sf_runoff.py\u001b[0m in \u001b[0;36mcreate_it_matrix\u001b[1;34m(daily_input, t_length)\u001b[0m\n\u001b[0;32m     37\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     38\u001b[0m     \u001b[1;31m# Compute the t_unit days average runoff\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m---> 39\u001b[1;33m     \u001b[0mrunoff_t_unit\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mrunoff\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mrolling\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mt_unit\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mmin_periods\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mt_unit\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mmean\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m     40\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     41\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;31mNameError\u001b[0m: name 't_unit' is not defined"
     ]
    }
   ],
   "source": [
    "C,eps,n=SVR_PCA_nested_CV_gridsearch(daily_input, C_range,  epsilon_range, components_range, t_range,n_splits,test_size)\n",
    "print(f'C={C}')\n",
    "print(f'eps={eps}')\n",
    "print(f'n={n}')"
   ]
  },
Marco Mazzolini's avatar
Marco Mazzolini committed
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "945f9297",
   "metadata": {},
   "outputs": [],
   "source": [
    "C=0.5994842503189409\n",
    "eps=5.994842503189409e-05\n",
    "n=180"
   ]
  },
Marco Mazzolini's avatar
Marco Mazzolini committed
  {
   "cell_type": "markdown",
   "id": "e9644d4d",
   "metadata": {},
   "source": [
    "C=0.21544346900318834\n",
    "epsilon=0.003593813663804626\n",
    "n=180\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c4d61ea5",
   "metadata": {},
   "source": [
    "### PREDICT RUNOFF ON TEST SET AND QUANTIFY THE PERFORMANCE"
   ]
  },
Loading
Loading full blame...