03-model-training.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Model Training\n",
    "Model training can be a never ending step in machine learning, so knowing when to stop experimenting is inmportant. In this notebook we will look at how different models compare and how to finetune them. We will use a series of techniques to improve the accuracy of a model.\n",
    "\n",
    "First we will load in our cleaned data from the last notebook."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 98,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "df = pd.read_parquet(\"solar_cleaned.parquet\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Split Data\n",
    "Now we will split the data into a testing set and a training set. Creating this test set is called hold-out validation and will be used as a sanity check for overfitting. Overfitting is when the model learns from the training set too well and preforms poorly on new data. There are other more in-depth validation techniques such as [cross-validation](https://en.wikipedia.org/wiki/Cross-validation_(statistics)), but we wont cover that in this demo."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 99,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Training shape: (30251, 39)\n",
      "Testing shape: (12965, 39)\n"
     ]
    }
   ],
   "source": [
    "middle = int(len(df) * 0.7)\n",
    "\n",
    "train = df.loc[:middle].copy()\n",
    "test = df.loc[middle:]\n",
    "\n",
    "print(\"Training shape: \" + str(train.shape))\n",
    "print(\"Testing shape: \" + str(test.shape))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will also load in the feature engineering pipeline from the previous notebook. Recall that this pipeline converts the raw data into normalized data that the model can use. We use this to seperate the inputs and outputs (x and y) of our training and testing data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 100,
   "metadata": {},
   "outputs": [],
   "source": [
    "import cloudpickle as cp\n",
    "#feature_pipeline = cp.load(open('feature_pipeline.sav', 'rb'))\n",
    "\n",
    "#X_train = feature_pipeline.fit_transform(train)\n",
    "X_train = train.drop(columns=['dni_efficiency', 'ghi_efficiency', 'dhi_efficiency', 'STATION','DATE',\n",
    "                                            'latitude','longitude'])\n",
    "y_train = train['dni_efficiency']\n",
    "\n",
    "#X_test = feature_pipeline.fit_transform(test)\n",
    "X_test = test.drop(columns=['dni_efficiency', 'ghi_efficiency', 'dhi_efficiency', 'STATION','DATE',\n",
    "                                            'latitude','longitude'])\n",
    "y_test = test['dni_efficiency']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Training\n",
    "We will train two different regression models: random forrest and xgboost. I removed the DHI and GHI efficiency outputs from the training set sense they were giving over 100% efficiency scores. This would cause the model to start predicting wrong values if they were used. Because we are only predicting efficiency, not training all of our outputs will help the model avoid more errors in the data.\n",
    "\n",
    "### Random Forest\n",
    "First we use sklearn's [Random Forest Regressor](https://en.wikipedia.org/wiki/Random_forest) for our first model. Random forest models work by creating many [decision trees](https://en.wikipedia.org/wiki/Decision_tree) and aggregating the average prediction of the combined trees. Each tree will split, for a desired amount of iterations, based on a random feature. There are many more parameters that can be tuned but to start we will just use the default."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 101,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "RandomForestRegressor()"
      ]
     },
     "execution_count": 101,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from sklearn.ensemble import RandomForestRegressor\n",
    "from sklearn import model_selection\n",
    "\n",
    "rfr = RandomForestRegressor()\n",
    "rfr.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### XGBoost\n",
    "In parallel, we will train a gradient boosted tree model called [XGBoost](https://en.wikipedia.org/wiki/XGBoost). These trees differ from random forest trees in that each tree builds off the last compared to aggregating each tree. We will train it on default parameters also."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 102,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "XGBRegressor(base_score=0.5, booster='gbtree', colsample_bylevel=1,\n",
       "             colsample_bynode=1, colsample_bytree=1, gamma=0, gpu_id=-1,\n",
       "             importance_type='gain', interaction_constraints='',\n",
       "             learning_rate=0.300000012, max_delta_step=0, max_depth=6,\n",
       "             min_child_weight=1, missing=nan, monotone_constraints='()',\n",
       "             n_estimators=100, n_jobs=0, num_parallel_tree=1, random_state=0,\n",
       "             reg_alpha=0, reg_lambda=1, scale_pos_weight=1, subsample=1,\n",
       "             tree_method='exact', validate_parameters=1, verbosity=None)"
      ]
     },
     "execution_count": 102,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import xgboost as xgb\n",
    "\n",
    "xg_reg = xgb.XGBRegressor()\n",
    "xg_reg.fit(X_train,y_train)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Model Validation\n",
    "Now we need to test to see how the models are preforming on data they have not seen. There are a number of metrics we can use, but to start we can use [R^2](https://en.wikipedia.org/wiki/Coefficient_of_determination). We use this metric to find the proportion of variance in our model which can show how well a model can predicts real values. The closer the R^2 value is to one, the more accurate our model is. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 103,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "R^2 Scores\n",
      "------------\n",
      "Random Forest: \t0.7169902009334733\n",
      "XGBoost: \t0.7151954636103839\n"
     ]
    }
   ],
   "source": [
    "from sklearn.metrics import r2_score, median_absolute_error\n",
    "\n",
    "predictions = rfr.predict(X_test)\n",
    "print(\"R^2 Scores\\n------------\\nRandom Forest: \\t\" + str(r2_score(y_test.values, predictions)))\n",
    "\n",
    "predictions = xg_reg.predict(X_test)\n",
    "print(\"XGBoost: \\t\" + str(r2_score(y_test.values, predictions)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Both scores are similer so we run another metric called the median absolute error. This metric is robust againts outliers which could help distingiush between the models. Smaller values equate to a better model preformance."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 104,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Median Absolute Error\n",
      "------------------------\n",
      "Random Forest: \t0.09300238718967224\n",
      "XGBoost: \t0.09418881670626102\n"
     ]
    }
   ],
   "source": [
    "predictions = rfr.predict(X_test)\n",
    "print(\"Median Absolute Error\\n------------------------\\nRandom Forest: \\t\" +\n",
    "      str(median_absolute_error(y_test.values, predictions)))\n",
    "\n",
    "predictions = xg_reg.predict(X_test)\n",
    "print(\"XGBoost: \\t\" + str(median_absolute_error(y_test.values, predictions)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Hyperparameter Tuning\n",
    "At the moment both models are preforming at an ok accuracy but it is not great. We can improve the models accuracy though a process called hyperparameter tuning. Here we will change the parameters for each model until a better accuracy is found. The method we will use to acomplish this is called random search. We will test a series of random parameters until we find the top preforming set.\n",
    "\n",
    "First we define a sampling distribution for each hyperparameter that we want to change.\n",
    "\n",
    "\n",
    "For the random forest model, we look at these main parameters: `n_estimators, max_depth, min_samples_split, min_samples_leaf, max_features`. Look at sklearns [documentation](https://scikit-learn.org/stable/modules/generated/sklearn.ensemble.RandomForestRegressor.html) on this model for what each parameter means.\n",
    "\n",
    "For the boosting model, we look at these main parameters: `n_estimators, max_depth, learning_rate, colsample_bynode, colsample_bytree, subsample`. Look at XGBoost's [documentation](https://xgboost.readthedocs.io/en/latest/python/python_api.html#module-xgboost.sklearn) on this model for what each parameter means."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 105,
   "metadata": {},
   "outputs": [],
   "source": [
    "from scipy.stats import randint as sp_randint\n",
    "from scipy.stats import uniform\n",
    "\n",
    "# both\n",
    "n_estimators = sp_randint(10, 1000)\n",
    "max_depth = sp_randint(1, 40)\n",
    "\n",
    "# rf\n",
    "min_samples_split = sp_randint(2, 5)\n",
    "min_samples_leaf = sp_randint(1, 5)\n",
    "max_features = sp_randint(1, 30)\n",
    "\n",
    "# boost\n",
    "learning_rate = uniform(loc=0, scale=0.05)\n",
    "colsample_bynode = uniform(loc=0, scale=1)\n",
    "colsample_bytree = uniform(loc=0, scale=1)\n",
    "subsample = uniform(loc=0, scale=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next we will create a function that will train each random model and output a Pandas DataFrame of its relevant values."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 106,
   "metadata": {},
   "outputs": [],
   "source": [
    "def record_model(model, model_type):\n",
    "    model.fit(X_train, y_train)\n",
    "    output = pd.DataFrame({\n",
    "        'model_type': model_type,\n",
    "        'model': [model],\n",
    "        'r2_score': r2_score(y_test.values, model.predict(X_test))\n",
    "    })\n",
    "    return output"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally we run though `itterations` amount of random tests."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 107,
   "metadata": {},
   "outputs": [],
   "source": [
    "itterations = 1\n",
    "\n",
    "models = pd.DataFrame(columns=['model_type', 'model', 'r2_score'])\n",
    "for i in range(itterations):\n",
    "    model = RandomForestRegressor(n_estimators=n_estimators.rvs(),\n",
    "                               max_depth=max_depth.rvs(),\n",
    "                               min_samples_split=min_samples_split.rvs(),\n",
    "                               min_samples_leaf=min_samples_leaf.rvs())\n",
    "    models = models.append(record_model(model, 'random_forest'))\n",
    "    \n",
    "    model = xgb.XGBRegressor(n_estimators=n_estimators.rvs(),\n",
    "                             max_depth=max_depth.rvs(),\n",
    "                            learning_rate=learning_rate.rvs(),\n",
    "                            colsample_bynode=colsample_bynode.rvs(),\n",
    "                            colsample_bytree=colsample_bytree.rvs(),\n",
    "                            subsample=subsample.rvs())\n",
    "    models = models.append(record_model(model, 'boost'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can repeat this process with smaller random ranges to further tune the model, but for the sake of this demo we will take the top models with the highest R^2 value for each model type.\n",
    "\n",
    "We plot the scores to see which sets of parameters preformed the best."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 108,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "\n",
       "<div id=\"altair-viz-bd7911a8162c463eba39f8eea257128f\"></div>\n",
       "<script type=\"text/javascript\">\n",
       "  (function(spec, embedOpt){\n",
       "    let outputDiv = document.currentScript.previousElementSibling;\n",
       "    if (outputDiv.id !== \"altair-viz-bd7911a8162c463eba39f8eea257128f\") {\n",
       "      outputDiv = document.getElementById(\"altair-viz-bd7911a8162c463eba39f8eea257128f\");\n",
       "    }\n",
       "    const paths = {\n",
       "      \"vega\": \"https://cdn.jsdelivr.net/npm//vega@5?noext\",\n",
       "      \"vega-lib\": \"https://cdn.jsdelivr.net/npm//vega-lib?noext\",\n",
       "      \"vega-lite\": \"https://cdn.jsdelivr.net/npm//vega-lite@4.8.1?noext\",\n",
       "      \"vega-embed\": \"https://cdn.jsdelivr.net/npm//vega-embed@6?noext\",\n",
       "    };\n",
       "\n",
       "    function loadScript(lib) {\n",
       "      return new Promise(function(resolve, reject) {\n",
       "        var s = document.createElement('script');\n",
       "        s.src = paths[lib];\n",
       "        s.async = true;\n",
       "        s.onload = () => resolve(paths[lib]);\n",
       "        s.onerror = () => reject(`Error loading script: ${paths[lib]}`);\n",
       "        document.getElementsByTagName(\"head\")[0].appendChild(s);\n",
       "      });\n",
       "    }\n",
       "\n",
       "    function showError(err) {\n",
       "      outputDiv.innerHTML = `<div class=\"error\" style=\"color:red;\">${err}</div>`;\n",
       "      throw err;\n",
       "    }\n",
       "\n",
       "    function displayChart(vegaEmbed) {\n",
       "      vegaEmbed(outputDiv, spec, embedOpt)\n",
       "        .catch(err => showError(`Javascript Error: ${err.message}<br>This usually means there's a typo in your chart specification. See the javascript console for the full traceback.`));\n",
       "    }\n",
       "\n",
       "    if(typeof define === \"function\" && define.amd) {\n",
       "      requirejs.config({paths});\n",
       "      require([\"vega-embed\"], displayChart, err => showError(`Error loading script: ${err.message}`));\n",
       "    } else if (typeof vegaEmbed === \"function\") {\n",
       "      displayChart(vegaEmbed);\n",
       "    } else {\n",
       "      loadScript(\"vega\")\n",
       "        .then(() => loadScript(\"vega-lite\"))\n",
       "        .then(() => loadScript(\"vega-embed\"))\n",
       "        .catch(showError)\n",
       "        .then(() => displayChart(vegaEmbed));\n",
       "    }\n",
       "  })({\"config\": {\"view\": {\"continuousWidth\": 400, \"continuousHeight\": 300}}, \"layer\": [{\"mark\": {\"type\": \"bar\", \"size\": 15}, \"encoding\": {\"color\": {\"type\": \"nominal\", \"field\": \"model_type\"}, \"x\": {\"type\": \"quantitative\", \"field\": \"index\", \"title\": \"Model\"}, \"y\": {\"type\": \"quantitative\", \"field\": \"r2_score\", \"scale\": {\"zero\": false}}}}, {\"mark\": {\"type\": \"text\", \"baseline\": \"top\", \"dy\": -15}, \"encoding\": {\"color\": {\"type\": \"nominal\", \"field\": \"model_type\"}, \"text\": {\"type\": \"quantitative\", \"field\": \"r2_score\", \"format\": \",.3r\"}, \"x\": {\"type\": \"quantitative\", \"field\": \"index\", \"title\": \"Model\"}, \"y\": {\"type\": \"quantitative\", \"field\": \"r2_score\", \"scale\": {\"zero\": false}}}}], \"data\": {\"name\": \"data-ebfb35a617582714f103857aebbbac0a\"}, \"width\": 700, \"$schema\": \"https://vega.github.io/schema/vega-lite/v4.8.1.json\", \"datasets\": {\"data-ebfb35a617582714f103857aebbbac0a\": [{\"model_type\": \"random_forest\", \"r2_score\": 0.722325999992195, \"index\": 0}, {\"model_type\": \"boost\", \"r2_score\": 0.6907924005233776, \"index\": 1}, {\"model_type\": \"random_forest\", \"r2_score\": 0.7162826354479085, \"index\": 2}, {\"model_type\": \"boost\", \"r2_score\": 0.5974194208097385, \"index\": 3}, {\"model_type\": \"random_forest\", \"r2_score\": 0.7064409918660777, \"index\": 4}, {\"model_type\": \"boost\", \"r2_score\": 0.6944168544342824, \"index\": 5}, {\"model_type\": \"random_forest\", \"r2_score\": 0.6782089227595868, \"index\": 6}, {\"model_type\": \"boost\", \"r2_score\": 0.7293830024216887, \"index\": 7}, {\"model_type\": \"random_forest\", \"r2_score\": 0.7221368292921438, \"index\": 8}, {\"model_type\": \"boost\", \"r2_score\": 0.46263451454019944, \"index\": 9}]}}, {\"mode\": \"vega-lite\"});\n",
       "</script>"
      ],
      "text/plain": [
       "alt.LayerChart(...)"
      ]
     },
     "execution_count": 108,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import altair as alt\n",
    "\n",
    "models_2 = models.reset_index(drop=True)\n",
    "models_2 = models_2.reset_index()\n",
    "\n",
    "base = alt.Chart(models_2[['model_type', 'r2_score', 'index']])\n",
    "\n",
    "bars = base.mark_bar(size=15).encode(\n",
    "x=alt.X('index:Q', title=\"Model\"),\n",
    "y= alt.Y('r2_score:Q', scale=alt.Scale(zero=False)),\n",
    "color='model_type:N')\n",
    "\n",
    "text = bars.mark_text(\n",
    "    baseline='top',\n",
    "    dy=-15\n",
    ").encode(\n",
    "    text=alt.Text('r2_score:Q', format=',.3r')\n",
    ")\n",
    "\n",
    "(bars + text).properties(width=700)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here we are extracting the best model for both model types and displaying the parameters that they used."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 109,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>n_estimators</th>\n",
       "      <th>max_depth</th>\n",
       "      <th>min_samples_split</th>\n",
       "      <th>min_samples_leaf</th>\n",
       "      <th>max_features</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>962</td>\n",
       "      <td>35</td>\n",
       "      <td>4</td>\n",
       "      <td>3</td>\n",
       "      <td>auto</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   n_estimators  max_depth  min_samples_split  min_samples_leaf max_features\n",
       "0           962         35                  4                 3         auto"
      ]
     },
     "execution_count": 109,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# top random forest\n",
    "rfr = models[models['model_type']=='random_forest']\n",
    "rfr = rfr.iloc[rfr['r2_score'].values.argmax()].values[1]\n",
    "pd.DataFrame(rfr.get_params(), index=[0])[['n_estimators', 'max_depth',\n",
    "                                           'min_samples_split', 'min_samples_leaf',\n",
    "                                           'max_features']]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 110,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>n_estimators</th>\n",
       "      <th>max_depth</th>\n",
       "      <th>learning_rate</th>\n",
       "      <th>colsample_bynode</th>\n",
       "      <th>colsample_bytree</th>\n",
       "      <th>subsample</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>948</td>\n",
       "      <td>12</td>\n",
       "      <td>0.006897</td>\n",
       "      <td>0.78998</td>\n",
       "      <td>0.471958</td>\n",
       "      <td>0.359216</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   n_estimators  max_depth  learning_rate  colsample_bynode  colsample_bytree  \\\n",
       "0           948         12       0.006897           0.78998          0.471958   \n",
       "\n",
       "   subsample  \n",
       "0   0.359216  "
      ]
     },
     "execution_count": 110,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# top boost\n",
    "xg_reg = models[models['model_type']=='boost']\n",
    "xg_reg = xg_reg.iloc[xg_reg['r2_score'].values.argmax()].values[1]\n",
    "\n",
    "# save this model for later export\n",
    "from copy import copy\n",
    "saved_model = copy(xg_reg)\n",
    "\n",
    "pd.DataFrame(xg_reg.get_params(), index=[0])[['n_estimators', 'max_depth',\n",
    "                                              'learning_rate', 'colsample_bynode',\n",
    "                                              'colsample_bytree', 'subsample']]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Feature Importances\n",
    "With hypertuning we were able to get about a 3 to 7 percent increase in accuracy. The model is still not preforming optimally, but we were able to rule out that the model was not the only issue. The issue might be attributed to the features we used to train. Initially we threw all the features at the model, but that may not be the best option. Luckily, an interesting option for tree models is that they provide feature importances. This allows us to see which features are influencing the model the most and which are doing nothing.\n",
    "\n",
    "We can grab the top five features from our random forest model and sort them from most influential to least influential."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 111,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[(7, 0.5560598822450813),\n",
       " (3, 0.08794821060257854),\n",
       " (11, 0.07518801192606869),\n",
       " (2, 0.053176975224588016),\n",
       " (4, 0.04076222646361437)]"
      ]
     },
     "execution_count": 111,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "l = list(enumerate(rfr.feature_importances_))\n",
    "l.sort(key=lambda x: -x[1])\n",
    "l[:5 ]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To visulaize this better, we use altair to graph the importances with their respective column names. We find that they are similer in distribution and both rely on cloud data for the majority of the model."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 112,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "\n",
       "<div id=\"altair-viz-61babdc4fd1343f585f7d4031d414348\"></div>\n",
       "<script type=\"text/javascript\">\n",
       "  (function(spec, embedOpt){\n",
       "    let outputDiv = document.currentScript.previousElementSibling;\n",
       "    if (outputDiv.id !== \"altair-viz-61babdc4fd1343f585f7d4031d414348\") {\n",
       "      outputDiv = document.getElementById(\"altair-viz-61babdc4fd1343f585f7d4031d414348\");\n",
       "    }\n",
       "    const paths = {\n",
       "      \"vega\": \"https://cdn.jsdelivr.net/npm//vega@5?noext\",\n",
       "      \"vega-lib\": \"https://cdn.jsdelivr.net/npm//vega-lib?noext\",\n",
       "      \"vega-lite\": \"https://cdn.jsdelivr.net/npm//vega-lite@4.8.1?noext\",\n",
       "      \"vega-embed\": \"https://cdn.jsdelivr.net/npm//vega-embed@6?noext\",\n",
       "    };\n",
       "\n",
       "    function loadScript(lib) {\n",
       "      return new Promise(function(resolve, reject) {\n",
       "        var s = document.createElement('script');\n",
       "        s.src = paths[lib];\n",
       "        s.async = true;\n",
       "        s.onload = () => resolve(paths[lib]);\n",
       "        s.onerror = () => reject(`Error loading script: ${paths[lib]}`);\n",
       "        document.getElementsByTagName(\"head\")[0].appendChild(s);\n",
       "      });\n",
       "    }\n",
       "\n",
       "    function showError(err) {\n",
       "      outputDiv.innerHTML = `<div class=\"error\" style=\"color:red;\">${err}</div>`;\n",
       "      throw err;\n",
       "    }\n",
       "\n",
       "    function displayChart(vegaEmbed) {\n",
       "      vegaEmbed(outputDiv, spec, embedOpt)\n",
       "        .catch(err => showError(`Javascript Error: ${err.message}<br>This usually means there's a typo in your chart specification. See the javascript console for the full traceback.`));\n",
       "    }\n",
       "\n",
       "    if(typeof define === \"function\" && define.amd) {\n",
       "      requirejs.config({paths});\n",
       "      require([\"vega-embed\"], displayChart, err => showError(`Error loading script: ${err.message}`));\n",
       "    } else if (typeof vegaEmbed === \"function\") {\n",
       "      displayChart(vegaEmbed);\n",
       "    } else {\n",
       "      loadScript(\"vega\")\n",
       "        .then(() => loadScript(\"vega-lite\"))\n",
       "        .then(() => loadScript(\"vega-embed\"))\n",
       "        .catch(showError)\n",
       "        .then(() => displayChart(vegaEmbed));\n",
       "    }\n",
       "  })({\"config\": {\"view\": {\"continuousWidth\": 400, \"continuousHeight\": 300}}, \"hconcat\": [{\"layer\": [{\"mark\": \"bar\", \"encoding\": {\"x\": {\"type\": \"quantitative\", \"field\": \"Importance\"}, \"y\": {\"type\": \"nominal\", \"field\": \"Feature\"}}, \"height\": 500, \"title\": \"Random Forest\", \"width\": 275}, {\"mark\": {\"type\": \"text\", \"align\": \"left\", \"baseline\": \"middle\", \"dx\": 3}, \"encoding\": {\"text\": {\"type\": \"quantitative\", \"field\": \"Importance\", \"format\": \",.2r\"}, \"x\": {\"type\": \"quantitative\", \"field\": \"Importance\"}, \"y\": {\"type\": \"nominal\", \"field\": \"Feature\"}}, \"height\": 500, \"title\": \"Random Forest\", \"width\": 275}], \"data\": {\"name\": \"data-8f7743bea9d32fe3e88e198a453dd700\"}}, {\"layer\": [{\"mark\": \"bar\", \"encoding\": {\"x\": {\"type\": \"quantitative\", \"field\": \"Importance\"}, \"y\": {\"type\": \"nominal\", \"field\": \"Feature\"}}, \"height\": 500, \"title\": \"XGBoost\", \"width\": 275}, {\"mark\": {\"type\": \"text\", \"align\": \"left\", \"baseline\": \"middle\", \"dx\": 3}, \"encoding\": {\"text\": {\"type\": \"quantitative\", \"field\": \"Importance\", \"format\": \",.2r\"}, \"x\": {\"type\": \"quantitative\", \"field\": \"Importance\"}, \"y\": {\"type\": \"nominal\", \"field\": \"Feature\"}}, \"height\": 500, \"title\": \"XGBoost\", \"width\": 275}], \"data\": {\"name\": \"data-6faa816ecd374c1d70e8125fdac14ac7\"}}], \"$schema\": \"https://vega.github.io/schema/vega-lite/v4.8.1.json\", \"datasets\": {\"data-8f7743bea9d32fe3e88e198a453dd700\": [{\"Feature\": \"temperature\", \"Importance\": 0.03562577280217239}, {\"Feature\": \"dew_point\", \"Importance\": 0.027152109307475164}, {\"Feature\": \"relative_humidity\", \"Importance\": 0.053176975224588016}, {\"Feature\": \"daily_precipitation\", \"Importance\": 0.08794821060257854}, {\"Feature\": \"station_pressure\", \"Importance\": 0.04076222646361437}, {\"Feature\": \"wind_speed\", \"Importance\": 0.03412401013482356}, {\"Feature\": \"hourly_visibility\", \"Importance\": 0.016641235395384806}, {\"Feature\": \"cloud_cover\", \"Importance\": 0.5560598822450813}, {\"Feature\": \"mostly_cloudy\", \"Importance\": 0.028800748213060876}, {\"Feature\": \"mostly_clear\", \"Importance\": 0.008065274596744334}, {\"Feature\": \"clear\", \"Importance\": 0.009030515708763495}, {\"Feature\": \"cloudy\", \"Importance\": 0.07518801192606869}, {\"Feature\": \"partly_cloudy\", \"Importance\": 0.0065103826811178786}, {\"Feature\": \"overcast\", \"Importance\": 0.0002649421605649014}, {\"Feature\": \"rain_light\", \"Importance\": 0.0076509382303586}, {\"Feature\": \"tstorm\", \"Importance\": 0.0069992045381454095}, {\"Feature\": \"drizzle\", \"Importance\": 0.002300154969895826}, {\"Feature\": \"rain_heavy\", \"Importance\": 0.0001050150131641543}, {\"Feature\": \"rain\", \"Importance\": 0.002430491168037719}, {\"Feature\": \"fog\", \"Importance\": 0.0010771293092536394}, {\"Feature\": \"snow_light\", \"Importance\": 6.134964202673294e-05}, {\"Feature\": \"snow\", \"Importance\": 9.74744293850254e-06}, {\"Feature\": \"snow_heavy\", \"Importance\": 6.686257895173773e-07}, {\"Feature\": \"freezing_rain\", \"Importance\": 1.4831345071227339e-05}, {\"Feature\": \"freezing_drizzle\", \"Importance\": 1.7225328037915427e-07}, {\"Feature\": \"ice_pellets\", \"Importance\": 0.0}, {\"Feature\": \"ice_pellets_light\", \"Importance\": 0.0}, {\"Feature\": \"ice_pellets_heavy\", \"Importance\": 0.0}, {\"Feature\": \"flurries\", \"Importance\": 0.0}, {\"Feature\": \"freezing_rain_heavy\", \"Importance\": 0.0}, {\"Feature\": \"freezing_rain_light\", \"Importance\": 0.0}, {\"Feature\": \"fog_light\", \"Importance\": 0.0}], \"data-6faa816ecd374c1d70e8125fdac14ac7\": [{\"Feature\": \"temperature\", \"Importance\": 0.015144347213208675}, {\"Feature\": \"dew_point\", \"Importance\": 0.015291542746126652}, {\"Feature\": \"relative_humidity\", \"Importance\": 0.029787147417664528}, {\"Feature\": \"daily_precipitation\", \"Importance\": 0.06699273735284805}, {\"Feature\": \"station_pressure\", \"Importance\": 0.01836101897060871}, {\"Feature\": \"wind_speed\", \"Importance\": 0.01620914787054062}, {\"Feature\": \"hourly_visibility\", \"Importance\": 0.017597422003746033}, {\"Feature\": \"cloud_cover\", \"Importance\": 0.13625745475292206}, {\"Feature\": \"mostly_cloudy\", \"Importance\": 0.03385981172323227}, {\"Feature\": \"mostly_clear\", \"Importance\": 0.030731279402971268}, {\"Feature\": \"clear\", \"Importance\": 0.05857749283313751}, {\"Feature\": \"cloudy\", \"Importance\": 0.22341768443584442}, {\"Feature\": \"partly_cloudy\", \"Importance\": 0.02034672163426876}, {\"Feature\": \"overcast\", \"Importance\": 0.013990417122840881}, {\"Feature\": \"rain_light\", \"Importance\": 0.09073188900947571}, {\"Feature\": \"tstorm\", \"Importance\": 0.04207774996757507}, {\"Feature\": \"drizzle\", \"Importance\": 0.017578262835741043}, {\"Feature\": \"rain_heavy\", \"Importance\": 0.014722248539328575}, {\"Feature\": \"rain\", \"Importance\": 0.037273578345775604}, {\"Feature\": \"fog\", \"Importance\": 0.015842387452721596}, {\"Feature\": \"snow_light\", \"Importance\": 0.023061620071530342}, {\"Feature\": \"snow\", \"Importance\": 0.017901046201586723}, {\"Feature\": \"snow_heavy\", \"Importance\": 0.005807942245155573}, {\"Feature\": \"freezing_rain\", \"Importance\": 0.017675580456852913}, {\"Feature\": \"freezing_drizzle\", \"Importance\": 0.006903701461851597}, {\"Feature\": \"ice_pellets\", \"Importance\": 0.0}, {\"Feature\": \"ice_pellets_light\", \"Importance\": 0.004317869897931814}, {\"Feature\": \"ice_pellets_heavy\", \"Importance\": 0.009541946463286877}, {\"Feature\": \"flurries\", \"Importance\": 0.0}, {\"Feature\": \"freezing_rain_heavy\", \"Importance\": 0.0}, {\"Feature\": \"freezing_rain_light\", \"Importance\": 0.0}, {\"Feature\": \"fog_light\", \"Importance\": 0.0}]}}, {\"mode\": \"vega-lite\"});\n",
       "</script>"
      ],
      "text/plain": [
       "alt.HConcatChart(...)"
      ]
     },
     "execution_count": 112,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def plot_importance(tree_model, drop=[], title=\"Feature Importance\"):\n",
    "    drop.extend(['dni_efficiency', 'ghi_efficiency', 'dhi_efficiency',\n",
    "                 'STATION','DATE', 'latitude','longitude'])\n",
    "\n",
    "    source = pd.DataFrame({\n",
    "        'Feature': list(train.drop(columns=drop)),\n",
    "        'Importance': list(tree_model.feature_importances_)\n",
    "    })\n",
    "\n",
    "    bars = alt.Chart(source, height=500, width = 275, title=title).mark_bar().encode(\n",
    "        y='Feature:N',\n",
    "        x=\"Importance:Q\"\n",
    "    )\n",
    "\n",
    "    text = bars.mark_text(\n",
    "        align='left',\n",
    "        baseline='middle',\n",
    "        dx=3\n",
    "    ).encode(\n",
    "        text=alt.Text('Importance:Q', format=',.2r')\n",
    "    )\n",
    "\n",
    "    return (bars + text)\n",
    "    \n",
    "alt.hconcat(plot_importance(rfr, title=\"Random Forest\"), plot_importance(xg_reg, title=\"XGBoost\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next we will dive into changing around the features to see what changes effect the final accuracy. First we will remove one of the top preforming features, cloud_cover, and see what results it gives.\n",
    "\n",
    "Here we created a function that takes in the columns that the model will not use, train the model, and return a R^2 score and graph. We will use the parameters we found through hypertuning."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 113,
   "metadata": {},
   "outputs": [],
   "source": [
    "def train_plot_model(model, title, cols):    \n",
    "    #get col list\n",
    "    dropped_train = train.drop(columns=['dni_efficiency', 'ghi_efficiency', 'dhi_efficiency', 'STATION','DATE',\n",
    "                                            'latitude','longitude'])\n",
    "    #drop = [dropped_train.columns.get_loc(col) for col in cols]\n",
    "    \n",
    "    # remove cols from transformed data\n",
    "    #svecs = np.delete(X_train, drop, axis=1)\n",
    "    svecs = X_train.drop(columns=cols)\n",
    "\n",
    "    model.fit(svecs, y_train)\n",
    "    \n",
    "    #np.delete(X_test, drop,axis=1)\n",
    "    \n",
    "    predictions = model.predict(X_test.drop(columns=cols))\n",
    "    print(title + \": \\tR^2 Score: \" + str(r2_score(y_test.values, predictions)))\n",
    "    \n",
    "    return plot_importance(model, cols.copy(), title=title)\n",
    "    \n",
    "dropped_columns = ['cloud_cover']\n",
    "#alt.hconcat(train_plot_model(rfr, 'Random Forest', dropped_columns),\n",
    " #          train_plot_model(xg_reg, 'XGBoost', dropped_columns))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With cloud_cover removed, we are able to see that not much changed in the model. Both models are preforming at the same accuracy and are roughly the same as the preformance prior with cloud_cover included. From this we can conclude that cloud_cover was highly correlated with a combintion of the other columns.\n",
    "\n",
    "Let's do the opposite of what we tested and remove all the cloud observation features."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 114,
   "metadata": {},
   "outputs": [],
   "source": [
    "dropped_columns = ['overcast', 'cloudy', 'mostly_cloudy',\n",
    "                   'partly_cloudy', 'mostly_clear', 'clear']\n",
    "#alt.hconcat(train_plot_model(rfr, 'Random Forest', dropped_columns),\n",
    "#          train_plot_model(xg_reg, 'XGBoost', dropped_columns))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The model with the cloud_cover feature preformed slightly worse than the cloud observations feature model, so it is safe to say that we can drop cloud_cover as a feature in favor of a better preforming set of features. On the other hand, having a less complex (less features) model could be favorable in some situations as long as the accuracy is not reduced too much. In our case we will opt for the more accurate set of features.\n",
    "\n",
    "Finally we can try trimming the feature set by removing features that have less than a .02% impact. This will give us the same advantage of having a less complex model but hopefully with less of a drop in accuracy."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 115,
   "metadata": {},
   "outputs": [],
   "source": [
    "dropped_columns = ['station_pressure',\n",
    "       'wind_speed', 'cloud_cover', 'overcast', 'fog', 'snow',\n",
    "       'snow_heavy', 'freezing_rain', 'freezing_drizzle', 'ice_pellets',\n",
    "       'ice_pellets_light', 'ice_pellets_heavy', 'flurries',\n",
    "       'freezing_rain_heavy', 'freezing_rain_light', 'fog_light']\n",
    "#alt.hconcat(train_plot_model(rfr, 'Random Forest', dropped_columns),\n",
    "#          train_plot_model(xg_reg, 'XGBoost', dropped_columns))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this test, there was an accuracy drop but that is expected when a good number of features are lost. We could continue to fine tune the features we use along with more hyperparameter tuning as well. But from the tests we have done, it looks like there is not much we can do to substantially increase the model's accuracy.\n",
    "\n",
    "From here we save the original boosting model so that we can use it outside of this notebook. Becouse of the random hypertuning we did, each iteration of this notebook will return a different model that will preform slightly differently. In the next notebook we will look at where the error in this model is coming from, and how to use that information to increase our model's accuracy."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 116,
   "metadata": {},
   "outputs": [],
   "source": [
    "import cloudpickle as cp\n",
    "import os\n",
    "\n",
    "def serialize_to(obj, default_filename):\n",
    "    filename = os.getenv(\"S2I_PIPELINE_STAGE_SAVE_FILE\", default_filename)\n",
    "    if filename != default_filename:\n",
    "        cp.dump(obj, open(default_filename, \"wb\"))\n",
    "    cp.dump(obj, open(filename, \"wb\"))\n",
    "\n",
    "serialize_to(saved_model, \"trained_model.sav\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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
}