notebooks/model_selection_techniques.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "1wMwXl-eL6Sx"
      },
      "source": [
        "# Testing on the training data\n",
        "\n",
        "**DO NOT DO IT**, since it is methodologically wrong!"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 1,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "Rg4PPnvSLe5Q",
        "outputId": "768a324a-8e0d-493b-9138-62d0eec1b0f0"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Accuracy: 0.97\n"
          ]
        }
      ],
      "source": [
        "from sklearn.datasets import load_iris\n",
        "from sklearn.ensemble import RandomForestClassifier\n",
        "from sklearn.metrics import accuracy_score\n",
        "\n",
        "iris = load_iris()\n",
        "\n",
        "X = iris.data\n",
        "y = iris.target\n",
        "\n",
        "clf = RandomForestClassifier(n_estimators=2, random_state=0)\n",
        "# X is our training data\n",
        "clf.fit(X, y)\n",
        "\n",
        "# This is an overly optimistic estimation since we are using X again!\n",
        "y_pred = clf.predict(X)\n",
        "acc = accuracy_score(y, y_pred)\n",
        "\n",
        "print(f'Accuracy: {acc:.2f}')"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "AGluy49qLjeJ"
      },
      "source": [
        "## Two-way holdout"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 2,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "aj3TlkZWLkj9",
        "outputId": "ce08d511-e61f-41a4-d2db-623a64b51e07"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Accuracy: 0.91\n"
          ]
        }
      ],
      "source": [
        "from sklearn.model_selection import train_test_split\n",
        "\n",
        "# split in train and test\n",
        "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=0)\n",
        "\n",
        "clf = RandomForestClassifier(n_estimators=2, random_state=0)\n",
        "clf.fit(X_train, y_train)\n",
        "\n",
        "# test with unseen data\n",
        "y_pred = clf.predict(X_test)\n",
        "acc = accuracy_score(y_test, y_pred)\n",
        "\n",
        "print(f'Accuracy: {acc:.2f}')"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "z4hoss4eL4IF"
      },
      "source": [
        "# k-fold cross validation"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 3,
      "metadata": {
        "id": "-_Yf27B-L4pL"
      },
      "outputs": [],
      "source": [
        "from sklearn.model_selection import cross_validate\n",
        "import timeit\n",
        "\n",
        "def do_cross_validation(clf, print_model=False, print_duration=False):\n",
        "    start = timeit.default_timer()\n",
        "    cv = cross_validate(clf, X, y, scoring='accuracy', cv=3)\n",
        "    scores = ' + '.join(f'{s:.2f}' for s in cv[\"test_score\"])\n",
        "    mean_ = cv[\"test_score\"].mean()\n",
        "    msg = f'Cross-validated accuracy: ({scores}) / 3 = {mean_:.2f}'\n",
        "\n",
        "    if print_model:\n",
        "        msg = f'\\nClassifier: {clf}\\n{msg}\\n'\n",
        "\n",
        "    if print_duration:\n",
        "        msg = f'Duration: {timeit.default_timer() - start}{msg}\\n'\n",
        "\n",
        "    print(msg)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 4,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "XfrxHB93Y_PP",
        "outputId": "aadc11b0-d8d8-4ddb-c8a0-374a0ec6b4a0"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Duration: 0.04880734300240874\n",
            "Classifier: RandomForestClassifier(n_estimators=2, random_state=0)\n",
            "Cross-validated accuracy: (0.98 + 0.92 + 0.96) / 3 = 0.95\n",
            "\n",
            "\n"
          ]
        }
      ],
      "source": [
        "clf = RandomForestClassifier(n_estimators=2, random_state=0)\n",
        "do_cross_validation(clf, True, True)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "U_qkItbFOEPP"
      },
      "source": [
        "## Applying cross-validation for model selection"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 5,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "AfExmrWAOHMc",
        "outputId": "c9133b2a-429e-430e-88f5-6b4937fc984d"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Default value for kernel:  rbf\n",
            "Duration: 0.036153421009657905\n",
            "Classifier: SVC(random_state=0)\n",
            "Cross-validated accuracy: (0.96 + 0.98 + 0.94) / 3 = 0.96\n",
            "\n",
            "\n"
          ]
        }
      ],
      "source": [
        "from sklearn.svm import SVC\n",
        "\n",
        "start = timeit.default_timer()\n",
        "svc = SVC(random_state=0)\n",
        "print('Default value for kernel: ', svc.kernel)\n",
        "do_cross_validation(svc, True, True)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "gz7sqaKQOYFo"
      },
      "source": []
    },
    {
      "cell_type": "code",
      "execution_count": 6,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "CoxscNzYObI6",
        "outputId": "fbffe700-569b-4886-90ae-6325c203a863"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "\n",
            "Classifier: SVC(kernel='linear', random_state=0)\n",
            "Cross-validated accuracy: (1.00 + 1.00 + 0.98) / 3 = 0.99\n",
            "\n",
            "\n",
            "Classifier: SVC(kernel='poly', random_state=0)\n",
            "Cross-validated accuracy: (0.98 + 0.94 + 0.98) / 3 = 0.97\n",
            "\n",
            "\n",
            "Classifier: RandomForestClassifier(n_estimators=2, random_state=0)\n",
            "Cross-validated accuracy: (0.98 + 0.92 + 0.96) / 3 = 0.95\n",
            "\n",
            "\n",
            "Classifier: RandomForestClassifier(n_estimators=5, random_state=0)\n",
            "Cross-validated accuracy: (0.98 + 0.94 + 0.94) / 3 = 0.95\n",
            "\n"
          ]
        }
      ],
      "source": [
        "do_cross_validation(SVC(kernel='linear', random_state=0), print_model=True)\n",
        "do_cross_validation(SVC(kernel='poly', random_state=0), print_model=True)\n",
        "do_cross_validation(RandomForestClassifier(n_estimators=2, random_state=0), print_model=True)\n",
        "do_cross_validation(RandomForestClassifier(n_estimators=5, random_state=0), print_model=True)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "YasAqT0fPeKe"
      },
      "source": [
        "# Nested cross-validation"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 9,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "2qggaX-dPhsQ",
        "outputId": "913308ae-92da-43d9-d69c-13dc928a99f5"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Duration: 0.5532464570133016\n",
            "Classifier: GridSearchCV(estimator=RandomForestClassifier(random_state=0),\n",
            "             param_grid={'n_estimators': [2, 5]})\n",
            "Cross-validated accuracy: (0.98 + 0.92 + 0.96) / 3 = 0.95\n",
            "\n",
            "\n",
            "Duration: 0.14918377198046073\n",
            "Classifier: GridSearchCV(estimator=SVC(random_state=0),\n",
            "             param_grid={'kernel': ['linear', 'poly']})\n",
            "Cross-validated accuracy: (1.00 + 0.94 + 0.98) / 3 = 0.97\n",
            "\n",
            "\n"
          ]
        }
      ],
      "source": [
        "from sklearn.model_selection import GridSearchCV\n",
        "\n",
        "start = timeit.default_timer()\n",
        "# random forest inner loop\n",
        "clf_grid = GridSearchCV(RandomForestClassifier(random_state=0), param_grid={'n_estimators': [2, 5]})\n",
        "# random forest outer loop\n",
        "do_cross_validation(clf_grid, print_model=True, print_duration=True)\n",
        "\n",
        "start = timeit.default_timer()\n",
        "# svc inner loop\n",
        "svc_grid = GridSearchCV(SVC(random_state=0), param_grid={'kernel': ['linear', 'poly']})\n",
        "# svc outer loop\n",
        "do_cross_validation(svc_grid, print_model=True, print_duration=True)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "# Nested CV - getting the final model\n",
        "\n",
        "\n",
        "Nested cross-validation itself doesn't directly produce a final model. Rather, it is a technique to get an unbiased estimated of the generalization error. \n",
        "\n",
        "There are three alternative approaches to produce the final model **after** using nested CV. \n",
        "1. The final model is produced by training on the entire dataset, and using the best hyperparameters found during the inner loop.\n",
        "2. The final model is produced using the algorithm selected in the inner loop, but performing an additional hyperparameter setting on the whole dataset.\n",
        "3. (Ensemble Model) The final model is built as an ensemble model by combining predictions from the multiple models trained in the inner loop.\n",
        "\n",
        "Approaches 1 and 2 are the most common ones. Both involve using the entire dataset to refit a model AFTER the generalization error has been estimated.\n",
        "\n",
        "Notice that in all of the three approaches described above, the estimate of the generalization error to be reported is the one resulting from the nested CV procedure. \n",
        "\n",
        "The two code blocks below provide examples of using the second approach."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Classification example"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 15,
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Model: svc\n",
            "Accuracy in the 3 outer folds: [0.94011976 0.96696697 0.93393393].\n",
            "Average acc: 0.9470068871266476\n",
            "\n",
            "Model: rfc\n",
            "Accuracy in the 3 outer folds: [0.95808383 0.97597598 0.93993994].\n",
            "Average acc: 0.9579999160837485\n",
            "\n",
            "Average score across the outer folds:  {'svc': 0.9470068871266476, 'rfc': 0.9579999160837485}\n",
            "\n",
            "****************************************************************************************************\n",
            "Now we choose the best model and refit on the whole dataset\n",
            "****************************************************************************************************\n",
            "\n",
            "Best model: \n",
            "\tRandomForestClassifier()\n",
            "\n",
            "Estimation of its generalization error (accuracy):\n",
            "\t0.9579999160837485\n",
            "\n",
            "Best parameter choice for this model: \n",
            "\t{'max_depth': 50}\n",
            "(according to cross-validation `KFold(n_splits=3, random_state=None, shuffle=False)` on the whole dataset).\n"
          ]
        }
      ],
      "source": [
        "from sklearn.model_selection import KFold, cross_val_score, GridSearchCV\n",
        "from sklearn.datasets import make_classification\n",
        "from sklearn.ensemble import RandomForestClassifier\n",
        "from sklearn.svm import SVC\n",
        "import numpy as np\n",
        "\n",
        "# `outer_cv` creates 3 folds for estimating generalization error\n",
        "outer_cv = KFold(3)\n",
        "\n",
        "# when we train on a certain fold, we use a second cross-validation\n",
        "# split in order to choose hyperparameters\n",
        "inner_cv = KFold(3)\n",
        "\n",
        "# create some regression data\n",
        "X, y = make_classification(n_samples=1000, n_features=10)\n",
        "\n",
        "# give shorthand names to models and use those as dictionary keys mapping\n",
        "# to models and parameter grids for that model\n",
        "models_and_parameters = {\n",
        "    'svc': (SVC(),\n",
        "            {'C': [0.01, 0.05, 0.1, 1]}),\n",
        "    'rfc': (RandomForestClassifier(),\n",
        "           {'max_depth': [5, 10, 50, 100, 200, 500]})}\n",
        "\n",
        "# we will collect the average of the scores on the 3 outer folds in this dictionary\n",
        "# with keys given by the names of the models in `models_and_parameters`\n",
        "average_scores_across_outer_folds_for_each_model = dict()\n",
        "\n",
        "# find the model with the best generalization error\n",
        "for name, (model, params) in models_and_parameters.items():\n",
        "    # this object is a regressor that also happens to choose\n",
        "    # its hyperparameters automatically using `inner_cv`\n",
        "    classifier_that_optimizes_its_hyperparams = GridSearchCV(\n",
        "        estimator=model, param_grid=params,\n",
        "        cv=inner_cv, scoring='accuracy')\n",
        "\n",
        "    # estimate generalization error on the 3-fold splits of the data\n",
        "    scores_across_outer_folds = cross_val_score(\n",
        "        classifier_that_optimizes_its_hyperparams,\n",
        "        X, y, cv=outer_cv, scoring='accuracy')\n",
        "\n",
        "    # get the mean accuracy across each of outer_cv's 3 folds\n",
        "    average_scores_across_outer_folds_for_each_model[name] = np.mean(scores_across_outer_folds)\n",
        "    error_summary = 'Model: {name}\\nAccuracy in the 3 outer folds: {scores}.\\nAverage acc: {avg}'\n",
        "    print(error_summary.format(\n",
        "        name=name, scores=scores_across_outer_folds,\n",
        "        avg=np.mean(scores_across_outer_folds)))\n",
        "    print()\n",
        "\n",
        "print('Average score across the outer folds: ',\n",
        "      average_scores_across_outer_folds_for_each_model)\n",
        "\n",
        "many_stars = '\\n' + '*' * 100 + '\\n'\n",
        "print(many_stars + 'Now we choose the best model and refit on the whole dataset' + many_stars)\n",
        "\n",
        "best_model_name, best_model_avg_score = max(\n",
        "    average_scores_across_outer_folds_for_each_model.items(),\n",
        "    key=(lambda name_averagescore: name_averagescore[1]))\n",
        "\n",
        "# get the best model and its associated parameter grid\n",
        "best_model, best_model_params = models_and_parameters[best_model_name]\n",
        "\n",
        "# Approach 1\n",
        "#best_model.fit(X, y)\n",
        "\n",
        "# Approach 2: # refit this best model on the whole dataset so that we can start\n",
        "# making predictions on other data, and now we have a reliable estimate of\n",
        "# this model's generalization error and we are confident this is the best model\n",
        "# among the ones we have tried\n",
        "final_classifier = GridSearchCV(best_model, best_model_params, cv=inner_cv)\n",
        "final_classifier.fit(X, y)\n",
        "\n",
        "print('Best model: \\n\\t{}'.format(best_model), end='\\n\\n')\n",
        "print('Estimation of its generalization error (accuracy):\\n\\t{}'.format(\n",
        "    best_model_avg_score), end='\\n\\n')\n",
        "print('Best parameter choice for this model: \\n\\t{params}'\n",
        "      '\\n(according to cross-validation `{cv}` on the whole dataset).'.format(\n",
        "      params=final_classifier.best_params_, cv=inner_cv))"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Regression example"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 12,
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Model: svr\n",
            "MSE in the 3 outer folds: [-35614.04476457 -35178.32923572 -32425.78493914].\n",
            "Average error: -34406.05297980857\n",
            "\n",
            "Model: rf\n",
            "MSE in the 3 outer folds: [-8666.61697967 -8442.3388464  -8122.1276498 ].\n",
            "Average error: -8410.361158623991\n",
            "\n",
            "Average score across the outer folds:  {'svr': -34406.05297980857, 'rf': -8410.361158623991}\n",
            "\n",
            "****************************************************************************************************\n",
            "Now we choose the best model and refit on the whole dataset\n",
            "****************************************************************************************************\n",
            "\n",
            "Best model: \n",
            "\tRandomForestRegressor()\n",
            "\n",
            "Estimation of its generalization error (negative mean squared error):\n",
            "\t-8410.361158623991\n",
            "\n",
            "Best parameter choice for this model: \n",
            "\t{'max_depth': 50}\n",
            "(according to cross-validation `KFold(n_splits=3, random_state=None, shuffle=False)` on the whole dataset).\n"
          ]
        }
      ],
      "source": [
        "from sklearn.model_selection import KFold, cross_val_score, GridSearchCV\n",
        "from sklearn.datasets import make_regression\n",
        "from sklearn.ensemble import RandomForestRegressor\n",
        "from sklearn.svm import SVR\n",
        "import numpy as np\n",
        "\n",
        "# `outer_cv` creates 3 folds for estimating generalization error\n",
        "outer_cv = KFold(3)\n",
        "\n",
        "# when we train on a certain fold, we use a second cross-validation\n",
        "# split in order to choose hyperparameters\n",
        "inner_cv = KFold(3)\n",
        "\n",
        "# create some regression data\n",
        "X, y = make_regression(n_samples=1000, n_features=10)\n",
        "\n",
        "# give shorthand names to models and use those as dictionary keys mapping\n",
        "# to models and parameter grids for that model\n",
        "models_and_parameters = {\n",
        "    'svr': (SVR(),\n",
        "            {'C': [0.01, 0.05, 0.1, 1]}),\n",
        "    'rfr': (RandomForestRegressor(),\n",
        "           {'max_depth': [5, 10, 50, 100, 200, 500]})}\n",
        "\n",
        "# we will collect the average of the scores on the 3 outer folds in this dictionary\n",
        "# with keys given by the names of the models in `models_and_parameters`\n",
        "average_scores_across_outer_folds_for_each_model = dict()\n",
        "\n",
        "# find the model with the best generalization error\n",
        "for name, (model, params) in models_and_parameters.items():\n",
        "    # this object is a regressor that also happens to choose\n",
        "    # its hyperparameters automatically using `inner_cv`\n",
        "    regressor_that_optimizes_its_hyperparams = GridSearchCV(\n",
        "        estimator=model, param_grid=params,\n",
        "        cv=inner_cv, scoring='neg_mean_squared_error')\n",
        "\n",
        "    # estimate generalization error on the 3-fold splits of the data\n",
        "    scores_across_outer_folds = cross_val_score(\n",
        "        regressor_that_optimizes_its_hyperparams,\n",
        "        X, y, cv=outer_cv, scoring='neg_mean_squared_error')\n",
        "\n",
        "    # get the mean MSE across each of outer_cv's 3 folds\n",
        "    average_scores_across_outer_folds_for_each_model[name] = np.mean(scores_across_outer_folds)\n",
        "    error_summary = 'Model: {name}\\nMSE in the 3 outer folds: {scores}.\\nAverage error: {avg}'\n",
        "    print(error_summary.format(\n",
        "        name=name, scores=scores_across_outer_folds,\n",
        "        avg=np.mean(scores_across_outer_folds)))\n",
        "    print()\n",
        "\n",
        "print('Average score across the outer folds: ',\n",
        "      average_scores_across_outer_folds_for_each_model)\n",
        "\n",
        "many_stars = '\\n' + '*' * 100 + '\\n'\n",
        "print(many_stars + 'Now we choose the best model and refit on the whole dataset' + many_stars)\n",
        "\n",
        "best_model_name, best_model_avg_score = max(\n",
        "    average_scores_across_outer_folds_for_each_model.items(),\n",
        "    key=(lambda name_averagescore: name_averagescore[1]))\n",
        "\n",
        "# get the best model and its associated parameter grid\n",
        "best_model, best_model_params = models_and_parameters[best_model_name]\n",
        "\n",
        "# now we refit this best model on the whole dataset so that we can start\n",
        "# making predictions on other data, and now we have a reliable estimate of\n",
        "# this model's generalization error and we are confident this is the best model\n",
        "# among the ones we have tried\n",
        "final_regressor = GridSearchCV(best_model, best_model_params, cv=inner_cv)\n",
        "final_regressor.fit(X, y)\n",
        "\n",
        "print('Best model: \\n\\t{}'.format(best_model), end='\\n\\n')\n",
        "print('Estimation of its generalization error (negative mean squared error):\\n\\t{}'.format(\n",
        "    best_model_avg_score), end='\\n\\n')\n",
        "print('Best parameter choice for this model: \\n\\t{params}'\n",
        "      '\\n(according to cross-validation `{cv}` on the whole dataset).'.format(\n",
        "      params=final_regressor.best_params_, cv=inner_cv))"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "# References\n",
        "\n",
        "1. [Model selection done right: A gentle introduction to nested cross-validation](https://ploomber.io/blog/nested-cv/).\n",
        "\n",
        "2. [Which is the final model from Nested Cross Validation: Accuracy or Frequency?](https://datascience.stackexchange.com/questions/116311/which-is-the-final-model-from-nested-cross-validation-accuracy-or-frequency)\n",
        "\n",
        "3. [What is the correct procedure for nested cross-validation?](https://stackoverflow.com/questions/64238730/what-is-the-correct-procedure-for-nested-cross-validation)\n",
        "\n",
        "4. [Nested Cross Validation (Cynthia Rudin)](https://youtu.be/az60jS7MQhU?list=PLNeXFnYrCJneoY_rKtWJy833YiMrCRi5f)\n",
        "\n",
        "5. [Nested cross-validation and selecting the best regression model - is this the right SKLearn process?](https://datascience.stackexchange.com/questions/13185/nested-cross-validation-and-selecting-the-best-regression-model-is-this-the-ri)\n",
        "\n",
        "6. [Model evaluation, model selection, and algorithm selection in machine learning](https://sebastianraschka.com/blog/2016/model-evaluation-selection-part1.html)\n",
        "\n",
        "7. [Nested cross validation for model selection](https://stats.stackexchange.com/questions/65128/nested-cross-validation-for-model-selection/65158#65158)\n",
        "\n",
        "8. [Training on the full dataset after cross-validation?](https://stats.stackexchange.com/questions/11602/training-on-the-full-dataset-after-cross-validation)\n",
        "\n",
        "9. [How to choose a predictive model after k-fold cross-validation?](https://stats.stackexchange.com/questions/52274/how-to-choose-a-predictive-model-after-k-fold-cross-validation)\n",
        "\n",
        "10. [How to Train a Final Machine Learning Model](https://machinelearningmastery.com/train-final-machine-learning-model/)\n",
        "\n",
        "11. [How to get from evaluation to final model](https://mindfulmodeler.substack.com/p/how-to-get-from-evaluation-to-final)"
      ]
    }
  ],
  "metadata": {
    "colab": {
      "provenance": [],
      "toc_visible": true
    },
    "kernelspec": {
      "display_name": "Python 3",
      "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.12.0"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}