{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "You can run this notebook in a [live session](https://binder.pangeo.io/v2/gh/pangeo-data/climpred/main?urlpath=lab/tree/docs/source/examples/subseasonal/weekly-subx-example.ipynb) [binder badge](https://binder.pangeo.io/v2/gh/pangeo-data/climpred/main?urlpath=lab/tree/docs/source/examples/subseasonal/weekly-subx-example.ipynb) or view it [on Github](https://github.com/pangeo-data/climpred/blob/main/docs/source/examples/subseasonal/weekly-subx-example.ipynb)." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Calculate skill of a MJO Index of SubX model GEOS_V2p1 as function of weekly lead time" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In this example, we demonstrate: \n", "\n", "1. How to remotely access data from the Subseasonal Experiment (SubX) hindcast database and set it up to be used in `climpred`. \n", "2. How to calculate the Anomaly Correlation Coefficient (ACC) using weekly data with `climpred`\n", "3. How to calculate and plot historical forecast skill of the real-time multivariate MJO (RMM) indices as function of lead time." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### The Subseasonal Experiment (SubX)\n", "\n", "Further information on SubX is available from [Pegion et al. 2019](https://journals.ametsoc.org/doi/full/10.1175/BAMS-D-18-0270.1) and the [SubX project website](http://cola.gmu.edu/subx/).\n", "\n", "The SubX public database is hosted on the International Research Institute for Climate and Society (IRI) data server [http://iridl.ldeo.columbia.edu/SOURCES/.Models/.SubX/](http://iridl.ldeo.columbia.edu/SOURCES/.Models/.SubX/).\n", "\n", "Since the SubX data server is accessed via this notebook, the time for the notebook to run may is several minutes and will vary depending on the speed that data can be downloaded. This is a large dataset, so please be patient. If you prefer to download SubX data locally, scripts are available from [https://github.com/kpegion/SubX](https://github.com/kpegion/SubX)." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Definitions\n", "\n", "**RMM**\n", ": Two indices (RMM1 and RMM2) are used to represent the MJO. Together they define the MJO based on 8 phases and can represent both the phase and amplitude of the MJO {cite:p}`Wheeler2004`. This example uses the observed RMM1 provided by Matthew Wheeler at the Center for Australian Weather and Climate Research. It is the version of the indices in which interannual variability has not been removed.\n", "\n", "**Skill of RMM**\n", ": Traditionally, the skill of the RMM is calculated as a bivariate correlation encompassing the skill of the two indices together {cite:p}`Rashid2011,Gottschalck2010`. Currently, `climpred` does not have the functionality to calculate the bivariate correlation, thus the anomaly correlation coefficient for RMM1 index is calculated here as a demonstration. The bivariate correlation metric will be added in a future version of `climpred`." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "data": { "application/javascript": [ "\n", " setTimeout(function() {\n", " var nbb_cell_id = 1;\n", " var nbb_unformatted_code = \"# linting\\n%load_ext nb_black\\n%load_ext lab_black\";\n", " var nbb_formatted_code = \"# linting\\n%load_ext nb_black\\n%load_ext lab_black\";\n", " var nbb_cells = Jupyter.notebook.get_cells();\n", " for (var i = 0; i < nbb_cells.length; ++i) {\n", " if (nbb_cells[i].input_prompt_number == nbb_cell_id) {\n", " if (nbb_cells[i].get_text() == nbb_unformatted_code) {\n", " nbb_cells[i].set_text(nbb_formatted_code);\n", " }\n", " break;\n", " }\n", " }\n", " }, 500);\n", " " ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# linting\n", "%load_ext nb_black\n", "%load_ext lab_black" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "data": { "application/javascript": [ "\n", " setTimeout(function() {\n", " var nbb_cell_id = 2;\n", " var nbb_unformatted_code = \"import warnings\\n\\nimport matplotlib.pyplot as plt\\n\\nimport xarray as xr\\nimport pandas as pd\\nimport numpy as np\\n\\nfrom climpred import HindcastEnsemble\\nimport climpred\\n\\nwarnings.filterwarnings(\\\"ignore\\\")\";\n", " var nbb_formatted_code = \"import warnings\\n\\nimport matplotlib.pyplot as plt\\n\\nimport xarray as xr\\nimport pandas as pd\\nimport numpy as np\\n\\nfrom climpred import HindcastEnsemble\\nimport climpred\\n\\nwarnings.filterwarnings(\\\"ignore\\\")\";\n", " var nbb_cells = Jupyter.notebook.get_cells();\n", " for (var i = 0; i < nbb_cells.length; ++i) {\n", " if (nbb_cells[i].input_prompt_number == nbb_cell_id) {\n", " if (nbb_cells[i].get_text() == nbb_unformatted_code) {\n", " nbb_cells[i].set_text(nbb_formatted_code);\n", " }\n", " break;\n", " }\n", " }\n", " }, 500);\n", " " ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "import warnings\n", "\n", "import matplotlib.pyplot as plt\n", "\n", "import xarray as xr\n", "import pandas as pd\n", "import numpy as np\n", "\n", "from climpred import HindcastEnsemble\n", "import climpred\n", "\n", "warnings.filterwarnings(\"ignore\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Read the observed RMM Indices" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.Dataset>\n",
       "Dimensions:  (time: 15613)\n",
       "Coordinates:\n",
       "  * time     (time) datetime64[ns] 1974-06-03 1974-06-04 ... 2017-07-24\n",
       "Data variables:\n",
       "    rmm1     (time) float64 1.867 1.868 1.924 1.574 ... 0.8918 0.9023 1.028
" ], "text/plain": [ "\n", "Dimensions: (time: 15613)\n", "Coordinates:\n", " * time (time) datetime64[ns] 1974-06-03 1974-06-04 ... 2017-07-24\n", "Data variables:\n", " rmm1 (time) float64 ..." ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" }, { "data": { "application/javascript": [ "\n", " setTimeout(function() {\n", " var nbb_cell_id = 3;\n", " var nbb_unformatted_code = \"obsds = climpred.tutorial.load_dataset(\\\"RMM-INTERANN-OBS\\\")[[\\\"rmm1\\\"]]\\nobsds\";\n", " var nbb_formatted_code = \"obsds = climpred.tutorial.load_dataset(\\\"RMM-INTERANN-OBS\\\")[[\\\"rmm1\\\"]]\\nobsds\";\n", " var nbb_cells = Jupyter.notebook.get_cells();\n", " for (var i = 0; i < nbb_cells.length; ++i) {\n", " if (nbb_cells[i].input_prompt_number == nbb_cell_id) {\n", " if (nbb_cells[i].get_text() == nbb_unformatted_code) {\n", " nbb_cells[i].set_text(nbb_formatted_code);\n", " }\n", " break;\n", " }\n", " }\n", " }, 500);\n", " " ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "obsds = climpred.tutorial.load_dataset(\"RMM-INTERANN-OBS\")[[\"rmm1\"]]\n", "obsds" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Read the SubX RMM1 data for the GMAO-GEOS_V2p1 model from the SubX data server. It is important to note that the SubX data contains weekly initialized forecasts where the `init` day varies by model. SubX data may have all NaNs for initial dates in which a model does not make a forecast, thus we apply `dropna` over the `S=init` dimension when `how=all` data for a given `S=init` is missing. This can be slow, but allows the rest of the calculations to go more quickly. \n", "\n", "Note that we ran the `dropna` operation offline and then uploaded the post-processed SubX dataset to the `climpred-data` repo for the purposes of this demo. This is how you can do this manually:\n", "\n", "```python\n", "url = 'http://iridl.ldeo.columbia.edu/SOURCES/.Models/.SubX/.GMAO/.GEOS_V2p1/.hindcast/.RMM/.RMM1/dods/'\n", "fcstds = xr.open_dataset(url, decode_times=False, chunks={'S': 1, 'L': 45}).dropna(dim='S',how='all')\n", "```" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.Dataset>\n",
       "Dimensions:  (S: 510, M: 4, L: 45)\n",
       "Coordinates:\n",
       "  * S        (S) float32 1.424e+04 1.425e+04 1.426e+04 ... 2.044e+04 2.045e+04\n",
       "  * M        (M) float32 1.0 2.0 3.0 4.0\n",
       "  * L        (L) float32 0.5 1.5 2.5 3.5 4.5 5.5 ... 40.5 41.5 42.5 43.5 44.5\n",
       "Data variables:\n",
       "    RMM1     (S, M, L) float32 -0.03237 0.1051 0.1946 ... -0.1499 0.08326 0.1389\n",
       "Attributes:\n",
       "    Conventions:  IRIDL
" ], "text/plain": [ "\n", "Dimensions: (S: 510, M: 4, L: 45)\n", "Coordinates:\n", " * S (S) float32 1.424e+04 1.425e+04 1.426e+04 ... 2.044e+04 2.045e+04\n", " * M (M) float32 1.0 2.0 3.0 4.0\n", " * L (L) float32 0.5 1.5 2.5 3.5 4.5 5.5 ... 40.5 41.5 42.5 43.5 44.5\n", "Data variables:\n", " RMM1 (S, M, L) float32 ...\n", "Attributes:\n", " Conventions: IRIDL" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" }, { "data": { "application/javascript": [ "\n", " setTimeout(function() {\n", " var nbb_cell_id = 4;\n", " var nbb_unformatted_code = \"fcstds = climpred.tutorial.load_dataset(\\\"GMAO-GEOS-RMM1\\\", decode_times=False)\\nfcstds\";\n", " var nbb_formatted_code = \"fcstds = climpred.tutorial.load_dataset(\\\"GMAO-GEOS-RMM1\\\", decode_times=False)\\nfcstds\";\n", " var nbb_cells = Jupyter.notebook.get_cells();\n", " for (var i = 0; i < nbb_cells.length; ++i) {\n", " if (nbb_cells[i].input_prompt_number == nbb_cell_id) {\n", " if (nbb_cells[i].get_text() == nbb_unformatted_code) {\n", " nbb_cells[i].set_text(nbb_formatted_code);\n", " }\n", " break;\n", " }\n", " }\n", " }, 500);\n", " " ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fcstds = climpred.tutorial.load_dataset(\"GMAO-GEOS-RMM1\", decode_times=False)\n", "fcstds" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The SubX data dimensions correspond to the following `climpred` dimension definitions: `X=lon`,`L=lead`,`Y=lat`,`M=member`, `S=init`. We will rename the dimensions to their `climpred` names." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "application/javascript": [ "\n", " setTimeout(function() {\n", " var nbb_cell_id = 5;\n", " var nbb_unformatted_code = \"fcstds = fcstds.rename({\\\"S\\\": \\\"init\\\", \\\"L\\\": \\\"lead\\\", \\\"M\\\": \\\"member\\\", \\\"RMM1\\\": \\\"rmm1\\\"})\";\n", " var nbb_formatted_code = \"fcstds = fcstds.rename({\\\"S\\\": \\\"init\\\", \\\"L\\\": \\\"lead\\\", \\\"M\\\": \\\"member\\\", \\\"RMM1\\\": \\\"rmm1\\\"})\";\n", " var nbb_cells = Jupyter.notebook.get_cells();\n", " for (var i = 0; i < nbb_cells.length; ++i) {\n", " if (nbb_cells[i].input_prompt_number == nbb_cell_id) {\n", " if (nbb_cells[i].get_text() == nbb_unformatted_code) {\n", " nbb_cells[i].set_text(nbb_formatted_code);\n", " }\n", " break;\n", " }\n", " }\n", " }, 500);\n", " " ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fcstds = fcstds.rename({\"S\": \"init\", \"L\": \"lead\", \"M\": \"member\", \"RMM1\": \"rmm1\"})" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's make sure that the `lead` dimension is set properly for `climpred`. SubX data stores `leads` as 0.5, 1.5, 2.5, etc, which correspond to 0, 1, 2, ... days since initialization. We will change the `lead` to be integers starting with zero. " ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "application/javascript": [ "\n", " setTimeout(function() {\n", " var nbb_cell_id = 6;\n", " var nbb_unformatted_code = \"fcstds[\\\"lead\\\"] = (fcstds[\\\"lead\\\"] - 0.5).astype(\\\"int\\\")\";\n", " var nbb_formatted_code = \"fcstds[\\\"lead\\\"] = (fcstds[\\\"lead\\\"] - 0.5).astype(\\\"int\\\")\";\n", " var nbb_cells = Jupyter.notebook.get_cells();\n", " for (var i = 0; i < nbb_cells.length; ++i) {\n", " if (nbb_cells[i].input_prompt_number == nbb_cell_id) {\n", " if (nbb_cells[i].get_text() == nbb_unformatted_code) {\n", " nbb_cells[i].set_text(nbb_formatted_code);\n", " }\n", " break;\n", " }\n", " }\n", " }, 500);\n", " " ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fcstds[\"lead\"] = (fcstds[\"lead\"] - 0.5).astype(\"int\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we need to make sure that the `init` dimension is set properly for `climpred`. We use an `xarray` convenience function to decode it into the proper CFTime calendar. It can detect that this is on a 360 day calendar." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "application/javascript": [ "\n", " setTimeout(function() {\n", " var nbb_cell_id = 7;\n", " var nbb_unformatted_code = \"fcstds = xr.decode_cf(fcstds, decode_times=True)\";\n", " var nbb_formatted_code = \"fcstds = xr.decode_cf(fcstds, decode_times=True)\";\n", " var nbb_cells = Jupyter.notebook.get_cells();\n", " for (var i = 0; i < nbb_cells.length; ++i) {\n", " if (nbb_cells[i].input_prompt_number == nbb_cell_id) {\n", " if (nbb_cells[i].get_text() == nbb_unformatted_code) {\n", " nbb_cells[i].set_text(nbb_formatted_code);\n", " }\n", " break;\n", " }\n", " }\n", " }, 500);\n", " " ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fcstds = xr.decode_cf(fcstds, decode_times=True)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Make Weekly/7-day Averages" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "Dimensions: (init: 510, member: 4, lead: 39)\n", "Coordinates:\n", " * init (init) datetime64[ns] 1999-01-01 1999-01-06 ... 2015-12-27\n", " * member (member) float32 1.0 2.0 3.0 4.0\n", " * lead (lead) int64 6 7 8 9 10 11 12 13 14 ... 36 37 38 39 40 41 42 43 44\n", "Data variables:\n", " rmm1 (init, member, lead) float32 0.1347 0.1095 ... 0.1404 0.07606\n", "Attributes:\n", " Conventions: IRIDL\n", "\n", "\n", "Dimensions: (time: 15456)\n", "Coordinates:\n", " * time (time) datetime64[ns] 1974-06-09 1974-06-10 ... 2017-07-24\n", "Data variables:\n", " rmm1 (time) float64 1.336 1.107 0.9046 0.695 ... 0.6265 0.673 0.7352\n" ] }, { "data": { "application/javascript": [ "\n", " setTimeout(function() {\n", " var nbb_cell_id = 8;\n", " var nbb_unformatted_code = \"fcstweekly = fcstds.rolling(lead=7, center=False).mean().dropna(dim=\\\"lead\\\")\\nobsweekly = obsds.rolling(time=7, center=False).mean().dropna(dim=\\\"time\\\")\\nprint(fcstweekly)\\nprint()\\nprint(obsweekly)\";\n", " var nbb_formatted_code = \"fcstweekly = fcstds.rolling(lead=7, center=False).mean().dropna(dim=\\\"lead\\\")\\nobsweekly = obsds.rolling(time=7, center=False).mean().dropna(dim=\\\"time\\\")\\nprint(fcstweekly)\\nprint()\\nprint(obsweekly)\";\n", " var nbb_cells = Jupyter.notebook.get_cells();\n", " for (var i = 0; i < nbb_cells.length; ++i) {\n", " if (nbb_cells[i].input_prompt_number == nbb_cell_id) {\n", " if (nbb_cells[i].get_text() == nbb_unformatted_code) {\n", " nbb_cells[i].set_text(nbb_formatted_code);\n", " }\n", " break;\n", " }\n", " }\n", " }, 500);\n", " " ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fcstweekly = fcstds.rolling(lead=7, center=False).mean().dropna(dim=\"lead\")\n", "obsweekly = obsds.rolling(time=7, center=False).mean().dropna(dim=\"time\")\n", "print(fcstweekly)\n", "print()\n", "print(obsweekly)" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "application/javascript": [ "\n", " setTimeout(function() {\n", " var nbb_cell_id = 9;\n", " var nbb_unformatted_code = \"fcstweekly = fcstweekly.isel(lead=slice(None, None, 7)).assign_coords(\\n lead=np.arange(1, 1 + fcstweekly.lead.size // 7 + 1)\\n)\";\n", " var nbb_formatted_code = \"fcstweekly = fcstweekly.isel(lead=slice(None, None, 7)).assign_coords(\\n lead=np.arange(1, 1 + fcstweekly.lead.size // 7 + 1)\\n)\";\n", " var nbb_cells = Jupyter.notebook.get_cells();\n", " for (var i = 0; i < nbb_cells.length; ++i) {\n", " if (nbb_cells[i].input_prompt_number == nbb_cell_id) {\n", " if (nbb_cells[i].get_text() == nbb_unformatted_code) {\n", " nbb_cells[i].set_text(nbb_formatted_code);\n", " }\n", " break;\n", " }\n", " }\n", " }, 500);\n", " " ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fcstweekly = fcstweekly.isel(lead=slice(None, None, 7)).assign_coords(\n", " lead=np.arange(1, 1 + fcstweekly.lead.size // 7 + 1)\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`climpred` requires that `lead` dimension has an attribute called `units` indicating what time units the `lead` is assocated with. Options are: `years`, `seasons`, `months`, `weeks`, `pentads`, `days`, `hours`, `minutes` or `seconds`. The `lead` `units` are `weeks`. " ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "application/javascript": [ "\n", " setTimeout(function() {\n", " var nbb_cell_id = 10;\n", " var nbb_unformatted_code = \"fcstweekly[\\\"lead\\\"].attrs = {\\\"units\\\": \\\"weeks\\\"}\";\n", " var nbb_formatted_code = \"fcstweekly[\\\"lead\\\"].attrs = {\\\"units\\\": \\\"weeks\\\"}\";\n", " var nbb_cells = Jupyter.notebook.get_cells();\n", " for (var i = 0; i < nbb_cells.length; ++i) {\n", " if (nbb_cells[i].input_prompt_number == nbb_cell_id) {\n", " if (nbb_cells[i].get_text() == nbb_unformatted_code) {\n", " nbb_cells[i].set_text(nbb_formatted_code);\n", " }\n", " break;\n", " }\n", " }\n", " }, 500);\n", " " ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fcstweekly[\"lead\"].attrs = {\"units\": \"weeks\"}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Create the {py:class}`.HindcastEnsemble` object and {py:meth}`.HindcastEnsemble.add_observations`." ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/html": [ "

climpred.HindcastEnsemble

" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<Initialized Ensemble>\n",
       "Dimensions:     (init: 510, member: 4, lead: 6)\n",
       "Coordinates:\n",
       "  * init        (init) object 1999-01-01 00:00:00 ... 2015-12-27 00:00:00\n",
       "  * member      (member) float32 1.0 2.0 3.0 4.0\n",
       "  * lead        (lead) int64 1 2 3 4 5 6\n",
       "    valid_time  (lead, init) object 1999-01-08 00:00:00 ... 2016-02-07 00:00:00\n",
       "Data variables:\n",
       "    rmm1        (init, member, lead) float32 0.1347 -0.7031 ... 0.2992 0.4478\n",
       "Attributes:\n",
       "    Conventions:  IRIDL
" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<Observations>\n",
       "Dimensions:  (time: 15456)\n",
       "Coordinates:\n",
       "  * time     (time) object 1974-06-09 00:00:00 ... 2017-07-24 00:00:00\n",
       "Data variables:\n",
       "    rmm1     (time) float64 1.336 1.107 0.9046 0.695 ... 0.6265 0.673 0.7352
" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" }, { "data": { "application/javascript": [ "\n", " setTimeout(function() {\n", " var nbb_cell_id = 11;\n", " var nbb_unformatted_code = \"hindcast = HindcastEnsemble(fcstweekly)\\nhindcast = hindcast.add_observations(obsweekly)\\nhindcast\";\n", " var nbb_formatted_code = \"hindcast = HindcastEnsemble(fcstweekly)\\nhindcast = hindcast.add_observations(obsweekly)\\nhindcast\";\n", " var nbb_cells = Jupyter.notebook.get_cells();\n", " for (var i = 0; i < nbb_cells.length; ++i) {\n", " if (nbb_cells[i].input_prompt_number == nbb_cell_id) {\n", " if (nbb_cells[i].get_text() == nbb_unformatted_code) {\n", " nbb_cells[i].set_text(nbb_formatted_code);\n", " }\n", " break;\n", " }\n", " }\n", " }, 500);\n", " " ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "hindcast = HindcastEnsemble(fcstweekly)\n", "hindcast = hindcast.add_observations(obsweekly)\n", "hindcast" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Calculate the Anomaly Correlation Coefficient (ACC) {py:func}`.climpred.metrics._pearson_r`." ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.Dataset>\n",
       "Dimensions:  (lead: 6)\n",
       "Coordinates:\n",
       "  * lead     (lead) int64 1 2 3 4 5 6\n",
       "    skill    <U11 'initialized'\n",
       "Data variables:\n",
       "    rmm1     (lead) float64 0.948 0.8719 0.7265 0.5658 0.4504 0.3607\n",
       "Attributes:\n",
       "    Conventions:                   IRIDL\n",
       "    prediction_skill_software:     climpred https://climpred.readthedocs.io/\n",
       "    skill_calculated_by_function:  HindcastEnsemble.verify()\n",
       "    number_of_initializations:     510\n",
       "    number_of_members:             4\n",
       "    alignment:                     maximize\n",
       "    metric:                        pearson_r\n",
       "    comparison:                    e2o\n",
       "    dim:                           init\n",
       "    reference:                     []
" ], "text/plain": [ "\n", "Dimensions: (lead: 6)\n", "Coordinates:\n", " * lead (lead) int64 1 2 3 4 5 6\n", " skill " ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "skill = hindcast.verify(\n", " metric=\"acc\", comparison=\"e2o\", alignment=\"maximize\", dim=\"init\"\n", ")\n", "skill" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Plot the skill as a function of lead time" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAETCAYAAAA7wAFvAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8/fFQqAAAACXBIWXMAAAsTAAALEwEAmpwYAAAaYUlEQVR4nO3debhdZX328e+dxDDPOSKSxGBliogMh4CACDLIJOkL6JtQgfiiUd8GtVBqWikCYovFSwq9oBInhMqMQiihoQqKDIEEgkBAaIQIQYYwSsIQgV//eJ6jKzt775OQrL3Jee7Pde0rez3rWWv91sk5695r3IoIzMysXIO6XYCZmXWXg8DMrHAOAjOzwjkIzMwK5yAwMyucg8DMrHAOAjPrOEnnSzqtxbjvSPrH/H5PSfMr4+ZJ2qdTdZbCQVAISeMk3S5pkaSn8/v/L0l5/PmSQtLYhunOzO0TGtr3zO1fabKs4ZJ+LOnZvLw7JB28jHVK0iRJ90h6WdKTkn4haVxDv4PzfBfl5fxY0vDK+AmS3pC0sOH17jx+d0m3SnpR0nOSbpG0U5u6dsnLWrvJuNm55i0kXS1pQZ7ndElbLst65/lMkjRL0muSzu+nb3X9/iDp19WfsaRR+f9ndsN0wyQtljSv0jYvtw1rsl4haVQe3kvSjflnNo9+SDpG0m8kvSTpKUnTJK3T33QR8fmI+Hp//WzlcRAUQNLxwFnAGcC7gI2BzwO7AUMrXR8CjqpMNwT4JPDbJrM9Gniu2j9PsyFwM7AYeD8wDDgTuEjS4ctQ7tnAl4HjgY2ATYETgf0ryzgcuAj41zz/9wOvATdL2qAyr9siYu2G1+8lrQv8J/BvwIZ5GafkeTQVETOA+cAS6yBpG2A0cDGwPjAV2JL0M74DuHoZ1rnP74HTgB8sY//bImLtvNxzgUskrd/QZ81cY58jgEeazOsRYHzfgKQPAGs29FmUazuhv8IkfQT4J2B8RKwDbA1c2t901iUR4dcAfgHrkf6AD+un3/nAt4CngA1y28HAdaQN+4RK37WAl4BxpA1+b2Xc14H7gEEN8/8K8DtAbWrYAnijOr8mfZTn83cN7YPyck/NwxOAm1vMoxd44S38LP8BuKGh7V+An7bovyEQwEZ5+GTgCtIG8SXgLuCDTaY7DTi/n1qWWD/SRjuAnfLwqDx8InBGpd8s4KvAvErbvNxvZqXtW7lfAKMalr1PdfoW9f0tcFU/v2+n5ffrADeSPgSoYdyewPyGWvfpxt/SQH55j2Dg+xCwGsv2yfTV3K/vMMxRwAVN+h0KLAQuB6aT9g767AtcGRFvNkxzGTCStLFv5aPAYxExq02fLfN8Lq825uVdmZffn4eANyT9SNIBDXsR7VwI7CFpBICkQaRP2D9q0X8P4MmIeLbSNjbXviFpr+YqSe9YxuU3JWkw8Gngj6SQrPoPYJykwZJGA2sDtzeZzQxgXUlb5/mNy9O+VbcDH5N0iqTdJK3WovaNgJ8Dt0TEFyNv7a2zHAQD3zDgmYh4va8hHxt/QdIrkvZo6H8BcFQ+xPAR4Kom8zwauDQi3iBtzMZVNmbDgCeaTPNEZXy7Wp+sNkian2t9VdJ7KtO3WkZ1/rvkaftevwWIiD8Au5M+7X4XWCBpqqSN29RGRDwG/AI4MjftTQrZaxv75vMV5wDHNYy6MyKuiIg/At8GVgd2abfcNnaR9AIpwL8FfCoinm7oMx94kPQp/ihSmLVyYe6zL/AA8PhbrIuI+BXpA8MOpJ/Ps5K+nUOmz7uBXwKXR8SJb3VZtuIcBAPfs8CwfLwfgIjYNSLWz+OW+B2IiJuBHtJhgf+MiFeq4/On4b2AH+emq0kbs4Py8DPAJk3q2KQyvl2tS0wbEcNJG/fVSIcN+qZvtYzq/GdExPqV119U5vtAREzI89+GtFH61za19fkRfw6CI4FL8kb9TyT1ANcD50bExQ3TP1ap4U3Shvrdy7DcZmbk/8cNSOcmPtyi3wWkQ0nj6T8Ijsh9m+0JLpeIuC4iPk7a+xmb5/uZSpeDgDWA76zosmzFOAgGvttIJ0HH9tex4j9IJ2ubbQyOJP3eXCPpSeBhUhD0HR76GXBoPmxS9UnSRvChNsu9ARguqbdNnwdJG89PVBvz8g4jHWZYLhHxG9Jx6W366Qrwk1zjXqRPvEscFsqHma4HpkbEN5pMP6Kh5uGkk8RvWUQsBL4AHClp+yZdriRtdB+OiEfbzOd3pJPGB5LWc6WIiDcj4uek/9/qz/i7wH8B0ySttbKWZ8vPQTDARcQLpCtizpV0uKR1JA2StB3ppG8zZ5MOD9zUZNzReX7bVV6HAQfm471nkk5Qf1/SuyStLmk8aQ/jhHbHgCPiQeA80tUv+0paIx9K2LXSJ0gnIk+UdESe/7uA7wHr5uW3JWkrScfnwzd9eznjScfJ24qIRaQTvj8Eflc9n5GvRppOOt49ucUsdpR0aN5D+zIppGfk6YdIWh0YDAzO6zakxXwa63qO9DM4qUXNH2XJT+OtHAN8NE+zhPx7szrwjjSo1SUNXWoOaeRYpUuWN1AyhnSosfFnPIkU7tdIWmMZ6rM6dPtstV+deQF/Rbqc8WVgAelk3kRgaB5/PvlKjSbT3kzard+FdDy6p0mfOcCk/H4k6XLK50hXLM0Exi5jnQK+CNwLvEI67v9L0h7FoEq/sXm+i/JyLgZGVMZPIF2BtLDhtRPpctHLSMfAF+V/zwPWXcYa9ySdX/hKQ/vRuX1RwzJH5vEns+RVQ7OBHSrTn5ynr75OblHDBBquiiLtXbwGbMufrxoa0mTaJa76ocWVOMAQKlcNVda7+vpFi/r2IO2dPZPX9SEqV3qx5JVBg0h7n9eT9i6r4/bEVw3V/lL+4ZpZzSSdDLwvIj7V7VrMqnxoyMyscLUFgaQfKD3K4L4W4yXpbElzlR4nsENdtdjbh6QPa+nHPiyUtLDbtcGfnnPTrD5f2WIDVm2HhvL16QuBCyJiqasxJB0IHEu6QmFn4KyI2LmWYszMrKXa9ggi4ibSSbxWxpJCIiI9x2V9Sc2uDTczsxp18xzBplRuriFdG75pl2oxMyvWMl2j3G2SJpIudWSttdbacauttupyRWZmq5Y777zzmYjoaTaum0HwOJW7LEnXQDd9tklETAGmAPT29sasWe2eSWZmZo0kNT6Q8E+6eWhoKunhZpK0C/BiRDR7kJiZmdWotj0CSReT7gocpvRVc18j3ZpORHwHmEa6Ymgu6W7XT9dVi5mZtVZbEETE+H7GB/DXdS3fzMyWje8sNjMrnIPAzKxwDgIzs8I5CMzMCucgMDMrnIPAzKxwq8QjJlaGUZOv7ejy5p1+UP+dzMzeBrxHYGZWOAeBmVnhHARmZoVzEJiZFc5BYGZWOAeBmVnhHARmZoVzEJiZFc5BYGZWOAeBmVnhHARmZoVzEJiZFc5BYGZWOAeBmVnhHARmZoVzEJiZFa6YL6YZ6PzFO2b2VnmPwMyscA4CM7PCOQjMzArnIDAzK5yDwMyscA4CM7PCOQjMzArnIDAzK5yDwMyscA4CM7PCOQjMzApXaxBI2l/Sg5LmSprcZPxISTdKmi3pHkkH1lmPmZktrbYgkDQYOAc4ABgNjJc0uqHbicBlEbE9MA44t656zMysuTr3CMYAcyPi4YhYDFwCjG3oE8C6+f16wO9rrMfMzJqoMwg2BR6rDM/PbVUnA5+SNB+YBhzbbEaSJkqaJWnWggUL6qjVzKxY3T5ZPB44PyKGAwcCF0paqqaImBIRvRHR29PT0/EizcwGsjqD4HFgRGV4eG6rOga4DCAibgNWB4bVWJOZmTWoMwhmAptL2kzSUNLJ4KkNfR4F9gaQtDUpCHzsx8ysg2oLgoh4HZgETAceIF0dNEfSqZIOyd2OBz4r6dfAxcCEiIi6ajIzs6XV+p3FETGNdBK42nZS5f39wG511mBmZu11+2SxmZl1mYPAzKxwDgIzs8I5CMzMCucgMDMrnIPAzKxwDgIzs8I5CMzMCucgMDMrnIPAzKxwtT5iwmxlGTX52o4ub97pB3V0eWbd5D0CM7PCOQjMzArnIDAzK5yDwMyscA4CM7PCOQjMzArnIDAzK5yDwMyscA4CM7PCOQjMzArnIDAzK5yDwMyscA4CM7PCOQjMzArnIDAzK5yDwMyscA4CM7PCOQjMzArnIDAzK5yDwMyscA4CM7PCOQjMzArnIDAzK1ytQSBpf0kPSporaXKLPp+UdL+kOZIuqrMeMzNb2pC6ZixpMHAOsC8wH5gpaWpE3F/psznw98BuEfG8pHfWVY+ZmTVX5x7BGGBuRDwcEYuBS4CxDX0+C5wTEc8DRMTTNdZjZmZN1BkEmwKPVYbn57aqLYAtJN0iaYak/Wusx8zMmqjt0NByLH9zYE9gOHCTpA9ExAvVTpImAhMBRo4c2eESzcwGtjr3CB4HRlSGh+e2qvnA1Ij4Y0Q8AjxECoYlRMSUiOiNiN6enp7aCjYzK1GdQTAT2FzSZpKGAuOAqQ19riLtDSBpGOlQ0cM11mRmZg1qC4KIeB2YBEwHHgAui4g5kk6VdEjuNh14VtL9wI3ACRHxbF01mZnZ0mo9RxAR04BpDW0nVd4HcFx+mZlZF/jOYjOzwjkIzMwK1+3LR82KN2rytR1d3rzTD+ro8uztz3sEZmaFcxCYmRXOQWBmVjgHgZlZ4RwEZmaFcxCYmRWuZRBIOkPS55q0f07S6fWWZWZmndJuj+CjwJQm7d8FDq6nHDMz67R2QbBafhbQEiLiTUD1lWRmZp3ULgheyd8pvITc9kp9JZmZWSe1e8TEScB1kk4D7sxtvaQvm/9yzXWZmVmHtAyCiLhO0l8CJwDH5ub7gMMi4t4O1GZmZh3QMggkrQ48FRFHN7T3SFo9Il6tvTozM6tdu3MEZwMfbtK+O3BmPeWYmVmntQuCHSPiJ42NEfFTYI/6SjIzs05qFwRrvsXpzMxsFdJug/60pDGNjbltQX0lmZlZJ7W7fPQE4DJJ57Pk5aNHAeNqrsvMzDqk5R5BRNwB7Ey6i3gC0Hf10NGkMDAzswGg7XcWR8RTwNck7QCMJ4XAHsCVHajNzMw6oN19BFuQNv7jgWeASwFFxF4dqs3MzDqg3R7Bb4BfAQdHxFwASX/TkarMzKxj2l01dCjwBHCjpO9K2hs/ddTMbMBpd7L4qogYB2wF3Eh60Nw7Jf27pP06VJ+ZmdWs7cligIhYBFwEXCRpA+ATwFeA62uuzcwGgFGTr+3o8uadflBHlzcQLNcdwhHxfERMiYi96yrIzMw6y4+KMDMrnIPAzKxwDgIzs8I5CMzMCucgMDMrnIPAzKxwtQaBpP0lPShprqTJbfodJikk9dZZj5mZLa22IJA0GDgHOAAYDYyXNLpJv3WALwG311WLmZm1VucewRhgbkQ8HBGLgUuAsU36fR34JvBqjbWYmVkLdQbBpsBjleH5ue1P8vccjIiItvegS5ooaZakWQsW+FsyzcxWpq6dLJY0CPg2cHx/ffNjLXojorenp6f+4szMClJnEDwOjKgMD89tfdYBtgF+IWkesAsw1SeMzcw6q84gmAlsLmkzSUNJX3g/tW9kRLwYEcMiYlREjAJmAIdExKwaazIzswa1BUFEvA5MAqYDDwCXRcQcSadKOqSu5ZqZ2fLp9/sIVkRETAOmNbSd1KLvnnXWYmZmzfnOYjOzwjkIzMwK5yAwMytcrecIzMwGuoHwnczeIzAzK5yDwMyscA4CM7PCOQjMzArnIDAzK5yDwMyscA4CM7PCOQjMzArnIDAzK5yDwMyscA4CM7PCOQjMzArnIDAzK5yDwMyscA4CM7PCOQjMzArnIDAzK5yDwMyscA4CM7PCOQjMzArnIDAzK5yDwMyscA4CM7PCOQjMzArnIDAzK5yDwMyscA4CM7PCOQjMzArnIDAzK5yDwMyscLUGgaT9JT0oaa6kyU3GHyfpfkn3SPq5pPfUWY+ZmS2ttiCQNBg4BzgAGA2MlzS6odtsoDcitgWuAP6lrnrMzKy5OvcIxgBzI+LhiFgMXAKMrXaIiBsj4uU8OAMYXmM9ZmbWRJ1BsCnwWGV4fm5r5RjguhrrMTOzJoZ0uwAASZ8CeoGPtBg/EZgIMHLkyA5WZmY28NW5R/A4MKIyPDy3LUHSPsBXgUMi4rVmM4qIKRHRGxG9PT09tRRrZlaqOoNgJrC5pM0kDQXGAVOrHSRtD5xHCoGna6zFzMxaqC0IIuJ1YBIwHXgAuCwi5kg6VdIhudsZwNrA5ZLuljS1xezMzKwmtZ4jiIhpwLSGtpMq7/epc/lmZtY/31lsZlY4B4GZWeEcBGZmhXMQmJkVzkFgZlY4B4GZWeEcBGZmhXMQmJkVzkFgZlY4B4GZWeEcBGZmhXMQmJkVzkFgZlY4B4GZWeEcBGZmhXMQmJkVzkFgZlY4B4GZWeEcBGZmhXMQmJkVzkFgZlY4B4GZWeEcBGZmhXMQmJkVzkFgZlY4B4GZWeEcBGZmhXMQmJkVzkFgZlY4B4GZWeEcBGZmhXMQmJkVzkFgZlY4B4GZWeEcBGZmhas1CCTtL+lBSXMlTW4yfjVJl+bxt0saVWc9Zma2tNqCQNJg4BzgAGA0MF7S6IZuxwDPR8T7gDOBb9ZVj5mZNVfnHsEYYG5EPBwRi4FLgLENfcYCP8rvrwD2lqQaazIzswaKiHpmLB0O7B8Rn8nDRwI7R8SkSp/7cp/5efi3uc8zDfOaCEzMg1sCD9ZSdHPDgGf67bXq8vqtugbyuoHXb2V7T0T0NBsxpINFvGURMQWY0o1lS5oVEb3dWHYneP1WXQN53cDr10l1Hhp6HBhRGR6e25r2kTQEWA94tsaazMysQZ1BMBPYXNJmkoYC44CpDX2mAkfn94cDN0Rdx6rMzKyp2g4NRcTrkiYB04HBwA8iYo6kU4FZETEV+D5woaS5wHOksHi76cohqQ7y+q26BvK6gdevY2o7WWxmZqsG31lsZlY4B4GZWeEcBGZmhXMQFEbSVpL2lrR2Q/v+3appZZE0RtJO+f1oScdJOrDbddVF0gXdrqEuknbP/3/7dbuWlUHSzpLWze/XkHSKpGskfVPSel2vzyeLl42kT0fED7tdx4qQ9EXgr4EHgO2AL0XE1XncXRGxQxfLWyGSvkZ6rtUQ4L+BnYEbgX2B6RHxjS6Wt8IkNV56LWAv4AaAiDik40WtRJLuiIgx+f1nSb+nPwX2A66JiNO7Wd+KkjQH+GC+mnIK8DL5sTq5/dCu1ucgWDaSHo2Ikd2uY0VIuhf4UEQszE96vQK4MCLOkjQ7IrbvboVvXV637YDVgCeB4RHxB0lrALdHxLbdrG9FSboLuB/4HhCkILiYfMl1RPyye9WtuOrvn6SZwIERsUDSWsCMiPhAdytcMZIeiIit8/slPnRJujsitutacawij5joFEn3tBoFbNzJWmoyKCIWAkTEPEl7AldIeg9pHVdlr0fEG8DLkn4bEX8AiIhXJL3Z5dpWhl7gS8BXgRMi4m5Jr6zqAVAxSNIGpMPViogFABGxSNLr3S1tpbivclTh15J6I2KWpC2AP3a7OAfBkjYGPgY839Au4NbOl7PSPSVpu4i4GyDvGRwM/ABYpT9xAYslrRkRLwM79jXm46+rfBBExJvAmZIuz/8+xcD6+10PuJP0txaSNomIJ/K5rFX9QwrAZ4CzJJ1IetDcbZIeAx7L47rKh4YqJH0f+GFE3Nxk3EURcUQXylppJA0nfXJ+ssm43SLili6UtVJIWi0iXmvSPgzYJCLu7UJZtZF0ELBbRPxDt2upk6Q1gY0j4pFu17Iy5BPGm5FCfH5EPNXlkgAHgZlZ8Xz5qJlZ4RwEZmaFcxCYmRXOQWBvO5IW1jDPefnEcbXtdkl3S3pU0oL8/m5Ju0q6YmXXkJf5l5JOkrS+pGf7vqNb0ockRT6hj6T1JD0nabn+RiWNyl8Buyx9h0q6KX8plBXMQWDFioid8408JwGXRsR2+XVrRBxe02L/Djg3Il4AngC2zu27ArPzvwC7AHfky0ZrERGLgZ8D/7euZdiqwUFgqwRJfyHpvyTdKelXkrbK7R/Pn+xnS/qZpI1z+0aSrpc0R9L3WI5r0aufqiVNkHSVpP/OexWT8jNwZkuaIWnDdvU1zHcL4LWI6PvC8lv584Z/V+DMhuFbJA2WdIakmZLukfS5yvxOqLSf0mR578117iTp/ZLuyHs890jaPHe7CvirZf3Z2MDkILBVxRTg2IjYEfhb4NzcfjOwS348wSWkT9wAXwNujoj3k55ZsyKPB9kGOBTYCfgG8HJe3m3AUf3UV7UbcFdl+Bb+vOF/L3A56Q5icvutwDHAixGxU17+Z5W+/nU/YHNgDOnRGjtK2qNvxpK2BK4EJkTETODzwFl5D6gXmJ+73pfnawXzsUF728t3l+4KXJ4PqUN6phDAcOBSSZsAQ4G+G4/2IG28iYhrJTXeLb48boyIl4CXJL0IXJPb7wW27ae+qk2ABZXhW4G/l7QZMC8iXlWyNunu6NuBL+Rl9B2qWo8UAPvl1+zcvnZufxToAa4GDo2I+/P424Cv5nMQP4mI/wGIiDckLZa0Tl5HK5CDwFYFg4AXWjyY69+Ab0fEVKVnJ51cw/Krdyy/WRl+k/Q31K6+qldIG3IAIuJ/JK0PfJy0oYb0mIVPk4JhYT6ZfGxETK/OSNLHgH+OiPMa2kcBL5ICYXfSg+qIiIsk3Q4cBEyT9LmIuCFPthrwaj+12wDmQ0P2tpcfIPeIpE8A5E/NH8yj1wMez++Prkx2E3BE7n8AsEGX6qt6AHhfQ9sM0sPk+oLgNuDLpMNGANOBL0h6R573FkpP5JwO/L+894CkTSW9M0+zGPg/wFGS+n4G7wUejoizSXsL2+b2jYBnIqLrDz6z7nEQ2NvRmpLmV17HkU5oHiPp18AcYGzuezLpkMydpId59TkF2EPpOfCHkj4h16lVfVU3Adv3XTKa3QKMAGbl4dtI5wv6HnL4PdKn+rvyCezzgCERcT1wEenhZfeSHim+Tt9MI2IRcDDwN5IOAT5JegLm3aRzHn1farMXcO0KrLcNAH7WkFkHSTqL9EUrP+t2LQCSfgJMjoiHul2LdY/3CMw665+ANbtdBKQbyoCrHALmPQIzs8J5j8DMrHAOAjOzwjkIzMwK5yAwMyucg8DMrHD/C3s+ExkV3TtQAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" }, { "data": { "application/javascript": [ "\n", " setTimeout(function() {\n", " var nbb_cell_id = 13;\n", " var nbb_unformatted_code = \"skill.rmm1.to_series().plot(kind=\\\"bar\\\")\\nplt.title(\\\"GMAO_GOES_V2p1 RMM1 Skill\\\")\\nplt.xlabel(\\\"Lead Time (Weeks)\\\")\\nplt.ylabel(\\\"ACC\\\")\\nplt.ylim(0.0, 1.0)\\nplt.show()\";\n", " var nbb_formatted_code = \"skill.rmm1.to_series().plot(kind=\\\"bar\\\")\\nplt.title(\\\"GMAO_GOES_V2p1 RMM1 Skill\\\")\\nplt.xlabel(\\\"Lead Time (Weeks)\\\")\\nplt.ylabel(\\\"ACC\\\")\\nplt.ylim(0.0, 1.0)\\nplt.show()\";\n", " var nbb_cells = Jupyter.notebook.get_cells();\n", " for (var i = 0; i < nbb_cells.length; ++i) {\n", " if (nbb_cells[i].input_prompt_number == nbb_cell_id) {\n", " if (nbb_cells[i].get_text() == nbb_unformatted_code) {\n", " nbb_cells[i].set_text(nbb_formatted_code);\n", " }\n", " break;\n", " }\n", " }\n", " }, 500);\n", " " ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "skill.rmm1.to_series().plot(kind=\"bar\")\n", "plt.title(\"GMAO_GOES_V2p1 RMM1 Skill\")\n", "plt.xlabel(\"Lead Time (Weeks)\")\n", "plt.ylabel(\"ACC\")\n", "plt.ylim(0.0, 1.0)\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### References\n", "\n", "```{eval-rst}\n", ".. bibliography::\n", " :filter: docname in docnames\n", "```" ] } ], "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.7.10" }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": true, "sideBar": true, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": true, "toc_window_display": false }, "varInspector": { "cols": { "lenName": 16, "lenType": 16, "lenVar": 40 }, "kernels_config": { "python": { "delete_cmd_postfix": "", "delete_cmd_prefix": "del ", "library": "var_list.py", "varRefreshCmd": "print(var_dic_list())" }, "r": { "delete_cmd_postfix": ") ", "delete_cmd_prefix": "rm(", "library": "var_list.r", "varRefreshCmd": "cat(var_dic_list()) " } }, "types_to_exclude": [ "module", "function", "builtin_function_or_method", "instance", "_Feature" ], "window_display": false } }, "nbformat": 4, "nbformat_minor": 4 }