Shap_XGboost_future_prediction_US.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# # Shap for US - future prediction for June 2020."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import xgboost as xgb\n",
    "import shap\n",
    "# Initialize your Jupyter notebook with initjs(), otherwise you will get an error message.\n",
    "shap.initjs()\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "from sklearn.model_selection import TimeSeriesSplit\n",
    "from sklearn.model_selection import GridSearchCV\n",
    "import numpy as np"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define our search space for grid search\n",
    "random_grid = {\n",
    "    'objective' : ['reg:squarederror'],\n",
    "    'n_estimators': [50, 100, 150, 200],\n",
    "    'learning_rate': [0.01, 0.1, 0.2, 0.3],\n",
    "    'max_depth': range(3, 10),\n",
    "    'colsample_bytree': [i/10.0 for i in range(1, 3)],\n",
    "    }"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Used to fix the shap error\n",
    "def myfun(self=None):\n",
    "        return model_bytearray"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def non_zero_variables(df):\n",
    "    cols = df.columns\n",
    "    cols_non_0 = []\n",
    "    for col in cols:\n",
    "        if ((df[col].eq(0).sum(axis=0)) < 0.4 * (len(df[col]))): #if we have more than 60% zeros\n",
    "            cols_non_0.append(col)\n",
    "    return cols_non_0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "for country in ['US']: \n",
    "     \n",
    "    #Read all variables\n",
    "    all_var = pd.read_csv('../../all_variables_and_GPI_monthly_all_countries/all_variables_%s.csv' \n",
    "                          %country, index_col = 0)\n",
    "\n",
    "    #Delete all columns that have more than 60% of their values 0\n",
    "    variables_non_0 = non_zero_variables(all_var) #Filter the variables that have many zeros\n",
    "    df_country = all_var[variables_non_0]\n",
    "    \n",
    "    print(country, len(all_var))\n",
    "\n",
    "    df_country = all_var\n",
    "    \n",
    "    all_var_extra = pd.read_csv('../../all_variables_200803_202009/all_variables_%s.csv'%country, index_col = 0)\n",
    "    \n",
    "    all_var_extra = all_var_extra[all_var.columns[1:].tolist()]\n",
    "\n",
    "    df_future = all_var_extra.loc[all_var_extra.index >= 202004]\n",
    "    \n",
    "    \n",
    "    #Set the target variable\n",
    "    Y = df_country['GPI']\n",
    "\n",
    "    #Set the independent variables\n",
    "    X = df_country.loc[:, df_country.columns != 'GPI']\n",
    "\n",
    "    #Set the training sets:\n",
    "    X_train = X.tail(72) #corresponds to the length of the dataset\n",
    "    Y_train = Y.tail(72) #corresponds to the length of the dataset\n",
    "    \n",
    "    tscv = TimeSeriesSplit(n_splits=10).split(X_train)\n",
    "\n",
    "    #Train the model\n",
    "\n",
    "    xgb_reg = xgb.XGBRegressor() #model to tune\n",
    "\n",
    "    xgb_reg_random = GridSearchCV(estimator = xgb_reg, param_grid = random_grid,\n",
    "                                   cv = tscv,  n_jobs = -1)\n",
    "\n",
    "    #Best model\n",
    "    bmodel = xgb_reg_random.fit(X_train, Y_train)   \n",
    "    \n",
    "    #Pull the best estimated model from the gridsearch and send it to the TreeExplainer\n",
    "    model = bmodel.best_estimator_\n",
    "    \n",
    "    #Fix a shap error\n",
    "    mybooster = model.get_booster()\n",
    "    model_bytearray = mybooster.save_raw()[4:]\n",
    "    mybooster.save_raw = myfun\n",
    "    #Finish Fix a shap error"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Get the shap values\n",
    "shap_values = shap.TreeExplainer(mybooster).shap_values(X_train)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Global variable importance plot:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "shap.summary_plot(shap_values, X_train, alpha=1, plot_type=\"bar\", max_display=12, show=False)\n",
    "        \n",
    "\n",
    "plt.title('Global variable importance United States XGBoost model ', loc='right', fontsize=25, fontweight='bold' )\n",
    "\n",
    "\n",
    "plt.savefig('../../shap_summary_bar_%s.pdf' %country, bbox_inches='tight')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Individual Shap Value plot for June 2020:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def shap_plot(j):\n",
    "    explainerModel = shap.TreeExplainer(mybooster)\n",
    "    \n",
    "    shap_values_Model = explainerModel.shap_values(df_future)\n",
    "    p = shap.force_plot(explainerModel.expected_value, shap_values_Model[j], df_future.iloc[[j]].round(0).astype(object), show=False,matplotlib=True)\n",
    "    p.set_figwidth(40)\n",
    "    p.set_figheight(3)\n",
    "    \n",
    "    plt.tick_params(axis='x'\n",
    "                    , labelsize = 12)\n",
    "    \n",
    "    plt.suptitle('United States XGBoost model\\n Prediction for June 2020 ', fontsize=25, fontweight='bold', y=1.27, x=0.45)\n",
    "\n",
    "    p.savefig('../../shap_plot_US.pdf', bbox_inches='tight')\n",
    "    return(p)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "shap_plot(2)"
   ]
  }
 ],
 "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.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}