Predictive_Maintenance_Final_Project.ipynb

{
  "metadata": {
    "kernelspec": {
      "language": "python",
      "display_name": "Python 3",
      "name": "python3"
    },
    "language_info": {
      "pygments_lexer": "ipython3",
      "nbconvert_exporter": "python",
      "version": "3.6.4",
      "file_extension": ".py",
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "name": "python",
      "mimetype": "text/x-python"
    },
    "colab": {
      "provenance": [],
      "include_colab_link": true
    }
  },
  "nbformat_minor": 0,
  "nbformat": 4,
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/ramesh6762/Predictive-Maintainance-for-iiot/blob/main/Predictive_Maintenance_Final_Project.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "## **Introduction**\n",
        "\n",
        "This work is placed in the field of Data Science with application to the\n",
        "area of predictive maintenance. The need to have a way to determine\n",
        "whether or not a particular machine will fail, as well as the nature of the\n",
        "failure, is essential for generation 4.0 industries. The main reason lies behind\n",
        "the following consideration: the repair or replacement of a faulty machine\n",
        "generally requires costs that are much higher than those required for the\n",
        "replacement of a single component. Therefore, the installation of sensors that\n",
        "monitor the state of the machines, collecting the appropriate information,\n",
        "can lead to great savings for industries.\n",
        "\n",
        "Here we use the AI4I Predictive Maintenance Dataset from the UCI Repository\n",
        "to carry out an analysis that aims to respond to the needs just reported.\n",
        "In particular, the work is presented through a lineup that characterize a typical\n",
        "Machine Learning application. In the first place the dataset is explored to\n",
        "obtain a deeper knowledge that can guide in fully understanding the ground\n",
        "truth. Then, some preprocessing techniques are applied to prepare the data\n",
        "for the algorithms we will use to make our predictions. We consider two main\n",
        "tasks: the first consists in establishing whether a generic machine is about\n",
        "to suffer a failure while the second concerns the determination of the nature\n",
        "of the fault. Finally, a comparison is provided between the results obtained\n",
        "by the latter, evaluating both their performance through appropriate metrics,\n",
        "and their interpretability."
      ],
      "metadata": {
        "id": "CFoUAOIuCBu7"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "## Table of Contents\n",
        "\n",
        "1. [Task and Data Description](#description)\n",
        "2. [Exploratory Analysis](#EDA) <br>\n",
        "    2.1 [ID Columns](#ID) <br>\n",
        "    2.2 [Target anomalies](#target) <br>\n",
        "    2.3 [Outliers Inspection](#outliers) <br>\n",
        "    2.4 [Resampling with SMOTE](#resampling) <br>\n",
        "    2.5 [Comparison after resampling](#resample_comparison) <br>\n",
        "    2.6 [Features scaling and Encoding](#encoding) <br>\n",
        "    2.7 [PCA and Correlation Heatmap](#pca) <br>\n",
        "    2.8 [Metrics](#metrics) <br>\n",
        "3. [Binary task](#binary) <br>\n",
        "    3.1 [Preliminaries](#preliminaries) <br>\n",
        "    3.2 [Feature Selection attempts](#selection) <br>\n",
        "    3.3 [Logistic Regression Benchmark](#binary_benchmark) <br>\n",
        "    3.4 [Models](#binary_models) <br>\n",
        "4. [Multi-class task](#multi) <br>\n",
        "    4.1 [Logistic Regression Benchmark](#multi_benchmark) <br>\n",
        "    4.2 [Models](#multi_models) <br>\n",
        "5. [Decision Paths](#decisionpath) <br>\n",
        "6. [Conclusions](#conclusions) <br>"
      ],
      "metadata": {
        "id": "YInh4bsjCBu9"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "## 1) **Task and Data description** <a id=\"description\"></a>\n",
        "\n",
        "Since real predictive maintenance datasets are generally difficult to obtain and in particular\n",
        "difficult to publish, the data provided by the UCI repository is a synthetic dataset that reflects\n",
        "real predictive maintenance encountered in industry to the best of their knowledge.\n",
        "The dataset consists of 10 000 data points stored as rows with 14 features in columns:\n",
        "* UID: unique identifier ranging from 1 to 10000;\n",
        "* Product ID: consisting of a letter L, M, or H for low (60% of all products), medium (30%)\n",
        "and high (10%) as product quality variants and a variant-specific serial number;\n",
        "* Air temperature [K]: generated using a random walk process later normalized to a standard\n",
        "deviation of 2 K around 300 K;\n",
        "* Process temperature [K]: generated using a random walk process normalized to a standard\n",
        "deviation of 1 K, added to the air temperature plus 10 K;\n",
        "* Rotational speed [rpm]: calculated from a power of 2860 W, overlaid with a normally\n",
        "distributed noise;\n",
        "* Torque [Nm]: torque values are normally distributed around 40 Nm with a standard deviation\n",
        "of 10 Nm and no negative values;\n",
        "* Tool wear [min]: The quality variants H/M/L add 5/3/2 minutes of tool wear to the used\n",
        "tool in the process;\n",
        "* Machine failure: label that indicates, whether the machine has failed in this particular data\n",
        "point for any of the following failure modes are true.\n",
        "The machine failure consists of five independent failure modes:\n",
        "* tool wear failure (TWF): the tool will be replaced of fail at a randomly selected tool wear\n",
        "time between 200 - 240 mins;\n",
        "* heat dissipation failure (HDF): heat dissipation causes a process failure, if the difference\n",
        "between air- and process temperature is below 8.6 K and the tools rotational speed is below\n",
        "1380 rpm;\n",
        "\n",
        "* power failure (PWF):the product of torque and rotational speed (in rad/s) equals the power\n",
        "required for the process. If this power is below 3500 W or above 9000 W, the process fails;\n",
        "* overstrain failure (OSF): if the product of tool wear and torque exceeds 11,000 minNm for\n",
        "the L product variant (12,000 M, 13,000 H), the process fails due to overstrain;\n",
        "* random failures (RNF): each process has a chance of 0,1 % to fail regardless of its process\n",
        "parameters.\n",
        "If at least one of the above failure modes is true, the process fails and the ’machine failure’\n",
        "label is set to 1. It is therefore not transparent to the machine learning method, which of the\n",
        "failure modes has caused the process to fail."
      ],
      "metadata": {
        "id": "n02KY-2BCBu-"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "## 2) **Exploratory Analysis** <a id=\"EDA\"></a>\n",
        "Our data exploration starts by checking that each entry is unique and there are no duplicates;\n",
        "this is done by veryfing that the number of unique ProductID corresponds to the number of\n",
        "observations. Then we print a report to look for missing values and check the data type for each\n",
        "column."
      ],
      "metadata": {
        "id": "ctO6LEaOCBu-"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "import numpy as np\n",
        "import pandas as pd\n",
        "import matplotlib.pyplot as plt\n",
        "import seaborn as sns\n",
        "\n",
        "# Import data\n",
        "data_path = '../input/machine-predictive-maintenance-classification/predictive_maintenance.csv'\n",
        "data = pd.read_csv(data_path)\n",
        "n = data.shape[0]\n",
        "# First checks\n",
        "print('Features non-null values and data type:')\n",
        "data.info()\n",
        "print('Check for duplicate values:',\n",
        "      data['Product ID'].unique().shape[0]!=n)"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:26.380943Z",
          "iopub.execute_input": "2022-07-18T13:14:26.381358Z",
          "iopub.status.idle": "2022-07-18T13:14:26.426027Z",
          "shell.execute_reply.started": "2022-07-18T13:14:26.381325Z",
          "shell.execute_reply": "2022-07-18T13:14:26.424754Z"
        },
        "trusted": true,
        "id": "R8g-K5NdCBu-"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "To sum up even more:\n",
        "* There is no missing data;\n",
        "* There are no duplicate values;\n",
        "* Six columns are numerical features, including UDI;\n",
        "* Three columns are categorical features, including ProductID.\n",
        "\n",
        "To make this distinction more clear we set numeric columns to float type."
      ],
      "metadata": {
        "id": "b_OF3ylQCBu-"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Set numeric columns dtype to float\n",
        "data['Tool wear [min]'] = data['Tool wear [min]'].astype('float64')\n",
        "data['Rotational speed [rpm]'] = data['Rotational speed [rpm]'].astype('float64')\n",
        "# Rename features\n",
        "data.rename(mapper={'Air temperature [K]': 'Air temperature',\n",
        "                    'Process temperature [K]': 'Process temperature',\n",
        "                    'Rotational speed [rpm]': 'Rotational speed',\n",
        "                    'Torque [Nm]': 'Torque',\n",
        "                    'Tool wear [min]': 'Tool wear'}, axis=1, inplace=True)"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:27.009857Z",
          "iopub.execute_input": "2022-07-18T13:14:27.010263Z",
          "iopub.status.idle": "2022-07-18T13:14:27.019852Z",
          "shell.execute_reply.started": "2022-07-18T13:14:27.010229Z",
          "shell.execute_reply": "2022-07-18T13:14:27.018844Z"
        },
        "trusted": true,
        "id": "_WMWedybCBu-"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "### 2.1) ID Columns <a id=\"ID\"></a>\n",
        "Before going into more technical matters we deal with the two ID columns as the model we will\n",
        "use could get confused by them, since it is unrealistic to think that the failure of a machine\n",
        "depends on its identifier. However, while UDI results in being a copy of the dataframe index,\n",
        "the column Product ID is made up of an initial letter followed by five numbers; there is a small\n",
        "chance that an hidden pattern lies behind this structure. However, the initial letter corresponds\n",
        "to the machine Type and the number sequences define three intervals based on the same feature;\n",
        "this allows to confirm that the Product ID column does not actually carry any more information\n",
        "than the feature Type and it is legit to drop it.\n",
        "The following histogram shows the number sequences:"
      ],
      "metadata": {
        "id": "M9tgEM7BCBu-"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Remove first character and set to numeric dtype\n",
        "data['Product ID'] = data['Product ID'].apply(lambda x: x[1:])\n",
        "data['Product ID'] = pd.to_numeric(data['Product ID'])\n",
        "\n",
        "# Histogram of ProductID\n",
        "sns.histplot(data=data, x='Product ID', hue='Type')\n",
        "plt.show()"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:27.578844Z",
          "iopub.execute_input": "2022-07-18T13:14:27.579684Z",
          "iopub.status.idle": "2022-07-18T13:14:27.965131Z",
          "shell.execute_reply.started": "2022-07-18T13:14:27.579648Z",
          "shell.execute_reply": "2022-07-18T13:14:27.963957Z"
        },
        "trusted": true,
        "id": "q__n7zp1CBu-"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# Drop ID columns\n",
        "df = data.copy()\n",
        "df.drop(columns=['UDI','Product ID'], inplace=True)"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:28.138893Z",
          "iopub.execute_input": "2022-07-18T13:14:28.139322Z",
          "iopub.status.idle": "2022-07-18T13:14:28.148255Z",
          "shell.execute_reply.started": "2022-07-18T13:14:28.139286Z",
          "shell.execute_reply": "2022-07-18T13:14:28.146944Z"
        },
        "trusted": true,
        "id": "ACCAPD_sCBu-"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "The following pie chart shows the percentages of machines by Type:"
      ],
      "metadata": {
        "id": "WpIVfnkpCBu-"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Pie chart of Type percentage\n",
        "value = data['Type'].value_counts()\n",
        "Type_percentage = 100*value/data.Type.shape[0]\n",
        "labels = Type_percentage.index.array\n",
        "x = Type_percentage.array\n",
        "plt.pie(x, labels = labels, colors=sns.color_palette('tab10')[0:3], autopct='%.0f%%')\n",
        "plt.title('Machine Type percentage')\n",
        "plt.show()"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:28.800947Z",
          "iopub.execute_input": "2022-07-18T13:14:28.801343Z",
          "iopub.status.idle": "2022-07-18T13:14:28.913364Z",
          "shell.execute_reply.started": "2022-07-18T13:14:28.801311Z",
          "shell.execute_reply": "2022-07-18T13:14:28.911903Z"
        },
        "trusted": true,
        "id": "qQsCZlbTCBu_"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "### 2.2) Target anomalies <a id=\"target\"></a>\n",
        "\n",
        "In this section we observe the distribution of the target to find any imbalances and correct them\n",
        "before dividing the dataset.\n",
        "The first anomaly respect to dataset’s description is that when the failure is random (RNF), the\n",
        "Machine Failure feature is not set to 1."
      ],
      "metadata": {
        "id": "fO8FCYDCCBu_"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Create lists of features and target names\n",
        "features = [col for col in df.columns\n",
        "            if df[col].dtype=='float64' or col =='Type']\n",
        "target = ['Target','Failure Type']\n",
        "# Portion of data where RNF=1\n",
        "idx_RNF = df.loc[df['Failure Type']=='Random Failures'].index\n",
        "df.loc[idx_RNF,target]"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:29.487652Z",
          "iopub.execute_input": "2022-07-18T13:14:29.488103Z",
          "iopub.status.idle": "2022-07-18T13:14:29.506997Z",
          "shell.execute_reply.started": "2022-07-18T13:14:29.488069Z",
          "shell.execute_reply": "2022-07-18T13:14:29.506087Z"
        },
        "trusted": true,
        "id": "LNax-PfkCBu_"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "Fortunately the machine failure RNF occurs in only 18 observations and it has a random nature\n",
        "therefore not predictable so we decide to remove these rows."
      ],
      "metadata": {
        "id": "IyHjurhNCBu_"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "first_drop = df.loc[idx_RNF,target].shape[0]\n",
        "print('Number of observations where RNF=1 but Machine failure=0:',first_drop)\n",
        "# Drop corresponding observations and RNF column\n",
        "df.drop(index=idx_RNF, inplace=True)"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:30.086197Z",
          "iopub.execute_input": "2022-07-18T13:14:30.086637Z",
          "iopub.status.idle": "2022-07-18T13:14:30.097817Z",
          "shell.execute_reply.started": "2022-07-18T13:14:30.0866Z",
          "shell.execute_reply": "2022-07-18T13:14:30.096401Z"
        },
        "trusted": true,
        "id": "6phYBXiACBu_"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "Going forward we find out that in 9 observations Machine failure is set to 1 when all types\n",
        "of failures are set to 0. We cannot understand if there really was a failure or not so let’s remove\n",
        "these observations too."
      ],
      "metadata": {
        "id": "Jt0K6d4ECBu_"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Portion of data where Machine failure=1 but no failure cause is specified\n",
        "idx_ambiguous = df.loc[(df['Target']==1) &\n",
        "                       (df['Failure Type']=='No Failure')].index\n",
        "second_drop = df.loc[idx_ambiguous].shape[0]\n",
        "print('Number of ambiguous observations:', second_drop)\n",
        "display(df.loc[idx_ambiguous,target])\n",
        "df.drop(index=idx_ambiguous, inplace=True)"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:30.600254Z",
          "iopub.execute_input": "2022-07-18T13:14:30.603078Z",
          "iopub.status.idle": "2022-07-18T13:14:30.62458Z",
          "shell.execute_reply.started": "2022-07-18T13:14:30.603027Z",
          "shell.execute_reply": "2022-07-18T13:14:30.623624Z"
        },
        "trusted": true,
        "id": "hfU9ygjJCBu_"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# Global percentage of removed observations\n",
        "print('Global percentage of removed observations:',\n",
        "     (100*(first_drop+second_drop)/n))\n",
        "df.reset_index(drop=True, inplace=True)   # Reset index\n",
        "n = df.shape[0]"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:31.175433Z",
          "iopub.execute_input": "2022-07-18T13:14:31.175905Z",
          "iopub.status.idle": "2022-07-18T13:14:31.183578Z",
          "shell.execute_reply.started": "2022-07-18T13:14:31.175858Z",
          "shell.execute_reply": "2022-07-18T13:14:31.182057Z"
        },
        "trusted": true,
        "id": "d32ikcjoCBu_"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "Our actions did not alterate the original data very much."
      ],
      "metadata": {
        "id": "vqdZEyfeCBu_"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "### 2.3) Outliers inspection <a id=\"outliers\"></a>\n",
        "\n",
        "The goal of this section is to check if the dataset contains any outlier, which are usually misleading\n",
        "for machine learning algorithms. We begin by looking at a statistical report of the numerical\n",
        "features."
      ],
      "metadata": {
        "id": "E9MTTrL7CBu_"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "df.describe()"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:31.722857Z",
          "iopub.execute_input": "2022-07-18T13:14:31.723581Z",
          "iopub.status.idle": "2022-07-18T13:14:31.759341Z",
          "shell.execute_reply.started": "2022-07-18T13:14:31.723536Z",
          "shell.execute_reply": "2022-07-18T13:14:31.758129Z"
        },
        "trusted": true,
        "id": "0pknuqaOCBu_"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "We can guess the presence of outliers in Rotational Speed and Torque because the maximum is\n",
        "very different from the third quartile. To make this consideration more concrete we take a closer\n",
        "look at the situation with boxplots, using histograms to understand the distribution."
      ],
      "metadata": {
        "id": "EC10Ca5NCBu_"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "num_features = [feature for feature in features if df[feature].dtype=='float64']\n",
        "# Histograms of numeric features\n",
        "fig, axs = plt.subplots(nrows=2, ncols=3, figsize=(18,7))\n",
        "fig.suptitle('Numeric features histogram')\n",
        "for j, feature in enumerate(num_features):\n",
        "    sns.histplot(ax=axs[j//3, j-3*(j//3)], data=df, x=feature)\n",
        "plt.show()\n",
        "\n",
        "# boxplot of numeric features\n",
        "fig, axs = plt.subplots(nrows=2, ncols=3, figsize=(18,7))\n",
        "fig.suptitle('Numeric features boxplot')\n",
        "for j, feature in enumerate(num_features):\n",
        "    sns.boxplot(ax=axs[j//3, j-3*(j//3)], data=df, x=feature)\n",
        "plt.show()"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:32.199679Z",
          "iopub.execute_input": "2022-07-18T13:14:32.200123Z",
          "iopub.status.idle": "2022-07-18T13:14:34.214477Z",
          "shell.execute_reply.started": "2022-07-18T13:14:32.200092Z",
          "shell.execute_reply": "2022-07-18T13:14:34.21307Z"
        },
        "trusted": true,
        "id": "8XWK8ALqCBu_"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "The boxplots highlight possible outliers in the features mentioned above, however in the case\n",
        "of Torque these are probably traceable to the way outliers are detected using boxplots (since the\n",
        "distribution is Gaussian it would be more appropriate to use the 3σ rule instead of the IQR); in\n",
        "the case of Rotational Speed the Gaussian distribution is skewed and it is not unrealistic to think\n",
        "that the few observation with high Rotational Speed are going to fail. As a result we keep the\n",
        "outliers for now and we reserve the right to decide whether to act on them or not after considering\n",
        "other aspects."
      ],
      "metadata": {
        "id": "B3TUPxIhCBu_"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "### 2.4) Resampling with SMOTE <a id=\"resampling\"></a>\n",
        "\n",
        "Another important consideration regards the extremely low occurrence of machine failures among\n",
        "the entire dataset, which percentage is equal only to 3.31%. Moreover, a pie plot showing the\n",
        "occurrence of the causes involved for each failure reveals a further degree of imbalance."
      ],
      "metadata": {
        "id": "UCx4aWbgCBu_"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Portion of df where there is a failure and causes percentage\n",
        "idx_fail = df.loc[df['Failure Type'] != 'No Failure'].index\n",
        "df_fail = df.loc[idx_fail]\n",
        "df_fail_percentage = 100*df_fail['Failure Type'].value_counts()/df_fail['Failure Type'].shape[0]\n",
        "print('Failures percentage in data:',\n",
        "      round(100*df['Target'].sum()/n,2))\n",
        "# Pie plot\n",
        "plt.title('Causes involved in Machine failures')\n",
        "plt.pie(x=df_fail_percentage.array, labels=df_fail_percentage.index.array,\n",
        "        colors=sns.color_palette('tab10')[0:4], autopct='%.0f%%')\n",
        "plt.show()"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:34.217396Z",
          "iopub.execute_input": "2022-07-18T13:14:34.217914Z",
          "iopub.status.idle": "2022-07-18T13:14:34.349351Z",
          "shell.execute_reply.started": "2022-07-18T13:14:34.217865Z",
          "shell.execute_reply": "2022-07-18T13:14:34.347913Z"
        },
        "trusted": true,
        "id": "T1JeyS3RCBvA"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "When dealing with machine learning problems classes imbalance is a great concern, as it can\n",
        "mislead both the models training process and our ability to interpret their results. As instance,\n",
        "if we build a model on this dataset that predicts that machines never fail, it should be 97%\n",
        "accurate. In order to avoid such effects and limit the preferential behaviour of the models with\n",
        "respect to individual classes we perform a data augmentation, with the aim of obtaining a ratio\n",
        "of 80 to 20 between functioning and faulty observations and the same percentage of occurrence\n",
        "between the causes involved in the failures.\n",
        "\n",
        "Among the most common data augmentation techniques we identify:\n",
        "* Under-sampling by deleting some data points from the majority class.\n",
        "* Over-Sampling by copying rows of data resulting in the minority class.\n",
        "* Over-Sampling with SMOTE (Synthetic Minority Oversampling Technique).\n",
        "\n",
        "The first two choices however result in extremely simplistic approaches; in particular the first one\n",
        "has the disadvantage of decreasing the length of the dataset in a context in which the available\n",
        "data are already limited. Therefore we use the SMOTE procedure to generate new samples,\n",
        "which is very much like slightly moving the data point in the direction of its neighbors. This way,\n",
        "the synthetic data point is not an exact copy of an existing data point but we can also be sure\n",
        "that it is also not too different from the known observations in the minority class. To be more\n",
        "precise, the SMOTE procedure works as follows: it draws a random sample from the minority class and for the observations in this sample, identifies the k nearest neighbors. It will then take\n",
        "one of those neighbors and identify the vector between the current data point and the selected\n",
        "neighbor. The vector will be multiplied by a random number between 0 and 1 and the synthetic\n",
        "data point is obtained by adding this vector to the current data point."
      ],
      "metadata": {
        "id": "fiz4vvFSCBvA"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "from sklearn.model_selection import train_test_split\n",
        "from imblearn.over_sampling import SMOTENC\n",
        "# n_working must represent 80% of the desired length of resampled dataframe\n",
        "n_working = df['Failure Type'].value_counts()['No Failure']\n",
        "desired_length = round(n_working/0.8)\n",
        "spc = round((desired_length-n_working)/4)  #samples per class\n",
        "# Resampling\n",
        "balance_cause = {'No Failure':n_working,\n",
        "                 'Overstrain Failure':spc,\n",
        "                 'Heat Dissipation Failure':spc,\n",
        "                 'Power Failure':spc,\n",
        "                 'Tool Wear Failure':spc}\n",
        "sm = SMOTENC(categorical_features=[0,7], sampling_strategy=balance_cause, random_state=0)\n",
        "df_res, y_res = sm.fit_resample(df, df['Failure Type'])"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:34.351244Z",
          "iopub.execute_input": "2022-07-18T13:14:34.352047Z",
          "iopub.status.idle": "2022-07-18T13:14:34.779888Z",
          "shell.execute_reply.started": "2022-07-18T13:14:34.351989Z",
          "shell.execute_reply": "2022-07-18T13:14:34.77851Z"
        },
        "trusted": true,
        "id": "tWHGfPG2CBvA"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "### 2.5) Comparison after resampling <a id=\"resample_comparison\"></a>\n",
        "\n",
        "The result is described in the following pie charts."
      ],
      "metadata": {
        "id": "wWnF8Fq2CBvA"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Portion of df_res where there is a failure and causes percentage\n",
        "idx_fail_res = df_res.loc[df_res['Failure Type'] != 'No Failure'].index\n",
        "df_res_fail = df_res.loc[idx_fail_res]\n",
        "fail_res_percentage = 100*df_res_fail['Failure Type'].value_counts()/df_res_fail.shape[0]\n",
        "\n",
        "# Percentages\n",
        "print('Percentage increment of observations after oversampling:',\n",
        "      round((df_res.shape[0]-df.shape[0])*100/df.shape[0],2))\n",
        "print('SMOTE Resampled Failures percentage:',\n",
        "      round(df_res_fail.shape[0]*100/df_res.shape[0],2))\n",
        "\n",
        "# Pie plot\n",
        "fig, axs = plt.subplots(ncols=2, figsize=(12,4))\n",
        "fig.suptitle('Causes involved in Machine failures')\n",
        "axs[0].pie(x=df_fail_percentage.array, labels=df_fail_percentage.index.array,\n",
        "        colors=sns.color_palette('tab10')[0:4], autopct='%.0f%%')\n",
        "axs[1].pie(x=fail_res_percentage.array, labels=fail_res_percentage.index.array,\n",
        "        colors=sns.color_palette('tab10')[0:4], autopct='%.0f%%')\n",
        "axs[0].title.set_text('Original')\n",
        "axs[1].title.set_text('After Resampling')\n",
        "plt.show()"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:34.782302Z",
          "iopub.execute_input": "2022-07-18T13:14:34.782677Z",
          "iopub.status.idle": "2022-07-18T13:14:35.026403Z",
          "shell.execute_reply.started": "2022-07-18T13:14:34.782645Z",
          "shell.execute_reply": "2022-07-18T13:14:35.024857Z"
        },
        "trusted": true,
        "id": "NNpW407_CBvA"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "As one can expect, the cases of Machine Failure mainly concern low quality machines, then\n",
        "those of medium quality and only a few times those of high quality. This difference is accentuated\n",
        "when the number of observations of non-functioning machines is (artificially) increased. However,\n",
        "from the kdeplots below it can be seen that this is not widely correlated with the features since\n",
        "differentiating according to the quality shows that distribution of the features does not present\n",
        "big differences, except for the two side peaks in Tool Wear (which is consistent with the data\n",
        "description). This suggests that probably the fact that the majority of failures concern type\n",
        "L machines is due to the greater presence of this type in the dataset and therefore that the\n",
        "correlation with the failure of the machine is due to statistical reasons."
      ],
      "metadata": {
        "id": "0l2XlKO8CBvA"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Kdeplot of numeric features (After resampling) - hue=Type\n",
        "fig, axs = plt.subplots(nrows=2, ncols=3, figsize=(19,7))\n",
        "fig.suptitle('Features distribution (After resampling)')\n",
        "custom_palette = {'L':'tab:blue', 'M':'tab:orange', 'H':'tab:green'}\n",
        "for j, feature in enumerate(num_features):\n",
        "    sns.kdeplot(ax=axs[j//3, j-3*(j//3)], data=df_res, x=feature,\n",
        "              hue='Type', fill=True, palette=custom_palette)\n",
        "plt.show()"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:35.028537Z",
          "iopub.execute_input": "2022-07-18T13:14:35.030032Z",
          "iopub.status.idle": "2022-07-18T13:14:36.942121Z",
          "shell.execute_reply.started": "2022-07-18T13:14:35.029975Z",
          "shell.execute_reply": "2022-07-18T13:14:36.941024Z"
        },
        "trusted": true,
        "id": "KMzX5EkXCBvA"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "Finally, let’s look at how the distribution of features has changed."
      ],
      "metadata": {
        "id": "NSRUTTEACBvA"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Kdeplot of numeric features (Original)\n",
        "fig, axs = plt.subplots(nrows=2, ncols=3, figsize=(18,7))\n",
        "fig.suptitle('Original Features distribution')\n",
        "enumerate_features = enumerate(num_features)\n",
        "for j, feature in enumerate_features:\n",
        "    sns.kdeplot(ax=axs[j//3, j-3*(j//3)], data=df, x=feature,\n",
        "                hue='Target', fill=True, palette='tab10')\n",
        "plt.show()\n",
        "# Kdeplot of numeric features (After resampling)\n",
        "fig, axs = plt.subplots(nrows=2, ncols=3, figsize=(18,7))\n",
        "fig.suptitle('Features distribution after oversampling')\n",
        "enumerate_features = enumerate(num_features)\n",
        "for j, feature in enumerate_features:\n",
        "    sns.kdeplot(ax=axs[j//3, j-3*(j//3)], data=df_res, x=feature,\n",
        "                hue=df_res['Target'], fill=True, palette='tab10')\n",
        "plt.show()\n",
        "# Kdeplot of numeric features (After resampling) - Diving deeper\n",
        "fig, axs = plt.subplots(nrows=2, ncols=3, figsize=(18,7))\n",
        "fig.suptitle('Features distribution after oversampling - Diving deeper')\n",
        "enumerate_features = enumerate(num_features)\n",
        "for j, feature in enumerate_features:\n",
        "    sns.kdeplot(ax=axs[j//3, j-3*(j//3)], data=df_res, x=feature,\n",
        "                hue=df_res['Failure Type'], fill=True, palette='tab10')\n",
        "plt.show()"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:36.943801Z",
          "iopub.execute_input": "2022-07-18T13:14:36.944269Z",
          "iopub.status.idle": "2022-07-18T13:14:40.964783Z",
          "shell.execute_reply.started": "2022-07-18T13:14:36.944217Z",
          "shell.execute_reply": "2022-07-18T13:14:40.963936Z"
        },
        "trusted": true,
        "id": "Q9DcyouQCBvA"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "The first thing we can observe is that the data augmentation was performed succesfully, as the\n",
        "feature distribution for faulty instancies have not been significantly distorted. It should also be\n",
        "noted that in Rotational Speed, Torque and Tool Wear the observations relating to failures have\n",
        "a density peak in extreme zones of the distribution. This implies that the outliers we discussed\n",
        "in Section 2.3 are not to be imputed to mistakes in the dataset building but rather to the natural\n",
        "variance of the same. This becomes even clearer when observing the distributions relative to the\n",
        "single causes of failure: in particular, an almost symmetrical behavior is recognized in Rotational\n",
        "Speed and Torque while in Tool Wear a clear separation is observed between PWF and HDF\n",
        "failures on lower values, and the peaks that are found at higher values relative to TWF and OSF.\n",
        "This is perfectly consistent with the description of the targets reported in the \"Task and dataset\n",
        "description\" section."
      ],
      "metadata": {
        "id": "DZLIG4rZCBvA"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "### 2.6) Features scaling and Encoding <a id=\"encoding\"></a>\n",
        "\n",
        "In order to make data exploitable for the algorithms we will run, we apply two transformations:\n",
        "* First, we apply a label encoding to the categorical columns, since Type is an ordinal feature\n",
        "and Cause must be represented in one column. The mapping follows this scheme:\n",
        "Type: {L=0, M=1, H=2}\n",
        "Cause: {Working=0, PWF=1, OSF=2, HDF=3, TWF=4}\n",
        "* Secondly we perform the scaling of the columns with StandardScaler. This is particularly\n",
        "useful for the good working of methods that rely on the metric space, such as PCA and KNN.\n",
        "It has been also verified that using StandardScaler leads to slightly better performances than\n",
        "using MinMaxScaler."
      ],
      "metadata": {
        "id": "gek6BvbGCBvA"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "from sklearn.preprocessing import StandardScaler\n",
        "from sklearn.decomposition import PCA\n",
        "\n",
        "sc = StandardScaler()\n",
        "type_dict = {'L': 0, 'M': 1, 'H': 2}\n",
        "cause_dict = {'No Failure': 0,\n",
        "              'Power Failure': 1,\n",
        "              'Overstrain Failure': 2,\n",
        "              'Heat Dissipation Failure': 3,\n",
        "              'Tool Wear Failure': 4}\n",
        "df_pre = df_res.copy()\n",
        "# Encoding\n",
        "df_pre['Type'].replace(to_replace=type_dict, inplace=True)\n",
        "df_pre['Failure Type'].replace(to_replace=cause_dict, inplace=True)\n",
        "# Scaling\n",
        "df_pre[num_features] = sc.fit_transform(df_pre[num_features])"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:40.967288Z",
          "iopub.execute_input": "2022-07-18T13:14:40.967981Z",
          "iopub.status.idle": "2022-07-18T13:14:41.002841Z",
          "shell.execute_reply.started": "2022-07-18T13:14:40.96794Z",
          "shell.execute_reply": "2022-07-18T13:14:41.001903Z"
        },
        "trusted": true,
        "id": "NBC9Ua6tCBvA"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "### 2.7) PCA and Correlation Heatmap <a id=\"pca\"></a>\n",
        "\n",
        "We run PCA to have a further way of displaying the data instead of making feature selection."
      ],
      "metadata": {
        "id": "nw-UBJ20CBvE"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "pca = PCA(n_components=len(num_features))\n",
        "X_pca = pd.DataFrame(data=pca.fit_transform(df_pre[num_features]), columns=['PC'+str(i+1) for i in range(len(num_features))])\n",
        "var_exp = pd.Series(data=100*pca.explained_variance_ratio_, index=['PC'+str(i+1) for i in range(len(num_features))])\n",
        "print('Explained variance ratio per component:', round(var_exp,2), sep='\\n')\n",
        "print('Explained variance ratio with 3 components: '+str(round(var_exp.values[:3].sum(),2)))"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:41.004452Z",
          "iopub.execute_input": "2022-07-18T13:14:41.004823Z",
          "iopub.status.idle": "2022-07-18T13:14:41.024215Z",
          "shell.execute_reply.started": "2022-07-18T13:14:41.00479Z",
          "shell.execute_reply": "2022-07-18T13:14:41.022814Z"
        },
        "trusted": true,
        "id": "AUc87yUOCBvE"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "Since the first three components are enough to almost fully represent the variance of the data\n",
        "we will project them in a three dimensional space."
      ],
      "metadata": {
        "id": "_APAFfPMCBvE"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# PCA for Data visualization\n",
        "pca3 = PCA(n_components=3)\n",
        "X_pca3 = pd.DataFrame(data=pca3.fit_transform(df_pre[num_features]), columns=['PC1','PC2','PC3'])\n",
        "\n",
        "# Loadings Analysis\n",
        "fig, axs = plt.subplots(ncols=3, figsize=(18,4))\n",
        "fig.suptitle('Loadings magnitude')\n",
        "pca_loadings = pd.DataFrame(data=pca3.components_, columns=num_features)\n",
        "for j in range(3):\n",
        "    ax = axs[j]\n",
        "    sns.barplot(ax=ax, x=pca_loadings.columns, y=pca_loadings.values[j])\n",
        "    ax.tick_params(axis='x', rotation=90)\n",
        "    ax.title.set_text('PC'+str(j+1))\n",
        "plt.show()"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:41.026537Z",
          "iopub.execute_input": "2022-07-18T13:14:41.02796Z",
          "iopub.status.idle": "2022-07-18T13:14:41.54442Z",
          "shell.execute_reply.started": "2022-07-18T13:14:41.027913Z",
          "shell.execute_reply": "2022-07-18T13:14:41.543274Z"
        },
        "trusted": true,
        "id": "79z6p2qcCBvF"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "The bar plot of Principal Components weights makes easy to understand what they represent:\n",
        "* PC1 is closely related to the two temperature data;\n",
        "* PC2 can be identified with the machine power, which is the product of Rotational Speed\n",
        "and Torque;\n",
        "* PC3 is identifiable with Tool Wear."
      ],
      "metadata": {
        "id": "V7Z8Z3IcCBvF"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "X_pca3.rename(mapper={'PC1':'Temperature',\n",
        "                      'PC2':'Power',\n",
        "                      'PC3':'Tool Wear'}, axis=1, inplace=True)\n",
        "\n",
        "# PCA plot\n",
        "color = []\n",
        "col = df_pre['Failure Type'].map({0:'tab:blue',1:'tab:orange',2:'tab:green',3:'tab:red',4:'tab:purple'})\n",
        "color.append(col)\n",
        "idx_w = col[col == 'tab:blue'].index\n",
        "color.append(col.drop(idx_w))\n",
        "colors = ['tab:blue','tab:orange','tab:green','tab:red','tab:purple']\n",
        "labelTups = [('No Failure','tab:blue'),\n",
        "             ('Power Failure', 'tab:orange'),\n",
        "             ('Overstrain Failure','tab:green'),\n",
        "             ('Heat Dissipation Failure', 'tab:red'),\n",
        "             ('Tool Wear Failure','tab:purple')]\n",
        "\n",
        "fig = plt.figure(figsize=(18,6))\n",
        "fig.suptitle('Data in 3D PCA space')\n",
        "full_idx = X_pca3.index\n",
        "\n",
        "for j, idx in enumerate([full_idx,idx_fail_res]):\n",
        "    ax = fig.add_subplot(1, 2, j+1, projection='3d')\n",
        "\n",
        "    lg = ax.scatter(X_pca3.loc[idx,'Temperature'],\n",
        "                    X_pca3.loc[idx,'Power'],\n",
        "                    X_pca3.loc[idx,'Tool Wear'],\n",
        "                    c=color[j])\n",
        "    ax.set_xlabel('$Temperature$')\n",
        "    ax.set_ylabel('$Power$')\n",
        "    ax.set_zlabel('$Tool Wear$')\n",
        "    ax.title.set_text('With'+str(j*'out')+' \"No Failure\" class')\n",
        "    ax.view_init(35, -10)\n",
        "    custom_lines = [plt.Line2D([],[], ls=\"\", marker='.',\n",
        "                               mec='k', mfc=c, mew=.1, ms=20) for c in colors[j:]]\n",
        "    ax.legend(custom_lines, [lt[0] for lt in labelTups[j:]],\n",
        "              loc='center left', bbox_to_anchor=(1.0, .5))\n",
        "\n",
        "plt.show()"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:41.546339Z",
          "iopub.execute_input": "2022-07-18T13:14:41.546669Z",
          "iopub.status.idle": "2022-07-18T13:14:42.353979Z",
          "shell.execute_reply.started": "2022-07-18T13:14:41.546639Z",
          "shell.execute_reply": "2022-07-18T13:14:42.352735Z"
        },
        "trusted": true,
        "id": "Gg0H0d53CBvF"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "The projection into the space generated by these three axes highlights that:\n",
        "* TWF is the class of failures best separated from all the others and seems to depend almost\n",
        "entirely on PC3 (Tool Wear);\n",
        "* PWF occupies two extreme bands along the PC2 (Power), it is independent of the other\n",
        "two components;\n",
        "18\n",
        "* The OSF and HDF classes are less separated than the others even if it can be observed\n",
        "that the first is characterized by a high Tool Wear and low power while the second is\n",
        "characterized by a high temperature and a low power."
      ],
      "metadata": {
        "id": "VifbG5DmCBvF"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Correlation Heatmap\n",
        "plt.figure(figsize=(7,4))\n",
        "sns.heatmap(data=df_pre.corr(), mask=np.triu(df_pre.corr()), annot=True, cmap='BrBG')\n",
        "plt.title('Correlation Heatmap')\n",
        "plt.show()"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:42.355681Z",
          "iopub.execute_input": "2022-07-18T13:14:42.356133Z",
          "iopub.status.idle": "2022-07-18T13:14:42.75554Z",
          "shell.execute_reply.started": "2022-07-18T13:14:42.35609Z",
          "shell.execute_reply": "2022-07-18T13:14:42.754182Z"
        },
        "trusted": true,
        "id": "GWEvoMZgCBvF"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "Unsurprisingly, we observe that the features related to temperature, as well as those related to\n",
        "power, are widely correlated. Furthermore, Tool Wear correlates well with both of our targets,\n",
        "confirming what we have observed by studying PCA. Finally, a less strong correlation is also\n",
        "observed between the torsion and the two targets."
      ],
      "metadata": {
        "id": "377mzsaDCBvF"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "### 2.8) Metrics <a id=\"metrics\"></a>\n",
        "\n",
        "To evaluate the models we will use from a quantitative point of view, we resort to some metrics\n",
        "that summarize some characteristics of the classification results:\n",
        "* Accuracy: expresses the fraction of instances that are classified correctly, it is the most intuitive metric that is usually used in classification tasks.\n",
        "$$ Accuracy = \\frac{TP + TN}{TP + TN + FT + FN} $$\n",
        "* AUC: can be considered as a measure of the separation between True Positives and True Negatives, that is, the ability of the model to distinguish between classes. In detail, it represents the area below the ROC curve, given by the estimate of the True Positive Rate (Recall) for each possible value of the True Negative Rate).\n",
        "* F1: reports the classification capacity of the model to Precision and Recall, giving both the same weight.\n",
        "$$F1 = 2\\frac{Precision * Recall}{Precision + Recall}$$\n",
        "Although generally effective, AUC can be optimistic in the case of highly unbalanced\n",
        "classes, as happens in the binary task, while the F1 score is more reliable in this kind of scenario.\n",
        "We consider this last metric particularly significant as it is able to mediate the cases in which\n",
        "the machines that are about to fail are classified as functioning (Recall) and the one in which\n",
        "functioning machines are classified as about to suffer a failure (Precision). To be more specific\n",
        "we will give more importance to Recall than Precision, by evaluating also an \"adjusted\" version\n",
        "of the F1 through a β parameter:\n",
        "$$F_\\beta = (1 + \\beta^2)\\frac{Precision * Recall}{\\beta^2  Precision + Recall}$$\n",
        "With the choice $\\beta = 2$ (common in literature) a greater influence of the Recall is obtained.\n",
        "This choice is motivated by the fact that in order to optimize the costs for the maintenance of\n",
        "the machinery it is a good thing to limit the purchase of unnecessary replacement materials but\n",
        "it results far more important to avoid the possibility of having to replace a machinery after it is\n",
        "broken, since this second scenario generally has higher costs."
      ],
      "metadata": {
        "id": "SPOwkLdVCBvF"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "## 3) **Binary task** <a id=\"binary\"></a>"
      ],
      "metadata": {
        "id": "E9aqVQoHCBvF"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "### 3.1) Preliminaries <a id=\"preliminaries\"></a>\n",
        "\n",
        "The goal of this section is to find the best model for binary classification of the dataset to predict whether or not there will be Machine Failure. Classification algorithms are part of data mining and use supervised machine learning methods to make predictions about data. In particular, a\n",
        "set of data already divided (”labeled”) into two or more classes of belonging is provided as input\n",
        "thanks to which a classification model is created, which will than be used on new (”unlabeled”)\n",
        "data to assign them to the appropriate class. The starting dataset is usually divided into three\n",
        "groups: the training dataset, i.e. the sample of data used to fit the model, the validation dataset,\n",
        "i.e. the sample of data used to provide an evaluation of a model fit on the training dataset while\n",
        "tuning model hyperparameters and the test dataset, which has the purpose of testing the model.\n",
        "At the beginning of a project a data scientist must make this division and the common ratios\n",
        "used are:\n",
        "* 70% train, 15% val, 15% test.\n",
        "* 80% train, 10% val, 10% test.\n",
        "* 60% train, 20% val, 20% test.\n",
        "\n",
        "In this project we use the ratio (80/10/10) for the split because we test the model for all of these\n",
        "strategies and find that it is the best one.\n",
        "The classification techniques we choose to implement are the following:\n",
        "* Logistic Regression: it estimates the probability of a dependent variable as a function\n",
        "of independent variables. The dependent variable is the output that we are trying to\n",
        "predict while the independent variables or explanatory variables are the factors that we\n",
        "feel could influence the output. For its simplicity and interpretability, we decide to use\n",
        "Logistic Regression as a Benchmark model, a basic model that represents the starting point\n",
        "for comparing the results obtained from other models.\n",
        "* K-nearest neighbors (K-NN): algorithm based on the calculation of the distance between the elements of the dataset. Data is assigned to a certain class if close enough to the other data of the same class. Parameter K represents the number of neighboring data taken into account when assigning classes.\n",
        "* Support Vector Machine: its aim is to find a hyperplane in an N-dimensional space (N—the number of features) that distinctly classifies the data points while maximizing the margin distance, i.e. the distance between data points of both classes.\n",
        "* Random Forest: it uses ensemble learning, which is a technique that combines many classifiers to provide solutions to complex problems. Random Forest uses bagging technique: it constructs a multitude of decision trees in parallel, all with the same importance, and the output is the class selected by most trees.\n",
        "* XGBoost: is a gradient-boosted decision tree (GBDT) machine learning library. A Gradient Boosting Decision Tree (GBDT) is a decision tree ensemble learning algorithm similar to Random Forest, from which differs because it uses a boosting technique: it iteratively trains an ensemble of shallow decision trees, with each iteration using the error residuals of the previous model to fit the next model. The final prediction is a weighted sum of all of the tree predictions."
      ],
      "metadata": {
        "id": "c05zm0-WCBvF"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "from sklearn.metrics import accuracy_score, roc_auc_score, f1_score, fbeta_score\n",
        "from sklearn.metrics import confusion_matrix, make_scorer\n",
        "from sklearn.inspection import permutation_importance\n",
        "from sklearn.model_selection import GridSearchCV\n",
        "from sklearn.linear_model import LogisticRegression\n",
        "from sklearn.neighbors import KNeighborsClassifier\n",
        "from sklearn.ensemble import RandomForestClassifier\n",
        "from xgboost import XGBClassifier\n",
        "from sklearn.svm import SVC\n",
        "import time\n",
        "\n",
        "# train-validation-test split\n",
        "X, y = df_pre[features], df_pre[['Target','Failure Type']]\n",
        "X_trainval, X_test, y_trainval, y_test = train_test_split(X, y, test_size=0.1, stratify=df_pre['Failure Type'], random_state=0)\n",
        "X_train, X_val, y_train, y_val = train_test_split(X_trainval, y_trainval, test_size=0.11, stratify=y_trainval['Failure Type'], random_state=0)"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:42.757538Z",
          "iopub.execute_input": "2022-07-18T13:14:42.758068Z",
          "iopub.status.idle": "2022-07-18T13:14:42.786281Z",
          "shell.execute_reply.started": "2022-07-18T13:14:42.75801Z",
          "shell.execute_reply": "2022-07-18T13:14:42.78519Z"
        },
        "trusted": true,
        "id": "8qC5aYMTCBvF"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "We define some functions to make the following subsections easier to read. If not interested in the details we suggest to skip to Section 3.1"
      ],
      "metadata": {
        "id": "gIwOdEzPCBvG"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "\"\"\"User-defined function: Evaluate cm, accurcay, AUC, F1 for a given classifier\n",
        "- model, fitted estimator.\n",
        "- X, data used to estimate class probabilities (paired with y_true)\n",
        "- y_true, ground truth with two columns\n",
        "- y_pred, predictions\n",
        "- task = 'binary','multi_class'\n",
        "\"\"\"\n",
        "def eval_preds(model,X,y_true,y_pred,task):\n",
        "    if task == 'binary':\n",
        "        # Extract task target\n",
        "        y_true = y_true['Target']\n",
        "        cm = confusion_matrix(y_true, y_pred)\n",
        "        # Probability of the minority class\n",
        "        proba = model.predict_proba(X)[:,1]\n",
        "        # Metrics\n",
        "        acc = accuracy_score(y_true, y_pred)\n",
        "        auc = roc_auc_score(y_true, proba)\n",
        "        f1 = f1_score(y_true, y_pred, pos_label=1)\n",
        "        f2 = fbeta_score(y_true, y_pred, pos_label=1, beta=2)\n",
        "    elif task == 'multi_class':\n",
        "        y_true = y_true['Failure Type']\n",
        "        cm = confusion_matrix(y_true, y_pred)\n",
        "        proba = model.predict_proba(X)\n",
        "        # Metrics\n",
        "        acc = accuracy_score(y_true, y_pred)\n",
        "        auc = roc_auc_score(y_true, proba, multi_class='ovr', average='weighted')\n",
        "        f1 = f1_score(y_true, y_pred, average='weighted')\n",
        "        f2 = fbeta_score(y_true, y_pred, beta=2, average='weighted')\n",
        "    metrics = pd.Series(data={'ACC':acc, 'AUC':auc, 'F1':f1, 'F2':f2})\n",
        "    metrics = round(metrics,3)\n",
        "    return cm, metrics\n",
        "\n",
        "\n",
        "\n",
        "\"\"\"User-defined function: Fits one estimator using GridSearch to search for the best parameters\n",
        "- clf, estimator\n",
        "- X, y = X_train, y_train\n",
        "- params, parameters grid for GridSearch\n",
        "- task = 'binary','multi_class'\n",
        "\"\"\"\n",
        "def tune_and_fit(clf,X,y,params,task):\n",
        "    if task=='binary':\n",
        "        f2_scorer = make_scorer(fbeta_score, pos_label=1, beta=2)\n",
        "        start_time = time.time()\n",
        "        grid_model = GridSearchCV(clf, param_grid=params,\n",
        "                                cv=5, scoring=f2_scorer)\n",
        "        grid_model.fit(X, y['Target'])\n",
        "    elif task=='multi_class':\n",
        "        f2_scorer = make_scorer(fbeta_score, beta=2, average='weighted')\n",
        "        start_time = time.time()\n",
        "        grid_model = GridSearchCV(clf, param_grid=params,\n",
        "                              cv=5, scoring=f2_scorer)\n",
        "        grid_model.fit(X, y['Failure Type'])\n",
        "\n",
        "    print('Best params:', grid_model.best_params_)\n",
        "    # Print training times\n",
        "    train_time = time.time()-start_time\n",
        "    mins = int(train_time//60)\n",
        "    print('Training time: '+str(mins)+'m '+str(round(train_time-mins*60))+'s')\n",
        "    return grid_model\n",
        "\n",
        "\n",
        "\n",
        "\"\"\"User-defined function: Makes predictions using the tuned classifiers.\n",
        "Then uses eval_preds to compute the relative metrics. Returns:\n",
        "- y_pred, DataFrame containing the predictions of each model\n",
        "- cm_list, confusion matrix list\n",
        "- metrics, DataFrame containing the metrics\n",
        "Input:\n",
        "- fitted_models, fitted estimators\n",
        "- X, data used to make predictions\n",
        "- y_true, true values for target\n",
        "- clf_str, list containing estimators names\n",
        "- task = 'binary','multi_class'\n",
        "\"\"\"\n",
        "def predict_and_evaluate(fitted_models,X,y_true,clf_str,task):\n",
        "    cm_dict = {key: np.nan for key in clf_str}\n",
        "    metrics = pd.DataFrame(columns=clf_str)\n",
        "    y_pred = pd.DataFrame(columns=clf_str)\n",
        "    for fit_model, model_name in zip(fitted_models,clf_str):\n",
        "        # Update predictions\n",
        "        y_pred[model_name] = fit_model.predict(X)\n",
        "        # Metrics\n",
        "        if task == 'binary':\n",
        "            cm, scores = eval_preds(fit_model,X,y_true,\n",
        "                                     y_pred[model_name],task)\n",
        "        elif task == 'multi_class':\n",
        "            cm, scores = eval_preds(fit_model,X,y_true,\n",
        "                                     y_pred[model_name],task)\n",
        "        # Update Confusion matrix and metrics\n",
        "        cm_dict[model_name] = cm\n",
        "        metrics[model_name] = scores\n",
        "    return y_pred, cm_dict, metrics\n",
        "\n",
        "\n",
        "\n",
        "\"\"\"User-defined function: Fit the estimators on multiple classifiers\n",
        "- clf, estimators\n",
        "- clf_str, list containing estimators names\n",
        "- X_train,y_train, data used to fit models\n",
        "- X_val,y_val, data used to validate models\n",
        "\"\"\"\n",
        "\n",
        "def fit_models(clf,clf_str,X_train,X_val,y_train,y_val):\n",
        "    metrics = pd.DataFrame(columns=clf_str)\n",
        "    for model, model_name in zip(clf, clf_str):\n",
        "        model.fit(X_train,y_train['Target'])\n",
        "        y_val_pred = model.predict(X_val)\n",
        "        metrics[model_name] = eval_preds(model,X_val,y_val,y_val_pred,'binary')[1]\n",
        "    return metrics"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:42.787974Z",
          "iopub.execute_input": "2022-07-18T13:14:42.788318Z",
          "iopub.status.idle": "2022-07-18T13:14:42.813969Z",
          "shell.execute_reply.started": "2022-07-18T13:14:42.788288Z",
          "shell.execute_reply": "2022-07-18T13:14:42.812738Z"
        },
        "trusted": true,
        "id": "qm3NGMmfCBvG"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "### 3.2) Feature selection attempts <a id=\"selection\"></a>\n",
        "\n",
        "Before going into the training of the models just mentioned we try to perform feature selection,\n",
        "exploiting the considerations we made about the correlation heatmap and the exploratory data\n",
        "analysis: just to remind, we noticed that the features \"Process temperature\" and \"Air temperature\"\n",
        "are positively correlated, and \"Torque\" and \"Rotational speed\" are negatively correlated.\n",
        "From the dateset description we see that the PWF failure occurs if the product between \"Torque\"\n",
        "and \"Rotational speed\" is in a certain range of values and, similarly, HDF failure occurs when\n",
        "the difference between \"Air temperature\" and \"Process temperature\" exceeds a certain value.\n",
        "For these reasons, completely deleting these columns seems to be a bad choice because important\n",
        "information can be lost but at the same time it is reasonable to see what happens if we combine\n",
        "them, taken by pairs, to create new features that still preserve a physical meaning. Therefore we\n",
        "proceed to compare the results obtained by fitting the classification models without tuning any\n",
        "parameter on the following datasets:\n",
        "* the original one;\n",
        "* the one obtained by removing the \"Process temperature\" and \"Air temperature\" columns, replacing them with a column of their product;\n",
        "* the one obtained by removing \"Torque\" and \"Rotational speed\", replacing them with a column of their product;\n",
        "* a combine the previous operations."
      ],
      "metadata": {
        "id": "qClwVAyTCBvG"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Models\n",
        "lr = LogisticRegression()\n",
        "knn = KNeighborsClassifier()\n",
        "svc = SVC(probability=True)\n",
        "rfc = RandomForestClassifier()\n",
        "xgb = XGBClassifier()\n",
        "\n",
        "clf = [lr,knn,svc,rfc,xgb]\n",
        "clf_str = ['LR','KNN','SVC','RFC','XGB']\n",
        "\n",
        "# Fit on raw train\n",
        "metrics_0 = fit_models(clf,clf_str,X_train,X_val,y_train,y_val)\n",
        "\n",
        "# Fit on temperature product train\n",
        "XX_train = X_train.drop(columns=['Process temperature','Air temperature'])\n",
        "XX_val = X_val.drop(columns=['Process temperature','Air temperature'])\n",
        "XX_train['Temperature']= X_train['Process temperature']*X_train['Air temperature']\n",
        "XX_val['Temperature'] = X_val['Process temperature']*X_val['Air temperature']\n",
        "metrics_1 = fit_models(clf,clf_str,XX_train,XX_val,y_train,y_val)\n",
        "\n",
        "# Fit on power product train\n",
        "XX_train = X_train.drop(columns=['Rotational speed','Torque'])\n",
        "XX_val = X_val.drop(columns=['Rotational speed','Torque'])\n",
        "XX_train['Power'] = X_train['Rotational speed']*X_train['Torque']\n",
        "XX_val['Power'] = X_val['Rotational speed']*X_val['Torque']\n",
        "metrics_2 = fit_models(clf,clf_str,XX_train,XX_val,y_train,y_val)\n",
        "\n",
        "# Fit on both products train\n",
        "XX_train = X_train.drop(columns=['Process temperature','Air temperature','Rotational speed','Torque'])\n",
        "XX_val = X_val.drop(columns=['Process temperature','Air temperature','Rotational speed','Torque'])\n",
        "XX_train['Temperature']= X_train['Process temperature']*X_train['Air temperature']\n",
        "XX_val['Temperature']= X_val['Process temperature']*X_val['Air temperature']\n",
        "XX_train['Power'] = X_train['Rotational speed']*X_train['Torque']\n",
        "XX_val['Power'] = X_val['Rotational speed']*X_val['Torque']\n",
        "metrics_3 = fit_models(clf,clf_str,XX_train,XX_val,y_train,y_val)\n",
        "\n",
        "# classification metrics barplot\n",
        "fig, axs = plt.subplots(nrows=2, ncols=3, figsize=(18,8))\n",
        "fig.suptitle('Classification metrics')\n",
        "for j, model in enumerate(clf_str):\n",
        "    ax = axs[j//3,j-3*(j//3)]\n",
        "    model_metrics = pd.DataFrame(data=[metrics_0[model],metrics_1[model],metrics_2[model],metrics_3[model]])\n",
        "    model_metrics.index = ['Original','Temperature','Power','Both']\n",
        "    model_metrics.transpose().plot(ax=ax, kind='bar', rot=0, )\n",
        "    ax.title.set_text(model)\n",
        "    ax.get_legend().remove()\n",
        "fig.subplots_adjust(top=0.9, left=0.1, right=0.9, bottom=0.12)\n",
        "axs.flatten()[-2].legend(title='Dataset', loc='upper center',\n",
        "                         bbox_to_anchor=(0.5, -0.12), ncol=4, fontsize=12)\n",
        "plt.show()"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:14:42.817961Z",
          "iopub.execute_input": "2022-07-18T13:14:42.818635Z",
          "iopub.status.idle": "2022-07-18T13:15:22.968833Z",
          "shell.execute_reply.started": "2022-07-18T13:14:42.818596Z",
          "shell.execute_reply": "2022-07-18T13:15:22.967678Z"
        },
        "trusted": true,
        "id": "1fMM4AqrCBvG"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "From the results obtained, we observe that all the models applied to the entire dataset perform\n",
        "better than when they are applied to the ones created by reducing the number of features.\n",
        "The best performances and the modest number of features from which our dataset is composed\n",
        "encourage us to opt to avoid the feature selection step."
      ],
      "metadata": {
        "id": "3adldgUdCBvG"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "### 3.3) Logistic Regression Benchmark <a id=\"binary_benchmark\"></a>\n",
        "\n",
        "We decide to use Logistic Regression as a Benchmark for our task. It\n",
        "represents an intermediate step between the basic model referred to in Section 2.4 and the more\n",
        "complex models that we have described and we will explore in depth in the following sections.\n",
        "Now we look at the results obtained and at the interpretability of the model."
      ],
      "metadata": {
        "id": "dYnpXQSSCBvG"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Make predictions\n",
        "lr = LogisticRegression(random_state=0)\n",
        "lr.fit(X_train, y_train['Target'])\n",
        "y_val_lr = lr.predict(X_val)\n",
        "y_test_lr = lr.predict(X_test)\n",
        "\n",
        "# Metrics\n",
        "cm_val_lr, metrics_val_lr = eval_preds(lr,X_val,y_val,y_val_lr,'binary')\n",
        "cm_test_lr, metrics_test_lr = eval_preds(lr,X_test,y_test,y_test_lr,'binary')\n",
        "print('Validation set metrics:',metrics_val_lr, sep='\\n')\n",
        "print('Test set metrics:',metrics_test_lr, sep='\\n')\n",
        "\n",
        "cm_labels = ['Not Failure', 'Failure']\n",
        "cm_lr = [cm_val_lr, cm_test_lr]\n",
        "# Show Confusion Matrices\n",
        "fig, axs = plt.subplots(ncols=2, figsize=(8,4))\n",
        "fig.suptitle('LR Confusion Matrices')\n",
        "for j, title in enumerate(['Validation Set', 'Test Set']):\n",
        "    ax = axs[j]\n",
        "    sns.heatmap(ax=ax, data=cm_lr[j], annot=True,\n",
        "              fmt='d', cmap='Blues', cbar=False)\n",
        "    axs[j].title.set_text(title)\n",
        "    axs[j].set_xticklabels(cm_labels)\n",
        "    axs[j].set_yticklabels(cm_labels)\n",
        "plt.show()\n",
        "\n",
        "# Odds for interpretation\n",
        "d = {'feature': X_train.columns, 'odds': np.exp(lr.coef_[0])}\n",
        "odds_df = pd.DataFrame(data=d).sort_values(by='odds', ascending=False)\n",
        "odds_df"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:15:22.970425Z",
          "iopub.execute_input": "2022-07-18T13:15:22.970926Z",
          "iopub.status.idle": "2022-07-18T13:15:23.350153Z",
          "shell.execute_reply.started": "2022-07-18T13:15:22.970881Z",
          "shell.execute_reply": "2022-07-18T13:15:23.349215Z"
        },
        "trusted": true,
        "id": "VOrWj3N3CBvG"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "The odds of logistic regression allow us to understand how the model is working. In particular,\n",
        "an unrealistically high importance is given to Torque and Rotational Speed. This is mainly due\n",
        "to the natural variance in these features, which is especially high when looking only at the failure\n",
        "cases and tends to \"deviate\" the model. However it is reasonable to believe, on the basis of\n",
        "exploratory analysis, that the first four features have a significantly greater relevance than the\n",
        "last two. We also expect greater reliability of the odds values when we apply logistic regression to\n",
        "the multiclass task, since the effects that are spread here appears to be localized around certain\n",
        "types of failures."
      ],
      "metadata": {
        "id": "bE5MLPU6CBvG"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "### 3.4) Models <a id=\"binary_models\"></a>\n"
      ],
      "metadata": {
        "id": "5f1js3aDCBvG"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Models\n",
        "knn = KNeighborsClassifier()\n",
        "svc = SVC()\n",
        "rfc = RandomForestClassifier()\n",
        "xgb = XGBClassifier()\n",
        "clf = [knn,svc,rfc,xgb]\n",
        "clf_str = ['KNN','SVC','RFC','XGB']\n",
        "\n",
        "# Parameter grids for GridSearch\n",
        "knn_params = {'n_neighbors':[1,3,5,8,10]}\n",
        "svc_params = {'C': [1, 10, 100],\n",
        "              'gamma': [0.1,1],\n",
        "              'kernel': ['rbf'],\n",
        "              'probability':[True],\n",
        "              'random_state':[0]}\n",
        "rfc_params = {'n_estimators':[100,300,500,700],\n",
        "              'max_depth':[5,7,10],\n",
        "              'random_state':[0]}\n",
        "xgb_params = {'n_estimators':[300,500,700],\n",
        "              'max_depth':[5,7],\n",
        "              'learning_rate':[0.01,0.1],\n",
        "              'objective':['binary:logistic']}\n",
        "params = pd.Series(data=[knn_params,svc_params,rfc_params,xgb_params],\n",
        "                   index=clf)\n",
        "\n",
        "# Tune hyperparameters with GridSearch (estimated time 8m)\n",
        "print('GridSearch start')\n",
        "fitted_models_binary = []\n",
        "for model, model_name in zip(clf, clf_str):\n",
        "    print('Training '+str(model_name))\n",
        "    fit_model = tune_and_fit(model,X_train,y_train,params[model],'binary')\n",
        "    fitted_models_binary.append(fit_model)"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:15:23.351537Z",
          "iopub.execute_input": "2022-07-18T13:15:23.352153Z",
          "iopub.status.idle": "2022-07-18T13:23:35.785342Z",
          "shell.execute_reply.started": "2022-07-18T13:15:23.352119Z",
          "shell.execute_reply": "2022-07-18T13:23:35.783572Z"
        },
        "trusted": true,
        "id": "h8N-6svFCBvG"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# Create evaluation metrics\n",
        "task = 'binary'\n",
        "y_pred_val, cm_dict_val, metrics_val = predict_and_evaluate(\n",
        "    fitted_models_binary,X_val,y_val,clf_str,task)\n",
        "y_pred_test, cm_dict_test, metrics_test = predict_and_evaluate(\n",
        "    fitted_models_binary,X_test,y_test,clf_str,task)\n",
        "\n",
        "# Show Validation Confusion Matrices\n",
        "fig, axs = plt.subplots(ncols=4, figsize=(20,4))\n",
        "fig.suptitle('Validation Set Confusion Matrices')\n",
        "for j, model_name in enumerate(clf_str):\n",
        "    ax = axs[j]\n",
        "    sns.heatmap(ax=ax, data=cm_dict_val[model_name], annot=True,\n",
        "                fmt='d', cmap='Blues', cbar=False)\n",
        "    ax.title.set_text(model_name)\n",
        "    ax.set_xticklabels(cm_labels)\n",
        "    ax.set_yticklabels(cm_labels)\n",
        "plt.show()\n",
        "\n",
        "# Show Test Confusion Matrices\n",
        "fig, axs = plt.subplots(ncols=4, figsize=(20,4))\n",
        "fig.suptitle('Test Set Confusion Matrices')\n",
        "for j, model_name in enumerate(clf_str):\n",
        "    ax = axs[j]\n",
        "    sns.heatmap(ax=ax, data=cm_dict_test[model_name], annot=True,\n",
        "                fmt='d', cmap='Blues', cbar=False)\n",
        "    ax.title.set_text(model_name)\n",
        "    ax.set_xticklabels(cm_labels)\n",
        "    ax.set_yticklabels(cm_labels)\n",
        "plt.show()\n",
        "\n",
        "# Print scores\n",
        "print('')\n",
        "print('Validation scores:', metrics_val, sep='\\n')\n",
        "print('Test scores:', metrics_test, sep='\\n')"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:23:35.786928Z",
          "iopub.execute_input": "2022-07-18T13:23:35.787601Z",
          "iopub.status.idle": "2022-07-18T13:23:37.90642Z",
          "shell.execute_reply.started": "2022-07-18T13:23:35.787563Z",
          "shell.execute_reply": "2022-07-18T13:23:37.905259Z"
        },
        "trusted": true,
        "id": "J6rrK5hnCBvG"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "All the selected models obtain similar results on the validation set (except KNN which is a little\n",
        "worse) and it is diffj7icult to determine if one works better than another by looking only at these\n",
        "values. Performance did not significantly drop when passing the test set, showing that overfitting\n",
        "was avoided. We comment on the results of the models by looking at the confusion matrices\n",
        "and the metrics obtained on the test set: in this way the formation of a hierarchy between the\n",
        "models used is slightly clearer, as all the metrics relating to a single model are smaller or larger\n",
        "than to the others and the time needed to search for the parameters is comparable, with the only\n",
        "exception of KNN. In particular KNN obtains the worst performances and XGB the best ones; in\n",
        "the middle we find SVC and RFC which achieve extremely similar results.\n",
        "\n",
        "About the parameters:\n",
        "* A Gridsearch has been started on the parameters which, looking in the literature, appear to be preponderant for each specific model;\n",
        "* The grid values to search for have been defined on the basis of literature and various tests, trying to keep the computational cost of finding the best values moderate.\n",
        "\n",
        "It is interesting to observe that the optimal parameters for RFC and XGB are the polar\n",
        "opposite: the former prefers to use a few estimators and go into depth while the latter uses more\n",
        "estimators with fewer splits. Furthermore, it must be taken into account that although XGB\n",
        "is the best classifier from a quantitative point of view, this is not true for what concerns the\n",
        "qualitative side. Both SVC and XGB in fact lack clear ways to interpret the results, while on\n",
        "the contrary RFC allows to have a complete understanding of how the algorithm worked. In\n",
        "any case, to get an idea of which features had greater importance in making the predictions, we\n",
        "report the permutation feature importances in a bar plot."
      ],
      "metadata": {
        "id": "26OVrLP7CBvG"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Evaluate Permutation Feature Importances\n",
        "f2_scorer = make_scorer(fbeta_score, pos_label=1, beta=2)\n",
        "importances = pd.DataFrame()\n",
        "for clf in fitted_models_binary:\n",
        "    result = permutation_importance(clf, X_train,y_train['Target'],\n",
        "                                  scoring=f2_scorer,random_state=0)\n",
        "    result_mean = pd.Series(data=result.importances_mean, index=X.columns)\n",
        "    importances = pd.concat(objs=[importances,result_mean],axis=1)\n",
        "importances.columns = clf_str\n",
        "\n",
        "# Barplot of Feature Importances\n",
        "fig, axs = plt.subplots(ncols=4, figsize=(20,4))\n",
        "fig.suptitle('Permutation Feature Importances')\n",
        "for j, name in enumerate(importances.columns):\n",
        "    sns.barplot(ax=axs[j], x=importances.index, y=importances[name].values)\n",
        "    axs[j].tick_params('x',labelrotation=90)\n",
        "    axs[j].set_ylabel('Importances')\n",
        "    axs[j].title.set_text(str(name))\n",
        "plt.show()"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:33:50.90093Z",
          "iopub.execute_input": "2022-07-18T13:33:50.901322Z",
          "iopub.status.idle": "2022-07-18T13:34:37.469245Z",
          "shell.execute_reply.started": "2022-07-18T13:33:50.90129Z",
          "shell.execute_reply": "2022-07-18T13:34:37.467878Z"
        },
        "trusted": true,
        "id": "dosLTae3CBvG"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "Remarks on Feature importances:\n",
        "* Type is the feature with the lowest significance, in accordance with what was observed during the exploratory analysis. However, its importance remains strictly positive in each of the cases considered and therefore removing it completely would have led to a decline in prediction performance, not justified by a significant computational gain;\n",
        "* Unlike Logistic Regression, the models tested place great emphasis on Tool wear as well as Torque and Rotational Speed. Since the former alone is related to a specific category of failures and strongly distorts the kdeplot of Machine failure, we have a sign that our models still worked well."
      ],
      "metadata": {
        "id": "Q94C_6oUCBvH"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "## 4) **Multi-class task** <a id=\"multi\"></a>\n",
        "We now proceed to the second task of this project, that is predict not only if there will be a\n",
        "failure, but also the type of failure that will occur. So we are in the case of multiclass classification problems that make the assumption that each sample is assigned to one and only one label. This hypothesis is verified because in data preprocessing we removed all the ambiguous observations that belonged to more than one class.\n",
        "\n",
        "For multiclass targets, when we calculate the values of AUC, F1 and F2 scores, we need to set the\n",
        "parameter \"average\". We choose \"average=weighted\", in order to account for class imbalance:\n",
        "in fact, at the end of data preprocessing, we have 80% WORKING machine and 20% that fail.\n",
        "As for binary classification task, we choose Logistic Regression as baseline model and we look for\n",
        "models that get higher values for the chosen metrics. In particular, we adapt to the multiclass\n",
        "case the models developed in the previous section. While many classification algorithms (such\n",
        "as K-nearest neighbor, Random Forest and XGBoost) naturally permit the use of more than\n",
        "two classes, some (like Logistic Regression and Support Vector Machines) are by nature binary\n",
        "algorithms; these can, however, be turned into multiclass classifiers by a variety of strategies. For\n",
        "our project, we decide to use \"OnevsRest\" approach, who involves training a single classifier per\n",
        "class, with the samples of that class as positive samples and all other samples as negatives. We\n",
        "choose it because it is computationally more efficient than other types of approach."
      ],
      "metadata": {
        "id": "LkEAOtQdCBvH"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "### 4.1) Logistic Regression Benchmark <a id=\"multi_benchmark\"></a>\n",
        "First let’s look at how the Logistic Regression behaves:"
      ],
      "metadata": {
        "id": "MN3pUNiPCBvH"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# multiclass classification\n",
        "lr = LogisticRegression(random_state=0,multi_class='ovr')\n",
        "lr.fit(X_train, y_train['Failure Type'])\n",
        "y_val_lr = lr.predict(X_val)\n",
        "y_test_lr = lr.predict(X_test)\n",
        "\n",
        "# Validation metrics\n",
        "cm_val_lr, metrics_val_lr = eval_preds(lr,X_val,y_val,y_val_lr,'multi_class')\n",
        "cm_test_lr, metrics_test_lr = eval_preds(lr,X_test,y_test,y_test_lr,'multi_class')\n",
        "print('Validation set metrics:',metrics_val_lr, sep='\\n')\n",
        "print('Test set metrics:',metrics_test_lr, sep='\\n')\n",
        "\n",
        "cm_lr = [cm_val_lr, cm_test_lr]\n",
        "cm_labels = ['No Fail','PWF','OSF','HDF','TWF']\n",
        "# Show Confusion Matrices\n",
        "fig, axs = plt.subplots(ncols=2, figsize=(9,4))\n",
        "fig.suptitle('LR Confusion Matrices')\n",
        "for j, title in enumerate(['Validation Set', 'Test Set']):\n",
        "    ax = axs[j]\n",
        "    sns.heatmap(ax=ax, data=cm_lr[j], annot=True,\n",
        "              fmt='d', cmap='Blues', cbar=False)\n",
        "    axs[j].title.set_text(title)\n",
        "    axs[j].set_xticklabels(cm_labels)\n",
        "    axs[j].set_yticklabels(cm_labels)\n",
        "plt.show()\n",
        "\n",
        "# Odds for interpretation\n",
        "odds_df = pd.DataFrame(data = np.exp(lr.coef_), columns = X_train.columns,\n",
        "                       index = df_res['Failure Type'].unique())\n",
        "odds_df"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:38:25.122715Z",
          "iopub.execute_input": "2022-07-18T13:38:25.124274Z",
          "iopub.status.idle": "2022-07-18T13:38:25.981333Z",
          "shell.execute_reply.started": "2022-07-18T13:38:25.124194Z",
          "shell.execute_reply": "2022-07-18T13:38:25.979714Z"
        },
        "trusted": true,
        "id": "gX4mNF7XCBvH"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "In the table above there are, for every class, the Logistic Regression’s odds that\n",
        "explain the contribution of each feature in the prediction of belonging to a specific class. By\n",
        "comparing this table with the PCA scatter and the comments we made, we understand that there\n",
        "is a complete agreement about the features that most affect the type of failure. For example, if\n",
        "we look at odds’ values of PWF, we see that Rotational Speed and Torque are the ones that are\n",
        "most important for the forecast of belonging to this class. In the analysis of the PCA we stated\n",
        "that PWF seems to be dependent only on PC2, i.e. the Power that is the product of Rotational\n",
        "Speed and Torque. We can make similar considerations for other classes."
      ],
      "metadata": {
        "id": "7JErmoDNCBvH"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "### 4.2) Models <a id=\"multi_models\"></a>\n",
        "For each model we launch the Gridsearch for hyperparameter optimization, using as metric to\n",
        "evaluate the model the weighted average F2 score. Similarly to the binary case, the Gridsearch\n",
        "has been started on the parameters that, looking in the literature, are found to be preponderant\n",
        "for each specific model and the grid values to look for have been defined according to the literature\n",
        "and several tests carried out."
      ],
      "metadata": {
        "id": "PmKx4KF1CBvH"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Models\n",
        "knn = KNeighborsClassifier()\n",
        "svc = SVC(decision_function_shape='ovr')\n",
        "rfc = RandomForestClassifier()\n",
        "xgb = XGBClassifier()\n",
        "clf = [knn,svc,rfc,xgb]\n",
        "clf_str = ['KNN','SVC','RFC','XGB']\n",
        "\n",
        "knn_params = {'n_neighbors':[1,3,5,8,10]}\n",
        "svc_params = {'C': [1, 10, 100],\n",
        "              'gamma': [0.1,1],\n",
        "              'kernel': ['rbf'],\n",
        "              'probability':[True],\n",
        "              'random_state':[0]}\n",
        "rfc_params = {'n_estimators':[100,300,500,700],\n",
        "              'max_depth':[5,7,10],\n",
        "              'random_state':[0]}\n",
        "xgb_params = {'n_estimators':[100,300,500],\n",
        "              'max_depth':[5,7,10],\n",
        "              'learning_rate':[0.01,0.1],\n",
        "              'objective':['multi:softprob']}\n",
        "\n",
        "params = pd.Series(data=[knn_params,svc_params,rfc_params,xgb_params],\n",
        "                    index=clf)\n",
        "\n",
        "\n",
        "# Tune hyperparameters with GridSearch (estimated time 8-10m)\n",
        "print('GridSearch start')\n",
        "fitted_models_multi = []\n",
        "for model, model_name in zip(clf, clf_str):\n",
        "    print('Training '+str(model_name))\n",
        "    fit_model = tune_and_fit(model,X_train,y_train,params[model],'multi_class')\n",
        "    fitted_models_multi.append(fit_model)"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:38:30.241714Z",
          "iopub.execute_input": "2022-07-18T13:38:30.242318Z",
          "iopub.status.idle": "2022-07-18T13:52:16.495056Z",
          "shell.execute_reply.started": "2022-07-18T13:38:30.242284Z",
          "shell.execute_reply": "2022-07-18T13:52:16.493729Z"
        },
        "trusted": true,
        "id": "502ZxNLoCBvH"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# Create evaluation metrics\n",
        "\n",
        "task = 'multi_class'\n",
        "y_pred_val, cm_dict_val, metrics_val = predict_and_evaluate(\n",
        "    fitted_models_multi,X_val,y_val,clf_str,task)\n",
        "y_pred_test, cm_dict_test, metrics_test = predict_and_evaluate(\n",
        "    fitted_models_multi,X_test,y_test,clf_str,task)\n",
        "\n",
        "# Show Validation Confusion Matrices\n",
        "fig, axs = plt.subplots(ncols=4, figsize=(20,4))\n",
        "fig.suptitle('Validation Set Confusion Matrices')\n",
        "for j, model_name in enumerate(clf_str):\n",
        "    ax = axs[j]\n",
        "    sns.heatmap(ax=ax, data=cm_dict_val[model_name], annot=True,\n",
        "                fmt='d', cmap='Blues', cbar=False)\n",
        "    ax.title.set_text(model_name)\n",
        "    ax.set_xticklabels(cm_labels)\n",
        "    ax.set_yticklabels(cm_labels)\n",
        "plt.show()\n",
        "\n",
        "# Show Test Confusion Matrices\n",
        "fig, axs = plt.subplots(ncols=4, figsize=(20,4))\n",
        "fig.suptitle('Test Set Confusion Matrices')\n",
        "for j, model_name in enumerate(clf_str):\n",
        "    ax = axs[j]\n",
        "    sns.heatmap(ax=ax, data=cm_dict_test[model_name], annot=True,\n",
        "                fmt='d', cmap='Blues', cbar=False)\n",
        "    ax.title.set_text(model_name)\n",
        "    ax.set_xticklabels(cm_labels)\n",
        "    ax.set_yticklabels(cm_labels)\n",
        "plt.show()\n",
        "\n",
        "# Print scores\n",
        "print('')\n",
        "print('Validation scores:', metrics_val, sep='\\n')\n",
        "print('Test scores:', metrics_test, sep='\\n')"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:55:13.218209Z",
          "iopub.execute_input": "2022-07-18T13:55:13.218612Z",
          "iopub.status.idle": "2022-07-18T13:55:16.535486Z",
          "shell.execute_reply.started": "2022-07-18T13:55:13.218572Z",
          "shell.execute_reply": "2022-07-18T13:55:16.534194Z"
        },
        "trusted": true,
        "id": "0kDnHm4RCBvH"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "By comparing the results obtained, we see that K-NN is the model that performs the worst and\n",
        "its accuracy is a little lower than Logistic Regression’s one. Despite this, we cannot exclude it a priori, as it still reaches high values for the metrics and, moreover, gives an immediate response.\n",
        "So, we can use it whenever we need to get an idea quickly about the situation and, then apply\n",
        "other models when we have more time.\n",
        "\n",
        "All other models perform better than the benchmark and they obtain high values for the chosen\n",
        "metrics both for validation and test set.\n",
        "SVC and RFC’s performances are very similar each other and XGB performs better than them.\n",
        "If we look at the training phase, SVC and RFC take the same time, while XGB takes more than\n",
        "four times as much as them. So, since, the improvement obtained with XGB is only 1.5%, one\n",
        "can choose which model he prefers according to his needs.\n",
        "While the best parameters for multiclass K-NN and SVC are the same as binary classification, for\n",
        "XGB and RFC the Gridsearch for the two types of task returns different parameters. Moreover,\n",
        "in the transition from binary to multiclass problem, the estimated training time remains the same\n",
        "for all models, except for XGB that triples it.\n",
        "In order to understand how features contribute to predictions, let’s look at the Permutation\n",
        "Feature Importances for each model."
      ],
      "metadata": {
        "id": "ME7BQi-HCBvI"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Evaluate Permutation Feature Importances\n",
        "f2_scorer = make_scorer(fbeta_score, beta=2, average='weighted')\n",
        "importances = pd.DataFrame()\n",
        "for clf in fitted_models_multi:\n",
        "    result = permutation_importance(clf, X_train,y_train['Failure Type'],\n",
        "                                  scoring=f2_scorer,random_state=0)\n",
        "    result_mean = pd.Series(data=result.importances_mean, index=X.columns)\n",
        "    importances = pd.concat(objs=[importances,result_mean],axis=1)\n",
        "\n",
        "importances.columns = clf_str\n",
        "\n",
        "# Barplot of Feature Importances\n",
        "fig, axs = plt.subplots(ncols=4, figsize=(20,4))\n",
        "fig.suptitle('Permutation Feature Importances')\n",
        "for j, name in enumerate(importances.columns):\n",
        "    sns.barplot(ax=axs[j], x=importances.index, y=importances[name].values)\n",
        "    axs[j].tick_params('x',labelrotation=90)\n",
        "    axs[j].set_ylabel('Importances')\n",
        "    axs[j].title.set_text(str(name))\n",
        "plt.show()"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:56:59.364784Z",
          "iopub.execute_input": "2022-07-18T13:56:59.365178Z",
          "iopub.status.idle": "2022-07-18T13:58:03.337653Z",
          "shell.execute_reply.started": "2022-07-18T13:56:59.365149Z",
          "shell.execute_reply": "2022-07-18T13:58:03.336745Z"
        },
        "trusted": true,
        "id": "jlSkd653CBvI"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "From previous barplots we see that the models give more importance to Torque, Tool wear\n",
        "and Rotational Speed while the Type contribution is very low. This is in accordance with the\n",
        "observations made in the exploration of the dataset in Section 1-2 and it is consistent with the\n",
        "Permutation Feature Importances of binary task.\n",
        "K-NN is the one who gives more importance to Type, but, different from binary case, here we\n",
        "see that for every model the Type contribution is almost zero. So, we test the model on a new\n",
        "dataset, the old one from which we removed the column Type. For K-NN and SVC there is an\n",
        "insignificant improvement in the metrics’ values, which were already very good. For RFC and XGB we do not see any change on metrics’ values. Since the training time for the different models\n",
        "is approximately equal in both cases, we let users choose which dataset to use."
      ],
      "metadata": {
        "id": "5J2RzuJuCBvI"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "## 5) **Decision Paths** <a id=\"decisionpath\"></a>\n",
        "\n",
        "Here we show the decision paths of one of the trees that make up the Random Forest for both\n",
        "tasks, truncated at depth=4. However this depth is enough to verify that trees require to be\n",
        "deep because the decision boundary are complex themselves and they are not overfitting. This\n",
        "is evident if one looks at the multi-class tree, where some kinds of failure do not appear before\n",
        "depth four, but also in the binary classification tree by looking at the evolution of the gini score\n",
        "while following most of the paths. A further remark can be made about the feature Type being\n",
        "the origin node of both graphs and separating the majority class (Low quality) from the other\n",
        "two at the first step. It appears just one time more in the upper side of the trees and shows\n",
        "sporadically again at the lowest floors, where its impact is scarce."
      ],
      "metadata": {
        "id": "KxOyvEvkCBvI"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Random Forest Decision Path\n",
        "from sklearn import tree\n",
        "import graphviz\n",
        "\n",
        "tree_binary = fitted_models_binary[2].best_estimator_.estimators_[0]\n",
        "tree_multi = fitted_models_multi[2].best_estimator_.estimators_[0]\n",
        "trees = [tree_binary,tree_multi]\n",
        "targets = ['Target', 'Failure Type']\n",
        "for decision_tree, target in zip(trees, targets):\n",
        "    decision_tree.fit(X_train,y_train[target])\n",
        "    classes = list(map(str,df_res[target].unique()))\n",
        "\n",
        "    dot_data = tree.export_graphviz(decision_tree, out_file=None,\n",
        "                                  feature_names=X.columns,\n",
        "                                  class_names=classes,\n",
        "                                  filled=True, rounded=True,\n",
        "                                  special_characters=True,\n",
        "                                  max_depth=4)  # uncomment to see full tree\n",
        "    graph = graphviz.Source(dot_data)\n",
        "    graph.render(target+\" Classification tree\")\n",
        "    display(graph)"
      ],
      "metadata": {
        "execution": {
          "iopub.status.busy": "2022-07-18T13:59:28.152908Z",
          "iopub.execute_input": "2022-07-18T13:59:28.153275Z",
          "iopub.status.idle": "2022-07-18T13:59:29.596684Z",
          "shell.execute_reply.started": "2022-07-18T13:59:28.153246Z",
          "shell.execute_reply": "2022-07-18T13:59:29.595256Z"
        },
        "trusted": true,
        "id": "WCQYdYzuCBvI"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "## 6) **Conclusions** <a id=\"conclusions\"></a>\n",
        "\n",
        "According to the analyses carried out and the results obtained, it is possible to make some\n",
        "conclusive considerations related to this project.\n",
        "\n",
        "We decided to tackle two tasks: predict whether a machine will fail or not and predict the type\n",
        "of failure that will occur. Before developing the models we did data preprocessing to ensure the\n",
        "validity of the assumptions of applicability of the models and ensure the best performances.\n",
        "Briefly, in preprocessing phase we have deleted some ambiguous samples, we applied a label\n",
        "encoding to the categorical columns and then we performed the scaling of the columns with\n",
        "StandardScaler. We also noticed the presence of some data points which at first we referred\n",
        "as outliers but later turned out to be part of the natural variance of the data and played an\n",
        "important role in the classification task. Then we ran PCA and found that most of the variance\n",
        "is explained by the first three components, that can be represented as the following features:\n",
        "combination of the two Temperatures, Machine Power (product of Rotational Speed and Torque)\n",
        "and Tool Wear. In according to this, we found that these are the features that contribute the\n",
        "most in the predictions when apply the models. Contrary to logical predictions, we demonstrated\n",
        "that the machine’s type does not affect the presence of failure.\n",
        "\n",
        "At the end, we can conclude that for both task the chosen models perform very well. For both\n",
        "tasks the best model is XGBoost and the worst is KNN; however the response time of KNN\n",
        "is instant while XGBoost takes more time and this further increase when we proceed with the\n",
        "multi-class classification task. The choice of the model depends on the needs of the company: for\n",
        "faster application one can use KNN while if one cares more about accuracy one can use XGBoost."
      ],
      "metadata": {
        "id": "Tat73nCdCBvI"
      }
    }
  ]
}