analyse_shap_paper.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Imports"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import json\n",
    "import os\n",
    "import torch\n",
    "\n",
    "import ccc\n",
    "import utils\n",
    "import utils_plots\n",
    "import utils_shap\n",
    "\n",
    "import ipywidgets as widgets\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import seaborn as sns\n",
    "\n",
    "from matplotlib.patches import Patch"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Settings"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### General"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "TARGET_MODE = 1\n",
    "GRADIENT_EXPLAINER = False\n",
    "MODEL = \"2022_02_21__11-11__ALDIS_paper\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Preparations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tmsubdir = f'targetmode_{TARGET_MODE}'\n",
    "modelpath = os.path.join(ccc.MODEL_ROOT_PATH, f'targetmode_{TARGET_MODE}', MODEL)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "_, model_name = utils.load_model(os.path.join(f'targetmode_{TARGET_MODE}', MODEL), torch.device(\"cpu\"), \"18\")\n",
    "\n",
    "prefix = \"_model_00018\"\n",
    "\n",
    "shap_path = os.path.join(modelpath, prefix + \"_shap_parquet_bg_by_lon_lat_no_flash\")\n",
    "df_path = os.path.join(modelpath, prefix + \"_test_df.pickle\")\n",
    "\n",
    "dd = pd.read_pickle(df_path)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dd.columns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "used_threshold = 0.8708981871604919  # F1 score based threshold of epoch 18 on validation set"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "vc_threshold = utils.getVeryConfidentThreshold(used_threshold)\n",
    "dd_transf = dd\n",
    "\n",
    "dd_transf.loc[:, \"pred_class\"] = np.where(dd_transf[\"output\"] > used_threshold, \"pred_flash\", \"pred_no_flash\")\n",
    "dd_transf.loc[:, \"real_class\"] = np.where(dd_transf[\"target\"] > 0.5, \"real_flash\", \"real_no_flash\")  # target col only contains 0s and 1s.\n",
    "\n",
    "dd_transf.loc[:, 'cat'] = np.select(\n",
    "    [\n",
    "        (dd_transf['pred_class'] == \"pred_flash\") & (dd_transf['real_class'] == \"real_flash\") & (dd_transf[\"output\"] < vc_threshold), \n",
    "        (dd_transf['pred_class'] == \"pred_flash\") & (dd_transf['real_class'] == \"real_flash\") & (dd_transf[\"output\"] >= vc_threshold), \n",
    "        (dd_transf['pred_class'] == \"pred_no_flash\") & (dd_transf['real_class'] == \"real_flash\"), \n",
    "        (dd_transf['pred_class'] == \"pred_flash\") & (dd_transf['real_class'] == \"real_no_flash\"), \n",
    "        (dd_transf['pred_class'] == \"pred_no_flash\") & (dd_transf['real_class'] == \"real_no_flash\"), \n",
    "    ], \n",
    "    [\n",
    "        'TP_LC', \n",
    "        'TP_VC',\n",
    "        'FN',\n",
    "        'FP',\n",
    "        'TN',\n",
    "    ], \n",
    "    default='ERROR'\n",
    ")\n",
    "\n",
    "dd_transf.loc[:, 'cluster'] = np.select(\n",
    "    [dd_transf['cat'] == \"TP_LC\", dd_transf['cat'] == \"TP_VC\", dd_transf['cat'] == \"FN\", dd_transf['cat'] == \"FP\", dd_transf['cat'] == \"TN\",],\n",
    "    [0, 1, 2, 3, 4,],\n",
    "    default=-1\n",
    ")\n",
    "\n",
    "dd_transf = dd_transf.rename(columns={\"output\": \"pred_score\"})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sns.histplot(dd_transf.query(\"cat in ['TP_LC', 'TP_VC']\")['pred_score'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "with open(os.path.join(modelpath, 'data_cfg.json'), 'r') as f:\n",
    "    config_data = json.load(f)\n",
    "\n",
    "with open(os.path.join(modelpath, 'model_cfg.json'), 'r') as f:\n",
    "    config_model = json.load(f)\n",
    "\n",
    "traincols = ccc.TRAIN_COLS\n",
    "\n",
    "cols = []\n",
    "\n",
    "for col in ccc.LVL_TRAIN_COLS:\n",
    "    for lvl_idx in range(74):\n",
    "        lvl = 64 + lvl_idx\n",
    "        cols.append(f\"{col}_lvl{lvl}\")\n",
    "\n",
    "cols.extend(ccc.INDEX_COLS)\n",
    "cols.extend(traincols)\n",
    "\n",
    "cols = list(set(cols))\n",
    "print(\"Load test data into spark df\", flush=True)\n",
    "test_package = utils.get_testdf_spark(config_data, cols + [\"cbh\", \"cth\"], None)\n",
    "sparkdd_test = test_package.drop(\"features\").drop(\"label\")\n",
    "\n",
    "print(\"Convert test data (excluding TNs) into spark df\", flush=True)\n",
    "spark = utils.getsparksession()\n",
    "sparkdd = spark.createDataFrame(dd_transf)\n",
    "\n",
    "print(\"Join the two dfs\")\n",
    "sparkdd = utils.joinDataframes(sparkdd, sparkdd_test)\n",
    "\n",
    "print(\"Convert to pandas df\")\n",
    "dd_enriched = sparkdd.toPandas()\n",
    "\n",
    "print(\"Free memory\")\n",
    "del sparkdd\n",
    "del sparkdd_test"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Clustering"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mask = (dd_enriched[\"cat\"] == 'TP_VC') | (dd_enriched[\"cat\"] == 'TP_LC')  # we only cluster for true positives\n",
    "dd_tp = dd_enriched.loc[mask, :]\n",
    "dd_fx = dd_enriched.loc[~mask, :]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dd_tp.head(10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "[len(dd_tp), len(dd_fx)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dd_enriched['cluster'].value_counts()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sample_size = 61431 # dd_enriched['cluster'].value_counts().min()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Plotting Profiles"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geoh_cols = [col for col in dd_enriched.columns if col.startswith(\"geoh_\")]\n",
    "\n",
    "df_many_cases = dd_enriched[ccc.INDEX_COLS + geoh_cols + ['cluster', 'cat']]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "len(df_many_cases)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df_many_cases_sampled = df_many_cases.groupby('cluster')[df_many_cases.columns].apply(lambda x: x.sample(n=sample_size) if len(x) > sample_size else x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df_many_cases_sampled['cluster'].value_counts()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "len(df_many_cases_sampled)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dshap = pd.read_parquet(shap_path) ## serves as test file\n",
    "dshap.rename(columns={utils_shap.colname_meta_infix(col) : col for col in ccc.INDEX_COLS}, inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df_many_cases_shap = utils.joinDataframes(df_many_cases_sampled, dshap)\n",
    "df_many_cases_shap.drop([\"flash_meta\"], axis=1, inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def invsig(y):\n",
    "    return np.log(y / (1 - y))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Used for debugging to quickly reload utils_plots.py\n",
    "import importlib\n",
    "importlib.reload(utils_plots)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ptype = \"q50\"  # can be mult, q50, q95\n",
    "use_cache = True\n",
    "write_cache = False\n",
    "\n",
    "separate_clusters = False\n",
    "\n",
    "plot_group = \"vartype_grouped\"  # \"confmat\": TP, FP, TN, FN; \"vartype\": CLOUD_HIGH, MASS_HIGH, WIND_HIGH; \"vartype_ext\": \"CLOUD_MASS_HIGH\", \"CLOUD_WIND_HIGH\", \"MASS_HIGH\", \"WIND_HIGH\"\n",
    "\n",
    "only_show_cols = []  # [\"ciwc\", \"cswc\"]\n",
    "\n",
    "y_axis = \"geopotential_altitude\"  # level, geopotential_altitude\n",
    "\n",
    "vis_save_path = os.path.join(modelpath, f\"{prefix}_shap_plots_bg_no_flash_{plot_group}\")\n",
    "\n",
    "if vis_save_path != \"\":\n",
    "    if not os.path.isdir(vis_save_path):\n",
    "        os.makedirs(vis_save_path)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "used_threshold = 0.8708981871604919  # maximum f1score threshold on validation set\n",
    "\n",
    "df_shap_to_plot = df_many_cases_shap.copy()\n",
    "\n",
    "df_plot_shap_cols = [df_shap_to_plot[[f\"{varname}_shapval_lvl{lvl}\" for lvl in range(64, 138)]].sum(axis=1) for varname in ccc.LVL_TRAIN_COLS]\n",
    "df_plot_shap_tmp = pd.concat(df_plot_shap_cols, axis=1)\n",
    "df_plot_shap_tmp.columns = [f\"{c}_shapsum\" for c in ccc.LVL_TRAIN_COLS]\n",
    "\n",
    "df_shap_to_plot = pd.concat([df_shap_to_plot, df_plot_shap_tmp], axis=1)\n",
    "\n",
    "traincols = []\n",
    "\n",
    "for traincol in ccc.LVL_TRAIN_COLS:\n",
    "    for idx in range(74):\n",
    "        lvl = 64 + idx\n",
    "        traincols.append(f\"{traincol}_shapval_lvl{lvl}\")\n",
    "            \n",
    "for c in traincols + [f\"{c}_shapsum\" for c in ccc.LVL_TRAIN_COLS]:\n",
    "    df_shap_to_plot[c] = df_shap_to_plot[c] / (invsig(used_threshold) - df_shap_to_plot[\"shap_base_value\"])\n",
    "\n",
    "if plot_group == \"vartype\":\n",
    "    plot_clusters = {\n",
    "                    1: 'TP_CLOUD_HIGH',\n",
    "                    2: 'TP_MASS_HIGH',\n",
    "                    3: 'TP_WIND_HIGH',\n",
    "                    4: 'TN',\n",
    "    }\n",
    "\n",
    "    cloud_sum = df_shap_to_plot['cswc_shapsum'] + df_shap_to_plot['ciwc_shapsum'] + df_shap_to_plot['crwc_shapsum'] + df_shap_to_plot['clwc_shapsum']\n",
    "    mass_sum = df_shap_to_plot['q_shapsum'] + df_shap_to_plot['t_shapsum']\n",
    "    wind_sum = df_shap_to_plot['u_shapsum'] + df_shap_to_plot['v_shapsum'] + df_shap_to_plot['w_shapsum']\n",
    "    \n",
    "    df_tempshap_cloud = df_shap_to_plot[df_shap_to_plot['cat'].isin(['TP_LC', 'TP_VC']) & (cloud_sum > 0.5)]\n",
    "    df_tempshap_mass = df_shap_to_plot[df_shap_to_plot['cat'].isin(['TP_LC', 'TP_VC']) & (mass_sum > 0.5)]\n",
    "    df_tempshap_wind = df_shap_to_plot[df_shap_to_plot['cat'].isin(['TP_LC', 'TP_VC']) & (wind_sum > 0.5)]\n",
    "    df_tempshap_tn = df_shap_to_plot[df_shap_to_plot['cat'].isin(['TN'])]\n",
    "    \n",
    "    df_tempshap_tn = df_shap_to_plot[df_shap_to_plot['cat'].isin(['TN'])]\n",
    "\n",
    "    df_tempshap_cloud.loc[:, \"cluster\"] = 1\n",
    "    df_tempshap_mass.loc[:, \"cluster\"] = 2\n",
    "    df_tempshap_wind.loc[:, \"cluster\"] = 3\n",
    "    df_tempshap_tn.loc[:, \"cluster\"] = 4\n",
    "\n",
    "    df_shap_to_plot = pd.concat([df_tempshap_cloud, df_tempshap_mass, df_tempshap_wind, df_tempshap_tn], ignore_index=True)\n",
    "elif plot_group == \"vartype_grouped\":\n",
    "    plot_clusters = {\n",
    "                    1: 'TP_CLOUD_MASS_HIGH',\n",
    "                    2: 'TP_CLOUD_WIND_HIGH',\n",
    "#                    3: 'TP_MASS_HIGH',\n",
    "#                    4: 'TP_WIND_HIGH',\n",
    "#                    5: 'TN',\n",
    "    }\n",
    "\n",
    "    cloud_sum = df_shap_to_plot['cswc_shapsum'] + df_shap_to_plot['ciwc_shapsum'] + df_shap_to_plot['crwc_shapsum'] + df_shap_to_plot['clwc_shapsum']\n",
    "    mass_sum = df_shap_to_plot['q_shapsum'] + df_shap_to_plot['t_shapsum']\n",
    "    wind_sum = df_shap_to_plot['u_shapsum'] + df_shap_to_plot['v_shapsum'] + df_shap_to_plot['w_shapsum']\n",
    "    \n",
    "    df_tempshap_cloudmass = df_shap_to_plot[df_shap_to_plot['cat'].isin(['TP_LC', 'TP_VC']) & (cloud_sum > 0.5) & (mass_sum > wind_sum)]\n",
    "    df_tempshap_cloudwind = df_shap_to_plot[df_shap_to_plot['cat'].isin(['TP_LC', 'TP_VC']) & (cloud_sum > 0.5) & (mass_sum <= wind_sum)]\n",
    "    df_tempshap_mass = df_shap_to_plot[df_shap_to_plot['cat'].isin(['TP_LC', 'TP_VC']) & (mass_sum > 0.5)]\n",
    "    df_tempshap_wind = df_shap_to_plot[df_shap_to_plot['cat'].isin(['TP_LC', 'TP_VC']) & (wind_sum > 0.5)]\n",
    "    df_tempshap_tn = df_shap_to_plot[df_shap_to_plot['cat'].isin(['TN'])]\n",
    "\n",
    "    df_tempshap_cloudmass.loc[:, \"cluster\"] = 1\n",
    "    df_tempshap_cloudwind.loc[:, \"cluster\"] = 2\n",
    "    df_tempshap_mass.loc[:, \"cluster\"] = 3\n",
    "    df_tempshap_wind.loc[:, \"cluster\"] = 4\n",
    "    df_tempshap_tn.loc[:, \"cluster\"] = 5\n",
    "\n",
    "    df_shap_to_plot = pd.concat([df_tempshap_cloudmass, df_tempshap_cloudwind, df_tempshap_mass, df_tempshap_wind, df_tempshap_tn], ignore_index=True)\n",
    "elif plot_group == \"confmat\":\n",
    "    plot_clusters = {\n",
    "                        0: 'TP less confident',\n",
    "#                        1: 'TP very confident',\n",
    "                        2: 'FN',\n",
    "#                        3: 'FP',\n",
    "                        4: 'TN',\n",
    "    }\n",
    "    \n",
    "    df_shap_to_plot.loc[:, 'cluster'] = np.select(\n",
    "        [df_shap_to_plot['cat'] == \"TP_LC\", df_shap_to_plot['cat'] == \"TP_VC\", df_shap_to_plot['cat'] == \"FN\", df_shap_to_plot['cat'] == \"FP\", df_shap_to_plot['cat'] == \"TN\",],\n",
    "        [0, 1, 2, 3, 4,],\n",
    "        default=-1\n",
    "    )\n",
    "\n",
    "palette = None\n",
    "\n",
    "df_shap_to_plot.query(\"cluster != -1\", inplace=True)\n",
    "df_shap_to_plot.drop(['cat'], axis=1, inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df_shap_to_plot[\"cluster\"].value_counts()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df_shap_to_plot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "utils_plots.plot_many_profiles(df_shap_to_plot, \"feature\", ptype=ptype, y_axis=y_axis, separate_clusters=separate_clusters, save_path=vis_save_path, use_cache=use_cache, plot_clusters=plot_clusters, only_show_cols=only_show_cols, write_cache=write_cache)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "utils_plots.plot_many_profiles(df_shap_to_plot, \"shap\", ptype=ptype, y_axis=y_axis, separate_clusters=separate_clusters, save_path=vis_save_path, use_cache=use_cache, plot_clusters=plot_clusters, only_show_cols=only_show_cols, write_cache=write_cache)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ucols = [\"longitude\", \"latitude\", \"year\", \"month\", \"day\", \"hour\", \"cbh\", \"cth\", \"cluster\"]\n",
    "df_cbh_cth_grouped = df_many_cases.reset_index()[ucols].groupby(ucols)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cbh_ls = []\n",
    "cth_ls = []\n",
    "cluster_ls = []\n",
    "\n",
    "for group, _ in df_cbh_cth_grouped:\n",
    "    cbh_ls.append(group[ucols.index(\"cbh\")])\n",
    "    cth_ls.append(group[ucols.index(\"cth\")])\n",
    "    cluster_ls.append(group[ucols.index(\"cluster\")])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df_cbh_cth = pd.DataFrame({\"cbh\" : cbh_ls, \"cth\" : cth_ls, \"cluster\" : cluster_ls})\n",
    "df_cbh_cth[\"ch\"] = df_cbh_cth[\"cth\"] - df_cbh_cth[\"cbh\"]\n",
    "df_cbh_cth = df_cbh_cth[df_cbh_cth[\"cluster\"].isin(plot_clusters.keys())]\n",
    "df_cbh_cth[\"cluster_labels\"] = df_cbh_cth[\"cluster\"].replace(plot_clusters)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "palette = sns.color_palette([utils_plots.CLUSTER_COLORS[col] for col in plot_clusters] if len(plot_clusters) > 0 else utils_plots.CLUSTER_COLORS)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "legend_elements = [Patch(color=utils_plots.CLUSTER_COLORS[plot_clusters[key]], label=plot_clusters[key]) for key in plot_clusters]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cloud_medians = df_cbh_cth.groupby(\"cluster\").median()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cloud_medians"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "g = plt.figure(figsize=(8, 15))\n",
    "graph = sns.violinplot(data=df_cbh_cth, x=\"cluster\", y=\"cth\", fliersize=3, palette=palette, cut=0)\n",
    "g.legend(handles=legend_elements, loc='center left', bbox_to_anchor=(0.9, 0.5))\n",
    "\n",
    "for idx in plot_clusters.keys():\n",
    "    graph.axhline(cloud_medians.query(f\"cluster == {idx}\")[\"cth\"].values[0], color=utils_plots.CLUSTER_COLORS[plot_clusters[idx]])\n",
    "\n",
    "ofile = os.path.join(vis_save_path, f\"cth_violinplot\")\n",
    "g.savefig(f\"{ofile}.pdf\")\n",
    "g.savefig(f\"{ofile}.png\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "g = plt.figure(figsize=(8, 15))\n",
    "graph = sns.violinplot(data=df_cbh_cth, x=\"cluster\", y=\"cbh\", fliersize=3, palette=palette, cut=0)\n",
    "g.legend(handles=legend_elements, loc='center left', bbox_to_anchor=(0.9, 0.5))\n",
    "\n",
    "for idx in plot_clusters.keys():\n",
    "    graph.axhline(cloud_medians.query(f\"cluster == {idx}\")[\"cbh\"].values[0], color=utils_plots.CLUSTER_COLORS[plot_clusters[idx]])\n",
    "\n",
    "ofile = os.path.join(vis_save_path, f\"violinplot_cbh\")\n",
    "g.savefig(f\"{ofile}.pdf\")\n",
    "g.savefig(f\"{ofile}.png\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "g = plt.figure(figsize=(8, 15))\n",
    "graph = sns.violinplot(data=df_cbh_cth, x=\"cluster\", y=\"ch\", fliersize=3, palette=palette, cut=0)\n",
    "g.legend(handles=legend_elements, loc='center left', bbox_to_anchor=(0.9, 0.5))\n",
    "\n",
    "for idx in plot_clusters.keys():\n",
    "    graph.axhline(cloud_medians.query(f\"cluster == {idx}\")[\"ch\"].values[0], color=utils_plots.CLUSTER_COLORS[plot_clusters[idx]])\n",
    "\n",
    "ofile = os.path.join(vis_save_path, f\"violinplot_ch\")\n",
    "g.savefig(f\"{ofile}.pdf\")\n",
    "g.savefig(f\"{ofile}.png\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "f77099b03e26f5f2a39d729b0e5ff7cce62b3e0b846f69c8115517d984996de7"
  },
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.10.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}