notebooks/evaluate_the_baseline.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![header](../figures/DC_SSH_QG_mapping-banner.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "***\n",
    "**Authors:**  Datlas & IGE <br>\n",
    "**Copyright:** 2022 Datlas & IGE <br>\n",
    "**License:** MIT\n",
    "\n",
    "<div class=\"alert alert-block alert-success\">\n",
    "    <h1><center>Evaluate the baseline</center></h1>\n",
    "    <h5><center>\n",
    "    <h5><center>This notebook shows how to address the data challenge by performing the baseline: an optimal.</h5> \n",
    "</div> \n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "The notebook is structured as follow:\n",
    "\n",
    "     1) The available scores \n",
    "     2) Download the reference and reconstructed fields\n",
    "     3) Regrid reconstructed fields to the reference fields\n",
    "     4) Statistical scores\n",
    "     5) Spectral scores\n",
    "     6) Save results\n",
    "     7) Display leaderboard metrics"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import xarray as xr\n",
    "import numpy\n",
    "import warnings\n",
    "import xrft\n",
    "import logging\n",
    "import os\n",
    "import sys\n",
    "import pandas as pd\n",
    "warnings.filterwarnings('ignore')\n",
    "\n",
    "logger = logging.getLogger()\n",
    "logger.setLevel(logging.INFO)\n",
    "\n",
    "sys.path.append('..')\n",
    "\n",
    "from src.mod_oi import *\n",
    "from src.mod_inout import *\n",
    "from src.mod_regrid import *\n",
    "from src.mod_eval import *\n",
    "from src.mod_plot import *"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-info\" role=\"alert\">\n",
    "\n",
    "<h2> The available scores </h2>\n",
    "\n",
    "</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "- the statistical comparison is based on the RMSE-based score $RMSE_{S}$ defined as:\n",
    "\n",
    "$$RMSE_{S}(t) = 1 - \\frac{RMSE(t)}{RMS(SSH_{true})}$$\n",
    "\n",
    "\n",
    "where RMS is the root mean square function, and with:\n",
    "\n",
    "\n",
    "$$RMSE(t) = \\sqrt{ \\frac{1}{N} \\sum_{i=1}^N (SSH_{reconstruction}(t,i) - SSH_{true}(t,i))^2   }$$\n",
    "\n",
    "\n",
    "N is the number of pixels included in the study domain. \n",
    "\n",
    "A score of 1 means a perfect reconstruction, a zero score means strong incoherency between the reconstruction & reference fields \n",
    "\n",
    "- the spectral analysis is based on the wavenumber-frequency power spectral density score $PSD_{S}^{wf}$ defined as:\n",
    "\n",
    "$$PSD_{S}^{wf} = 1 - \\frac{PSD^{wf}(SSH_{reconstruction} - SSH_{true})}{PSD^{wf}(SSH_{true})}$$\n",
    "\n",
    "A score of 1 means a perfect reconstruction, a negative score means strong incoherency between the reconstruction & reference fields \n",
    "\n",
    "--- \n",
    "\n",
    "The **Leaderboard** summarizes the following key metrics:\n",
    "\n",
    "   1) $\\overline{RMSE_{S}(t)}$ : mean RMSE score (normalized RMSE) \n",
    " \n",
    "   2) $\\sigma(RMSE_{S}(t))$ : standard deviation of the RMSE(t) => give an insight on the temporal stability of the reconstruction\n",
    " \n",
    "   3) $\\lambda_{x}$ : the minimum spatial scale resolved (wavelength in degree) \n",
    " \n",
    "   4) $\\lambda_{t}$ : the minimum temporal scale resolved (wavelength in days)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-info\" role=\"alert\">\n",
    "\n",
    "<h2> Download the reference and reconstructed fields </h2>\n",
    "\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "time_min = numpy.datetime64('2012-10-22')                # domain min time\n",
    "time_max = numpy.datetime64('2012-12-02')                # domain max time"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# reconstrcuted SSH field\n",
    "input_file = '../results/ssh_reconstruction_2012-10-22-2012-12-02_jason1.nc'\n",
    "ds_oi1_grid = xr.open_dataset(input_file)\n",
    "ds_oi1_grid"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# reference SSH field\n",
    "# Note: dc_ref is used for reggriding step\n",
    "dc_ref = xr.open_mfdataset('../dc_ref/*.nc', combine='nested', concat_dim='time', parallel=True)\n",
    "dc_ref"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Note: dc_ref_sample is used for reggriding step (daily mean is enougth !!!!!)\n",
    "dc_ref_sample = dc_ref.sel(time=slice(time_min, time_max)).resample(time='1D').mean()\n",
    "del dc_ref\n",
    "dc_ref_sample"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-info\" role=\"alert\">\n",
    "\n",
    "<h2> Regrid reconstructed fields to the reference fields </h2>\n",
    "\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Regrid fields to same spatio-temporal grid\n",
    "ds_oi1_regrid = oi_regrid(ds_oi1_grid, dc_ref_sample)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-info\" role=\"alert\">\n",
    "\n",
    "<h2> Statistical scores </h2>\n",
    "\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Eval\n",
    "rmse_t_oi1, rmse_xy_oi1, leaderboard_nrmse, leaderboard_nrmse_std = rmse_based_scores(ds_oi1_regrid, dc_ref_sample)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-info\" role=\"alert\">\n",
    "\n",
    "<h2> Spectral scores </h2>\n",
    "\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "psd_oi1, leaderboard_psds_score, leaderboard_psdt_score  = psd_based_scores(ds_oi1_regrid, dc_ref_sample)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-info\" role=\"alert\">\n",
    "\n",
    "<h2> Save results </h2>\n",
    "\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define outputs\n",
    "output_directory = '../results/'\n",
    "if not os.path.exists(output_directory):\n",
    "    os.mkdir(output_directory)\n",
    "filename_rmse_t = output_directory + 'rmse_t_ssh_reconstruction_2012-10-22-2012-12-02_jason1.nc'\n",
    "filename_rmse_xy = output_directory + 'rmse_xy_ssh_reconstruction_2012-10-22-2012-12-02_jason1.nc'\n",
    "filename_psd = output_directory + 'psd_ssh_reconstruction_2012-10-22-2012-12-02_jason1.nc'\n",
    "filename_dc_ref_sample = output_directory + 'dc_ref_2012-10-22-2012-12-02_sample.nc'\n",
    "filename_oi_regrid = output_directory + 'ssh_reconstruction_regridded_2012-10-22-2012-12-02_jason1.nc'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Save to netcdf files\n",
    "rmse_t_oi1.to_netcdf(filename_rmse_t)\n",
    "rmse_xy_oi1.to_netcdf(filename_rmse_xy)\n",
    "psd_oi1.name = 'psd_score'\n",
    "psd_oi1.to_netcdf(filename_psd)\n",
    "dc_ref_sample.to_netcdf(filename_dc_ref_sample)\n",
    "ds_oi1_regrid.to_netcdf(filename_oi_regrid)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-info\" role=\"alert\">\n",
    "\n",
    "<h2> Display leaderboard metrics </h2>\n",
    "\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = [['demo 1 nadir', \n",
    "         leaderboard_nrmse, \n",
    "         leaderboard_nrmse_std, \n",
    "         leaderboard_psds_score, \n",
    "         leaderboard_psdt_score,\n",
    "        'Covariances not optimized',\n",
    "        'example_data_eval.ipynb']]\n",
    "Leaderboard = pd.DataFrame(data, \n",
    "                           columns=['Method', \n",
    "                                    \"µ(RMSE) \", \n",
    "                                    \"σ(RMSE)\", \n",
    "                                    'λx (degree)', \n",
    "                                    'λt (days)', \n",
    "                                    'Notes',\n",
    "                                    'Reference'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"Summary of the leaderboard metrics:\")\n",
    "Leaderboard"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(Leaderboard.to_markdown())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.6.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}