docs/notebooks/time_series_analysis.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Analyzing time series\n",
    "\n",
    "We will use the 'xclim' package and it's powerful time-series analysis tools to analyze the streamflow observations of the Salmon River basin. We will compute a few indicators, but you can refer to the xclim documentation to see how you can best make use of it for your specific needs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "\n",
    "import xarray as xr\n",
    "import xclim\n",
    "from pandas.plotting import register_matplotlib_converters\n",
    "\n",
    "from ravenpy.utilities.testdata import get_file, open_dataset\n",
    "\n",
    "register_matplotlib_converters()\n",
    "\n",
    "# Get the file we will use to analyze flows\n",
    "file = \"hydro_simulations/raven-gr4j-cemaneige-sim_hmets-0_Hydrographs.nc\"\n",
    "ds = open_dataset(file)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Base flow index\n",
    "\n",
    "The base flow index is the minimum 7-day average flow divided by the mean flow."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "help(xclim.land.base_flow_index)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The base flow index needs as input arguments a DataArray storing the stream flow time series, and the frequency at which the index is computed (`YS`: yearly, `QS-DEC`: seasonally)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "out = xclim.land.base_flow_index(ds.q_sim)\n",
    "out.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To compute generic statistics of a time series, use the `stats` process."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "help(xclim.generic.stats)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Here we compute the annual summer (JJA) minimum\n",
    "out = xclim.generic.stats(ds.q_sim, op=\"min\", season=\"JJA\")\n",
    "out.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Frequency analysis\n",
    "\n",
    "The process `freq_analysis` is similar to the previous stat in that it fits a series of annual maxima or minima to a statistical distribution, and returns the values corresponding to different return periods."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "help(xclim.generic.return_level)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For example, computing the Q(2,7), the minimum 7-days streamflow with a two-year reoccurrence, can be done using the following."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "out = xclim.generic.return_level(ds.q_sim, mode=\"min\", t=2, dist=\"gumbel_r\", window=7)\n",
    "out"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "An array of return periods can be passed."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "out = xclim.generic.return_level(\n",
    "    ds.q_sim, mode=\"max\", t=(2, 5, 10, 25, 50, 100), dist=\"gumbel_r\"\n",
    ")\n",
    "out.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Getting the parameters of the distribution and comparing the fit\n",
    "\n",
    "It's sometimes more useful to store the fitted parameters of the distribution rather than storing only the quantiles. In the example below, we're first computing the annual maxima of the simulated time series, then fitting them to a gumbel distribution using the `fit` process."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import json\n",
    "\n",
    "with xclim.set_options(\n",
    "    check_missing=\"pct\", missing_options={\"pct\": {\"tolerance\": 0.05}}\n",
    "):\n",
    "    ts = xclim.generic.stats(ds.q_sim, op=\"max\")\n",
    "\n",
    "ts"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "with xclim.set_options(check_missing=\"skip\"):\n",
    "    pa = xclim.generic.fit(ts.isel(nbasins=0), dist=\"gumbel_r\")\n",
    "pa"
   ]
  }
 ],
 "metadata": {
  "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.10.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}