bias boxplot

bopen · Mar 7, 2024 · 1b3a61d · 1b3a61d
1 parent 91808c4
commit 1b3a61d
Showing 1 changed file with 39 additions and 36 deletions.
diff --git a/notebooks/wp4/extreme_temperature_indices.ipynb b/notebooks/wp4/extreme_temperature_indices.ipynb
@@ -85,7 +85,7 @@
     "assert timeseries in (\"annual\", \"DJF\", \"MAM\", \"JJA\", \"SON\")\n",
     "\n",
     "# Choose CORDEX or CMIP6\n",
-    "collection_id = \"CORDEX\"\n",
+    "collection_id = \"CMIP6\"\n",
     "assert collection_id in (\"CORDEX\", \"CMIP6\")\n",
     "\n",
     "# Define region for analysis\n",
@@ -789,48 +789,51 @@
    "metadata": {},
    "outputs": [],
    "source": [
+    "weights = collection_id == \"CMIP6\"\n",
     "mean_datasets = [\n",
-    "    diagnostics.spatial_weighted_mean(\n",
-    "        ds.expand_dims(model=[model]), weights=collection_id == \"CMIP6\"\n",
-    "    )\n",
+    "    diagnostics.spatial_weighted_mean(ds.expand_dims(model=[model]), weights=weights)\n",
     "    for model, ds in model_datasets.items()\n",
     "]\n",
     "mean_ds = xr.concat(mean_datasets, \"model\")\n",
-    "for index, da in mean_ds[\"trend\"].groupby(\"index\"):\n",
-    "    df_slope = da.to_dataframe()[[\"trend\"]]\n",
-    "    ax = df_slope.boxplot()\n",
-    "    ax.scatter(\n",
-    "        x=[1] * len(df_slope),\n",
-    "        y=df_slope,\n",
-    "        color=\"grey\",\n",
-    "        marker=\".\",\n",
-    "        label=\"models\",\n",
-    "    )\n",
+    "mean_bias_ds = diagnostics.spatial_weighted_mean(bias, weights=weights)\n",
+    "for is_bias, ds in zip((False, True), (mean_ds, mean_bias_ds)):\n",
+    "    for index, da in ds[\"trend\"].groupby(\"index\"):\n",
+    "        df_slope = da.to_dataframe()[[\"trend\"]]\n",
+    "        ax = df_slope.boxplot()\n",
+    "        ax.scatter(\n",
+    "            x=[1] * len(df_slope),\n",
+    "            y=df_slope,\n",
+    "            color=\"grey\",\n",
+    "            marker=\".\",\n",
+    "            label=\"models\",\n",
+    "        )\n",
     "\n",
-    "    # Ensemble mean\n",
-    "    ax.scatter(\n",
-    "        x=1,\n",
-    "        y=da.mean(\"model\"),\n",
-    "        marker=\"o\",\n",
-    "        label=f\"{collection_id} Ensemble Mean\",\n",
-    "    )\n",
+    "        # Ensemble mean\n",
+    "        ax.scatter(\n",
+    "            x=1,\n",
+    "            y=da.mean(\"model\"),\n",
+    "            marker=\"o\",\n",
+    "            label=f\"{collection_id} Ensemble Mean\",\n",
+    "        )\n",
     "\n",
-    "    # ERA5\n",
-    "    da = ds_era5[\"trend\"].sel(index=index)\n",
-    "    da = diagnostics.spatial_weighted_mean(da)\n",
-    "    ax.scatter(\n",
-    "        x=2,\n",
-    "        y=da.values,\n",
-    "        marker=\"o\",\n",
-    "        label=\"ERA5\",\n",
-    "    )\n",
+    "        # ERA5\n",
+    "        labels = [f\"{collection_id} Ensemble\"]\n",
+    "        if not is_bias:\n",
+    "            da = ds_era5[\"trend\"].sel(index=index)\n",
+    "            da = diagnostics.spatial_weighted_mean(da)\n",
+    "            ax.scatter(\n",
+    "                x=2,\n",
+    "                y=da.values,\n",
+    "                marker=\"o\",\n",
+    "                label=\"ERA5\",\n",
+    "            )\n",
+    "            labels.append(\"ERA5\")\n",
     "\n",
-    "    # Figure settings\n",
-    "    ax.set_xticks([1, 2], [f\"{collection_id} Ensemble\", \"ERA5\"])\n",
-    "    ax.set_ylabel(da.attrs[\"units\"])\n",
-    "    plt.suptitle(f\"Trend of {index}\")\n",
-    "    plt.legend()\n",
-    "    plt.show()"
+    "        ax.set_xticks(range(1, len(labels) + 1), labels)\n",
+    "        ax.set_ylabel(da.attrs[\"units\"])\n",
+    "        plt.suptitle(f\"Trend{' bias ' if is_bias else ' '}of {index}\")\n",
+    "        plt.legend()\n",
+    "        plt.show()"
    ]
   }
  ],