{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "0de7782b",
   "metadata": {},
   "source": [
    "# Comparing Different Ranking Methods\n",
    "\n",
    "There are several options with skexplain to rank features. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "8571d857",
   "metadata": {},
   "outputs": [],
   "source": [
    "import sys, os \n",
    "current_dir = os.getcwd()\n",
    "path = os.path.dirname(current_dir)\n",
    "sys.path.append(path)\n",
    "import skexplain\n",
    "import plotting_config\n",
    "import shap\n",
    "from skexplain.common.utils import shap_values_to_importance, coefficients_to_importance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "5fd55d65",
   "metadata": {},
   "outputs": [],
   "source": [
    "estimators = skexplain.load_models()\n",
    "X,y = skexplain.load_data()\n",
    "X = X.astype({'urban': 'category', 'rural':'category'})"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b727eef3",
   "metadata": {},
   "source": [
    "For demonstration purposes, we will evaluate how important predictors vary for the logistic regression model. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "61749bbd",
   "metadata": {},
   "outputs": [],
   "source": "explainer = skexplain.InterpretToolkit(estimators=estimators[-1], X=X, y=y,)\nexplainer.set_plotting_config(\n    display_feature_names=plotting_config.display_feature_names,\n    display_units=plotting_config.display_units,\n    feature_colors=plotting_config.color_dict,\n)"
  },
  {
   "cell_type": "markdown",
   "id": "4e2895ad",
   "metadata": {},
   "source": [
    "## Using ALE for feature rankings \n",
    "\n",
    "Though feature importance is defined with respect to model performance, we can also rank features by their magnitude of their first-order effect. Inspired by Greenwell et al. (2018, https://arxiv.org/abs/1805.04755), we can compute the standard deviation of the ALE (rather than the partial dependence) and those features with the highest standard deviation have the biggest range of contributions to a model's performance. This method has limitations though: \n",
    "* The standard deviation is affected by outliers and may skew the results. \n",
    "* When feature interactions are strong, the 1-D expected contribution may not fully capture the strength of the first-order effect.\n",
    "* Only considers first-order effects\n",
    "* Does not explicitly measure the contribution of features to the model's accuracy\n",
    "* Unclear when the importance score is significantly/practically different from zero"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "ee4cc11c",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "3dc2317ec85e4a4a92482f4d75853395",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "  0%|          | 0/28 [00:00<?, ?it/s]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "fab95171e30747599b8d08bd1fbe6f60",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "  0%|          | 0/2 [00:00<?, ?it/s]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Arrays of bytes/strings is being converted to decimal numbers if dtype='numeric'. This behavior is deprecated in 0.24 and will be removed in 1.1 (renaming of 0.26). Please convert your data to numeric values explicitly instead.\n",
      "Arrays of bytes/strings is being converted to decimal numbers if dtype='numeric'. This behavior is deprecated in 0.24 and will be removed in 1.1 (renaming of 0.26). Please convert your data to numeric values explicitly instead.\n",
      "Arrays of bytes/strings is being converted to decimal numbers if dtype='numeric'. This behavior is deprecated in 0.24 and will be removed in 1.1 (renaming of 0.26). Please convert your data to numeric values explicitly instead.\n",
      "Arrays of bytes/strings is being converted to decimal numbers if dtype='numeric'. This behavior is deprecated in 0.24 and will be removed in 1.1 (renaming of 0.26). Please convert your data to numeric values explicitly instead.\n",
      "Arrays of bytes/strings is being converted to decimal numbers if dtype='numeric'. This behavior is deprecated in 0.24 and will be removed in 1.1 (renaming of 0.26). Please convert your data to numeric values explicitly instead.\n",
      "Arrays of bytes/strings is being converted to decimal numbers if dtype='numeric'. This behavior is deprecated in 0.24 and will be removed in 1.1 (renaming of 0.26). Please convert your data to numeric values explicitly instead.\n"
     ]
    }
   ],
   "source": [
    "ale_1d_ds = explainer.ale(features='all',n_bootstrap=1,subsample=1000,n_jobs=8, n_bins=20)\n",
    "ale_var_1d = explainer.ale_variance(ale=ale_1d_ds)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8b10f0cd",
   "metadata": {},
   "source": [
    "## SHAP-based importance \n",
    "\n",
    "We can compile SHAP values over several examples as a measure of importance. In skexplain, we have two ways of summarizing the SHAP values. First, `method=\"sum\"` in `shap_values_to_importance` will take the sum of the absolute value of SHAP values per feature. The highest total is considered the most important feature and so on. Second, `method=\"std\"`, will take the standard deviation of the SHAP values per feature similar to the ALE variance. The higher the variance, the more important the feature is considered. \n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c4d2f3ce",
   "metadata": {},
   "outputs": [],
   "source": "X_subset = shap.sample(X, 10, random_state=22)\nexplainer = skexplain.InterpretToolkit(estimators[-1], X=X_subset,)\nexplainer.set_plotting_config(\n    display_feature_names=plotting_config.display_feature_names,\n    display_units=plotting_config.display_units,\n    feature_colors=plotting_config.color_dict,\n)\nbackground_dataset = shap.sample(X, 100)\nresults = explainer.local_attributions(method='shap', shap_kws={'masker': background_dataset, 'algorithm': 'auto'})\nshap_values = results['shap_values__Logistic Regression'].values\nshap_rank = shap_values_to_importance(shap_values, estimator_name='Logistic Regression', feature_names=X.columns)"
  },
  {
   "cell_type": "markdown",
   "id": "61b26bf1",
   "metadata": {},
   "source": [
    "## Multi- and Single-pass Permutation Importance \n",
    "\n",
    "In the [permutatation importance notebook](https://github.com/monte-flora/mintpy/blob/master/tutorial_notebooks/permutation_importance_tutorial.ipynb), we already compute the multi-pass and single-pass permutation importance. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "e28140eb",
   "metadata": {},
   "outputs": [],
   "source": [
    "perm_results = explainer.load(fnames='multipass_importance_naupdc.nc')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fb6e9362",
   "metadata": {},
   "source": [
    "## Model Coefficients \n",
    "\n",
    "For a logistic regression model in particular, we have access to the model coefficients which can be interpreted as importance. There are several factors to be aware of when interpreting model coefficients as importance\n",
    "\n",
    "* Does not consider model performance \n",
    "* Assumes the variance between each feature is identical \n",
    "* Assumes independent features "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "504c2c12",
   "metadata": {},
   "outputs": [],
   "source": [
    "coef_rank = coefficients_to_importance(coefficients=estimators[-1][1].coef_[0], \n",
    "                                       estimator_name='Logistic Regression', \n",
    "                                       feature_names=X.columns)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a40bdf3d",
   "metadata": {},
   "outputs": [],
   "source": "data = [perm_results, perm_results, ale_var_1d,  shap_rank, coef_rank]\npanels = [('multipass', 'Logistic Regression'), \n          ('singlepass', 'Logistic Regression'), \n          ('ale_variance', 'Logistic Regression'), \n          ('shap', 'Logistic Regression'), \n          ('coefs', 'Logistic Regression'), \n          ]\n\nfig = explainer.plot_importance(data=data, \n                                panels=panels, \n                                num_vars_to_plot=10, \n                               )"
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "dc5b64ef",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "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.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}